using ode45 for coupled equations

I am attempting to simplify my forward euler solution to an ODE solution, but I don't know how to "call" another script to solve.

Here is my forward euler:

dt=.01;
tmin=5;
tmax=10;
t=tmin:dt:tmax;
n=length(t);
z=zeros(1,n);
x=zeros(1,n);
x(1)=1;
z(1)=1;

for i=1:n-1
    x(i+1)=x(i)+dt*(t(i)+x(i)+z(i));
    z(i+1)=z(i)+dt*(t(i)+x(i)-(3*x(i)));
end

And my current script for ode45:

 tspan = [5 10];
 t0 = 5;
 dt=.01;
 tmin=5;
 tmax=10;
 t=tmin:dt:tmax;
 n=length(t);
 z=zeros(1,n);
 x=zeros(1,n);
 x(1)=1;
 z(1)=1;
 
 [t,x] = ode45(@(t,x)myfun, tspan, t0);

Answer 1

[James Tursa]

You are mixing syntax here:

@(t,x)myfun

Either use the syntax @(t,x)expression, where expression is the derivative

or use the syntax @myfun, where myfun is a function that returns the derivative.

E.g., using y as the input state where we define y(1) = x and y(2) = z

@(t,y)[t+y(1)+y(2);t+y(1)-(3*y(1))]

But if that last x(i) in your z derivative was supposed to be z(i), then it would be this:

@(t,y)[t+y(1)+y(2);t+y(1)-(3*y(2))]

Or have a separate function for this in a file called myfun.m and use the @myfun syntax:

function dy = myfun(t,y);
dy = [t+y(1)+y(2);t+y(1)-(3*y(2));  % either y(1) or y(2) for that last term, you need to check this
end

Finally, that last input argument for ode45 needs to be the initial values of x and z, or in this case y(1) and y(2). So that last argument should be [1;1], not 5.

[James Tursa]

Kolleggerm1's Answer moved here:

Thank you for clearing that up. I didn't realize that I didn't need to call another script. Now that I have adjusted my code to reflect that, the only issue is in line 14 (y(1)=x).

The error says that the indices on the left are not compatible with the right. Any thoughts on clearing that up?

 tspan = [5 10];
 t0 = 5;
 dt=.01;
 tmin=5;
 tmax=10;
 t=tmin:dt:tmax;
 n=length(t);
 z=zeros(1,n);
 x=zeros(1,n);
 x(1)=1;
 z(1)=1;
 y(1)=x;
 y(2)=z;
 
[t,x]=ode45(@(t,y)[t+y(1)+y(2);t+y(1)-3*y(1)],tspan,t0);

[James Tursa]

The ode45( ) function will create the outputs ... these are not something that you pre-allocate ahead of time. Also, you still don't have that last input argument correct. Instead of t0, it should be the initial value of y. So something like this:

 tspan = [5 10];
 t0 = 5;
 dt=.01;
 tmin=5;
 tmax=10;
 t=tmin:dt:tmax;
 n=length(t);
% z=zeros(1,n); ode45 will create this for you
% x=zeros(1,n); ode45 will create this for you
 x0=1; % changed nomenclature
 z0=1; % changed nomenclature
 y0(1)=x0; % changed nomenclature
 y0(2)=z0; % changed nomenclature
 
[t,y]=ode45(@(t,y)[t+y(1)+y(2);t+y(1)-3*y(1)],tspan,y0); % using y and y0

Then the y output variable will contain both the x and z values.

And, I would still advise you to double check that z derivative. It looks strange to have y(1) appearing twice and I strongly suspect that second y(1) should be y(2) instead.

Answer 2

[Kolleggerm1]

Thank you for your help with this. I understand the confusion about y(1) appearing twice and will rechck my derivatives. This was my first try ever deriving and working with coupled equations so I don't doubt that I may need to practice.

Thank you!