 Open Access
 Total Downloads : 1705
 Authors : Muzaffar Ali Quazi , Ponnusamy Baskar
 Paper ID : IJERTV1IS5100
 Volume & Issue : Volume 01, Issue 05 (July 2012)
 Published (First Online): 02082012
 ISSN (Online) : 22780181
 Publisher Name : IJERT
 License: This work is licensed under a Creative Commons Attribution 4.0 International License
Computer Simulation Of Pneumatic Engine Operation
Muzaffar Ali Quazi , Ponnusamy BASKAR *
A mathematical model for the three cylinder pneumatic engine is proposed, which allows calculating both the dynamic characteristics of piston motion and flowing gas parameters without using any fitting procedures. The corresponding computer code in MATLABSIMULINK software is developed and numerical simulation of the engine operation has been accomplished With the help of MATLAB SIMULINK software and AUTOMOTIVE STUDIO software. The approach proposed allows calculation of a wide set of thermodynamic and operational parameters for various pneumatic cylinders and can be used for development of the highly efficient pneumatic engine intended for vehicle propulsion.

At the present time, a new direction in designing automobile using the compressed air technologies and pneumatic power plants is being developed[13].compressed air engine having the high efficiency as compared to gasoline engine and correspondingly low consumption of compressed air are necessary for development of nonpolluting pneumatic automobile that run on compressed air. The piston expansion machine based on pneumatic cylinder most closely correspond to this criterion[4].recent developments in pneumatic servo system and innovative pneumatic components[5,6] show important advances, which are expected in vehicle applications also.
Design optimization for a pneumatic engine of a given set of characteristics is possible as a result of mathematical modeling of working cycle. Therefore, the development of a adequate mathematical model is a reasonable scientific and technical approach.
The purpose of making the new design of the pneumatic engine and mathematical model developed is to determine the basic dynamic parameter and the effectiveness of new design. The dynamic parameters such as air pressure in the cylinder, position and speed of the piston in time, cycle frequency and the operational characteristics such as power, efficiency, specific work, air consumption, etc. of the pneumatic engine are being considered in this analysis of three cylinder pneumatic engine.

Each component in engine frame and cylinder has design limit. To ensure that these are not exceeded in operation, each frame and each cylinder has a design rating above which it may not be used. The loads used to rate pneumatic engine are discussed below. All three cylinders have a maximum allowable inlet pressure. All engine components are subjected to alternating loads and the rated pressure of a cylinder will be based on fatigue considerations.
Every cylinder has a minimum clearance it can be built with. This controls the volumetric efficiency of the cylinder and hence the capacity for a given pressure ratio and air consumption. The clearance of a cylinder can usually be increased if the maximum capacity is not needed for a given application.
Cylinders has a fixed number of valves and valve size. A cylinder with a few or small valves for its size will have losses and will give poor efficiency.
Each cylinder exerts a rod load on the running gear components, and a frame load on the stationary components. These can be evaluated by considering the forces acting on the various components.
Frame Load=PHE AP
Where
PCE ( AP
AROD )
PHE
and PCE
= Pressure in the head end and crank end of the cylinder
AP and AROD =The area of the piston and piston rod
The frame load will vary through the cycle as the pressures in the head end and crank end of the cylinder vary. The maximum tensile and maximum compressive stresses are calculated. These are the loads the stationary components and bolting must be designed to resist. The rod load, the force exerted on the piston rod, connecting rod and crankshaft is different for each component. It is the frame load plus the inertia of all the parts outboard of the component that is of interest. The inertia is the mass times the piston acceleration, and varies through the cycle.Once the suction and discharge pressure, the suction air temperature, the required flow rate and the air composition are determined, a engine can be selected to design.
The engine speed and the stroke will then be determined by the horsepower requirement. Alow horsepower application will require a light, low stroke, high speed engine. A high horse power application will require a heavy, long stroke, low speed engine.

The schematic diagram of a pneumatic engine is represented on fig. In this we consider only one cylinder in spite of three, as the working of the all three cylinder are same and the mathematical model have no effect on, which cylinder we used under observation.
In the above cylinder, we will first consider a stroke of the piston from left to right. The movement of the piston from point A to B corresponds to the process of filling of the cylinder by
compressed air from inlet manifold
Pm1
and exhaust of simple air from the exhaust manifold Pm2. At
the point B the supply from the inlet manifold is cut off by closing of the inlet valve, as the axis of the
cylinder and supporting assembly orifice are not coinside due to cylinder harmonic motion. The movement of the piston between B and D corresponding to the process of expansion of the working air in the cylinder, as the compressed air is expand. At the point D the exhaust valve is open, as the axis of the cylinder and supporting assemblys exhaust orifice are coinside, due to the cylinder harmonic motion. All three cylinder which are placed at 1200 to each other, done the same work vice versa.
For a description of the dynamics of the piston movement between the point A and D, it is necessary to determine the parameters of working gas state. The equation of the piston motion generally can be written as:
d 2 x
dt 2

p1s1 F
(1)
Where p1, is the pressure in the working cylinder,s1is the useful areas of the piston for cylinder is the resistance force; which consist of the force of friction Ffr, and loading force FL; M is the mass of piston with all moving parts(roads, crank shaft mechanism, etc.)
In the above equation of the piston motion. The quantity p1 is unknown. Now we have to derive the equation of pressure at every stage of motion of piston measured from the right edge of the piston.
Let we suppose that all the thermal energy dQm1, which is admitted with the air in cylinder is changed to internal energy dU1, and the work of the gas expansion dW1, and according to the first law of thermodynamics, the equation of the energy balance is written as
dQm1
dU1
dW1
(2)
Assuming, that pressure in the system of receiver manifold does not vary during filling of the
working cylinder, we use the relation dQm1 dHm1. In this case the quantity of thermal energy, which
has arrived from the inlet( Pm1) to the cylinder, is equal to the product of the mass of gas dmm1and the
specific enthalpy( dQm1
hm1dmm1) and the change of the gas internal energy dU1and work
dW1
made by it are equal accordingly the following form
dU1
d (u1m1 )and
dW1
p1dV. Therefore, eq.2 can be written in
hm1dm1
u1dm1
m1du1
p1dV1 , (3)
Where u1is specific internal energy of gas in the cylinder, V1is a volume of cylinder, m1is amass of
gas entering in the cylinder, quantities with index mrelate to the manifold of pipeline.
We can express in equation 3 the values of enthalpy and internal energy through the product of temperature and het capacity at constant pressure cpand volume cv, according to equation
cpTm1dm1
cvT1dm1
cv m1du1
p1dV1
(4)
The equation of state of the real gas in the working cavity is written as
p1V1
zm1RT1
(5)
Where R is the gas constant and compressibility factor z determines the extent of nonideality of the working fluid. Numerical calculation of the compressibility factor z for air or nitrogen, accomplished in the approach proposed in ref.[7] for the real gases, showed that in the pressure range considered (P=0.1..3.5MPa) and temperature range T 210k the assumption that z = 1 is highly accurate, and the gas can be treated as ideal.
Substituting in equation 4 the value m1dT1,taken from equation 5 with the approximation z=1, and
using the common notation cp
cv
kand
cp cv R where k is the adiabatic exponent, after simple
transformations one can obtain the following expression
kRTm1dm1
V1dp1
kp1dV1
(6)
We can then replace the mass of air dmm1, entering the volume V1during time dt in equation 6, by the
corresponding value of the consumption function G1 (define as dmm1
equation relative to pressure
G1dt ) and then express the
dp1
kG1RTm1dt V
kp1
dV1 V
(7)
1 1
We can also determine the function G1 that describes the gas consumption. For this purpose we can find the parameters from Bernoullis equation. According to Bernoullis equation;
dh vdv 0
For isentropic flow this can be written as
dp1
1
vdv 0
(1a)
If the flow is incompressible than constant
Then equation (1a) on integration fields
1
dp1
1
1 dv2
2
constant
p1 1 v2
1 2
constant
(2a)
As we know when the flow is isentropically decelerated to zero velocity, the resultant pressure is the stagnation pressure. Therefore, when v=0,
p1 pm1,
1 m1,
Thus, constant in equation (2a) , pm1
m1
constant
Therefore, p1
1
1 v2
2
pm1 ,
m1
In incompressible flow since,
m1
1 constant,
v
p
1 2
1 2 1
pm1,
But for adiabatic,
h c T k RT , (3a)
0 p 0 k 1
RT pm1 , h
0
k
k pm1 , h
1 v2
0
1
h , h c T
k RT , kRT
2And RT
p1 ,
0 0
m1
2 k p
1 m1 2
p k 1
h 1 , (4a)
k 1 k 1 1
From equation (4a) & (3a),
2 1
v2
k 1 2
constant, (5a)
At T
0, h
0, v vmax, From equation (3a),
h 1 v2
, (6a)
0 2 max
At v 0,
0, From equation (5a),
2
constant h 0 , (7a)
0 k 1
By equation (5a),(6a)and (7a),
2 1 1 2
v2 v2 0 h ,
k 1 2 2
max k 1 0
Therefore,
k p1
1 v2
k pm1 , (8a)
1
k 1 1 2 k 1 m1
Then from (8a) we have,
k p 1 v2
constant,
k 1 2
Of the adiabatic braked stream(parameters of braking).in this case we will put the velocity of a stream equal to zero over the area of inlet valve into Bernoullis, then we have
k p1
1 v2
k pm1
(8)
1
k 1 1 2 k 1 m1
Now one can find the velocity of air entering the first cylinder from Eq. 8
v 2k
( pm1 p1 )
(9)
1 k 1
m1 1
From the adiabatic equation
p 1
1 ( 1 )k
m1 pm1
(10)
One can find the value of air density 1 in the working cylinder,
p 1
1 m1 ( 1 )k ,
pm1
By substituting expression (10) into formula (9) we have,
2k p k 1
v
1 k 1
RTm1[1 ( 1 ) k
pm1
], (11)
The gas consumption function can be defined as G1 v1 1 1 1, where a1 is the area of the cross section
of the inlet valve; 1is a coefficient of gas consumption through the inlet valve. Let us substitute the velocity of gas, determined by formula 11, into the expression for the gas consumption. Then we obtain
2k p k 1
G1 1 1 1 k
RTm1[1 ( 1 ) k ]
1 pm1
(12)
If we substitute the gas density in the cylinder from equation (10) to (12) we can derive
p 1 2k
p k 1
G1 m1 ( 1 )k 1 1
pm1 k
RTm1[1 ( 1 ) k ] =
1 pm1
2k p 2
p k 1
m1
G 2 RT [( 1 )k
( 1 ) k ] =
1 1 1 k
1 m1
pm1 pm1
2k 1
p 2 p k 1
G1 1 1 pm1 k
1 RTm1
[( 1 )kpm1
( 1 ) k ]
pm1
As a result, the function G1 of the gas consumption from the unlimited volume (manifold) can be determined by the formula [811]
2k 1
p 2 p k 1
G1 1 1 pm1 k
1 RTm1
[( 1 )kpm1
( 1 ) k ]
pm1
(13)
Where pm1 and Tm1 are the gas pressure and temperature in manifold pm1..
Note that the losses of gas pressure in the pipeline and local resistances are taken into account by introducing the coefficient of consumption Âµ [9, 11], which besides takes into account the compression of flow stream during exhaust, the speed of the gas as it approach the aperture and other factors. More often, this coefficient of consumption is determined experimentally or with the help of approximations. When the flow of gas occurs within a short span of a pipe at high velocity, the exhaust process is considered to be adiabatic, and it is possible to use formula 13. The process of compression of the gas in a working cavity is described by Eq. 6, which we can solve simultaneously with Eq. 13. As a result it is possible to determine the pressure p1 in cylinder as a function of time. We must note that this process cannot be described by means of elementary polytropic process with a constant parameter, n. In reality, it follows from formula 13, that the consumption G1 is a function of the pressure ratio, in which the numerator is the pressure of the medium into which gas flows, and the denominator is the pressure of the medium from which this gas moves.
We will present formula 13 for gas flow from a pipe in a more convenient form [9]
1
G 1 1kpm1 1
RTm1
2 k 1 p 2k
(14)
Where,
( )k
( ) k ; 1
1 ; K
pm1 k 1
In order to find the maximum of the gas consumption factor ,let us set its derivative to zero, which one can be written as
2 2 k 1 1
k 1 k 0,
k k
From which one can obtain the critical ratio of pressures
2 k
* [ ]k 1
(15)
k 1
It is necessary to distinguish between two regimes of the flow; subcritical, when the consumption function G1 is determined by the formula 13, and supercritical, at which the maximum critical consumption of gas G* is obtained after substitution of the critical ratio of pressures from Equation 15 into Equation 14[9]
*
G 1 1
pm1
* 2k
(k 1)RTm1
(16)
Where * 0.5282and for k=1.4
If we substitute the value of the critical gas consumption G* into Equation 7 instead of G1, we can obtain the equation that describes the process during the supercritical mode in a cavity of changing
volume
kG*RT dV
dp m1 dt kp 1 (17)
1 V 1 V
1 1
The analysis of Equation. 7, 13, 16 and 17 shows, that the process of change of gas state in a cylinder being filled, both for variable and for constant volume does not coincide with any one of the elementary thermodynamic processes, which occur with a constant polytropic exponent. The process can be described using the variable polytropic exponent n, which, in the beginning of the process, equals the adiabatic exponent, and then monotonically decreases,
n 1 [
0 (k
1) ]
(18)
Where, 0
p0 pm1
; At = 1, i. e. at the end of process, the value of polytropic exponent asymptotically
approaches the isothermal exponent n =1.
Change of the gas state by polytropic process (with a constat polytropic exponent) is possible in a working cylinder when there is a constant mass of gas. Some examples include the internal combustion engines after closure of the inlet valve, in flapvalve pneumatic motors, where the chamber is completely isolated by plates from the inlet and outlet ports, and in compressiondriven devices and accumulators, when there is expansion of the compressed gas. In the case of variable gas quantity in a cylinder, it is necessary to investigate the process by the energy balance stated in equ (6).
Thus, we can obtain the differential equation for determination of the gas pressure during the filling of cavity 1 in a general form, by substituting the value of gas consumption G1 from expression 14 into Equation (7)
dp1
1 1kp1 1
RTm1
kp1 dV1
(19)
dt V1 V1 dt
Lets consider the expansion stage of the pneumatic cylinder operation; the process of gas expansion during the movement of piston from point B to point D. We will describe this polytropic process with a parameter 1 n 1.4 , that allows us to take into account, to a first approximation, the possible process of heat exchange in the pneumatic cylinder.
Then the pressure in the cylinder can be determined by the following expression,
p1 p1B
( xB ts )n
x ts
(20)
Where xB is the corresponding distance traveled by the piston to arrive at position B (see Fig. 1); x is the current piston position; ts is a thickness of the piston; p1B is the gas pressure at the beginning of the expansion stage.
As a result it is necessary to solve the next equation for determining the dynamic characteristics of the piston motion in the pneumatic cylinder
0, x xA
d 2 x p s p s F
1 1 2 2 , xA x xD
(21)
dt 2 M
0, x xD
Thus the total set of the equation, describing the dynamics of the pneumatic cylinder can be presented as
dp1
1 1kp1 1
RTm1
kp1 dV1
dt V1
if x xB;
V1 dt
p1 p1B
0, x xA
( xB ts )n if x xB; x ts
d 2 x p s p s F
1 1 2 2 , xA x xD
dt 2 M
0, x xD
For the numerical solution of the derived differential parameters, describe the physical sizes and operating conditions of the pneumatic cylinder.


A computer program in MATLABSIMULINK and modeling in AUTOMOTIVE STUDIO(PNEUMATIC SIMULATION) has been developed for calculation of the dynamic parameters of the pneumatic cylinder by using the derived equations.
The following parameters and initial conditions were chosen for the simulation.
working fluidair/nitrogen, R=298.8 J/kgk ; adiabatic k=1.4 ; polytropic n=1.25exponents; initial pressure in the cylinder=.942bar ; pressure in the manifold=7bar ; gas temperature in the manifold=300k ; loading force=2N ; coefficient of friction=.15 ; gas consumption coefficient in cylinder=.7 ; diameter of input and output valve=6mm ; diameter of the piston=23mm ; initial piston speed=0 ; mass of piston=300gm ; piston stroke=50mm ; distance xA=2mm,xB=4mm,xC=47mm,xD=48mm ;
The result of numerical modeling of the gas pressure in the cylinder as the function of the piston
position,pressure,speed and time are presented in fig.5 ,6 and 7.
Now, lets analyze the physical process, which occur within the cylinder, in accordance with the numerical results obtained.
Until the force arising from the pressure in cylinder exceeds the force of friction and force of useful resistance. During this period there is filling of the constant volume of cylinder, and the piston remain motionless at point A.On the graph of pressure as a function of the piston motion this process is displayed by a practically vertical segment near the value of x=2mm.
The following stage of the pneumatic cylinder operation is the continuing process of filling of cylinder, for an already moving piston from point A up to B, the pressure in cylinder, which was establish after filling the cylinder while the piston was motionless, falls to 6bar.
After the piston passes point B, the outlet valve of exhaust is open, the inlet valve of the working cylinder is closed, and the piston continues to move to the right under the action of pressure force in cylinder, during this period pressure in cylinder decreases.Further movement of the piston there corresponds to the gas expansion with the closet inlet valve and open outlet valve. The final pressure in the cylinder (.924bar) accurately coincides with the pressure in the exhaust pipeline.
Note that the mathematical model developed here allows us to calculate the time dependence of the piston and velocity of the piston for pneumatic cylinder. The PV diagram is presented for the process under consideration, using the pressure dependences obtained for piston motion in the cylinder.

Data for each piston was entered into MATLAB Simulink computer simulation and AUTOMOTIVE STUDIO software, of pneumatic engine operation to generate flow velocities .These inputs include: the mass of the piston and all moving parts (piston, rod, rack, gear, etc.), useful areas of each piston, input and output pressures, gas temperature at inlet, and required piston velocity. The general equation of the piston motion can be written as follows:
d 2 x
M dt 2
p1s1
p2 s2 F
The MFR required for cylinder can be calculated with the following equation .
.
M VA
Power produced by the engine:
n p n 1
Power =

Q1 p2V1
i ( 1 ) n
n 1 p2
n p n 1
Work per cycle =
n
V 1
1 p2V1
i [( 1 ) n 1]
p2
1 c [(r )n 1]
V
vol p
s


Power(b.p.):0.683hp
Power(I.p.):0.853hp
Volumetric efficiency(vol ):0.72
Isothermal efficiency(iso ):0.89
A Good Approximation
Typical Compressor produces 4 CFM per 1 Hp 1 Hp = 0.746/0.9 = 0.829kW
Therefore, 1 CFM = 0.207kW 4rs/kwhr, 1 cfm = .828rs/hr It takes 0.69rs/km.
The three cylinder engine is fabricated, and mathematical model developed in this paper allows one to accomplish the numerical simulation of the working process and to determine the main dynamic characteristics of the pneumatic cylinder. By using the results of the numerical calculations, the analysis of the particulars of changes of gas pressure in cylinder can be accomplished. The subsequent stages of pneumatic cylinder operation can also be studied to calculate the main operational characteristics.
NOMENCLATURE: 

Asp 
specific useful work, kJ/kg; 
a1, a2 
Cross sections of inlet and outlet valves, m2; 
B 
Coeffient of friction; 
Cv 
Heat capacity at constant volume J / Kg. K ; 
Cp 
Heat capacity at constant pressure J / Kg. K; 
D 
Diameter of piston, m; 
Dd1,Dd2 
Rod diameters in cylinder, m; 
d1, d2 
Diameters of input and output valves, m; 
F,FL,FFR 
Resistance, loading and friction forces, N; 
Gas consumption factor; 

G1, G 2 
Gas consumption functions for cylinder, kg/s; 
GN2 
Specific gas consumption Kg / Kwhr 
Hm1 
Enthalpy of gas entering the cylinder, J; 
hm1 
Specific enthalpy, J/kg; 
K 
Adiabatic exponent; 
M 
Mass of piston, kg; 
/tr>
1 , 2 
Coefficients for gas consumption in cylinder; 
m1 
Mass of gas entering the cylinder, kg; 
N 
Working power of cylinder, kW; n polytropic exponent; 
F 
Frequency of operation, rot/min; p1, p 2 pressures in cylinder, MPa; 
pm1, pm2 
Pressures in 1st and 2nd manifolds, MPa; 
p10, p20 
Initial pressures in cylinder, MPa; 
R 
Gas constant J / Kg.K 
Qm1 
Thermal energy entering the cylinder, J; 
Qm2 
Thermal energy, which is removed from cylinder, J; 
1, 2 
Gas densities in cylinder, kg/m3; 
m1, m2 
Gas densities in 1st and 2nd pipelines, kg/m3; 
S1, S2 
Useful areas of piston for cylinder, m2; 
0, 1, 2, * 
Ratios of pressures; 
Tm1, Tm2 
Gas temperatures in 1st and 2nd manifolds, K; 
T1, T2 
Gas temperatures in cylinder, K; 
U1, U2 
Internal energies of gas in cylinder, J; 
u1 
Specific internal energy of gas in the cylinder, J/kg; 
V1, V2 
Volumes of cylinder, m3; 
1 
Velocity of air entering the 1st cylinder ,m/s; 
m2 
Velocity of air exhaust from thecylinder, m/s; 
W1, W2 
Works for gas expansion in cylinder, J; 
X 
Current piston position, m; 
xA,xB,xC,xD 
Distances (see Fig. 1), m; xE piston stroke, m; 
Z 
Compressibility factor 

Cryogenic engineering book, second Edison by Thomas m.flym,cryocee, inc., Louisville ,
Colorado.USA

Air Engine Design for Machining Class, Hugh Currin,April 11, 2007.

Sasa Trajkovic The Pneumatic Hybrid Vehicle, A New Concept for Fuel Consumption Reduction, Doctoral thesis, Division of Combustion Engines, Department of Energy Sciences, Faculty of Engineering, Lund university.

Michael Beeman; Design and Evaluation of an Advanced Adiabatic Compressed Air Energy Storage System at the MichiganUtah Mine ,A thesis submitted to the faculty of The University of Utah in partial fulfillment of the requirements for the degree of Master of Science Department of Mechanical Engineering, The University of Utah, August 2010

Gas dynamics maurite j.zucrow,jo D Hoffman, school of mechanical engineering Purdue university

Compressor hand book,Paul C. Hanlon Editor,McGrawHill.

Thermodynamics fundamentals for application john p.oconnell, j.m.haile,university of Virginia, Cambridge university.

Fluid power with application,6th ed., by Anthony Esposito,professor emeritus, department of manufacturing Engg.Miami University, oxford, Ohio.

J.m.tressler, T.clement, H.kazerooni, M.lim, Dynamic behavior of pneumatic system for lower extremity extenders, university of California at Berkeley,(2002 IEEE).

Y.Y. LinChen, J. Wang, Q.H. Wu Department of Electrical Engineering and Electronics, A software tool development for pneumatic actuator system simulation and design, The University of Liverpool, Brownlow Hill, Liverpool L69 3GJ, UK . 2003.

Tony Wong, Pascal Bigras, Department of Automated Manufacturing Engineering, A software application for visualizing and understanding hydraulic and pneumatic networks, Ã‰cole de technologies supÃ©rieure University of Quebec, Daniel Cervera ,Department of Mechanical Engineering Technologies, College de Valleyfield

Simulink, Getting Started Guide,R2011b (Matlab&Simulink)

Reciprocating Compressor Performance and Sizing Fundamentals, A Practical Guide to Compressor Technology, Second Edition, By Heinz P. Bloch Copyright Â© 2006 John Wiley & Sons, Inc.
P.G.Scholar ,
School of Mechanical & Building Sciences VIT University, Vellore14,
Tamilnadu, India
School of Mechanical & Building Sciences VIT University, Vellore14,
Tamilnadu, India