close all
clear all

gamma = 1.4;

V0 = 1.5*231;  %in^3

V1 = .9*V0;  

Aa = pi*3^2;  %6 Inch diameter

lo = V0/Aa;
l1 = Aa*lo/V1;

Ap = (pi/4)*(2.0^2-1^2);

z1 = 0;  %Min stroke
z2 = 75;  %Max stroke
z = z1:(z2-z1)/100:z2;

F1 = 200;
P1 = F1/Ap;  %Precharge (psi)


l = Ap*(z-z1)/Aa+l1;
V = Aa*(lo-l);

CR = V1./V;


Piso = P1*(lo-l1)./(lo-l1-Ap*(z-z1)/Aa);
Fiso = Piso*Ap;

Padiabatic= P1*((lo-l1)./(lo-l1-Ap*(z-z1)/Aa)).^gamma;
Fadiabatic = Padiabatic*Ap;









figure
subplot(3,1,1);
plot(z,Piso,z,Padiabatic);
ylabel('Pressure (psi)');
xlabel('Rod position (in)');
legend('isothermal','adiabatic','Location','NorthWest');

subplot(3,1,2);
plot(z,Fiso,z,Fadiabatic);
ylabel('Force (lbs)');
xlabel('Rod position (in)');
legend('isothermal','adiabatic','Location','NorthWest');

subplot(3,1,3);
plot(z,Fiso,z,Fadiabatic);
ylabel('Force (lbs)');
xlabel('Rod position (in)');
legend('isothermal','adiabatic','Location','NorthWest');




figure
subplot(2,1,1);
plot(z,l);
ylabel('Accumulator Piston Position (in)');
xlabel('Rod position (in)');

subplot(2,1,2);
plot(z,CR);
ylabel('Compression Ratio');
xlabel('Rod position (in)');