
Trying to plot specific values from a vector
11 views (last 30 days)
Show older comments
Caroline Kenemer
on 3 May 2021
Commented: Mathieu NOE
on 4 May 2021
Hi, I am trying to plot the TSFC values and Pcr values for specfic values of To4 which are 1500,1600, and 1700. So that I will end up with three lines. How could I do this?
0 Comments
Accepted Answer
Mathieu NOE
on 3 May 2021
hello
this is my suggestion - I put the equations in a sub function to make it easy after to do a for loop for multiple To4 values
also changes the Pcr from linspace to logspace vector as it seems to better render in the semilogx plot (my opinion)

clc
clearvars
%
%Turbojet - Rolls Royce Olympus 593
%Range = 6440 km
%Engine Weight = 3386 kg;
%Ramjet is based off the Blackbird SR-71
%% Takeoff TJ
Pa = 101325; % in Pa
Pe = Pa;
rho = 1.225; %kg/m^3
Ai = 0.4902; %inlet area in meters
ma = 64.4; %mass flow rate in kg/s
u = ma/(rho*Ai);
h = 0; %km
T = 169000; % Thrust in N
R = 287; %J/kg*K
Ta = 288.15; % ambient temperature at Sea Level
gamma = 1.4;
Pcr = logspace(0,2,25);
To4 = [1500;1600;1700];
for ci = 1:length(To4)
TSFC(ci,:) = mysubfn1(Pcr,u,gamma,R,Ta,Pa,ma,T,To4(ci))
leg_str{ci} = ['To4 = ' num2str(To4(ci))];
end
semilogx(Pcr,TSFC)
legend(leg_str);
function [TSFC] = mysubfn1(Pcr,u,gamma,R,Ta,Pa,ma,T,To4)
M = u/sqrt(gamma*R*Ta);
%compressor inlet conditions
To2 = Ta*(1+((gamma-1)/2)*M^2);
yd = 1.4;
nd = 0.97;
Po2 = Pa*(1+nd*((To2/Ta)-1))^(yd/(yd-1));
% Pcr = linspace(0,100,10);
% To4 = linspace(1500,1700,10);
Po3 = Po2.*Pcr;
yc = 1.37;
nc = 0.85;
To3 = To2.*(1+(1./nc).*(Pcr.^((yc-1)./yc))-1);
Q = 45000000; %heating value in K/kg
Cp = 1000; %J/kgK
f = ((To4./To3)-1)./(((Q.*To3)./Cp)-To4./To3);
Po4 = Po3;
To5 = To4-(To3-To2);
nt = 0.90;
yt =1.33;
Po5 = Po4.*(1-(1/nt).*(1-(To5./To4))).^(yt/(yt-1));
%nozzle conditions with no afterburner
To6 = To5;
Po6 = Po5;
yn = 1.36;
nn = 0.98;
ue = sqrt(2.*nn.*(yn./(yn-1)).*R.*To6.*(1-(Pa./Po6).^(yn-1)./yn));
mf = f.*ma;
TSFC = mf./T;
end
4 Comments
Mathieu NOE
on 4 May 2021
hello again
1/ if you have multiple functions in one code , the functions must always be at the very end of the code, all "main" code must be written in the upper section; you can have multiple functions , simply make sure they have individual names
2/ xlabel, ylabel, and title : it's exactly this way to code it :
semilogx(Pcr,TSFC)
title('TSFC vs. Pcr vs. To4')
legend(leg_str);
xlabel('Pcr')
ylabel('TSFC')
3/ sure, you can always get back to linspace
or this way : Pcr = (0:10:100);
did you notice that Pcr = 0 will issue a NaN result because you make a division by zero in this line :
f = ((To4./To3)-1)./(((Q.*To3)./Cp)-To4./To3);
I send back the full code :
clc
clearvars
%
%Turbojet - Rolls Royce Olympus 593
%Range = 6440 km
%Engine Weight = 3386 kg;
%Ramjet is based off the Blackbird SR-71
%% Takeoff TJ
Pa = 101325; % in Pa
Pe = Pa;
rho = 1.225; %kg/m^3
Ai = 0.4902; %inlet area in meters
ma = 64.4; %mass flow rate in kg/s
u = ma/(rho*Ai);
h = 0; %km
T = 169000; % Thrust in N
R = 287; %J/kg*K
Ta = 288.15; % ambient temperature at Sea Level
gamma = 1.4;
% Pcr = logspace(0,2,25);
Pcr = (0:10:100);
To4 = [1500;1600;1700];
for ci = 1:length(To4)
TSFC(ci,:) = mysubfn1(Pcr,u,gamma,R,Ta,Pa,ma,T,To4(ci))
leg_str{ci} = ['To4 = ' num2str(To4(ci))];
end
figure(1),
plot(Pcr,TSFC)
title('TSFC vs. Pcr vs. To4')
legend(leg_str);
xlabel('Pcr')
ylabel('TSFC')
function [TSFC] = mysubfn1(Pcr,u,gamma,R,Ta,Pa,ma,T,To4)
M = u/sqrt(gamma*R*Ta);
%compressor inlet conditions
To2 = Ta*(1+((gamma-1)/2)*M^2);
yd = 1.4;
nd = 0.97;
Po2 = Pa*(1+nd*((To2/Ta)-1))^(yd/(yd-1));
% Pcr = linspace(0,100,10);
% To4 = linspace(1500,1700,10);
Po3 = Po2.*Pcr;
yc = 1.37;
nc = 0.85;
To3 = To2.*(1+(1./nc).*(Pcr.^((yc-1)./yc))-1);
Q = 45000000; %heating value in K/kg
Cp = 1000; %J/kgK
f = ((To4./To3)-1)./(((Q.*To3)./Cp)-To4./To3);
Po4 = Po3;
To5 = To4-(To3-To2);
nt = 0.90;
yt =1.33;
Po5 = Po4.*(1-(1/nt).*(1-(To5./To4))).^(yt/(yt-1));
%nozzle conditions with no afterburner
To6 = To5;
Po6 = Po5;
yn = 1.36;
nn = 0.98;
ue = sqrt(2.*nn.*(yn./(yn-1)).*R.*To6.*(1-(Pa./Po6).^(yn-1)./yn));
mf = f.*ma;
TSFC = mf./T;
end
More Answers (0)
See Also
Categories
Find more on Gas Dynamics in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!