Paper Final

Proceedings of the Third International Conference on Modeling, Simulation and Applied Optimization Sharjah,U.A.E January 20-22 2009

THE SIMULATION OF TRANSIENT HEAT TRANSFER WITHIN PISTON IN SPARK IGNITION COMBUSTION ENGINES Dr. Mahmood Farzaneh Gord

Hamid Hajializadeh

Shahrood University of Technology The Faculty of Mechanical Engineering, Shahrood, Iran [email protected]

Shahrood University of Technology The Faculty of Mechanical Engineering, Shahrood, Iran [email protected]

ABSTRACT

The equation of state for an ideal gas is

A combined thermodynamic simulation and transient conduction heat transfer through piston has been carried out in this study. In the thermodynamic simulation, a computer program has been developed to study the full operation cycle of a four stroke internal combustion engine. The simulation used to calculate the pressure and temperature field existing in realistic engine combustion chambers for various engine parameters. The axisymmetric conduction heat transfer equation has been solved numerically to study transient temperature distribution within piston. The results show a good agreement with previous studies where applicable. INTRODUCTION

PV = mRT

(1)

Taking the logarithm of both sides and differentiating with respect to crank angle gives

1 ( PdV + VdP) = mdT + Tdm R

(2)

It can be assumed that the internal energy is only a function of the temperature and the combustion process is modeled as a single zone. The first law of thermodynamics can be written in differential form for an open thermodynamic system. If changes in potential energy are neglected, then

mC CV (dT ) + CV TC (dm) = ∂Q − PdV + (dm)hi

It is well known that heat transfer influences internal combustion engine (ICE) performance, efficiency and emissions significantly. Cooling the cylinder wall and piston head is desired because of problems such as thermal stresses, deterioration of the lubricating oil film and knock and pre ignition. On the other hand, an increase of heat transfer to the combustion chamber walls will lower the gas temperature and pressure within the cylinder, and this reduces the work per cycle transferred to the piston. The piston is one of critical places in term of thermal load as no effective way of cooling exists in reality. In order to obtain the temperature distribution in a typical ICE piston, a combined thermodynamic simulation and transient conduction heat transfer through piston has been done in this study. The thermodynamic simulation was based on first law of thermodynamics in which the combustion chamber was assumed an open system. The simulation has been carried out for full cycle (720 degree of crank angle) so the intake and exhaust stroke has been modeled. Axisymmetric transient conduction heat transfer equation was solved numerically in order to obtain temperature distribution through piston. The engine simulator and heat Conduction solver has been coupled to predict in cylinder properties and piston temperature. This approach allows more accurate simulations of engine combustion and heat transfer.

+ mi C P (dT ) − (dm)he + m e C P (dT )

(3)

With combining equation (2) and (4) and arranging it we get

CV C ) PdV = ∂Q − V VdP + (dm)hi + mi C P (dT ) (4) R R − (dm)he + m e C P (dT )

(1 +

Introducing

CV 1 and C P = CV K then rearranging = R K −1

equation (4) we have

dP K − 1 ∂Q K dV dm dT = − P + ( i ) RKTi + mi RK ( i ) dθ V ∂θ V dθ dθ dθ dme dTe −( ) RKTe − me RK ( ) dθ dθ (5) This equation has to be solved iteratively, and to do so, we should calculate its terms.

1.1.Engine & Heat Transfer Modeling As analytic functions cannot be used to describe engine processes, it is necessary to solve the governing equations on a step-wise basis.

ICMSA0’09-1

Proceedings of the Third International Conference on Modeling, Simulation and Applied Optimization Sharjah,U.A.E January 20-22 2009

  θ − θ  0 x b (θ ) = 1 − exp − a   ∆ θ   

m +1 

  

(10)

Where, θ is the crank angle, θ 0 is the start of combustion,

∆θ is the total combustion duration ( x b = 0 to x b = 1 ). a and m are adjustable parameters. Varying a and m changes the shape of the mass fraction burned curves. Actual mass fraction burned curves have fitted with a = 5 and m = 2. Ql , is the lost heat due to convection in combustion chamber. It can be calculated from ∂Ql =

hA π (T g − TW ) ω 180

(11)

Where, h is the heat transfer coefficient, A is the surface area in contact with the gases, T g is the mean charge temperature and TW is the wall temperature. Figure1.Engine geometric properties

To calculate convective heat transfer coefficient we use woschni correlation.

According to figure1, the volume of the engine cylinder and the rate of its change can be expressed, respectively, by

1 πB 2  2 a  V = s( ) ( ) + 1 − cos θ + (1 − cos 2θ ) 2 4  Cr − 1 4r 

(6)

dV 1 πB 2  a  = s( ) sin θ + sin 2θ  (7) dθ 2 4  2r 

h = 3.26b −0.2 P 0.8 T −0.55 v 0.8

Where, p is instantaneous cylinder pressure, b is bore diameter, T is instantaneous gas temperature and v is the characteristic velocity of gases. According to woschni correlation, the characteristic velocity of in cylinder gas for a four stroke engine without swirl is

Where, B is the cylinder internal diameter, S is the engine stroke, Cr is the engine compression ratio, a is the crank radius of the engine and r is the length of the connecting rod.

v = c1 S p + c 2

Equation (5) is valid whether or not dQ is interpreted as the

−180 ≤ θ ≤ θ 0

heat addition due to combustion or the heat lost by the gases in the cylinder because of convection. To admit both possibilities, we can write

dQ = Qin dx − dQl

V d Tr ( p − pm ) Pr V r

(13)

The constants are:

c1 = 2.28

c2 = 0

c1 = 2.28

c 2 = 3.24 × 10 −3

(14)

θ 0 ≤ θ ≤ +180

(15)

(8)

p m , is the motoring pressure that can be calculated from

Where Qin is the total value of the heat that can be released from the combustion of the quantity, m f , of fuel.

Qin = η C m f C fl

(12)

[( pm =

(9)

Where, C fl is the lower calorific value of the fuel and η c is the combustion efficiency. We consider the combustion efficiency, the value of 0.9. The heat addition for spark ignition engines may be a prescribed function of crank angle. The function that is generally used for calculating mass fraction burned during combustion is weib function.

CrV d γ ) pa ] Cr − 1 Vγ

(16)

The mass flow rate through a poppet valve is usually described by the equation for the compressible flow through a flow restriction. This equation is derived from a one dimensional isentropic flow analysis and the real gas flow effects are included by means of an experimentally determined discharge coefficient, CD . The mass flow rate is related to the upstream stagnation pressure p 0 and stagnation temperature T0 , static pressure just downstream the flow restriction p r and a reference area AR characteristic of the valve design.

ICMSA0’09-2

Proceedings of the Third International Conference on Modeling, Simulation and Applied Optimization Sharjah,U.A.E January 20-22 2009 γ −1     C A P P P  2γ   m = D R 0.50 ( T )  1 − ( T ) γ    P γ − 1 P ( RT0 ) 0 0     1 γ

.

0.5

(17)

chamber was a specified heat flux, which was obtained through engine simulation in first part of analysis. Computational model for this analysis is shown below.

γ

When the flow is choked, PT ≤  2  γ −1 the appropriate  

P0

 γ + 1

equation is .

m=

C D AR P0 ( RT0 ) 0.5

1 (γ ) 2

γ +1

 2  2(γ −1)    γ + 1

(18)

The value of C D for both intake and exhaust valves is taken 0.7 [3] and the reference area AR is the valve head area. The approach to determine the piston temperature distributions is deduced from the energy conservation law. For simplicity it is assumed that the engine geometry is axisymmetric. In this paper we calculate two dimensional temperature variations in piston. For two dimensional temperature distribution, the unsteady heat flow equation is

1 ∂ ∂T ∂ ∂T 1 ∂T (r )+ ( )= ( ) r ∂r ∂r ∂z ∂z α ∂t Where r

(21)

and z are the coordinates in the radial and axial

directions, respectively, T is the wall temperature, α =

k is ρc p

the thermal diffusivity and t is time. Using the explicit Saul’yev finite difference method, the explicit finite difference expression can be obtained. Ti n, j+1 − Ti n, j ∆t

=

h h 2α [h1 (1 + 2 )(Ti n+1, j − Ti n, j ) − h2 (1 − 1 )(Ti n, j+1 − Ti n−+1,1j )] + h1h2 h 2ri 2ri

2α [l1 (Ti n, j +1 − Ti n, j ) − l 2 (Ti n, j+1 − Ti n, j+−11 )] l1l 2l

(22)

Where

h1 = ri − ri −1 h2 = ri +1 − ri h = h1 + h2 l1 = z j − z j −1

(23)

l2 = z j +1 − z j l = l1 + l2 The equation then can be solved using iterative method. As the boundary conditions, the temperatures on the outer surfaces of the piston were treated as constant surface temperature conditions, in which the temperature equals either the measured coolant temperature or the oil temperature in the crankshaft case [4]. This is reasonable because the variation of the temperature on these surfaces is much smaller than that on the inner surfaces of the combustion chamber. The boundary condition on the inner gas-side surface of the combustion

Figure2.The grid structure and boundary condition used in heat transfer computation The two dimensional model has been considered as three blocks and temperature distribution on the joint grids of blocks boundary, preserved as an average temperature of up and down grids. 1.2.Results & Discussion Now the equation (5) can be solved iteratively using governed functions and equations of thermal parameters. Throughout the present work, a prototype single cylinder, four stroke engine is used to illustrate all of the thermodynamic points to be made. It has a bore of 100mm, a stroke of 111.1mm and a swept volume of 872.5cm3 and it works with speed of 3000rpm. The fuel is octane with a declared calorific value of 44.3 MJ/kg. Initial temperature of piston was considered 400K and engine simulation started with this value of temperature. In this case the piston was stationary and the material was an aluminum alloy with conductivity k=150W/mK, density p=2790kg/m3 and specific heat cp = 883J/kgK. The temperature on the outer boundaries was considered 375K [4]. After first simulation the temperature of piston will be modified and the program re simulates engine cycle till results get a good precision. In figures 3 and 4 the simulated pressure and work of the engine has been shown for various compression ratio of engine. The increasing cylinder pressure with greater compression ratio is easily seen. The higher compression ratio means compression into a smaller volume at TDC, raising the pressure and temperature at the end of compression.

ICMSA0’09-3

Proceedings of the Third International Conference on Modeling, Simulation and Applied Optimization Sharjah,U.A.E January 20-22 2009 In figure 5 the convective heat transfer coefficient that has been calculated using equation (12) is drawn for various compression ratios. Because of the direct dependency of convective heat transfer coefficient to pressure of the cylinder, at higher compression ratios there is an obvious increment in heat transfer coefficient. In cylinder temperature for engine cycle is also drawn in figure 6. Increasing pressure of the cylinder causes increment of in-cylinder temperature. But increasing the convective heat transfer coefficient, transfer this energy of fuel to the coolant and causes lower in-cylinder temperature at higher compression ratio.

Figure3. Simulated pressure of engine for compression ratio 5 to 11

Figure6. In cylinder temperature variation with crank angle for compression ratio 5 to 11

Figure4. Simulated pressure of engine for compression ratio 5 to 11 versus volume

Using the in-cylinder thermal parameters, the boundary condition for upper boundary of piston is achieved. The two dimensional temperature distribution of piston is calculated for the full cycle of engine. In figures 7 and 8 temperature distribution of piston are shown for compression ratio 5 and 11 at the moment that piston position is at TDC after ignition.

Figure5. Heat transfer coefficient variation with crank angle for compression ratio 5 to 11

Figure7. Temperature distribution of piston at TDC for compression ratio 5

ICMSA0’09-4

Proceedings of the Third International Conference on Modeling, Simulation and Applied Optimization Sharjah,U.A.E January 20-22 2009

Figure8. Temperature distribution of piston at TDC for compression ratio 11

Figure10. Temperature distribution of piston at TDC for fuel METHANOL

From the figures, it is obvious that mean temperature of piston is higher in higher compression ratio. Simulation repeated for various fuels and temperature profile calculated in the piston when compression ratio is 8. At table 1 some properties of these fuels have shown. Fuel

Formula

Molecular weight

LHV(Kj/Kg)

(A/F)stoi

Methane

CH4

16

50010

17.2

Propane

C3H8

44

46360

15.6

Methanol

CH4O

32

19910

6.5

Ethanol

C2H6O

46

26820

9

Table1. Alternative fuels with properties Figures 9 to 12 show the temperature distribution in the piston when piston is at TDC after ignition.

Figure11. Temperature distribution of piston at TDC for fuel ETHANOL

Figure9. Temperature distribution of piston at TDC for fuel METHANE

Figure12. Temperature distribution of piston at TDC for fuel PROPANE

ICMSA0’09-5

Proceedings of the Third International Conference on Modeling, Simulation and Applied Optimization Sharjah,U.A.E January 20-22 2009

CONCLUSIONS A combined thermodynamic simulation and transient conduction heat transfer through piston has been carried out in this study. In thermodynamic simulation, a computer program has been developed to study the full operation cycle of a four stroke internal combustion engine. The simulation used to calculate the pressure and temperature field existing in realistic engine combustion chambers for various engine parameters. The results of the thermodynamic simulation are agreed well with experimental and analytical results of previous studies by other researchers. The one dimensional conduction heat transfer equation has been solved numerically using an explicit finite difference method. The thermodynamic simulation and heat transfer solver has been combined to study transient temperature profile within piston. The computed temperature swings at different depths into the metal below the piston surface demonstrate that for a depth exceeding about 1.2mm, the temperature profile exhibits pseudo steady state characteristics. These results also have a good agreement with the results of Liu and Reitz [4]. REFERENCES Abd Alla, G.H. Computer simulation of a four stroke spark ignition engine, Energy Conversion and Management 43 (2002) 1043–1061 Farzaneh-Gord M., Maghrebi, M. J., Hajializadeh, H.,Optimizing Four Stroke Spark Ignition Engine Performance, The second International conference on Modeling, Simulation, And Applied optimization, Abu Dhabi, UAE, March 24-27 2007 Alkidas, A.C., Cole, R.M., Transient heat-flux measurements in a divided-chamber diesel engine, Trans. ASME, J. Heat Transfer 107 (1985) 439–444. Liu, Yong, Reitz, R. D., Modeling of heat conduction within chamber walls for multidimensional internal combustion engine simulations, Int J. Heat Mass Transfer. Vol. 41. Nos 6-7. pp. 859869. 1998 Rakopoulos, C. D.,. Rakopoulos, D. C., Mavropoulos, G. C., Giakoumis, E. G., Experimental and theoretical study of the short term response temperature transients in the cylinder walls of a diesel engine at various operating conditions, Applied Thermal Engineering 24 (2004), 679–702 Stone, Richard, ”Introduction to internal combustion engines”, Department of Engineering Science, University of Oxford, 1999 Heywood, John B.,” Internal Combustion Engine Fundamentals”, McGraw-Hill, Inc. 1988 M. Weclas, A. Melling and F. Durst, ” FLOW SEPARATION IN THE INLET VALVE GAP OF PISTON ENGINES”, Prog. Energy Combust. Sci. Vol 24, pp 165-195, 1998 Woschni G, Fieger J. Determination of local heat transfer coefficients at the piston of a high speed diesel engine by evaluation of measured temperature distribution. SAE 790834, 1979. Zeng Pin, Assanis, Dennis N, “Cylinder Pressure Reconstruction and its Application to Heat Transfer Analysis”, 2004 SAE World Congress Detroit, Michigan, March 8-11, 2004, 2004-01-0922

ICMSA0’09-6