using ode45 for coupled equations

1 次查看（过去 30 天）
Kolleggerm1
Kolleggerm1 2019-10-23
编辑： James Tursa 2019-10-24
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);

请先登录，再回答此问题。

采纳的回答

James Tursa
James Tursa 2019-10-24
编辑：James Tursa 2019-10-24
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.
  2 个评论
James Tursa
James Tursa 2019-10-24
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
James Tursa 2019-10-24
编辑：James Tursa 2019-10-24
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.

请先登录，再进行评论。

更多回答（1 个）

Kolleggerm1
Kolleggerm1 2019-10-24
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!

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by