Plotting Air Density as a function of altitude (in SI Units)

hello, i am trying to plot air density as a function of altitude, but i am having trouble with the plotting itself.

when i run the code, there are no calculation errors, however the graph i get out of the plot is odd and doesnt look right. Can someone explain what my errors are?

Any help would be greatly appreciated, thanks in advance.

All units are in SI units too btw

% ---------------------------------------------------------
%  AIR DENSITY AS A FUNCTION OF ALTITUDE
% ---------------------------------------------------------
clc
clear

X = 40000;                                                  %defining the altitude limit
rho = zeros(1,X);                                           %creating empty rho vector

for z = 1:X 
    
    if z < 10999                                            %inside the troposphere
        T = 15.04 - 0.00649 * z;                        
        P = (101.29) * ((T + 273.15)/288.08).^(5.256);
        rho(z) =  P./(0.2869*(T + 273.15));
    end
   
    if (11000 < z) && ( z < 24999)                          %inside the the lower stratosphere
        
        T = -56.46;
        P = 22.65 * exp(1.73 - 0.000157 * z);
        rho(z) = P./(0.2869*(T + 273.15));
    end 
    
    if z > 25000                                            %upper stratosphere and beyond
       
        T = -131.21 + 0.00299 * z;
        P = 2.488 * ((T + 273.15)/288.08).^(-11.388);
        rho(z) = P./(0.2869*(T + 273.15));
    end 
    rho = rho.';                                            %transpose of rho vector
    L = length(rho);
    z = 1:L;                                                %re-defining altitude vector to match rho vector dimensions
                                                            
                                                            
end

plot(z,rho);                                                %plotting rho as a function of z

[Image Analyst]

If not this,

then what do you expect?

[Jared Guzman]

Its just most graphs i have seen of air density changing with height dont include the vertical legs, or appear periodic...? I am not very experienced with these calculations so I could be wrong to assume that my graph has an error?

[dpb]

Just because you don't get an actual error doesn't mean calculations are correct -- you've got a discontinuity in the correlations between the 2nd and third sections; there's something wrong there in your implementation.

You also have a missed calculation point between each of the three section breakpoints resulting in one returned zero value --

if z < 10999  
     ...
if (11000 < z) && ( z < 24999)                          %inside the the lower stratosphere
  ...    
if z > 25000                                            %upper stratosphere and beyond
  ...

You didn't calculate anything for z==11000 nor 25000...

There's a missing normaliztion or offset or something in the last section to match up to the previous results; would have to see the source of the correlations to determine what.

But, just looking at what the correlations as written generate,

z=25000;                          % look at breakpoint
T = -131.21 + 0.00299 * z         % upper stratosphere correlation T
  -56.4600                       
% matches the lower strosphere so that's ok
P = 2.488 * ((T + 273.15)/288.08).^(-11.388)
P =
   63.7166
>>  P = 22.65 * exp(1.73 - 0.000157 * z)
P =
    2.5223
>> 

But the two P values aren't even close; something's not right in that correlation it would seem.

[Jared Guzman]

Okay i see what youre saying, thank you.

https://www.grc.nasa.gov/www/k-12/airplane/atmosmet.html

the above link is where i get the equations from, so one would expect the pressure to be different > 25000 m than at 24999 m technically, but youre right that is a huge difference between the two pressures.

I greatly appreciate your help with this!

[Jared Guzman]

I found my error in the calculation, changing the P equation for the upper strat to P = 2.488 * ((T + 273.15)/216.6).^(-11.388) fixes it

Thanks again! I would upvote/accept your answer but these are comments so i cant. This helped a great deal

Answer 1

[dpb]

Just because you don't get an actual error doesn't mean calculations are correct -- you've got a discontinuity in the correlations between the 2nd and third sections; there's something wrong there in your implementation.

You also have a missed calculation point between each of the three section breakpoints resulting in one returned zero value at each breakpoint--two, overall.

if z < 10999  
     ...
if (11000 < z) && ( z < 24999)                          %inside the the lower stratosphere
  ...    
if z > 25000                                            %upper stratosphere and beyond
  ...

You didn't calculate anything for z==11000 nor 25000...

There's a missing normaliztion or offset or something in the last section to match up to the previous results; would have to see the source of the correlations to determine what.

But, just looking at what the correlations as written generate,

z=25000;                          % look at breakpoint
T = -131.21 + 0.00299 * z         % upper stratosphere correlation T
  -56.4600                       
% matches the lower strosphere so that's ok
P = 2.488 * ((T + 273.15)/288.08).^(-11.388)
P =
   63.7166
>>  P = 22.65 * exp(1.73 - 0.000157 * z)
P =
    2.5223
>> 

But the two P values aren't even close; something's not right in that correlation it would seem.