error in ode45 - must return a column vector

1 view (last 30 days)
Hossein
Hossein on 16 Aug 2014
Edited: Star Strider on 18 Aug 2014
Hi everyone im trying to solve a system of coupled ordinary differential equations using ODE45 and i keep getting this error:
Error using odearguments (line 91) THERMAL_BUBBLE must return a column vector.
Error in ode45 (line 114) [neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
i put my code below for you guys to help me find my error thanks
function dy=thermal_bubble(t,y)
dy=zeros(1,2,3);
d=998;
c=1500;
r0=10^-5;
s=0.0725;
e=3/4;
k=0.016;
u=0.000001;
p0=10^5;
f=25000;
pa=1000;
cv=1.5*10^3;
a=2.338*10^-5;
Tinf=300;
g=1.31;
dy(1)=y(2);
dy(2)=((1/d)*(1+(y(2))/c)*(p0+2*s/r0)*(((y(1))/r0)^(3*e))+((p0+2*s/r0)*(((y(1))/r0)^(3*e))*g*((-3*y(2))/y(1))+(g-1)*((3*k/y(1))*(y(3)-Tinf)*sqrt(3*y(2)*(g-1)/(a*y(1))))-4*u*(((y(2))/(y(1)))^2)+2*s*y(2)/((y(1))^2)+2*pi*f*pa*cos(2*pi*f*t))*((y(1))/(d*c))-(3/2)*((y(2))^2)*(1-(y(2))/c))*(((y(1))*(1-(y(2))/c))^-1)+(1+(4*u/(d*c))*(y(1)*((1-(y(2))/c))^-1))^-1;
dy(3)=((4*pi*(y(1))^2)/((cv*4*pi*r0^3)/3))*(k*sqrt(abs(y(2))/(a*y(1))))*(y(3)-Tinf)-y(2)*(p0+2*s/r0)*((y(1))/r0)^(3*e);
end
---------------------------------------------------------------------------------------------------------------------------------------------- and then i run this in command window and error appears
>>[t,y]=ode45(@thermal_bubble,[0 20/25000],[0.00001 0 25]);

Answers (1)

Star Strider
Star Strider on 16 Aug 2014
Edited: Star Strider on 16 Aug 2014
This is one problem:
dy=zeros(1,2,3);
Replace it with:
dy=zeros(3,1);
It will run then.
Also, this is a stiff problem, so change the solver:
[t,y]=ode15s(@thermal_bubble,[0 20/25000],[0.00001 0 25]);
  2 Comments
Star Strider
Star Strider on 18 Aug 2014
Edited: Star Strider on 18 Aug 2014
My pleasure!
(The sincerest expression of appreciation here on MATLAB Answers is to Accept the Answer that most closely solves your problem.)

Sign in to comment.

Categories

Find more on Programming in Help Center and File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!