ode 45 y componer un vector con valores de una matriz
Publicado por Ana (1 intervención) el 17/03/2014 11:28:57
Buenas a todos,
Estoy usando ode45 para resolver una ecuación diferencial ordinaria. En este caso lo que quiero es obtener el valor de la Presión (Pr) en función del ángulo theta (theta3). Como ya sabe, ode45 da como resultados dos columnas, que en este caso serán [theta, Pr]
Lo que yo quiero es que los valores numéricos obtenidos en la columna Pr, se me almacenen como componentes de un vector "P" definido anteriormente en el código, y cuyas componentes i=1:280 son conocidas. Quiero que estos valores se almacenen en i=281:340
Además, no sé cómo definir 'theta3' para que se considere como variable en el problema
phi=0.8 ;
rc=15;
lambda=0.25;
b=0.086;
s=0.086;
r=0.043;
R=287;
cv=717;
x2=280;
x3=340;
%-------------------------------------------------------------------
k1= 2 ;
k2= 5000 ;
k3= 14.2./phi.^0.644 ;
k4= 0.79.*k3.^0.25 ;
beta= 0.25 ;
%------------------------------------------------------------------
syms theta2;
fpre=1-(1-((theta2-x2)./b).^k1).^k2 ;
fdif=1-exp(-k3.*((theta2-x2)./b).^k4) ;
Xb= beta.*fpre+(1-beta).*fdif ;
dXb=diff(Xb);
%---------------------------------------------------------------------
DPt=@(theta3,Pr)((dQ-Pr.*dVt).*R./cv-Pr.*dVt)./V ;
V=(pi.*b^2.*s./4).*(1./(rc-1)+1./2.*(1./lambda+1-cosd(theta3)-(1./lambda.^2-sind(theta3.^2).^0.5))) ;
dVt=pi.*b.^2.*r./4.*sind(theta3).*(1+lambda.*cosd(theta3)./(sqrt(1-lambda.^2.*sind(theta3).^2))) ;
dXbsust=subs(dXb, theta2, theta3);
dQ=Qt.*dXbsust ;
tspan=[x2+1 x3];
%P(i)=dsolve('DP(i)=((dQ(i)-P(i).*dVt(i)).*R./Cv-P(i).*dVt(i))./V(i)', 'P(0)=100' , 'theta')
[theta3,Pr]=ode45(DPt,tspan, 100);
for i=281:340
P(i)=[theta3,Pr]
V=(pi.*b^2.*s./4).*(1./(rc-1)+1./2.*(1./lambda+1-cosd(theta(i))-(1./lambda.^2-sind(theta(i).^2).^0.5))) ;
%Derivada de V respecto de theta
dVt=pi.*b.^2.*r./4.*sind(theta(i)).*(1+lambda.*cosd(theta(i))./(sqrt(1-lambda.^2.*sind(theta(i)).^2))) ;
T(i)=P(i).*V(i)./(m.*R) ;
end
figure(1)
title ('Presión')
hold all
grid on
plot (theta,P(i))
Muchas gracias de antemano.
Estoy usando ode45 para resolver una ecuación diferencial ordinaria. En este caso lo que quiero es obtener el valor de la Presión (Pr) en función del ángulo theta (theta3). Como ya sabe, ode45 da como resultados dos columnas, que en este caso serán [theta, Pr]
Lo que yo quiero es que los valores numéricos obtenidos en la columna Pr, se me almacenen como componentes de un vector "P" definido anteriormente en el código, y cuyas componentes i=1:280 son conocidas. Quiero que estos valores se almacenen en i=281:340
Además, no sé cómo definir 'theta3' para que se considere como variable en el problema
phi=0.8 ;
rc=15;
lambda=0.25;
b=0.086;
s=0.086;
r=0.043;
R=287;
cv=717;
x2=280;
x3=340;
%-------------------------------------------------------------------
k1= 2 ;
k2= 5000 ;
k3= 14.2./phi.^0.644 ;
k4= 0.79.*k3.^0.25 ;
beta= 0.25 ;
%------------------------------------------------------------------
syms theta2;
fpre=1-(1-((theta2-x2)./b).^k1).^k2 ;
fdif=1-exp(-k3.*((theta2-x2)./b).^k4) ;
Xb= beta.*fpre+(1-beta).*fdif ;
dXb=diff(Xb);
%---------------------------------------------------------------------
DPt=@(theta3,Pr)((dQ-Pr.*dVt).*R./cv-Pr.*dVt)./V ;
V=(pi.*b^2.*s./4).*(1./(rc-1)+1./2.*(1./lambda+1-cosd(theta3)-(1./lambda.^2-sind(theta3.^2).^0.5))) ;
dVt=pi.*b.^2.*r./4.*sind(theta3).*(1+lambda.*cosd(theta3)./(sqrt(1-lambda.^2.*sind(theta3).^2))) ;
dXbsust=subs(dXb, theta2, theta3);
dQ=Qt.*dXbsust ;
tspan=[x2+1 x3];
%P(i)=dsolve('DP(i)=((dQ(i)-P(i).*dVt(i)).*R./Cv-P(i).*dVt(i))./V(i)', 'P(0)=100' , 'theta')
[theta3,Pr]=ode45(DPt,tspan, 100);
for i=281:340
P(i)=[theta3,Pr]
V=(pi.*b^2.*s./4).*(1./(rc-1)+1./2.*(1./lambda+1-cosd(theta(i))-(1./lambda.^2-sind(theta(i).^2).^0.5))) ;
%Derivada de V respecto de theta
dVt=pi.*b.^2.*r./4.*sind(theta(i)).*(1+lambda.*cosd(theta(i))./(sqrt(1-lambda.^2.*sind(theta(i)).^2))) ;
T(i)=P(i).*V(i)./(m.*R) ;
end
figure(1)
title ('Presión')
hold all
grid on
plot (theta,P(i))
Muchas gracias de antemano.
Valora esta pregunta


0