Solve differential equation Systems with ODE45
2 views (last 30 days)
Show older comments
Brendan Görres
on 22 Oct 2018
Commented: Brendan Görres
on 23 Oct 2018
I want to solve a differental equation system with ODE45.
close all;clc;
%basic planet/moon Parametres%
re=6378e3;
g0=9.81;
rho0=1.225;
sh=7500;
%parametres related to the vehicle%
m0=68e3;
Cdrag=0.5;
diameter=5;
area= pi*((diameter/2)^2);
Adrag= area*Cdrag;
Isp=390;
thrust= 933910;
R=15;
%initial calculations%
ceff= Isp*g0;
propflow= thrust/ceff;
mf=m0/R;
mp= m0-mf;
tburn=mp/propflow;
%Functions%
syms h;
rhoh=rho0*exp(-h/sh);
gh=g0*((re/re+h)^2);
r(h)=re+h;
%Pitchover altitude mass and velocity%
syms x
y0=89.85;
hpo=-(g0/2)*(((m0-x)/propflow)^2)+(ceff/propflow)*(x*log(x/m0)+m0-x)==130;
mpo=vpasolve(hpo,x);
vpo=-g0*((m0-mpo)/propflow)-ceff*log(mpo/m0);
tpo=(m0-mpo)/propflow;
%inital conditions%
tRange=[0 tburn];
Y0=[y0;vpo;hpo;0;mpo];
%Solve ODE§
[tSol,YSol]=ode45(@AscentProblem,tRange,Y0);
%Define Functions%
function dYdt= AscentProblem(t,Y)
y=Y(1);
v=Y(2);
h=Y(3);
x=Y(4);
m=Y(5);
dydt=((v(t)/(re+h(t)))-(g0/v(t)).*((re/re+h(t)).^2)).*cos(y(t));
dvdt=(propflow*ceff/m(t))-(Adrag* rhoh.*(v(t).^2))/2*m(t)- gh .*sin(y(t));
dhdt=v(t).*sin(y(t));
dxdt=(re/(re+h(t))).*v(t).*cos(y(t));
dmdt=-propflow;
%Define dYdt%
dYdt=[dydt;dvdt;dhdt;dxdt;dmdt];
end
Array indices must be positive integers or logical values.
Error in sym/subsref (line 870) R_tilde = builtin('subsref',L_tilde,Idx);
Error in APode45>AscentProblem (line 50) dydt=((v(t)/(re+h(t)))-(g0/v(t)).*((re/re+h(t)).^2)).*cos(y(t));
Error in odearguments (line 90) f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 115) odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in APode45 (line 40) [tSol,YSol]=ode45(@AscentProblem,tRange,Y0);
PLease help, I cannot find my mistake.
Accepted Answer
Torsten
on 22 Oct 2018
Y0=[y0;vpo;hpo;0;mpo];
must be a vector of doubles.
Check the type of every entry of this vector and modify your code accordingly.
3 Comments
Torsten
on 22 Oct 2018
Replace v(t), h(t) y(t) and m(t) by v,h,y and m.
Further make all the constants you use in AscentProblem accessible there (g0, re, propflow,...)
More Answers (1)
Steven Lord
on 22 Oct 2018
The variables you define inside your ODE function as each containing one element of the Y input with which ode45 calls your ODE function (y, v, h, x, and m) are scalars. You shouldn't be trying to index into them using the variable t. t is unlikely to be the exact value 1, meaning it is unlikely to be a valid index into a scalar. Just use y, v, h, x, and m.
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!