Question: local stability of ODE system

my system [u(t),v(t),w(t)] is highly non-linear. I want to understand the local behavior at [u=0,v=0,w=0]. The boundary conditions are as follows: u(0) is given but v(0) and w(0) are not. I want to know if there is a suitable choice of v(0) and w(0) such that the system will converge to the critical point [0,0,0].

Can this system be approximated by a linear system? Or by a two-dimensional system?

[
diff(u(t),t) = (1-u(t))*(u(t)*(1/40+v(t))-1/3*v(t)-1/3*w(t)), 
diff(v(t),t) = -3/2*u(t)*2^(1/2)/(3^(1/2)*2^(1/2)
              *(1/(3/40+v(t)))^(1/2))^(3/2)*(3/2*2^(1/2)*(3^(1/2)*2^(1/2)
              *(1/(3/40+v(t)))^(1/2))^(1/2)-1/40*3^(1/2)*2^(1/2)*(1/(3/40+v(t)))^(1/2)), 
diff(w(t),t) = -(w(t)-1/200)*(w(t)+1/200)+1/3*w(t)*(v(t)+w(t))/u(t)*(1-u(t))
];

Below I show what the 2-D fieldplots look like ...

Any help will be much appreciated.


u=0:

 

 

v=0:

 

w=0:

 

Can I steer my system towards the origin?

 


Below I show what I've done so far. The code may be copy-pasted into a Maple session.

# Parameter and Functional Assignments

params := [ A=3, d=1/20, s=1/3, a=1/2, b1=2/100, b2=3/100, w1=1/2, w2=1/2, b=1/40 ];
f   := x -> A*x^(a)-d*x:
fx  := x -> A*a*x^(a-1)-d:
fxx := x -> A*a*(a-1)*x^(a-2):
ifx := y -> (A*a/(d+y))^(1-a):
h   := y -> (fxx(ifx(y)))*(f(ifx(y))): h(y):

# Differential Equation System
udot := diff(u(t),t) = (1-u(t))*(u(t)*(b+v(t))-s*(v(t)+w(t))):
vdot := diff(v(t),t) = u(t)*h(b+v(t)):
wdot := diff(w(t),t) = -(w(t)+b1-b)*(w(t)+b2-b)+s*w(t)*(v(t)+w(t))/u(t)*(1-u(t)):
sys:=[udot,vdot,wdot]: eval(sys,params);

# Transformation
newsys := subs(u(t)=u,v(t)=v,w(t)=w,map(rhs,sys *~ u(t))):
newsys := simplify(newsys) assuming positive:
eval(newsys,params);

# Does the fieldplot provide clues?
plots:-fieldplot3d(
[eval(newsys,params)[1],eval(newsys,params)[2],eval(newsys,params)[3]],
u=-1e-2..1e-2, v = -1e-2..1e-2, w=-1e-2..1e-2,
fieldstrength=fixed,
axes=box
);

# Does the 2-D fieldplot provide clues?
plots:-fieldplot( [eval(eval(newsys,params),u=0)[2],eval(eval(newsys,params),u=0)[3]] ,v=-1e-2..1e-2, w = -1e-2..1e-2, fieldstrength=fixed);
plots:-fieldplot( [eval(eval(newsys,params),v=0)[1],eval(eval(newsys,params),v=0)[3]] ,u=-1e-2..1e-2, w = -1e-2..1e-2, fieldstrength=fixed);
plots:-fieldplot( [eval(eval(newsys,params),w=0)[1],eval(eval(newsys,params),w=0)[2]] ,u=-1e-2..1e-2, v = -1e-2..1e-2, fieldstrength=fixed);

# I can't even decide how to let [u,v,w] tend to [0,0,0], as the limit is clearly path-dependent:
Jac := VectorCalculus:-Jacobian(newsys,[u,v,w]):
e := LinearAlgebra:-Eigenvalues(Jac):
series('leadterm'( series('leadterm'( series('leadterm'( e[1]), u=0 )), v=0 )), w=0);
series('leadterm'( series('leadterm'( series('leadterm'( e[1]), v=0 )), u=0 )), w=0);
series('leadterm'( series('leadterm'( series('leadterm'( e[1]), v=0 )), w=0 )), u=0);
series('leadterm'( series('leadterm'( series('leadterm'( e[1]), w=0 )), v=0 )), u=0);
series('leadterm'( series('leadterm'( series('leadterm'( e[1]), u=0 )), w=0 )), v=0);
series('leadterm'( series('leadterm'( series('leadterm'( e[1]), w=0 )), u=0 )), v=0);
 

Please Wait...