how to find equlibrium point of 5 non linear system with numerical method

Question

Latifah 2023-4-2

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1939779-how-to-find-equlibrium-point-of-5-non-linear-system-with-numerical-method

评论： Latifah 2023-4-2

I have a system of 5 non linear ordinary differential equations with variable coefficients . I am trying to find the equilibrium points by hand but it seems like it is not possible without the help of a numerical method. What would be a good method to calculate equilibrium points of the system?

Another question (somehow related to the problem above): Would it be possible to check the stability of the equilibrium points and then draw a bifurcation diagram? If so, please suggest some way out!

 e1 = 13;
    g  = 0.0125;
    h  = 0.284253;
    f  = 0.05;
    q1 = 0;
    k1 = 10;
    d1 = 0.0412;
    e2 = 0.0188;
    j  = 0.0082;
    q2 = 0;
    k2 = 10;
    d2 = 0.0288;
    b  = 2;
    d4 = 0.1152;
    e3 = 0.166667;
    a  = 1.7;
    q3 = 0;
    k3 = 10;
    d3 = 0.1152;
    r  = 0.5;
    m  = 1.02;
    q4 = 0;
    k4 = 10;
    dydt(1) = e1 + (g*y(3)*y(1)/(h + y(3))) + (f*y(3)*y(1)) - (y(1)*(1 + ((q1/k1)*y(1)))) - (d1*y(1));
    dydt(2) = e2*y(2) + (f*y(1)*y(3)) - (j*y(2)) - (y(2)*(1 + ((q2/k2)*y(2)))) - (d2*y(2));
    dydt(3) = (b*e2*y(2)) - (d4*y(3));
    dydt(4) = e3*y(4) + (j*y(2)) - (a*y(4)) - (y(4)*(1 + ((q3/k3)*y(4)))) - (d3*y(4));
    dydt(5) = (r*y(5)*(1 - (m*y(5)))) + (a*y(4)) - y(5)*(1 + ((q4/k4)*y(5)));

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Torsten 2023-4-2

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1939779-how-to-find-equlibrium-point-of-5-non-linear-system-with-numerical-method#answer_1206954

在 MATLAB Online 中打开

I don't know if the result is as expected.

format long
yequi = fsolve(@fun,rand(5,1))
Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.
yequi = 5×1
  12.485593545908634
   0.000000000000008
   0.000000000000003
   0.000000000000000
   0.000000000000000
fun(yequi)
ans = 1×5
1.0e-13 *

  -0.657252030578093  -0.063585506384309  -0.000000000537866   0.000000000086870  -0.000243217417049
function dydt = fun(y)
e1 = 13;
g  = 0.0125;
h  = 0.284253;
f  = 0.05;
q1 = 0;
k1 = 10;
d1 = 0.0412;
e2 = 0.0188;
j  = 0.0082;
q2 = 0;
k2 = 10;
d2 = 0.0288;
b  = 2;
d4 = 0.1152;
e3 = 0.166667;
a  = 1.7;
q3 = 0;
k3 = 10;
d3 = 0.1152;
r  = 0.5;
m  = 1.02;
q4 = 0;
k4 = 10;
dydt(1) = e1 + (g*y(3)*y(1)/(h + y(3))) + (f*y(3)*y(1)) - (y(1)*(1 + ((q1/k1)*y(1)))) - (d1*y(1));
dydt(2) = e2*y(2) + (f*y(1)*y(3)) - (j*y(2)) - (y(2)*(1 + ((q2/k2)*y(2)))) - (d2*y(2));
dydt(3) = (b*e2*y(2)) - (d4*y(3));
dydt(4) = e3*y(4) + (j*y(2)) - (a*y(4)) - (y(4)*(1 + ((q3/k3)*y(4)))) - (d3*y(4));
dydt(5) = (r*y(5)*(1 - (m*y(5)))) + (a*y(4)) - y(5)*(1 + ((q4/k4)*y(5)));
end

4 个评论
显示 2更早的评论隐藏 2更早的评论

Torsten 2023-4-2

编辑：Torsten 2023-4-2

在 MATLAB Online 中打开

I left it to you to glue the solutions together.

syms e1 g h f q1 k1 d1 e2 f j q2 k2 d2 b e2 d4 e3 a q3 k3 d3 r m q4 k4
y = sym('y',[5 1]);
dydt(1) = e1 + (g*y(3)*y(1)/(h + y(3))) + (f*y(3)*y(1)) - (y(1)*(1 + ((q1/k1)*y(1)))) - (d1*y(1));
dydt(2) = e2*y(2) + (f*y(1)*y(3)) - (j*y(2)) - (y(2)*(1 + ((q2/k2)*y(2)))) - (d2*y(2));
dydt(3) = (b*e2*y(2)) - (d4*y(3));
dydt(4) = e3*y(4) + (j*y(2)) - (a*y(4)) - (y(4)*(1 + ((q3/k3)*y(4)))) - (d3*y(4));
dydt(5) = (r*y(5)*(1 - (m*y(5)))) + (a*y(4)) - y(5)*(1 + ((q4/k4)*y(5)));
y2 = solve(dydt(3)==0,y(2))
y1 = simplify(solve(subs(dydt(2),y(2),y2)==0,y(1)))
y3 = solve(subs(dydt(1),y(1),y1)==0,y(3),'MaxDegree',3)
y31 = y3(1);
y32 = y3(2);
y33 = y3(3);
y11 = subs(y1,y(3),y31);
y12 = subs(y1,y(3),y32);
y13 = subs(y1,y(3),y33);
y21 = subs(y2,y(3),y31);
y22 = subs(y2,y(3),y32);
y23 = subs(y2,y(3),y33);
y41 = solve(subs(dydt(4),y(2),y21)==0,y(4))
y411 = y41(1);
y412 = y41(2);
y42 = solve(subs(dydt(4),y(2),y22)==0,y(4))
y421 = y42(1);
y422 = y42(2);
y43 = solve(subs(dydt(4),y(2),y23)==0,y(4))
y431 = y43(1);
y432 = y43(2);
y511 = solve(subs(dydt(5),y(4),y411)==0,y(5))
y512 = solve(subs(dydt(5),y(4),y412)==0,y(5))
y521 = solve(subs(dydt(5),y(4),y421)==0,y(5))
y522 = solve(subs(dydt(5),y(4),y422)==0,y(5))
y531 = solve(subs(dydt(5),y(4),y431)==0,y(5))
y532 = solve(subs(dydt(5),y(4),y432)==0,y(5))

Torsten 2023-4-2

Is it possible for the equilibrium point to be negative?

It's the numerical error in the equation fun(yequi) that is negative, not the solution itself (yequi) (if this is what you mean with your question).

Latifah 2023-4-2

thank you , it helps me a lot

请先登录，再进行评论。

how to find equlibrium point of 5 non linear system with numerical method

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

4 个评论
显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

how to find equlibrium point of 5 non linear system with numerical method

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

4 个评论 显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论