Plotting a p-V diagram

154 views (last 30 days)
Jonathan Bird
Jonathan Bird on 4 Nov 2019
Commented: Jonathan Bird on 3 May 2020
I need to plot a p-V diagram for an ideal cycle using the ideal gas law: pV=mRT
Ive chosen arbitrary values of 373K and 273K for the max and min temperatures, 5kg for the mass, 0.2871 for R and [p1, V1] is [7.8378, 50]. Using this I can calculate the values for [p2 V2], [p3 V3] and [p4 V4]. However I need [p1 V1] and [p2 V2] to be connected by an isothermal and the same for [p3 V3] and [p4 V4]. Does anyone know how I can do this?? The plot should look like this:
p-V.png
Many thanks!

Accepted Answer

Star Strider
Star Strider on 4 Nov 2019
I’ve not considered this since senior Thermodynamics (back in the phlogiston era).
This should get you started:
Tv1 = linspace(273, 323, 10);
Tv2 = linspace(323, 373, 10);
Vv1 = linspace(150, 200, 10);
[T1m,V1m] = ndgrid(Tv1,Vv1);
[T2m,V1m] = ndgrid(Tv2,Vv1);
R = 0.2871;
m = 5;
% p1 = 7.8378;
% V1 = 50;
p = @(T,V) m*R*T./V;
figure
surf(T1m, V1m, p(T1m,V1m))
hold on
surf(T2m, V1m, p(T2m,V1m))
hold off
grid on
xlabel('T')
ylabel('V')
zlabel('p')
view(90,0)
figure
plot(Vv1, p(Tv1,Vv1), 'LineWidth',2)
hold on
plot(Vv1, p(Tv2,Vv1), 'LineWidth',2)
plot(Vv1(1)*[1 1], p([Tv1(1) Tv1(end)],[1 1]*Vv1(1)), '-g', 'LineWidth',2)
plot(Vv1(end)*[1 1], p([Tv2(1) Tv2(end)],[1 1]*Vv1(end)), '-g', 'LineWidth',2)
hold off
grid
xlabel('V')
ylabel('p')
axis([125 225 2.0 3.3])
legend('T(273 - 323)', 'T(323 - 373)')
If you want to plot the arrows, use the quiver function. You may need to change the plot arguments in figure(2) to get the correct directions.
Make any necessary changes to get the result you want.
Plotting a p-V diagram - 2019 11 04.png
  2 Comments
Jonathan Bird
Jonathan Bird on 5 Nov 2019
Thanks for your help.
Im just a bit confused about where the 323 comes from? Surely temperature is constant along the red and blue lines?
Star Strider
Star Strider on 5 Nov 2019
It only works correctly if the temperature vector is split into two separate ranges of equal differences. That’s just the nature of the adiabatic relation.
Feel free to experiment with it, for example changing the limits of ‘Tv3’, ‘Tv4’, and ‘Vv2’ (the two temperature vectors and the volume vector respectively) in figure(3):
VL = 50; % Vector Lengths
Tv3 = linspace(173, 273, VL);
Tv4 = linspace(273, 373, VL);
Vv2 = linspace(50, 150, VL);
figure
plot(Vv2, p(Tv3,Vv2), 'LineWidth',2)
hold on
plot(Vv2, p(Tv4,Vv2), 'LineWidth',2)
plot(Vv2(1)*[1 1], p([Tv3(1) Tv3(end)],[1 1]*Vv2(1)), '-g', 'LineWidth',2)
plot(Vv2(end)*[1 1], p([Tv4(1) Tv4(end)],[1 1]*Vv2(end)), '-g', 'LineWidth',2)
hold off
grid
xlabel('V')
ylabel('p')
axis([25 175 1.0 9])
legend(sprintf('T(%d - %d)',min(Tv3),max(Tv3)), sprintf('T(%d - %d)',min(Tv4),max(Tv4)))
Just be sure that all the vector lengths are always the same (set by ‘VL’ here). (The axis call makes the plot more readable. If the plot seems cut off or too small, comment it out and then restate its limits so it looks acceptable.)

Sign in to comment.

More Answers (1)

Fernanda Silverio
Fernanda Silverio on 9 Apr 2020
Hi,
This is correct for the Carnot cycle where you have the adiabatic transformation. For Stirling you should have isotherms, therefore create an array with VL elements with the same temperature.
Look at the difference
:
R = 0.2871;
m = 5;
VL = 10;
%Define Temperatures ( It only works correctly
% if the temperature vector is split into two separate
% ranges of equal differences. )
Th = 373; %K
Tc = 273; %K
Tm = (Th + Tc)/2;
%Define Volume Delta
V1 = 200;
V2 = 150;
Tv1 = linspace(Tc, Tm, VL); % Temp low -> Adiabatic transformation
Tv2 = linspace(Tm, Th, VL); % Temp high -> Adiabatic transformation
Vv1 = linspace(V2, V1, VL);
Thigh = zeros(1, 10) + Th;
Tlow = zeros(1, 10) + Tc;
%Ideal gases law
p = @(T,V) m*R*T./V;
%P-v Diagram
figure('DefaultAxesFontSize',12)
set(0,'DefaultTextFontSize',12)
plot(Vv1, p(Tv1,Vv1), '-k', 'LineWidth',2) %adiabatic
hold on
plot(Vv1, p(Thigh, Vv1), '-r', 'LineWidth',2); %isotherm
plot(Vv1, p(Tv2,Vv1), '-k', 'LineWidth',2) %adiabatic
plot(Vv1, p(Tlow, Vv1), '-r', 'LineWidth',2);%isotherm
%adiabatic
plot(Vv1(1)*[1 1], p([Tv1(1) Tv1(end)],[1 1]*Vv1(1)), '-k', 'LineWidth',2)
plot(Vv1(end)*[1 1], p([Tv2(1) Tv2(end)],[1 1]*Vv1(end)), '-k', 'LineWidth',2)
%isotherm
plot(Vv1(1)*[1 1], [p(Tc, Vv1(1)), p(Th, Vv1(1))], '-r', 'LineWidth',2)
text(mean(Vv1), p(Tc, mean(Vv1))+0.05, 'T cold');
plot(Vv1(end)*[1 1], [p(Tc, Vv1(end)), p(Th, Vv1(end))], '-r', 'LineWidth',2)
text(mean(Vv1), p(Th, mean(Vv1))+0.05, 'T hot');
%adiabatic
plot(V1, min(p(Tv1,Vv1)), 'k.', 'MarkerSize', 23) %POINT 1
text(V1+2, min(p(Tv1,Vv1)), '1');
plot(V2, max(p(Tv1,Vv1)), 'k.', 'MarkerSize', 23) %POINT 2
text(V2-3, max(p(Tv1,Vv1)), '2');
plot(V2, max(p(Tv2,Vv1)), 'k.', 'MarkerSize', 23) %POINT 3
text(V2-3, max(p(Tv2,Vv1)), '3');
plot(V1, min(p(Tv2,Vv1)), 'k.', 'MarkerSize', 23) %POINT 4
text(V1+2, min(p(Tv2,Vv1)), '4', 'FontSize',12);
%isotherm
plot(V1, min(p(Tlow,Vv1)), 'r.', 'MarkerSize', 23) %POINT 1
text(V1+2, min(p(Tlow,Vv1)), '1');
plot(V2, max(p(Tlow,Vv1)), 'r.', 'MarkerSize', 23) %POINT 2
text(V2-3, max(p(Tlow,Vv1)), '2');
plot(V2, max(p(Thigh,Vv1)), 'r.', 'MarkerSize', 23) %POINT 3
text(V2-3, max(p(Thigh,Vv1)), '3');
plot(V1, min(p(Thigh,Vv1)), 'r.', 'MarkerSize', 23) %POINT 4
text(V1+2, min(p(Thigh,Vv1)), '4');
set(gca,'XTick',[], 'YTick', []) % REMOVE VALUE FROM AXIS
hold off
grid
xlabel('V')
ylabel('p')
axis([135 215 1.7 3.8])
legend('adiabatic', 'isotherm')
  1 Comment
Jonathan Bird
Jonathan Bird on 3 May 2020
Many thanks for your help, do you by any chance know if this adiabatic transformation is consistent with that suggested by Finkelstein?

Sign in to comment.

Categories

Find more on Thermodynamics and Heat Transfer 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!