## 3 Control Method
### Basic Equations for the HEV Powertrain  
The control target, i.e., the HEV drivetrain system, was configured to have a torsional damper, a planetary gear set, and two electric motors. The engine torque was input as a system disturbance. The equation of the system's motion was derived using the Lagrangian equation. Considering the characteristics of the planetary gears, the relational expression for the rotational speed can be defined as follows:  
\[\omega_{m2}=(\lambda+1)\omega_{c}-\lambda\omega_{m1},\ \ \ \omega_{p}=-\frac{ \lambda+1}{\lambda-1}\omega_{c}+\frac{2\lambda}{\lambda-1}\omega_{m1} \tag{1}\]  
where \(\lambda\) and \(\omega\) denote the planetary gear ratio and rotation speed, respectively. The subscripts \(c\), \(m1\), \(m2\), and \(p\) correspond to the carrier, motor 1, motor 2, and pinion gear, respectively. The relationship between each torque can be expressed as follows:  
\[T_{r}=\frac{T_{m1}}{\lambda}+T_{m2},\ \ \ T_{c}=(\lambda+1)T_{r}. \tag{2}\]  
where \(T_{r}\), \(T_{c}\), \(T_{m1}\), and \(T_{m2}\) indicate the torque of the ring gear, carrier, and motors 1 and 2, respectively. The kinetic energy of the entire system \(E_{k}\), potential energy owing to the damper mechanism \(E_{u}\), and dissipation energy \(E_{d}\) can be expressed as follows:  
\[E_{k}=\frac{1}{2}\{l_{e}\omega_{e}^{2}+l_{e}\omega_{e}^{2}+l_{m1}\omega_{m1}^{2 }+l_{m2}\omega_{m2}^{2}+l_{p}\omega_{p}^{2}+n_{p}m_{p}(R_{s}+R_{p})^{2}\omega_{ c}^{2}\},\ E_{u}=\frac{1}{2}k_{dp}(\theta_{e}-\theta_{c})^{2},\ \ E_{d}=\frac{1}{2}c_{dp}(\omega_{e}-\omega_{c})^{2} \tag{3}\]  
where \(\theta\) is the rotation angle; \(I\) is the moment of inertia; \(k_{dp}\) and \(c_{dp}\) denote the torsional stiffness and damping coefficient of the damper, respectively; \(R_{s}\) and \(R_{p}\) denote the sun gear and pinion gear radius, respectively; \(m_{p}\) is the pinion gear mass; and \(n_{p}\) denotes the number of pinion gears. The subscript \(e\) used with \(\theta\), \(I\), and \(\omega\) indicates the engine crankshaft. In the Lagrangian equation, the Lagrangian is expressed as \(L=E_{k}-E_{u}\), the dissipation energy is denoted as \(D=E_{d}\), and the external force in this system is represented as \(f_{i}\). The rotational speed of each element is as follows: \(\omega_{e}=\dot{\theta}_{e}\), \(\omega_{c}=\dot{\theta}_{c}\), \(\omega_{m1}=\dot{\theta}_{m1}\), \(\omega_{m2}=\dot{\theta}_{m2}\), and \(\omega_{p}=\dot{\theta}_{p}\). The Lagrangian equation can be expressed as follows:  
\[\frac{d}{dt}\bigg{(}\frac{\partial L}{\partial\dot{\theta}_{i}}\bigg{)}-\frac {\partial L}{\partial\dot{\theta}_{i}}+\frac{\partial D}{\partial\dot{\theta}_ {i}}=f_{i}\ \ \ \ (i=1,2,3). \tag{4}\]  
In the planetary gear system, \(\omega_{m1}\) and \(\omega_{p}\) were eliminated using the relationship from Equation (1). The basic equations were implemented as a linear system with three degrees of freedom. The rotational speed and angle of each element \(i\) in Equation (4)  
Figure 2: The controller structure of the ESN-based MPC. The solid and dashed blue lines indicate an observed past disturbance and a predicted future disturbance using ESN, respectively. The control horizon \(H_{u}\) and prediction horizon \(H_{p}\) were applied to obtain the MPC optimization to evaluate the cost function. are set as follows: \(\dot{\theta}_{1}=\omega_{e}\), \(\dot{\theta}_{2}=\omega_{c}\), \(\dot{\theta}_{3}=\omega_{m2}\), \(\theta_{1}=\theta_{e}\), \(\theta_{2}=\theta_{c}\), and \(\theta_{3}=\theta_{m2}\). By solving Equations (1)-(4) for \(i=1\)-3, we obtain Equation (5), which is the equation of the system's motion:  
\[\left[ \begin{array}{c}l_{e}\dot{\omega}_{e}\\ l_{1}\dot{\omega}_{c}\\ l_{1}\dot{\omega}_{m2}\end{array} \right]=\left[ \begin{array}{cc}k_{dp}&c_{ dp}\\ l_{2}k_{dp}&l_{2}c_{dp}\\ l_{3}k_{dp}&l_{3}c_{dp}\end{array} \right]\left[ \begin{array}{c}\theta_{e}- \theta_{c}\\ \omega_{e}-\omega_{c}\end{array} \right]+\left[ \begin{array}{c}f_{1}\\ f_{2}\\ f_{3}\end{array}