How can I make ode45 work for me?

Trying unsuccessfully to make this work:
main test.m file
clear all;
clc;
Kc=1;
Dab=1;
Dl=1;
N=1000;
Cp0=ones(1,N);
Cf0=zeros(1,N);
t0=1;
tf=300;
[T,C]=ode45(@(t,c) dcdt(t,c,Kc,Dab,Dl),[t0 tf],[Cp0 Cf0]);
function dcdt.m file
function dcdt=dcdt(t,c,Kc,Dab,Dl)
%Constants
N=1000;
r=2.5*10^-3;
delr=r/N;
z=5;
delz=z/N;
ez=0.4;
u=2.24*10^-2;
dcdt=zeros(4,N);
%Node R=0
dcdt(1,1)=6*Dab/2*(c(1,2)-c(1,1));
%Nodes R=2 to N-1
for i=2:(N-1)
dcdt(2,i)=(Dab/delr^2+Dab/(r*delr))*c(1,i+1)...
-2*(Dab/delr^2)*c(1,i)...
+(Dab/delr^2-Dab/(r*delr))*c(1,i-1);
end
%Node R=N
dcdt(3,N)=((-2*Dab/(delr^2))*(4*r^2/delr^2-4*r/delr+1)*c(1,N-1)+...
(2*Dab/(delr^2)*(4*r^2/delr^2-4*r/delr+1)-Kc*(4*r^2/delr^2))*c(1,N)+...
Kc*(4*r^2/delr^2)*c(2,1))/(7*r^3/(6*delr^2)-r^2/(2*delr)+r/2-delr/6);
%Node R>N
for i=1:N
dcdt(4,i)=(Dl/delz^2-u/2*delz)*c(2,i+1)...
-(2/delz^2+(1-e/e))*(3*kc/r)*c(2,i)+...
(Dl/delz^2+u/2*delz)*c(2,i-1)...
+(1-ez/ez)*(3*Kc/r)*c(1,N);
end
Running the main file produces this:
??? Attempted to access c(1,2); index out of bounds because size(c)=[1000,1].
Error in ==> dcdt at 14 dcdt(1,1)=6*Dab/2*(c(1,2)-c(1,1));
Error in ==> @(t,c)dcdt(t,c,Kc,Dab,Dl)
Error in ==> odearguments at 110 f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ==> ode45 at 173 [neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in ==> test at 15 [T,C]=ode45(@(t,c) dcdt(t,c,Kc,Dab,Dl),[t0 tf],Cp0)
I have two variables C(1,) and C(2,). My main goal is obtaining C in two dimensions: space and time, with c(1,2) being the 1st variable in the second position in whatever time it simulated. I have scoured the web and still cannot find a solution to this. Don't see anything wrong in the ode45 syntax. Not sure how can I make this work.
Thanks!

2 Comments

Hello Red
I will recommend to start using breakpoints into the function dcdt.m, to identify problems. The variable "c" is a column vector of 2000 row and 1 column. In this case the line 12 in your code is not correct.
dcdt(1,1)=6*Dab/2*(c(1,2)-c(1,1));
In that case c(1,2) doesn't exist. If c(1,2) doesn't exist, there is no value for dcdt(1,1) (second error).
And the error continue from here. Errors are related with code implementation, not ode45 function.
Hoper it helps. Best regards
Javier
Red
Red on 21 Sep 2012
Edited: Red on 21 Sep 2012
Hello Javier, thanks for your help and suggestions.
I am new to working with ode45 but have some experience in coding.
As much as I get what you mean, my aim is to make c a (2,1000) matrix, not a (2000,1) matrix. I have gone through several examples of ode45 implementation which results in the function code.
I know that the initial error will lead to more errors but I can't really work out where the code is implemented wrongly which brings me to seek help here. However my final aim is written in the original question which is to use ode45 to obtain 2 C variables for each point in space and time. I believe the desired outcome would be a three dimensional matrix of C(2,N,t) or something similar.
Hope that this could highlight my plight as I am at my wits end towards solving this.
Once again, thanks and regards, Red

Sign in to comment.

Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Tags

Asked:

Red
on 18 Sep 2012

Community Treasure Hunt

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

Start Hunting!