solving differential equation for diffrent intial condition

9 ビュー (過去 30 日間)
samee baig
samee baig 2016 年 7 月 14 日
コメント済み: samee baig 2016 年 7 月 14 日
I am trying to solve differential equation for different initial conditions and plotting the result for each initial condition.i need help on writing a code using a loop to do the same and plotting the result against each initial condition. my current code is as below function xdot=fifteenstate2(t,x)
T=1;
v(1) =-x(1)-x(2)*x(3)-x(3)*x(4)-x(5)*x(6)-x(7)*x(8)-x(9)*x(10)-x(11)*x(12)-x(13)*x(14);
v(2) = -x(2);
v(3) =-x(5);
v(4) =-x(5)*x(6)-x(7);
v(5) =-x(9);
v(6) =-x(11);
v(7) =-x(13);
v(8) =x(3);
v(9) =x(6);
v(10) =x(8);
v(11) =x(10);
v(12) =x(12);
v(13) =x(14);
v(14) =x(4);
v=[v(1);v(2);v(3);v(4);v(5);v(6);v(7);v(8);v(9);v(10);v(11);v(12);v(13);v(14)];
%uu=[1;1;1;1;1;1;1];
d(1)=3.54491*sign(v(8))*sqrt(abs(v(8))/T);
d(2)=3.54491*sign(v(9))*sqrt(abs(v(9))/T);
d(3)=3.54491*sign(v(10))*sqrt(abs(v(10))/T);
d(4)=3.54491*sign(v(11))*sqrt(abs(v(11))/T);
d(5)=3.54491*sign(v(12))*sqrt(abs(v(12))/T);
d(6)=3.54491*sign(v(13))*sqrt(abs(v(13))/T);
d(7)=3.54491*sign(v(14))*sqrt(abs(v(14))/T);
d=[d(1);d(2);d(3);d(4);d(5);d(6);d(7)];
y(1)=sin((2*pi*t)/T); y(2)=cos((2*pi*t)/T);
y=[y(1);y(2)];
uu(1)=v(1)+d(1)*y(2)+d(2)*y(2)-d(3)*y(2)-d(4)*y(1)-d(5)*y(1)+d(6)*y(2);%+d(7)*y(1)+d(1)*d(7)*y(2)*y(2);
uu(2)=v(2)+d(1)*y(1);%-d(1)*d(7)*y(1)*y(2);
uu(3)=v(3)+d(2)*y(1);
uu(4)=v(4)+d(3)*y(1);
uu(5)=v(5)+d(4)*y(2);
uu(6)=v(6)+d(5)*y(2);
uu(7)=v(7)+d(6)*y(1);
uu=[uu(1);uu(2);uu(3);uu(4);uu(5);uu(6);uu(7)];
xdot(1)=uu(1); xdot(2)=uu(2); xdot(3)=x(2)*uu(1); xdot(4)=x(3)*uu(1); xdot(5)=uu(3); xdot(6)=x(5)*uu(1); xdot(7)=uu(4); xdot(8)=x(7)*uu(1); xdot(9)=uu(5); xdot(10)=x(9)*uu(1); xdot(11)=uu(6); xdot(12)=x(11)*uu(1); xdot(13)=uu(7); xdot(14)=x(13)*uu(1); xdot=xdot';
clear all close all clc
Tspan =linspace(0,100,500);
IC=[1;-3;2;-2;3;-4;1;-2;5;3;-2;-2;2;-3]; %options = odeset('RelTol',1e-2,'AbsTol',[1e-2 1e-2 1e-2 1e-2 1e-2 1e-2]); [T,Y] = ode45(@fifteenstate2, Tspan,IC)
x1 = Y(:,1); x2 = Y(:,2); x3 = Y(:,3); x4 = Y(:,4); x5 = Y(:,5); x6 = Y(:,6); x7 = Y(:,7); x8 = Y(:,8); x9 = Y(:,9); x10 = Y(:,10); x11 = Y(:,11); x12 = Y(:,12); x13 = Y(:,13); x14 = Y(:,14);
figure(1) plot(T,x1,T,x3,T,x4,T,x6,T,x8,T,x10,T,x12,T,x14,'linewidth',2) legend('x1','x3','x4','x6','x8','x10','x12','x14',8); xlabel('t(s)') ylabel('x_1,x_3,x_4,x_6,x_8,x_10,x_12,x_14') grid on
figure(2) plot(T,x2,T,x5,T,x7,T,x9,T,x11,T,x13,'linewidth',2) legend('x2','x5','x7','x9','x11','x13',6); xlabel('t(s)') ylabel('x_2,x_5,x_7,x_9,x_11,x_13') grid on
  2 件のコメント
KSSV
KSSV 2016 年 7 月 14 日
ok.....go ahead, do it.......is some thing objecting you?
samee baig
samee baig 2016 年 7 月 14 日
yes,cannot seem to figure out how to write the code

サインインしてコメントする。

採用された回答

Mischa Kim
Mischa Kim 2016 年 7 月 14 日
Samee, check out this answer. It does exactly what you are looking for.
  1 件のコメント
samee baig
samee baig 2016 年 7 月 14 日
its not working,is it because the initial condition i wrote are in vector form ?

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeOrdinary Differential Equations についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by