Consider the 2-dimensional optimal control problem of the LQR kind $$ \min_u \int_0^\infty (x^T Q x + u^TRu) \, dt \quad\text{such that}\quad \begin{cases}\dot x(t) = Ax(t)+Bu(t) \\ x(0) = \begin{pmatrix}1\\-1\end{pmatrix}\end{cases} $$ with $x=\begin{pmatrix}x_1 \\ x_2\end{pmatrix}$, $u=\begin{pmatrix}u_1\\u_2\end{pmatrix}$, $Q=B=I$ (identity), $R=\gamma I$, $A=\begin{pmatrix}-\alpha & \alpha\\\beta & -\beta\end{pmatrix}$, $\alpha>\beta>0$ and $\gamma>0$.

Solve it using the stationary Riccati equation $$ 0 = Q+A^TS+SA-SBR^{-1}BS,\quad \text{with}\quad S=\begin{pmatrix}s_1 & s_2 \\ s_2 & s_3\end{pmatrix}. $$

In order to find $x$ and $u$ we have to combine the equation for $\dot x$ with the equation for the optimal control $u = -R^{-1}B^TSx = -\frac1\gamma Sx$. So, firstly, we have to find the expression for $S$.

The Riccati equation is simplified as $$ 0 = I+A^TS+SA-\frac1\gamma S^2 $$ which is equivalent to the following nonlinear system $$ \begin{cases} -\dfrac{s_1^2}{\gamma}-2\alpha s_1 -\dfrac{s_2^2}{\gamma}+2\beta s_2+1=0\\ \alpha s_1 - \alpha s_2 - \beta s_2 + \beta s_3 - \dfrac{s_1s_2}{\gamma}-\dfrac{s_2s_3}{\gamma}=0\\ -\dfrac{s_2^2}{\gamma}+2\alpha s_2 - \dfrac{s_3^2}{\gamma} - 2\beta s_3+1 = 0 \end{cases} $$

How to solve such a non-linear system to find expressions for $s_1, s_2$ and $s_3$? I think there is a fast way to solve either the system or directly the Riccati equation in matrix form, but I don't know how.

The solutions provided by Matlab, using the code below, are so long that it let me think is not the correct way to solve the problem

syms x y z a b g
eqn1 = 0 == -x^2/g-2*a*x-y^2/g+2*b*y+1;
eqn2 = 0 == a*x-a*y-b*y+b*z-x*y/g-y*z/g;
eqn3 = 0 == -y^2/g+2*a*y-z^2/g-2*b*z+1;
[x,y,z] = solve([eqn1, eqn2, eqn3], [x, y, z])

The solution of a Riccati equation can be found by determining the eigenvectors of the Hamiltonian; see e.g. here on Wikipedia. This is your best hope for a closed-form symbolic solution, in my view. In your case, the Hamiltonian is [ \begin{bmatrix} -\alpha & \alpha & -\gamma^{-1} & 0\\ \beta & -\beta & 0 & -\gamma^{-1}\\ -1 & 0 & \alpha & -\beta\\ 0 & -1 & -\alpha & \beta \end{bmatrix}. ] I don't see immediately how to diagonalize it, but maybe someone here has a more trained eye than me.

  • $\begingroup$ Thanks, I did not know about this method. Calling H the Hamiltonian matrix, the diagonalization can be performed in matlab with [V,D] = eig(H), where V is the matrix of eigenvectors and D the diagonal matrix of eigenvalues. However, also in this case the values computed are very long. $\endgroup$ – sound wave Mar 23 '20 at 18:21
  • $\begingroup$ Using V and D, how can I construct the matrices $U_1$ and $U_2$ cited in the page you linked? In particular, see i.imgur.com/WbfKw7I.png $\endgroup$ – sound wave Mar 23 '20 at 22:31
  • 1
    $\begingroup$ $\begin{bmatrix}U_1 \\ U_2\end{bmatrix}\in\mathbb{C}^{2n\times n}$ ($n=2$ in your case) is obtained by concatenating horizontally the $n$ eigenvectors that are associated to eigenvalues in the left half-plane. In (succinct) Matlab, U = V(:, real(diag(D)) < 0). If you have eigenvalues with real part zero, then it is trickier (but that shouldn't happen in your case). $\endgroup$ – Federico Poloni Mar 24 '20 at 7:40
  • 1
    $\begingroup$ If the eigenvalues are complex, $U$ (as we constructed it) is complex. If the eigenvalues are real, $U$ is real. That wiki page is assuming that $H$ may have complex eigenvalues, but not imaginary eigenvalues (that's two different things). $\endgroup$ – Federico Poloni Mar 25 '20 at 7:30
  • 1
    $\begingroup$ Depending what you mean by "solve the problem", yes, you might have to do it. $\endgroup$ – Federico Poloni Mar 25 '20 at 16:01

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.