1
$\begingroup$

I have a four dimensional dynamical system where I have a Hopf bifurcation. My aim is to figure out if limit cycles are stable or not.

The calibration in the system is

paramFinal = {ρ -> 0.85, u -> 0.95, h -> 60, a -> 1, b -> 0.1, δ -> 0.001, γ -> 0.005, ψ -> 60, θ -> 0.04586401037966653, β ->  0.8772372953238555};

KK4 = 33.5081
SS4 = 0.212842
lambdaSS4 = 0.19493
muSS4 = 4.02085
Subscript[a, 1] = 151.467

I have the following system :

kdot = (a - δ) k - (h (2 - u) - λ Subscript[a, 1]);
sdot = s (1 - s) - γ a k;
lambdadot = ((ρ + θ) - (a - δ)) λ + μ γ a;
mudot = (ρ + θ) μ - μ (1 - 2 s) - θ ψ b s^-1 ;

To be able to have Jacobian matrix, I write

Subscript[f, 11] = D[kdot, k];
Subscript[f, 12] = D[kdot, s];
Subscript[f, 13] = D[kdot, λ];
Subscript[f, 14] = D[kdot, μ]; 
Subscript[f, 21] = D[sdot, k];
Subscript[f, 22] = D[sdot, s];
Subscript[f, 23] = D[sdot, λ];
Subscript[f, 24] = D[sdot, μ];
Subscript[f, 31] = D[lambdadot, k];
Subscript[f, 32] = D[lambdadot, s];
Subscript[f, 33] = D[lambdadot, λ];
Subscript[f, 34] = D[lambdadot, μ];
Subscript[f, 41] = D[mudot, k];
Subscript[f, 42] = D[mudot, s];
Subscript[f, 43] = D[mudot, λ];
Subscript[f, 44] = D[mudot, μ];

To write the fourth order characteristic polynominal equation, I write down the the sum of minors of matrix and determinant.

Subscript[b, 1] = (Subscript[f, 11] + Subscript[f, 22] + Subscript[f, 33] + Subscript[f, 44]) /. paramFinal /. k -> KK4 /. s ->  SS4 /. μ ->  muSS4 /. λ ->   lambdaSS4

Subscript[b, 2] = Det[{{Subscript[f, 11], Subscript[f, 12]}, {Subscript[f, 21], Subscript[f, 22]}}] + Det[{{Subscript[f, 11], Subscript[f, 13]}, {Subscript[f, 31], Subscript[f, 33]}}] + Det[{{Subscript[f, 11], Subscript[f, 14]}, {Subscript[f, 41], Subscript[f, 44]}}] + Det[{{Subscript[f, 22], Subscript[f, 23]}, {Subscript[f, 32], Subscript[f, 33]}}] + Det[{{Subscript[f, 22], Subscript[f, 24]}, {Subscript[f, 42], Subscript[f, 44]}}] + Det[{{Subscript[f, 33], Subscript[f, 34]}, {Subscript[f, 43], Subscript[f,44]}}]/. paramFinal /. k -> KK4 /. s ->  SS4 /. μ ->  muSS4 /. [Lambda] ->   lambdaSS4;

Subscript[b, 3] = (Det[{{Subscript[f, 22], Subscript[f, 23], Subscript[f, 24]}, {Subscript[f, 32], Subscript[f, 33], Subscript[f, 34]}, {Subscript[f, 42], Subscript[f, 43], Subscript[f, 44]}}] + Det[{{Subscript[f, 11], Subscript[f, 13], Subscript[f, 14]}, {Subscript[f, 31], Subscript[f, 33], Subscript[f, 34]}, {Subscript[f, 41], Subscript[f, 43], Subscript[f, 44]}}] + Det[{{Subscript[f, 11], Subscript[f, 12], Subscript[f, 14]}, {Subscript[f, 21], Subscript[f, 22], Subscript[f, 24]}, {Subscript[f, 41], Subscript[f, 42], Subscript[f, 44]}}] +Det[{{Subscript[f, 11], Subscript[f, 12], Subscript[f, 13]}, {Subscript[f, 21],Subscript[f, 22], Subscript[f, 23]}, {Subscript[f, 31], Subscript[f, 32],Subscript[f, 33]}}]) /. paramFinal /. k -> KK4 /. s ->  SS4 /. μ ->muSS4/. λ -> lambdaSS4


Subscript[b, 4] = Det[{{Subscript[f, 11], Subscript[f, 12], Subscript[f, 13], Subscript[f, 14]}, {Subscript[f, 21], Subscript[f, 22], Subscript[f, 23], Subscript[f, 24]}, {Subscript[f, 31], Subscript[f, 32], Subscript[f, 33], Subscript[f, 34]}, {Subscript[f, 41], Subscript[f, 42], Subscript[f, 43],Subscript[f, 44]}}] /. paramFinal /. k -> KK4 /. s ->  SS4 /. μ -> muSS4 /. λ ->   lambdaSS4

Then, I compute the eigenvalues.

Solve[λ^4 + Subscript[b, 1] λ^3 + Subscript[b, 2] λ^2 + Subscript[b, 3] λ + Subscript[b, 4] == 0, {λ}

There are four roots, where a pair of eigenvalues are purely imaginary and two others are conjugate with real and complex parts, which shows the existence of limit cycles. My problem is that I can not figure out how to find numerically Lyapunov number to know if cycles are stable or not. Any hints or suggestions ? or another way to know if cycles are stable or not ?

$\endgroup$
1
  • 1
    $\begingroup$ May I suggest that you use f[n] instead of Subscript[f, n] in the question. It is more compact and occasionally more reliable. Also, consider including a plot of the bifurcation diagram. $\endgroup$
    – bbgodfrey
    May 18, 2017 at 20:05

1 Answer 1

2
$\begingroup$

I recommend that you review the Liu criterion for the Hopf bifurcation, since it is the most suitable for detecting limit cycles in systems with three or more equations.

It is more convenient to write your system and the Jacobian matrix as follows:

f1 = (a - δ) k - (h (2 - u) - λ a1);
f2 = s (1 - s) - γ a k;
f3 = ((ρ + θ) - (a - δ)) λ + μ γ a;
f4 = (ρ + θ) μ - μ (1 - 2 s) - θ ψ b s^-1;
F = {f1, f2, f3, f4};
V = {k, s, λ, μ};

J = FullSimplify[D[F, {V}]]
$\endgroup$
0

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.