ODE to solve for density

7 visualizaciones (últimos 30 días)
Christian Boaitey
Christian Boaitey el 2 de Mzo. de 2021
Comentada: Bjorn Gustavsson el 2 de Mzo. de 2021
Hi,
I'm currently trying to calculate air density which involves solving an ODE for altitude. I have plotted a graph for it, I get the right shape but the axis and numbers are wrong. The first figure is what I get for my output, whereas the second is the correct one (using the line vector). Can anybody help me? Much appreciated
I get this for my output
[h,p]=ode45(@p_prime,[0 7.5],0);
plot(p,h)
ylabel('Altitude (km)')
xlabel('Air density (kg/m^3)')
function dp=p_prime(h,p)
p0=1.225
hscale=7.5
dp=(p0*exp(-h/hscale))
end

Respuesta aceptada

Mathieu NOE
Mathieu NOE el 2 de Mzo. de 2021
hello
as far as i can tell , the second plot is simply the derivative given by you and not the solution of the ODE
hx = linspace(0,7.5,100);
p0=1.225;
hscale=7.5;
dp=(p0*exp(-hx/hscale)); % derivative
figure(1), plot(dp,hx);
The curve altitude /density should have a negative slope , so at least there is a sign error in the derivative equation, drho/dh must be negative - beside that I don't see consistence with what I read here :
  4 comentarios
Christian Boaitey
Christian Boaitey el 2 de Mzo. de 2021
Editada: Christian Boaitey el 2 de Mzo. de 2021
The equation we were given was po*e^(h/hscale), this is an approximation of density as a function of altitude. These come from the curved earth equations. Other than that, I'm still confused as to why the plot is reversed?
Bjorn Gustavsson
Bjorn Gustavsson el 2 de Mzo. de 2021
You get the plot you get because you interpret things wrong. That expression fo the pressure-variation are the solution to the ODE I gave you below. You can derive that ODE from the equation-of-motion for a slice of atmosphere under steady-state conditions (air-parcel at rest so, no acceleration) for which the difference bewteen the pressure-force upward from below and downward from above is balanced by the force of gravity on the air-parcel. That leads to the ODE below. That ode you can implement in a matlab-function, then with some initial condition, like ground-pressure, you can call ode45 and obtain a solution for the pressure-variation.
Where you go wrong is that you plug in the solution of the barometric ODE in a function for an ODE and integrate that one, starting with an initial condition for the ground-pressure that is zero (which in itself is wrong, we have approximately 1 atm pressure at ground-level...).
The dp/dz should be strictly negative, since the weight of each atmospheric layer is kept steady by a correspondinly larger pressure from below than from above. Your pressure-gradient is set to be larger than zero a all altitudes which is unphysical. This also does not come from any "curved earth equation" (other than that gravity in itself makes planetary-sized bodies rather curved...)

Iniciar sesión para comentar.

Más respuestas (1)

Bjorn Gustavsson
Bjorn Gustavsson el 2 de Mzo. de 2021
Your ODE is not correct for the pressure-variation. What you have in the p_prime function is closer to the solution of the pressure-balance. If you go back to the book you'll see something like this:
That is what you should have in the p_prime-function that you then integrate with altitude.

Categorías

Más información sobre Programming en Help Center y File Exchange.

Productos


Versión

R2020b

Community Treasure Hunt

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

Start Hunting!

Translated by