# I did not understand what the problem is ?

3 views (last 30 days)
Ous Chkiri on 3 Nov 2019
Answered: Cris LaPierre on 4 Nov 2019
R=8.314/32;
T0=2930;
a=(12)/((2.027*10^6)^0.45);
rhoP=1920;
Astar=pi*0.25^2;
k=1.35;
n=0.45;
P0=101325;
v0=pi*rp^2*8;
%syms P(t)
%at beginning of the integration set initial values for the persistent variables
t1=0; %initial time step
dP=@(t,P)Fun(t,P,R,T0,rp,a,n,Ab,P0,rhoP,Astar,k,v0);%@(t,P)(Ab.*a.*P.^n.*(rhoP-rhoO)-P.*Astar.*sqrt(k/(R.*T0)).*(2/(k+1)).^((k+1)/(2.*(k-1)))).*R.*T0./v0;
[t,P]=ode45(dP, [0,0.1], P0);
figure(1)
plot(t,y)
xlabel("Time (s)")
ylabel("Chamber Pressure (Pa)")
title("Chamber Pressure vs Time (Start-Up)")
dP=@(t1,P)Fun(t,P,R,T0,rp,a,n,t1,Ab,P0,rhoP,Astar,k,v0);%@(t,P)(Ab.*a.*P.^n.*(rhoP-rhoO)-P.*Astar.*sqrt(k/(R.*T0)).*(2/(k+1)).^((k+1)/(2.*(k-1)))).*R.*T0./v0;
[t,P]=ode45(dP, [0,60], P0);
figure(2)
plot(t,y)
xlabel("Time (s)")
ylabel("Chamber Pressure (Pa)")
title("Chamber Pressure vs Time ")
function dP = Fun(t,P,Ab,R,T0,rp,a,n,rhoP,Astar,k)
if t==0
rp=0.35;
end
Ab=2*pi*rp*8;
rhoO=P/(R*T0);
if rp>=0.7
Ab=0;
end
v0=pi*rp^2*8;
t1=t;
rp=min(rp+((a*P^n)*10^-3)*(t-t1),0.7);
Ab=2*pi*rp*8;
dP = (Ab.*a.*P.^n.*(rhoP-rhoO)-P.*Astar.*sqrt(k/(R.*T0)).*(2/(k+1)).^((k+1)/(2.*(k-1)))).*R.*T0./v0;
end
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in HW08A (line 15)
[t,P]=ode45(dP, [0,0.1], P0);

#### 1 Comment

darova on 3 Nov 2019
Try
[t,P]=ode45(@(t1,P)Fun(t,P,R,T0,rp,a,n,t1,Ab,P0,rhoP,Astar,k,v0), [0,60], P0);
Also you can use nested function to get rid of number of arguments
function main
% constants
[t,y] = ode45(@Fun,ts,p0);
% plot ...
function dP = Fun(t,P)
% do stuff
end
end

Cris LaPierre on 4 Nov 2019
Not sure there is enough info for us to know either. Some observations, though.
1. You use variable rp before you define it. That is causing an error.
2. You use the variable Ab as an input to the function, but only define it inside the function. You must consider variable scope.
3. You call the function Fun twice but each time you use different input arguments. You must call it with exactly the same inputs that you defined it with.
4. You plot t vs y, but y does not exist. Do you mean P?
It appears you are a little unsure of what you are doing. It might be good to find an example in the documentation that you can adapt to fit your needs. It would definitely be a good idea to get some more practice with creating and using functions.
Here's what I could do to get your code to run.
R=8.314/32;
T0=2930;
a=(12)/((2.027*10^6)^0.45);
rhoP=1920;
Astar=pi*0.25^2;
k=1.35;
n=0.45;
%at beginning of the integration set initial values for the persistent variables
P0=101325;
t1=0; %initial time step
dP=@(t,P)Fun(t,P,R,T0,rp,a,n,rhoP,Astar,k);
[t,P]=ode45(dP, [0,0.1], P0);
plot(t,P)
xlabel("Time (s)")
ylabel("Chamber Pressure (Pa)")
title("Chamber Pressure vs Time (Start-Up)")
[t,P]=ode45(dP, [0,60], P0);
figure
plot(t,P)
xlabel("Time (s)")
ylabel("Chamber Pressure (Pa)")
title("Chamber Pressure vs Time ")
function dP = Fun(t,P,R,T0,rp,a,n,rhoP,Astar,k)
if t==0
rp=0.35;
end
Ab=2*pi*rp*8;
rhoO=P/(R*T0);
if rp>=0.7
Ab=0;
end
v0=pi*rp^2*8;
t1=t;
rp=min(rp+((a*P^n)*10^-3)*(t-t1),0.7);
Ab=2*pi*rp*8;
dP = (Ab.*a.*P.^n.*(rhoP-rhoO)-P.*Astar.*sqrt(k/(R.*T0)).*(2/(k+1)).^((k+1)/(2.*(k-1)))).*R.*T0./v0;
end