The differential equation problem with variable solution by using ode45

2019 10 31
7 个回答
回答已采纳
更新时间：2019 11 1
6 次查看（30 天）

0 个投票

I have four coupled diffrential equation shown bellow :-
In which a ,b ,c ,d ,e,f are constant and x[t] we will get from the solution of this second order diffrential equation .how to write code for it in matlab.plz help

4 个评论
显示 2更早的评论 隐藏 2更早的评论

Walter Roberson
Walter Roberson 2019-10-31
If you have the symbolic toolbox then I recommend that you look at the flow of calls in the first example under the odeFunction documentation.
darova
darova 2019-10-31
Solve this equation first
image.png
Use results to solve the system of ODE
abhishek singh
abhishek singh 2019-10-31

I have solved individually the second order differential equation how to use it output as x(t) in coupled differential equation.

请先登录，再进行评论。

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

 采纳的回答

darova
darova 2019-10-31

1 个投票

Here is an idea
[t1,x1] = ode45(@F1,ts,x0); % solve second order
[t2,s] = ode45(@(t,s)F2(t,s,t1,x1(:,1)),ts,s0); % solve system
function ds = F2(t,s,t1,x1)
x = interp1(t1,x1,t); % extract x
ds(1,1) = a*x-b ...
ds(2,1) = c*x*s(1) ...
end

3 个评论
显示 1更早的评论 隐藏 1更早的评论

abhishek singh
abhishek singh 2019-10-31
This is the matlab file i have constructed what else i have to change in this.noscllator have second order diffrential equation and coupled equation is in xotss. what else i have to do.
darova
darova 2019-10-31
编辑：darova 2019-10-31
Here is how your code should look like
function main
%% your constants
[t1,x1] = ode45(@noscillator,[0 10],[0 1]);
[t2,s1] = ode45(@(t,s) xotss(t,s,t1,x1(:,1)), [0 10], [1 0 0 0]);
plot(t2,s1(:,1))
function xdot=noscillator(t,x)
xdot(1) = x(2);
xdot(2) = -(omega^2)*x(1)-3*(gamma/m)*x(1)^2 - 4*(beta/m)*x(1)^3 + (V/m)*cos(w *t);
xdot=xdot';
end
function dxdt = xotss(t,s,t1,x1)
x = interp1(t1,x1,t);
dxdt(1) = (a*x-b)*s(1) + c*x*s(2);
% ...
dxdt = -1i*dxdt'
end
end

请先登录，再进行评论。

更多回答（6 个）

abhishek singh
abhishek singh 2019-10-31

0 个投票

Sir ,I need to calculate rho11 at diffrent time but i am getting an error plz help.
function F1
%% your constants
[t1,x1] = ode45(@noscillator,[0 100],[0 1]);
[t2,s1] = ode45(@(t,s) xotss(t,s,t1,x1(:,1)), [0 100], [1 0 0 0]);
for ti = 0:1:100
rho11(ti)=s1(ti,1).*s1(ti,1)'-s1(ti,3).*s1(ti,3)';
rho12(ti)=s1(ti,1).*s1(ti,2)'+s1(ti,3).*s1(ti,4)';
rho21(ti)=s1(ti,2).*s1(ti,1)'+s1(ti,4).*s1(ti,3)';
rho22(ti)=s1(ti,2).*s1(ti,2)'-s1(ti,4).*s1(ti,4)';
end
s1(:,1)
plot(t2,rho11+rho22)
function xdot=noscillator(t,x)
m=1; omega=1; beta=0.01;gamma=0.01; V=.5; w=(3.0)/(2*pi);
xdot(1) = x(2);
xdot(2) = -(omega^2)*x(1)-3*(gamma/m)*x(1)^2 - 4*(beta/m)*x(1)^3 + (V/m)*cos(w *t);
xdot=xdot';
end
function dsdt = xotss(t,s,t1,x1)
x = interp1(t1,x1,t);
dsdt = zeros(4,1); % this creates a empty coloumn
%vector that you can fill with four derivatives:
a=1;b=0.05625;c=0.010825;d=0.5;e=0.5;f=0.01875;h=0.00625;% define constants
dsdt(1) = -1i* s(1) *(a+h*x-b) -1i* c*x* s(2);
dsdt(2) = -1i*c*x *s(1) -1i* (d-e-h*x)*s(2) +1i* f*s(3);
dsdt(3) = 1i* f*s(2) -1i*(e+h*x-d) *s(3)-1i* c*x* s(4);
dsdt(4) = -1i* c*x* s(3) -1i* (b-a-h*x)* s(4);
end
end
abhishek singh
abhishek singh 2019-10-31

0 个投票

Subscript indices must either be real positive integers or logicals.
Error in F1 (line 6)
rho11(ti)=s1(ti,1).*s1(ti,1)'-s1(ti,3).*s1(ti,3)';

1 个评论
显示 -1更早的评论 隐藏 -1更早的评论

darova
darova 2019-10-31
It means
सदस्यता सूचकांकों को वास्तविक धनात्मक पूर्णांक या तार्किक होना चाहिए।
In your language. Any ideas what the problem it might be?

请先登录，再进行评论。

Walter Roberson
Walter Roberson 2019-10-31

0 个投票

[t1,x1] = ode45(@noscillator,[0:100],[0 1]);
[t2,s1] = ode45(@(t,s) xotss(t,s,t1,x1(:,1)), [0:100], [1 0 0 0]);
for ti = 0:1:100
rho11(ti+1)=s1(ti+1,1).*s1(ti+1,1)'-s1(ti+1,3).*s1(ti+1,3)';
rho12(ti+1)=s1(ti+1,1).*s1(ti+1,2)'+s1(ti+1,3).*s1(ti+1,4)';
rho21(ti+1)=s1(ti+1,2).*s1(ti+1,1)'+s1(ti+1,4).*s1(ti+1,3)';
rho22(ti+1)=s1(ti+1,2).*s1(ti+1,2)'-s1(ti+1,4).*s1(ti+1,4)';
end

3 个评论
显示 1更早的评论 隐藏 1更早的评论

Walter Roberson
Walter Roberson 2019-11-1
Note that s1(ti+1,1)' means the conjugate complex transpose of s1(ti+1,1) . It is, however, a scalar, so transpose does not make any change. The You are also expecting real-valued results, so the conjugate is probably not makeing any changes. I suspect you are doing the equivalent of squaring the value.
I worry that you might have that that s1(ti+1,1)' is the derivative of s1(ti+1,1) .

请先登录，再进行评论。

abhishek singh
abhishek singh 2019-11-1

0 个投票

No , I am not using derivative ,it complex conjugate only.

1 个评论
显示 -1更早的评论 隐藏 -1更早的评论

Rik
Rik 2019-11-1
Please do not post your comments as answer. Their order can change, which makes it confusing.

请先登录，再进行评论。

类别

在 帮助中心 和 File Exchange 中查找有关 Programming 的更多信息

产品

标签

Community Treasure Hunt

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

Start Hunting!

Translated by