Area enclosed within a curve

I'm trying to find the enclosed area for two curves, shown in the plot titled "Ideal and Schmidt pV Diagrams". I think using the "trapz" function should work but im not sure how to discretise the x-axis, thank you for any suggestions.

V_se=0.00008747;
V_sc=0.00004106;
v=V_sc/V_se;
T_e=333;
T_c=293;
t=T_c/T_e;
phi=pi/2;
T_r=312.6;
%let V_de=1000, V_r=500, V_dc=50
V_de=0.000010;
V_r=0.000050;
V_dc=0.000010;
V_d=V_de+V_r+V_dc;
T_d=324.7;
s=(V_d/V_se)*(T_c/T_d);
B=t+1+v+2*s;

A=sqrt(t^2-2*t+1+2*(t-1)*v*cos(phi)+v^2);
c=A/B;
p_m=101325;
p_max=p_m*(sqrt(1+c)/sqrt(1-c));
p_min=p_m*(sqrt(1-c)/sqrt(1+c));
delta=atan((v*sin(phi))/(t-1+v*cos(phi)));
alpha=0:0.2:2*pi;
for k = 1:numel(alpha)
  p(k) = (p_m*sqrt(1-c^2))/(1-c*cos(alpha(k)-delta));
end
figure
plot (alpha,p)
xlabel('\alpha')
ylabel('p')
% min = min(p);



alpha=0:0.2:2*pi;
for k = 1:numel(alpha)
  V_total(k) = (V_se/2)*(1-cos(alpha(k)))+(V_se/2)*(1+cos(alpha(k)))+(V_sc/2)*(1-cos(alpha(k)-phi))+V_d;
end
figure
plot (alpha,V_total)
xlabel('\alpha')
ylabel('V_total')


alpha=0:0.0001:2*pi;
for k = 1:numel(alpha)
  p_1(k) = (p_m*sqrt(1-c^2))/(1-c*cos(alpha(k)+delta));
  V_1_total(k) = (V_se/2)*(1-cos(alpha(k)))+(V_se/2)*(1+cos(alpha(k)))+(V_sc/2)*(1-cos(alpha(k)-phi))+V_d;
end
p_start = (p_m*sqrt(1-c^2))/(1-c*cos(0-delta));
P_mid = (p_m*sqrt(1-c^2))/(1-c*cos(pi-delta));
figure
plot (V_1_total,p_1)
xlabel('V_total')
ylabel('p')

hold on
plot(Vv1, p_2(Tv1,Vv1), 'LineWidth',2)
plot(Vv1, p_2(Tv2,Vv1), 'LineWidth',2)
plot(Vv1(1)*[1 1], p_2([Tv1(1) Tv1(end)],[1 1]*Vv1(1)), '-k', 'LineWidth',2)
plot(Vv1(end)*[1 1], p_2([Tv2(1) Tv2(end)],[1 1]*Vv1(end)), '-k', 'LineWidth',2)
title('Ideal and Schmidt pV Diagrams')
hold off


% I = trapz(V_1_total,p)

[dpb]

I've not enough time to really look at the moment, but use MATLAB vector operations instead of loops as a style...

...
alpha=0:0.2:2*pi;
%for k = 1:numel(alpha)
%  p(k) = (p_m*sqrt(1-c^2))/(1-c*cos(alpha(k)-delta));
%end
p=(p_m*sqrt(1-c^2))/(1-c*cos(alpha-delta));
...
    
...
alpha=0:0.2:2*pi;
%for k = 1:numel(alpha)
%  V_total(k) = (V_se/2)*(1-cos(alpha(k)))+(V_se/2)*(1+cos(alpha(k)))+(V_sc/2)*(1-cos(alpha(k)-phi))+V_d;
%end
V_total=(V_se/2).*(1-cos(alpha))+(V_se/2).*(1+cos(alpha))+(V_sc/2).*(1-cos(alpha-phi))+V_d;
...

Similarly for the remaining loop(s)...

I used the "dot" multiply operator altho I think here it doesn't matter as I think the multipliers are all just constants. Mandatory for operating on two vectors element-wise, however.

See the documentation for the gory details...

Just a guess I'm thinking you'll probably want to define a grid of points and use polyarea but I've not run the code to look at the area in question so it is just that, a guess.

[David Wilson]

Your code above fails. What is variable Vv1 for example?

[Dhananjay Kumar]

Can you first clearly write the mathematical equations please? And then what you want to do with them?

[Jonathan Bird]

I'm trying to plot an ideal cycle and a practical cycle on the same pressure-volume plot. I then want to find the area enclosed within each of the 2 regions on the plot I have titled "Ideal and Schmidt pV diagrams". The ideal cycle is governed by the ideal gas law, pV=mRT and for the

practical cycle the pressure is calculated as function of the crank angle alpha.

The volumes are also calculated as function of crank angle and then summed with the dead volume V_d to give the total volume V_total.

Sorr that the code didnt work, I've adapted it as below, many thanks for your help.

% variables for practical cycle

V_se=0.00008747;

V_sc=0.00004106;

v=V_sc/V_se;

T_e=333;

T_c=293;

t=T_c/T_e;

phi=pi/2;

T_r=312.6;

V_de=0.000010;

V_r=0.000050;

V_dc=0.000010;

V_d=V_de+V_r+V_dc;

T_d=324.7;

s=(V_d/V_se)*(T_c/T_d);

B=t+1+v+2*s;

A=sqrt(t^2-2*t+1+2*(t-1)*v*cos(phi)+v^2);

c=A/B;

p_m=101325;

p_max=p_m*(sqrt(1+c)/sqrt(1-c));

p_min=p_m*(sqrt(1-c)/sqrt(1+c));

delta=atan((v*sin(phi))/(t-1+v*cos(phi)));

alpha=0:0.2:2*pi;

%calculating the pressure as a function of crank angle

for k = 1:numel(alpha)

p(k) = (p_m*sqrt(1-c^2))/(1-c*cos(alpha(k)-delta));

end

figure

plot (alpha,p)

xlabel('\alpha')

ylabel('p')

% min = min(p);

%calculating volume as a function of crank angle

alpha=0:0.2:2*pi;

for k = 1:numel(alpha)

V_total(k) = (V_se/2)*(1-cos(alpha(k)))+(V_se/2)*(1+cos(alpha(k)))+(V_sc/2)*(1-cos(alpha(k)-phi))+V_d;

end

figure

plot (alpha,V_total)

xlabel('\alpha')

ylabel('V_total')

% Plotting pressure as a function of volume for practical cycle

alpha=0:0.0001:2*pi;

for k = 1:numel(alpha)

p_1(k) = (p_m*sqrt(1-c^2))/(1-c*cos(alpha(k)+delta));

V_1_total(k) = (V_se/2)*(1-cos(alpha(k)))+(V_se/2)*(1+cos(alpha(k)))+(V_sc/2)*(1-cos(alpha(k)-phi))+V_d;

end

% p_start = (p_m*sqrt(1-c^2))/(1-c*cos(0-delta));

% P_mid = (p_m*sqrt(1-c^2))/(1-c*cos(pi-delta));

figure

plot (V_1_total,p_1)

xlabel('V_total')

ylabel('p')

%setting variables of the ideal cycle

Tv1 = linspace(253, 293, 10);

Tv2 = linspace(293, 333, 10);

Vv1 = linspace(0.0001575, 0.0001985, 10);

[T1m,V1m] = ndgrid(Tv1,Vv1);

[T2m,V1m] = ndgrid(Tv2,Vv1);

R = 287.1;

m = 0.0001897;

%ideal gas law

p_2 = @(T,V) m*R*T./V;

hold on

%plotting ideal cycle

plot(Vv1, p_2(Tv1,Vv1), 'LineWidth',2)

plot(Vv1, p_2(Tv2,Vv1), 'LineWidth',2)

plot(Vv1(1)*[1 1], p_2([Tv1(1) Tv1(end)],[1 1]*Vv1(1)), '-k', 'LineWidth',2)

plot(Vv1(end)*[1 1], p_2([Tv2(1) Tv2(end)],[1 1]*Vv1(end)), '-k', 'LineWidth',2)

title('Ideal and Schmidt pV Diagrams')

hold off

% I = trapz(V_1_total,p)

[dpb]

Update the code in the original; then use the 'Code' button to format it legibily/neatly...

Probably if you just generated the final data to create the plot and attached it would be simpler--how you generated the curves is really of not much import here I think...unless you want to try symbolic integration (and I don't have toolbox to try).