How do I use a for loop in my ode15s based code for shooting method and get multiple graphs?

2 ビュー (過去 30 日間)
naygarp
naygarp 2018 年 5 月 28 日
コメント済み: MINATI 2019 年 2 月 11 日
I am using ode15s solver to solve a set of odes by shooting technique and obtain graphs of the solutions.
while trying to vary some parameters in the equations I have used a for loop.
But I am not getting the proper graphs.
How do I use the 'for' loop properly to get the accurate graphs.
Here in the code below I need to vary for three values of the 'lambda' parameter, so that I could obtain 3 different graphs in one figure
function shooting_method
clc;clf;clear;
global lambda gama Pr Rd Lew Nb Nt Mn
gama=1;
Mn=6;
Rd=0.1;
Pr=10;
Nb=0.3;
Lew=10;
Nt=0.3;
pp = [0.5 1 1.5];
for i=1:numel(pp)
lambda = pp(i);
%lambda=0.5;
x=[1 1 1];
format long
options= optimset('Display','iter');
x1= fsolve(@solver,x);
end
end
function F=solver(x)
options= odeset('RelTol',1e-8,'AbsTol',[1e-8 1e-8 1e-8 1e-8 1e-8 1e-8 1e-8]);
[t,u] = ode15s(@equation,linspace(0,2),[0 1 x(1) 1 x(2) 1 x(3)],options)
%sol= [t,u];
s=length(t);
F= [u(s,2),u(s,4),u(s,6)];
%y1 = deval(u(0,:))
plot(t,u(:,2),'LineWidth',2);hold on %t,u(:,4));
end
function dy=equation(t,y)
global lambda gama Pr Rd Lew Nb Nt Mn
dy= zeros(7,1);
dy(1)= y(2);
dy(2)= y(3)*(y(3)^2+gama^2)/(y(3)^2+lambda*gama^2);
dy(3)= y(2)^2/3-(2*y(1)*y(3)*(y(3)^2+gama^2))/(3*(y(3)^2+lambda*gama^2))+Mn*y(2);
dy(4)= y(5);
dy(5)= -(2*Pr*y(1)*y(5))/(3*(1+Rd)) - (Nb*y(5)*y(7))/(1+Rd) - (Nt*y(5)^2)/(1+Rd);
dy(6)= y(7);
dy(7)= -((2*Lew*y(1)*y(7))/3)+ (Nt/Nb)*((2*Pr*y(1)*y(5))/(3*(1+Rd)) + (Nb*y(5)*y(7))/(1+Rd) + (Nt*y(5)^2)/(1+Rd));
end
  1 件のコメント
MINATI
MINATI 2019 年 2 月 2 日
Dear naygarp
Can you share the question or the paper of whose this is the code
minatipatra456@gmail.com

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

回答 (1 件)

Torsten
Torsten 2018 年 5 月 29 日
編集済み: Torsten 2018 年 5 月 30 日
function shooting_method
clc;clf;clear;
global lambda gama Pr Rd Lew Nb Nt Mn
gama=1;
Mn=6;
Rd=0.1;
Pr=10;
Nb=0.3;
Lew=10;
Nt=0.3;
pp = [0.5 1 1.5];
for i=1:numel(pp)
lambda = pp(i);
%lambda=0.5;
x=[1 1 1];
format long
options= optimset('Display','iter');
x1= fsolve(@solver,x);
options= odeset('RelTol',1e-8,'AbsTol',[1e-8 1e-8 1e-8 1e-8 1e-8 1e-8 1e-8]);
[t,u] = ode15s(@equation,linspace(0,2,20),[0 1 x1(1) 1 x1(2) 1 x1(3)],options)
U(i,:,:)=u;
T(i,:)=t;
end
plot(T(1,:),U(1,:,2),T(2,:),U(2,:,2),T(3,:),U(3,:,2))
end
  3 件のコメント
Torsten
Torsten 2018 年 6 月 1 日
Of course,you will have to run the above modified function "shooting_method" together with your two functions "solver" and "equation" from above.
MINATI
MINATI 2019 年 2 月 11 日
Dear naygarp
Can you share the question or the paper of whose this is the code
minatipatra456@gmail.com

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

カテゴリ

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

Community Treasure Hunt

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

Start Hunting!

Translated by