Problem solving a system of differential equations

6 views (last 30 days)
I am trying to solve a system of differential equations in Matlab.
dn/du=(-2*u*n-K*(n*u-(1+g)))/(1+u^2+K*u*(u-(1+g)/n))
dxi/du=(1-u^2)/(1+u^2+K*u*(u-(1+g)/n))
df/du=(2*u+K*u^2*(u-(1+g)/n))/(1+u^2+K*u*(u-(1+g)/n))
K and gamma are constants. and i have the initial conditions ( n(0)=1, xi(0)=0, f(0)=0 )
First i tried to use 'dsolve' as shown below.
g=0.1;
K=3;
syms n(u) u
n(u)=dsolve(diff(n,u)== (-2*u*n-K*(n*u-(1+g)))/(1+u^2+K*u*(u-(1+g)/n)),n(0)==1)
syms x(u) u n
x(u)= dsolve(diff(x,u)== (1-u^2)/(1+u^2+K*u*(u-(1+g)/n)),x(0)==0)
syms f(u) u n
f(u)=dsolve(diff(f,u)== (2*u+K*u^2*(u-(1+g)/n))/(1+u^2+K*u*(u-(1+g)/n)),f(0)==0)
from this i got that there is not explicit solution for the first equation and returned an [empty sym] as shown below.
Warning: Explicit solution could not be found.
In dsolve (line 201)
In Untitled (line 37)
n(u) = [ empty sym ]
x(u) = (1089*atan(33/(1600*n^2 - 1089)^(1/2) - (80*n*u)/(1600*n^2 - 1089)^(1/2)))/(160*n*(1600*n^2 - 1089)^(1/2)) - (35937*log(40*n*u^2 - 33*u + 10*n))/(2*(- 256000*n^3 + 174240*n)) - u/4 + (33*log(10*n)*(1600*n^2 - 1089)^(1/2) - 2178*atan(33/(1600*n^2 - 1089)^(1/2)) + 8000*n^2*atan(33/(1600*n^2 - 1089)^(1/2)))/(320*n*(1600*n^2 - 1089)^(1/2)) - (25*n*atan(33/(1600*n^2 - 1089)^(1/2) - (80*n*u)/(1600*n^2 - 1089)^(1/2)))/(1600*n^2 - 1089)^(1/2) + (26400*n^2*log(40*n*u^2 - 33*u + 10*n))/(- 256000*n^3 + 174240*n)
f(u) = (3*u^2)/8 - (33*u)/(160*n) - (log(40*n*u^2 - 33*u + 10*n)*(3200000*n^4 - 3920400*n^2 + 1185921))/(2*(- 10240000*n^4 + 6969600*n^2)) + (log(10*n)*(3200000*n^4 - 3920400*n^2 + 1185921))/(2*(- 10240000*n^4 + 6969600*n^2)) + (33*atan((1089*(2800*n^2 - 1089))/((1600*n^2 - 1089)^(1/2)*(92400*n^2 - 35937)))*(2800*n^2 - 1089))/(6400*n^2* (1600*n^2 - 1089)^(1/2)) - (33*atan((1089*(2800*n^2 - 1089))/((1600*n^2 - 1089)^(1/2)*(92400*n^2 - 35937)) - (2640*nu(2800*n^2 - 1089))/((1600*n^2 - 1089)^(1/2)*(92400*n^2 - 35937)))*(2800*n^2 - 1089))/(6400*n^2*(1600*n^2 - 1089)^(1/2))
Then i tried to use 'ode45'.
u=0:0.1:1;
K=3;
g=0.1;
dndu=@(u,n) (-2*u*n-K*(n*u-(1+g)))/(1+u^2+K*u*(u-(1+g)/n));
[u,n]=ode45(dndu, u, 1) % initial n=1
dxidu=@(u,xi) (1-u^2)/(1+u^2+K*u*(u-(1+g)/n));
[u,xi]=ode45(dxidu, u, 0); %initial xi=0
dfdu=@(u,f) (2*u+K*u^2*(u-(1+g)/n))/(1+u^2+K*u*(u-(1+g)/n));
[u,f]=ode45(dfdu, u, 0) %initial f=0
K and g are constants and u gets values for 0 to 1. When i run this the first equation is solved and got an answer for u and n but i get an error for the second equation that matrix dimension must agree as shown below.
u = 0 0.1000 0.2000 0.3000 0.4000 0.5000 0.6000 0.7000 0.8000 0.9000 1.0000
n = 1.0000 1.3409 1.6243 1.7746 1.7945 1.7263 1.6119 1.4800 1.3470 1.2208 1.1048
Error using /
Matrix dimensions must agree.
Error in @(u,xi)(1-u^2)/(1+u^2+Ku(u-(1+g)/n))
Error in odearguments (line 87)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 113)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in ode45method (line 24)
[u,xi]=ode45(dxidu, u, 0); %initial xi=0
Can anyone help me how to approach this type of problem? Thanks in advance!

Accepted Answer

Torsten
Torsten on 8 Jul 2015
uspan=0:0.1:1;
K=3;
g=0.1;
dydu=@(u,y)[(-2*u*y(1)-K*(y(1)*u-(1+g)))/(1+u^2+K*u*(u-(1+g)/y(1)))
(1-u^2)/(1+u^2+K*u*(u-(1+g)/y(1)))
(2*u+K*u^2*(u-(1+g)/y(1)))/(1+u^2+K*u*(u-(1+g)/y(1)))];
[U,Y]=ode45(dydu,uspan,[1 0 0]);
plot(U,Y(:,1),U,Y(:,2),U,Y(:,3));
Best wishes
Torsten.
  3 Comments
Torsten
Torsten on 8 Jul 2015
xi and f depend on n. This means that you have to solve the three equations simultaneously.
The reason for the error you encounter is that in your second and third equation you divide by the complete vector n which is wrong. You have to divide by the value of n at the corresponding value of u.
Best wishes
Torsten.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!