The multimodeling structure is defined by a linear dynamical system that has one slow and two fast subsystems: \[ \begin{aligned} \dot x_0(t)&= A_{00}x_0(t) + A_{01}x_1(t) + A_{02}x_2(t) +B_{01}u_1(t)+B_{02}u_2(t),\\ \varepsilon_1\dot x_1(t)&=A_{10}x_0(t) + A_{11}x_1(t) +\varepsilon_3A_{12}x_2(t) +B_{11}u_1(t)+\varepsilon_3B_{12}u_2(t),\\ \varepsilon_2\dot x_2(t)&=A_{20}x_0(t) +\varepsilon_3A_{21}x_1(t) + A_{22}x_2(t) +\varepsilon_3B_{21}u_1(t)+B_{22}u_2(t),\end{aligned} \] where \(x_0\in \mathbb R^{n_0}\) represents slow state variables; \(x_1\in \mathbb R^{n_1}\), \(x_2\in \mathbb R^{n_2}\) are fast state variables; \(u_1\in \mathbb R^{m_1}\), \(u_2\in \mathbb R^{m_2}\) are control inputs; \(\varepsilon_3\) is a small weak coupling parameter; \(\varepsilon_1\) and \(\varepsilon_2\) are small positive singular perturbation parameters of the same order of magnitude. The quadratic performance criterion has to be minimized by the proper choice of the control variable \(u_1(t)\) and \(u_2(t)\). The authors propose a method for the exact decomposition of the optimal control such that the optimal solution is obtained in terms of reduced-order nonsymmetric one pure-slow and two pure-fast algebraic Riccati equations for which the Newton method is perfectly suited. This proposed approach is conceptually simpler and numerically more efficient than the ones previously used.