MATLAB Answers

Coupled ODE with ode45

250 views (last 30 days)
Rick
Rick on 2 Nov 2014
Commented: santosh kumar on 5 Oct 2017
Hello, I am trying to solve these two coupled differential equations, but I can't seem to get it to work. I always have difficulty using ODE45...but why isn't the variable X being recognized? The image shows the differential equations I am trying to solve using MATLAB
alpha = 0.001; % L^-1
FA0 = 2.5; % mol/min
CA0 = 0.305; % mol/L
eps = 2;
k = 0.044;
CA = @(y,X) CA0*(1-X).*y./(1+eps*X);
rA = @(y,X) -k*CA(y,X);
vSpan = [0 500]; % L
dXdV = @(V,X) -rA(y,X)/FA0;
dydV = @(V,y) -alpha*(1+eps.*X)./(2*y);
[V, Y] = ode45(dydV,vSpan,1);
[V, X] = ode45(dXdV,vSpan,0);
But this is the error I encounter
Undefined function or variable 'X'.
Error in @(V,y)-alpha*(1+eps.*X)./(2*y)
Error in odearguments (line 88)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn,
...

Accepted Answer

Zoltán Csáti
Zoltán Csáti on 2 Nov 2014
I found the mistake. It was a sign. :) The correct function is
function dz = myode(v,z)
alpha = 0.001;
C0 = 0.3;
esp = 2;
k = 0.044;
f0 = 2.5;
dz = zeros(2,1);
dz(1) = k*C0/f0*(1-z(1)).*z(2)./(1+esp*z(1));
dz(2) = -alpha*(1+esp*z(1))./(2*z(2));
  5 Comments
santosh kumar
santosh kumar on 5 Oct 2017
to define a function, use "syms function" example: to define x y z syntax: "syms x y z"

Sign in to comment.

More Answers (2)

Zoltán Csáti
Zoltán Csáti on 2 Nov 2014
On line dydV = @(V,y) -alpha*(1+eps.*X)./(2*y); X is not included. I would solve it as a coupled system or solve it analytically by hand. I don't see the impact of C in the differential equations. Is C needed?
  1 Comment
Rick
Rick on 2 Nov 2014
the Concentration (CA) is needed for rA, which is in the dXdV differential equation

Sign in to comment.


Zoltán Csáti
Zoltán Csáti on 2 Nov 2014
Code your function as follows.
function dz = myode(v,z)
alpha = 0.001;
C0 = 0.3;
esp = 2;
k = 0.044;
f0 = 2.5;
dz = zeros(2,1);
dz(1) = k*C0/f0*(1-z(1)).*z(2)./(1-esp*z(1));
dz(2) = -alpha*(1+esp*z(1))./(2*z(2));
Then call it as
[v z] = ode45(@myode,[0 500],[0 1]);
  2 Comments
Zoltán Csáti
Zoltán Csáti on 2 Nov 2014
Yes, I plotted it. It is very sensitive to the parameters. Try different parameter values.

Sign in to comment.

Tags

Community Treasure Hunt

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

Start Hunting!