MATLAB: While Loop Nested in For Loop Not Transferring Values to Next Column (1 Column is Created, Every Other Column is Zeros)

errorfor loopmatrixmatrix manipulationnestednested loopwhile loop

I'm trying to use a while loop nested in a for loop to create a matrix of values from a preset list of 1×15 vectors. When I run the code, it does the first iteration perfectly well, and then when z=2, all the values in the following columns are 0 (except for the top row, which is extrapolated from the preset vectors). I've tried to use a preset matrix of zeros (1000,15) and have the nested loop replace the values until the while loop condition is satisfied (It should be no more than 956 elements based on the largest M2 value (which occurs repeatedly in the vector). I suspect that the reason that the values are not transferring over is because of the unevenness in the vector values and their interaction with the while loop. I'm not too concerned with getting a matrix with an even number of values for each initial M2 (M_max). I'm more concerned with getting the values I need to move forward with the calculations.
T0=300; %Stagnation temperature, K
R=287; %gas constant, J/kgK
gam=1.4; %specific heat ratio
gamneg=gam-1;
gampos=gam+1;
T_crit=T0*0.8333; %critical temperature, K
v_crit=1*sqrt(gam*R*T_crit); %velocity of the throat, m/s
rho=1.225; %static density of air, assuming this is also the stagnation density of air
rho_crit=rho*0.6339; %critical density of air
Pressure_Conv=1E+5; %conversion factor from bar to Pa
Liter_Conv=1.667E-5; %conversion factor from LPM to cubic meters/s
Patm=101325; %Atmospheric Pressure, Pa
%% Vectors and Pump Data
m=[15 .27 .195 .360 .210 .430 .480 .380 1.8 1.6 2.4 4.3 5.2 5.6 7.8]; %mass vecotr
W=9.81.*m; %weight vecotr
Q_LPM=[55 2.5 3.1 4.2 4.5 8.0 3 2.7 12 15 24 25 30 50 78]; %Volumetic flow rate vector, Q in LPM
Q=Liter_Conv.*Q_LPM; %Q in cubmic meters/s
m_dot=rho_crit.*Q; %mass flow rate determined by volumetric flow rate and critical density at the throat, kg/s
P0abs=[8 2.4 2 2.5 2.5 2.5 6.5 6.5 8 8 8 5 6 3 3]; %Total pressure in bar-abs
P0=Pressure_Conv.*P0abs; %Total Pressure in Pa-abs
A_crit=m_dot*sqrt(T0)./(0.04042.*P0); %area at the throat
d_crit=sqrt(pi.*A_crit)*2; %diameter at the throat
Prat_max=P0./Patm; %Pressure ratio bt. stagnation and atmospheric for use in the maximum Mach # calculation
M_max=sqrt((2/gamneg)*((Prat_max.^(gamneg/gam))-1)); %maximum Mach # attainable with pressure drop
 for z=1:1:15
    n=1; %counting variable for while loop
    M2(n,z)=M_max(1,z); %initialized M2
    while M2>1.05
        M2(n,z)=M_max(1,z)-n*0.001; %decreasing maximum Mach # for each pump
        Trat(n,z)=1+((gamneg/2)*M2(n,z)^2); %Temperature ratio 
        T2(n,z)=T0/Trat(n,z); %Exit Temperature 
        v2(n,z)=M2(n)*sqrt(gam*R*T2(n,z)); %Exit Velocity 
        Arat(n,z)=(1/M2(n,z))*((2/gampos)*(1+((gamneg/2)*M2(n,z)^2)))^(gampos/(2*gamneg)); %Area ratio based on Exit Mach #
        Ae(n,z)=A_crit(z)*Arat(n,z); %Calculating Exit Area 
        Prat(n,z)=(1+((gamneg/2)*M2(n,z)^2))^(gam/gamneg); %Pressure ratio as dictated by Mach #
        Pe(n,z)=P0(z)/Prat(n,z); %Exit pressure changing with Mach #
        Thrust(n,z)=(m_dot(z)*v2(n,z))+((Pe(n,z)-Patm)*Ae(n,z)); %Thrust
        n=n+1;
    end
end

Best Answer

M2=zeros(956,15);
for z=1:1:15
    n=1; %counting variable for while loop
    M2(:,z)=M_max(1,z); %initialized M2
    while M2(n,z)>1.05
        if n==956;break;end
        M2(n,z)=M_max(1,z)-n*0.001; %decreasing maximum Mach # for each pump
        Trat(n,z)=1+((gamneg/2)*M2(n,z)^2); %Temperature ratio 
        T2(n,z)=T0/Trat(n,z); %Exit Temperature 
        v2(n,z)=M2(n)*sqrt(gam*R*T2(n,z)); %Exit Velocity 
        Arat(n,z)=(1/M2(n,z))*((2/gampos)*(1+((gamneg/2)*M2(n,z)^2)))^(gampos/(2*gamneg)); %Area ratio based on Exit Mach #
        Ae(n,z)=A_crit(z)*Arat(n,z); %Calculating Exit Area 
        Prat(n,z)=(1+((gamneg/2)*M2(n,z)^2))^(gam/gamneg); %Pressure ratio as dictated by Mach #
        Pe(n,z)=P0(z)/Prat(n,z); %Exit pressure changing with Mach #
        Thrust(n,z)=(m_dot(z)*v2(n,z))+((Pe(n,z)-Patm)*Ae(n,z)); %Thrust
        n=n+1;
    end
end
Related Question