7
$\begingroup$

For didactic purpose I investigate the stability of the ODEs system {x'[t] == x[t] + y[t], y'[t] == -5 x[t] - 2 y[t]}.

Its stream plot

StreamPlot[{x + y, -5 x - 2 y}, {x, -10, 10}, {y, -10, 10}]

enter image description here

clearly shows a spiral sink at the origin.

How to prove it, making use of Lyapunov function? First of all, how to explicitly write down that Lyapunov function for the system under consideration with help of LyapunovSolve? I don't find it here. The first example in the Applications section of the documentation to LyapunovSolve is useful, but it produces the conclusion, not a Lyapunov function.

Improve this question
$\endgroup$

2 Answers 2

Reset to default
6
$\begingroup$

Maybe this?:

ClearAll[v, x, y];
a = CoefficientArrays[{x + y, -5 x - 2 y}, {x, y}] // Last;
p = LyapunovSolve[Transpose[a], -IdentityMatrix[2]];
v[x_, y_] =  p . {x, y} . {x, y} // Simplify (* Lyapunov function *)

(* Check "energy dissipation" dv/dt *)
Dt[v[x, y]] /. Thread[Dt[{x, y}] -> {x + y, -5 x - 2 y}] // Simplify

(*  1/6 (32 x^2 + 14 x y + 5 y^2)  -- Lyapunov function *)

(*  -x^2 - y^2  -- dv/dt is neg. def. *)
Improve this answer
$\endgroup$
4
  • 1
    $\begingroup$ As usually, very good. $\endgroup$
    – user64494
    Jun 2, 2021 at 19:05
  • $\begingroup$ @Michael E2 could you do this to find a Lyapunov function for "any" system? $\endgroup$
    – Math
    Mar 17, 2022 at 14:47
  • $\begingroup$ @Math I doubt it. It fails on this: {x'[t] == x[t]^3 + y[t], y'[t] == -5 x[t] - 2 y[t]^3}. But it succeeds here: ode = {x'[t] == Sin[x[t]] + y[t], y'[t] == -5 x[t] - 2 Sin[y[t]]}; vf = First@SolveValues[ode, {x'[t], y'[t]}]; a = D[vf, {{x[t], y[t]}}] /. {x[t] -> 0, y[t] -> 0}; p = LyapunovSolve[Transpose[a], -IdentityMatrix[2]]; v[x_, y_] = p . {x, y} . {x, y} // Simplify. I don't think there's any guarantee that there exists a quadratic Lyapunov function. $\endgroup$
    – Michael E2
    Mar 17, 2022 at 16:51
  • $\begingroup$ @MichaelE2 Thank you! Can you please have a look at this question: mathematica.stackexchange.com/questions/265288/… $\endgroup$
    – Math
    Mar 18, 2022 at 11:35
3
$\begingroup$

Stability of {x',y'}]== f[x,y] can be investigated by searching for critical points of f (see e.g: https://en.wikipedia.org/wiki/Stability_theory). In your case this is obvious {0,0}. Then in the general case f is linearized at the critical points, again in your case that is not needed, as f is already linear. Then the Eigensystem is

mat= m = {{1, 1}, {-5, -2}};
{x'[t], y'[t]} = mat . {x[t], y[t]}
{eigenval, eigenvec}=Eigensystem[{{1, 1}, {-5, -2}}] // N

{-0.5 + 1.65831 I, -0.5 - 1.65831 I}, {{-0.3 - 0.331662 I, 
   1.}, {-0.3 + 0.331662 I, 1.}}}

The stability depends on the real parts of the eigenvalues. In your case mat is real, what makes the 2 eigenvalues conjugate complex and there is only one real part of -0.5.

If the real part is negative, as in your case, we have an attractive critical point, called a sink. As the imaginary part is different from zero, it is a spiral sink. The trajectory will approach zero the faster, the larger the absolute value of the eigenvalue.

Improve this answer
$\endgroup$
2
  • 1
    $\begingroup$ Sorry, this is not it. The question is about the construction of Lyapunov function with help of Mathematica. Hope I am and was clear. $\endgroup$
    – user64494
    Jun 2, 2021 at 18:40
  • $\begingroup$ You have a very simple problem that does not need the Liapunov coefficient. This coefficient is needed if you want to investigate the rate of separation of infinitesimally close trajectories. .The Liapunov theory is more evolved than the simple stability analysis. $\endgroup$
    – Daniel Huber
    Jun 2, 2021 at 18:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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