Iteratively using 2 For Loops to create 4 data sets

Hi,

I am trying to use 2 for loops, one nested into the first one. I want to change a variable 4 times in the first for loop, so that i get 20 values for the output from the second for loop. I want to repeat this the 4 times in the first one. From this I want to save the data individually to then be plotted.

When it runs, it isnt executing properly so either has a size of 4 or 20.

load('Thermo CW.mat')
P1 = 24990;%alt 40000ft
V1 =  238.32; %Mach 0.9
T1 = 220.79;
%Changing Parameters
CompStage = 9;
CompOPR = linspace(15,30,4);
for n = 1:length(CompOPR)
    
    
    CompSPR = CompOPR.^(1/CompStage);
    %Parameters
    CompEff = Values(5);
    CPair = [1005 1005 1005 1005];
    kAir = Values(28);
    CombdP = 0.02;
    TOPR = Values(24);
    TStages = 3;
    TSPR = TOPR^(1/TStages);
    kGas = [1.369863013698630];
    TEff = Values(9);
    CPgas = [1100 1100 1100 1100];
    
    T2s = T1 + ((V1^2)./(2.*CPair));
    P2s = P1 * ((T2s/T1).^(kAir/(kAir-1)));
    T2 = T2s;
    P2 = P2s;
    CompT1 = T2;
    CompP1 = P2;
    CompT1s = T2;
    CompP2s = P2;
    CompT2s = CompT1.*((CompSPR).^((kAir-1)/kAir));
    CompT2 = (((CompT2s - CompT1)/CompEff)+CompT1);
    ComP2s = CompP1 * ((CompT2s/CompT1).^(kAir/(kAir-1)));
    CompP2 = CompP1 * ((CompT2/CompT1).^(kAir/(kAir-1)));
    CompT3s = CompT2.*((CompSPR).^((kAir-1)/kAir));
    CompT3 = (((CompT3s - CompT2)/CompEff)+CompT2);
    ComP3s = CompP2 * ((CompT3s/CompT2).^(kAir/(kAir-1)));
    CompP3 = CompP2 * ((CompT3/CompT2).^(kAir/(kAir-1)));
    CompT4s = CompT3.*((CompSPR).^((kAir-1)/kAir));
    CompT4 = (((CompT4s - CompT3)/CompEff)+CompT3);
    CompP4s = CompP3 * ((CompT4s/CompT3).^(kAir/(kAir-1)));
    CompP4 = CompP3 * ((CompT4/CompT3).^(kAir/(kAir-1)));
    CompT5s = CompT4.*((CompSPR).^((kAir-1)/kAir));
    CompT5 = (((CompT5s - CompT4)/CompEff)+CompT4);
    CompP5s = CompP4 * ((CompT5s/CompT4).^(kAir/(kAir-1)));
    CompP5 = CompP4 * ((CompT5/CompT4).^(kAir/(kAir-1)));
    CompT6s = CompT5.*((CompSPR).^((kAir-1)/kAir));
    CompT6= (((CompT6s - CompT5)/CompEff)+CompT5);
    CompP6s = CompP5 * ((CompT6s/CompT5).^(kAir/(kAir-1)));
    CompP6 = CompP5 * ((CompT6/CompT5).^(kAir/(kAir-1)));
    CompT7s = CompT6.*((CompSPR).^((kAir-1)/kAir));
    CompT7= (((CompT7s - CompT6)/CompEff)+CompT6);
    CompP7s = CompP6 * ((CompT7s/CompT6).^(kAir/(kAir-1)));
    CompP7 = CompP6 * ((CompT7/CompT6).^(kAir/(kAir-1)));
    CompT8s = CompT7.*((CompSPR).^((kAir-1)/kAir));
    CompT8= (((CompT8s - CompT7)/CompEff)+CompT7);
    CompP8s = CompP7 * ((CompT8s/CompT7).^(kAir/(kAir-1)));
    CompP8 = CompP7 * ((CompT8/CompT7).^(kAir/(kAir-1)));
    CompT9s = CompT8.*((CompSPR).^((kAir-1)/kAir));
    CompT9= (((CompT9s - CompT8)/CompEff)+CompT8);
    CompP9s = CompP8 * ((CompT9s/CompT8).^(kAir/(kAir-1)));
    CompP9 = CompP8 * ((CompT9/CompT8).^(kAir/(kAir-1)));
    CompT10s = CompT9.*((CompSPR).^((kAir-1)/kAir));
    CompT10= (((CompT10s - CompT9)/CompEff)+CompT9);
    CompP10s = CompP9 * ((CompT10s/CompT9).^(kAir/(kAir-1)));
    CompP10 = CompP9 * ((CompT10/CompT9).^(kAir/(kAir-1)));
    T3s = T2s.*((CompOPR).^((kAir-1)/kAir))
    T3 = CompT10;
    P3s = CompP10s;
    P3 = CompP10;
    CompWorks = CPair .*(T3s-T2s);
    CompWork = CPair .*(T3 - T2);
    T4 = linspace(1500,1900,20);
    for k = 1:length(T4)
        T4s = T4;
        P4s = P3;
        P4 = P4s * (1-CombdP);
        TT1 = T4;
        TT1s = T4s;
        TP1s = P4;
        TP1 = P4;
        TT2s = TT1 * ((TSPR^((kGas-1)/kGas)));
        TT2 = (((TT2s - TT1)*TEff)+TT1);
        TP2s = TP1 *((TT2s/TT1)^(kGas/(kGas-1)));
        TP2 = TP1 * ((TT2/TT1)^(kGas/(kGas-1)));
        TT3s = TT2 * ((TSPR^((kGas-1)/kGas)));
        TT3 = (((TT3s - TT2)*TEff)+TT2);
        TP3s = TP2 *((TT3s/TT2)^(kGas/(kGas-1)));
        TP3 = TP2 * ((TT3/TT2)^(kGas/(kGas-1)));
        TT4s = TT3 * ((TSPR^((kGas-1)/kGas)));
        TT4 = (((TT4s - TT3)*TEff)+TT3);
        TP4s = TP3 *((TT4s/TT3)^(kGas/(kGas-1)));
        TP4 = TP3 * ((TT4/TT3)^(kGas/(kGas-1)));
        AuxPwr = (0.08) * CompWork;
        TrueTurbOutT = AuxPwr + CompWork;
        TurbWorkMechLoss = TrueTurbOutT ./ (0.98);
        
        T51 = reshape(T5,1,20)
        T5s = (-(TurbWorkMechLoss / CPgas)) + T4;
       
        P5 = P4 .* ((T51/T4).^(kGas./(kGas - 1)));
        P5s = P5;
        P6s = P1;
        P6 = P1;
        
        T6s1 = size(((P6s ./ P5s) .^((kGas - 1) ./ kGas)))
        
        T5s1 = reshape(T5s,20,1)
        T5 = T5s1
        size(T5s1)
        
        T6s = (T5s1 .* ((P6s ./ P5s) .^((kGas - 1) ./ kGas)));
        T6 = T5 .* ((P6 ./ P5).^((kGas - 1)./ kGas));
        V6 = (0.98.*(((T5-T6).*(2.*CPgas)).^(1/2)).^(2)).^(1/2);
        mDotair = Values(31);
        CPavg = (CPair + CPgas)/2;
        
        T31 = reshape(T3,4,1)
        
        Qin = mDotair * CPavg *(T4-T31);
        LHV = Values(17)*(10^6);
        mDotfuel = Qin/LHV;
        
        
        mOut = (mDotair + mDotfuel)
        size(mOut)
        size(V6)
        momOut = (mOut)*V6
        momIn = (mDotair.*V1)
        
       
        
        
        size(mDotair)
        
        size(V1)
       
        Thrust = momOut-momIn;
        mDotexhaust = mDotair + mDotfuel;
        SpecificThrust = Thrust/(mDotair);
    end
end

Any Help would be brilliant, I have seen examples for similar questions but nothing that I can apply to this.

Thanks in advance

Answer 1

[Kevin Phung]

You dont seem to be doing any sort of indexing with CompOPR. so every loop iteration is the same.

for example, I think you meant to do this for this line:

    CompSPR = CompOPR(n)^(1/CompStage);

[Matthew Hole]

Yes, thats what I indeeded to put. However, when I was writing the code I did the second For Loop first (shown below) and it works perfectly. So I then tried to do the next For Loop being the one you referenced above. So now I dont even get the result of below in the one above. Ie the 20 values I would get in the code below for Thrust, specific thrust etc, I now only get 4 values from the one above. Any ideas?

load('Thermo CW.mat')

P1 = 24990;%alt 40000ft
V1 =  238.32; %Mach 0.9
T1 = 220.79;


%Changing Parameters
CompStage = 9;
CompOPR = 15;
CompSPR = CompOPR^(1/CompStage);

%Parameters
CompEff = Values(5);
CPair = Values(4);
kAir = Values(28);
CombdP = 0.02;
TOPR = Values(24);
TStages = 3;
TSPR = TOPR^(1/TStages);
kGas = Values(30);
TEff = Values(9);
CPgas = Values(8);


T2s = T1 + ((V1^2)/(2*CPair));
P2s = P1 * ((T2s/T1)^(kAir/(kAir-1)));
T2 = T2s;
P2 = P2s;

CompT1 = T2;
CompP1 = P2;
CompT1s = T2;
CompP2s = P2;

CompT2s = CompT1*((CompSPR)^((kAir-1)/kAir));
CompT2 = (((CompT2s - CompT1)/CompEff)+CompT1);
ComP2s = CompP1 * ((CompT2s/CompT1)^(kAir/(kAir-1)));
CompP2 = CompP1 * ((CompT2/CompT1)^(kAir/(kAir-1)));

CompT3s = CompT2*((CompSPR)^((kAir-1)/kAir));
CompT3 = (((CompT3s - CompT2)/CompEff)+CompT2);
ComP3s = CompP2 * ((CompT3s/CompT2)^(kAir/(kAir-1)));
CompP3 = CompP2 * ((CompT3/CompT2)^(kAir/(kAir-1)));

CompT4s = CompT3*((CompSPR)^((kAir-1)/kAir));
CompT4 = (((CompT4s - CompT3)/CompEff)+CompT3);
CompP4s = CompP3 * ((CompT4s/CompT3)^(kAir/(kAir-1)));
CompP4 = CompP3 * ((CompT4/CompT3)^(kAir/(kAir-1)));

CompT5s = CompT4*((CompSPR)^((kAir-1)/kAir));
CompT5 = (((CompT5s - CompT4)/CompEff)+CompT4);
CompP5s = CompP4 * ((CompT5s/CompT4)^(kAir/(kAir-1)));
CompP5 = CompP4 * ((CompT5/CompT4)^(kAir/(kAir-1)));

CompT6s = CompT5*((CompSPR)^((kAir-1)/kAir));
CompT6= (((CompT6s - CompT5)/CompEff)+CompT5);
CompP6s = CompP5 * ((CompT6s/CompT5)^(kAir/(kAir-1)));
CompP6 = CompP5 * ((CompT6/CompT5)^(kAir/(kAir-1)));

CompT7s = CompT6*((CompSPR)^((kAir-1)/kAir));
CompT7= (((CompT7s - CompT6)/CompEff)+CompT6);
CompP7s = CompP6 * ((CompT7s/CompT6)^(kAir/(kAir-1)));
CompP7 = CompP6 * ((CompT7/CompT6)^(kAir/(kAir-1)));

CompT8s = CompT7*((CompSPR)^((kAir-1)/kAir));
CompT8= (((CompT8s - CompT7)/CompEff)+CompT7);
CompP8s = CompP7 * ((CompT8s/CompT7)^(kAir/(kAir-1)));
CompP8 = CompP7 * ((CompT8/CompT7)^(kAir/(kAir-1)));

CompT9s = CompT8*((CompSPR)^((kAir-1)/kAir));
CompT9= (((CompT9s - CompT8)/CompEff)+CompT8);
CompP9s = CompP8 * ((CompT9s/CompT8)^(kAir/(kAir-1)));
CompP9 = CompP8 * ((CompT9/CompT8)^(kAir/(kAir-1)));

CompT10s = CompT9*((CompSPR)^((kAir-1)/kAir));
CompT10= (((CompT10s - CompT9)/CompEff)+CompT9);
CompP10s = CompP9 * ((CompT10s/CompT9)^(kAir/(kAir-1)));
CompP10 = CompP9 * ((CompT10/CompT9)^(kAir/(kAir-1)));

T3s = T2s*((CompOPR)^((kAir-1)/kAir))
T3 = CompT10;
P3s = CompP10s;
P3 = CompP10;

CompWorks = CPair *(T3s-T2s);
CompWork = CPair *(T3 - T2);

T4 = linspace(1500,1900,20);


for k = 1:length(T4)
    T4s = T4;
    P4s = P3;
    P4 = P4s * (1-CombdP);

    TT1 = T4;
    TT1s = T4s;
    TP1s = P4;
    TP1 = P4;

    TT2s = TT1 * ((TSPR^((kGas-1)/kGas)));
    TT2 = (((TT2s - TT1)*TEff)+TT1);
    TP2s = TP1 *((TT2s/TT1)^(kGas/(kGas-1)));
    TP2 = TP1 * ((TT2/TT1)^(kGas/(kGas-1)));

    TT3s = TT2 * ((TSPR^((kGas-1)/kGas)));
    TT3 = (((TT3s - TT2)*TEff)+TT2);
    TP3s = TP2 *((TT3s/TT2)^(kGas/(kGas-1)));
    TP3 = TP2 * ((TT3/TT2)^(kGas/(kGas-1)));

    TT4s = TT3 * ((TSPR^((kGas-1)/kGas)));
    TT4 = (((TT4s - TT3)*TEff)+TT3);
    TP4s = TP3 *((TT4s/TT3)^(kGas/(kGas-1)));
    TP4 = TP3 * ((TT4/TT3)^(kGas/(kGas-1)));

    AuxPwr = (0.08) * CompWork;

    TrueTurbOutT = AuxPwr + CompWork;

    TurbWorkMechLoss = TrueTurbOutT/(0.98);


    T5s = (-(TurbWorkMechLoss/CPgas)) + T4;
    T5 = T5s;
    P5 = P4 * ((T5/T4)^(kGas/(kGas - 1)));
    P5s = P5;

    P6s = P1;
    P6 = P1;
    T6s = T5s *((P6s/P5s)^((kGas-1)/kGas));
    T6 = T5 * (((P6/P5)^((kGas-1)/kGas)));

    V6 = (0.98*(((T5-T6)*(2*CPgas)).^(1/2)).^(2)).^(1/2);

    mDotair = Values(31);

    CPavg = (CPair + CPgas)/2;

    Qin = mDotair * CPavg *(T4-T3);

    LHV = Values(17)*(10^6);

    mDotfuel = Qin/LHV;

    Thrust = ((mDotair + mDotfuel).*V6)-(mDotair.*V1);
    
    mDotexhaust = mDotair + mDotfuel;

    SpecificThrust = Thrust/(mDotair);
end


plot(T4,SpecificThrust)
xlabel('Temperature, K')
ylabel('Specific Thrust, N s / kg')