01 Abstract
The ATM Milano tramway fleet includes the historic "Carelli 1928" vehicles — Series 1500 cars that have been in continuous urban service for nearly a century. Each tram is propelled by four separately excited DC motors, each rated at 21 kW, fed from a 600 V DC overhead line.
Despite their age, the electromechanical architecture of these drives presents a rich control problem: the decoupling of flux and torque requires coordinated management of the excitation and armature circuits, and the variable speed profile over a 10 km route — including uphill and downhill gradients — demands robust disturbance rejection combined with strict actuator constraint handling.
This project develops and validates a complete cascade control solution, covering:
- Parameter identification and KVL-based consistency verification of the machine constant
- Derivation of the armature, excitation and mechanical transfer functions
- Pole-zero cancellation tuning of three nested PI regulators
- Dynamic field-weakening strategy for above-rated-speed operation
- Anti-windup back-calculation to handle actuator saturation cleanly
- Full Simulink simulation over the prescribed 10 km kinematic profile
The Carelli trams still operate on ATM lines 1 and 9 in Milan. Their four DC motors with separate excitation allow independent control of torque and flux — a concept that predates modern vector control by half a century, yet maps directly onto contemporary cascade PI architectures.
02 System Identification
2.1 Motor Parameters
The motors are characterised by the datasheet parameters below. A critical step in model identification concerns the machine constant. Although the datasheet provides a torque constant \(K_t = 1.06\,\text{Nm/A}\), applying Kirchhoff's Voltage Law to the armature circuit at rated conditions reveals that this value must be interpreted as the total machine constant \(K_s\), i.e. it already incorporates the rated excitation current \(I_e = 5\,\text{A}\):
$$ E_n \big|_{\text{KVL}} = V_{DC} - R_a I_{tot} = 600 - 0.39 \times 156 = 539.16\,\text{V} $$
$$ E_n \big|_{K_s} = K_s \cdot I_e \cdot \Omega_n = 1.06 \times 5 \times 101.6 = 538.5\,\text{V} \approx 539\,\text{V} \;\checkmark $$
$$ T_n \big|_{\text{motor}} = K_s \cdot I_e \cdot I_{tot} = 1.06 \times 5 \times 156 = 826.8\,\text{Nm} $$
$$ T_n \big|_{\text{power}} = \frac{4 P_r}{\Omega_n} = \frac{84\,000}{101.6} = 827.2\,\text{Nm} \;\checkmark $$
Both checks confirm \(K_s = K_t = 1.06\,\text{Nm/A}^2\). Any other interpretation produces a voltage mismatch that would make rated performance physically impossible within the 600 V supply constraint.
2.2 Vehicle Parameters
2.3 Kinematic Profile — 10 km Urban Track
The longitudinal speed profile is defined by the following seven segments. Speed transitions are rate-limited to \(a_{max} = 1.0\,\text{m/s}^2\) for passenger comfort. Slope variations are modelled as measurable torque disturbances acting on the mechanical plant.
| Segment [km] | Slope | Speed | v [m/s] | Regime |
|---|---|---|---|---|
| 0 – 1 | 0 % | \(v_r / 2\) | 5.83 | Low speed — urban entry |
| 1 – 3 | 0 % | \(v_r\) | 11.67 | Rated speed — flat |
| 3 – 4 | +5 % | \(v_r\) | 11.67 | Uphill — extra load disturbance |
| 4 – 6 | 0 % | \(v_{max}\) | 11.67* | Max speed — field weakening active |
| 6 – 8 | 0 % | \(v_r\) | 11.67 | Back to rated — flat |
| 8 – 9 | −5 % | \(v_r\) | 11.67 | Downhill — reduced traction demand |
| 9 – 10 | 0 % | \(v_r / 2\) | 5.83 | Low speed — urban exit |
* vmax = 42 km/h = 11.67 m/s = vr in this case. The field-weakening region is activated when the motor shaft speed exceeds the base speed ωbase.
03 Modelling of the Physical System
3.1 Electrical Model — Decoupling Strategy
The separately excited DC motor has two independent electrical circuits: the armature circuit, which carries the traction current and generates electromagnetic torque, and the excitation (field) circuit, which establishes the air-gap flux. The two are coupled through the back-EMF \(e = K_s I_e \omega\), which acts as a feedback from the mechanical state into the armature voltage equation.
The control architecture is therefore built on three nested loops, ordered from fastest to slowest: excitation current → armature current → angular speed. To enable independent tuning of each loop, a feedforward back-EMF compensation \(\hat{e} = K_s \hat{I}_e \hat{\omega}\) is injected at the armature voltage reference, decoupling the electrical dynamics from the mechanical ones.
3.2 Mechanical Model
The tram is modelled as a single equivalent wheel of diameter \(d\) and total loaded mass \(M\), rolling on a slope angle \(\gamma\) under no-slip conditions. Resolving forces along the longitudinal axis and moments about the wheel-rail contact point yields the two governing equations:
$$ \sum F_x:\quad f - Mg\sin\gamma = M\ddot{x} \;\Rightarrow\; f = Mg\sin\gamma + M\ddot{x} $$
$$ \sum M_P:\quad m_e - \beta\omega - Mg\frac{d}{2}\sin\gamma\cdot\rho = J\dot{\omega} $$
The gravity-slope term \(Mg\,(d/2)\sin\gamma\cdot\rho\) is treated as a measurable disturbance torque on the motor shaft. With \(\rho = 13/74\) and \(d = 0.68\,\text{m}\), a 5 % grade produces a disturbance torque of approximately 118 Nm — about 14 % of rated torque \(T_n \approx 827\,\text{Nm}\).
3.3 Transfer Functions
All three plant transfer functions are asymptotically stable first-order systems, allowing straightforward pole-zero cancellation PI design:
$$ G_A(s) = \frac{1}{R_a + sL_a} \qquad L_a = R_a\cdot\tau_a = 0.39 \times 0.01 = 3.9\,\text{mH} $$
$$ G_E(s) = \frac{1}{R_e + sL_e} \qquad L_e = R_e\cdot\tau_e = 12 \times 0.1 = 1.2\,\text{H} $$
$$ G_M(s) = \frac{1}{\beta + sJ} \qquad J = M\!\left(\rho\frac{d}{2}\right)^{\!2} = 25400\times(0.1757\times0.34)^2 \approx 9.10\,\text{kg·m}^2 $$
| Parameter | Symbol | Value |
|---|---|---|
| Armature inductance | La | 3.9 mH |
| Excitation inductance | Le | 1.2 H |
| Nominal back-EMF | En | ≈ 539 V |
| Equivalent inertia | J | ≈ 9.10 kg·m² |
| Mechanical time constant | τM = J/β | ≈ 11.2 s |
| Base speed (FW onset) | ωbase = En/(KsIe) | ≈ 101.6 rad/s |
04 Control Design
The cascade architecture requires three PI controllers, ordered by increasing time constant. A PID structure is deliberately avoided: the derivative action would amplify measurement noise and introduce overshoot — unacceptable for a public transport drive. Each PI is designed via pole-zero cancellation, placing the controller zero at the plant pole \(-R/L\) and selecting the gain to achieve the target crossover frequency with 90° phase margin.
4.1 PI Controller Tuning
Design rule: for a first-order plant \(G(s) = 1/(R+sL)\), the PI gains that achieve crossover frequency \(\omega_c\) with 90° phase margin are \(K_p = \omega_c L\) and \(K_i = \omega_c R\). This zero-cancellation condition is necessary; without it, the resulting phase margin would be less than 90° and the transient would show overshoot.
| Loop | ωc [rad/s] | Phase Margin | Kp | Ki |
|---|---|---|---|---|
| Excitation current Ie | 40 | 90° | 48 | 480 |
| Armature current Ia | 20 | 90° | 0.078 | 7.8 |
| Angular speed ω | 2 | 90° | 181.2 | 1.62 |
Loop bandwidth separation is mandatory: each inner loop must settle at least 5–10 times faster than its outer loop. The chosen values (40 ≫ 20 ≫ 2 rad/s) respect this condition. Violating the separation would cause inter-loop coupling and destabilise the cascade.
4.2 Actuator Saturation & Anti-windup
Every loop output is bounded by physical actuator limits. During saturation, an unconstrained PI integrator continues accumulating the control error, so when the system exits saturation it responds with a large, delayed overshoot — the classic integrator windup problem. The chosen remedy is the back-calculation anti-windup scheme, which feeds back the saturation difference to discharge the integrator at a rate governed by \(K_b = 1/\tau\):
$$ \dot{x}_I(t) = K_i\,e(t) + K_b\!\left[u_{sat}(t) - u(t)\right] $$
$$ u(t) = K_p\,e(t) + x_I(t) \;,\quad u_{sat}(t) = \text{sat}_{[u_{min},\,u_{max}]}\!\left[u(t)\right] $$
When the output saturates, the feedback term \(K_b[u_{sat}-u]\) actively reduces the integrator state, preventing further windup. Back-calculation is preferred in industrial drives for its simplicity, smooth recovery transient and robustness to model uncertainty.
| Loop | τ | Kb = 1/τ | Saturation limits |
|---|---|---|---|
| Excitation current | 0.1 s | 10 | [0 V, 60 V] |
| Armature current | 0.01 s | 100 | [−600 V, +600 V] |
| Speed (torque out) | ≈ 11.2 s | 0.0089 | [−827 Nm, +827 Nm] |
4.3 Field-weakening Strategy
At rated speed the back-EMF already occupies most of the 600 V supply. To operate above the rated speed (segments km 4–6), the excitation current must be dynamically reduced to keep the back-EMF below the voltage limit, entering the constant-power (field-weakening) region:
$$ I_{e,ref} = \begin{cases} I_{e,rated} & \text{if } \omega \leq \omega_{base} \\[6pt] \dfrac{E_n}{K_s\,\omega} & \text{if } \omega > \omega_{base} \end{cases} $$
$$ \omega_{base} = \frac{E_n}{K_s I_e} = \frac{539}{1.06 \times 5} \approx 101.6\,\text{rad/s} \;(= \Omega_r \;\checkmark) $$
The \(1/\omega\) law is implemented in Simulink as a lookup table driven by the measured motor speed. The excitation PI then tracks the reduced reference, reducing \(I_e\) smoothly without transient spikes.
05 Simulation Results
The complete Simulink model is simulated over the full 10 km track. All three control loops perform within specification: zero overshoot and zero steady-state error on the speed reference, armature current clamped at 156 A during acceleration with smooth anti-windup release, and excitation current correctly reduced to ≈ 2.6 A in the field-weakening region (km 4–6).
Risultati completi nella pagina dedicata
I grafici Simulink originali — velocità v(t), corrente di armatura i_a(t), corrente di eccitazione i_e(t) — con osservazioni dettagliate e schemi di controllo sono disponibili nella pagina dei risultati.
| Performance Index | Result | Requirement / Cause |
|---|---|---|
| Steady-state speed error | ≈ 0 | PI integral action |
| Speed overshoot | 0 % | PM = 90° by design |
| Max armature current | 156 A | Saturation block active |
| Anti-windup transient | No spikes | Back-calculation |
| Field weakening | Active km 4–6 | \(1/\omega\) law |
| Slope disturbance rejection | Full recovery | Integral action |
06 Conclusions
The cascade PI architecture proves capable of controlling the Carelli 1928 tram drive across its full operating envelope — including rated-speed traction, above-rated field-weakening and slope disturbances — while respecting every electrical and mechanical constraint.
Three design decisions proved essential for this result. First, interpreting \(K_t\) as the total machine constant \(K_s\) was the only physically consistent choice, as demonstrated by the dual KVL/power-balance verification. Second, strict bandwidth separation (40 / 20 / 2 rad/s) was necessary for independent loop stability. Third, back-calculation anti-windup was required to prevent integrator windup during the frequent current-saturation events typical of urban traction cycles.
The study also highlights a broader engineering lesson: sophisticated control performance can be achieved with simple PI structures, provided that the plant is correctly modelled, the actuator constraints are explicitly handled, and the loops are designed with a clear hierarchy of time scales.
This project was submitted for the course Dynamics of Electrical Machines and Drives (10 CFU) of the MSc in Automation and Control Engineering, Politecnico di Milano, under the supervision of Prof. Francesco Castelli Dezza.
07 MATLAB Code — Key Excerpts
The full script performs parameter identification, KVL consistency checks, PI gain computation and anti-windup setup. The workspace variables are then loaded directly into the Simulink model.
Parameter identification
%% Motor parameters Ra = 0.39; ta = 10e-3; La = ta * Ra; % La = 3.9 mH Re = 12; te = 0.1; Le = te * Re; % Le = 1.2 H Ks = 1.06; Ie = 5; ome = 970*2*pi/60; % 101.6 rad/s %% KVL consistency check En_motor = Ks * Ie * ome; % 538.5 V ← from machine constant En_KVL = 600 - Ra * 156; % 539.2 V ← from KVL → match ✓ Tn_motor = Ks * Ie * 156; % 826.8 Nm ← motor equation Tn_mech = 4 * 21e3 / ome; % 827.2 Nm ← power balance → match ✓ %% Mechanical parameters mt = 15e3; mp = 80; np = 130; M = mt + mp * np; % 25 400 kg d = 680e-3; rho = 13/74; beta = 0.81; J = M * (rho * d/2)^2; % ≈ 9.10 kg·m²
PI design — pole-zero cancellation
%% Target bandwidths [rad/s] wi = 20; we = 40; wm = 2; %% PI gains: Kp = omega*L, Ki = omega*R KpA = wi*La; KiA = wi*Ra; % Armature: Kp=0.078, Ki=7.8 KpE = we*Le; KiE = we*Re; % Excitation: Kp=48, Ki=480 KpM = wm*J; KiM = wm*beta;% Speed: Kp=181.2, Ki=1.62
Anti-windup back-calculation
%% Kb = 1/tau for each loop Kb_a = 1/ta; % = 100 (armature, tau = 10 ms) Kb_e = 1/te; % = 10 (excitation, tau = 0.1 s) Kb_m = 1/(J/beta); % ≈ 0.0089 (speed, tau ≈ 11.2 s) %% Saturation limits Te_max = Ks * Ie * 156; % ±826.8 Nm (speed PI output) Va_max = 600; % ±600 V (armature voltage)
Field-weakening
%% Base speed and FW computation En = Ks * Ie * ome; % nominal back-EMF ≈ 539 V ome_base = En / (Ks * Ie); % ≈ 101.6 rad/s (= omega_rated ✓) ome_max = v_max_ms/(rho*d/2); % max shaft speed %% Simulink lookup: ie_ref = Ie if omega <= ome_base %% En/(Ks*omega) if omega > ome_base
08File di Progetto
Il progetto completo include il codice MATLAB per l'identificazione dei parametri e la sintesi dei controllori PI, il modello Simulink con tutti i blocchi di controllo in cascata, e il report tecnico completo. Tutti i file sono disponibili nel repository GitHub.