On the Parametrisation of Turbocharger Power and Heat Transfer Models

6th IFAC Symposium Advances in Automotive Control Munich, Germany, July 12-14, 2010

On the Parametrisation of Turbocharger Power and Heat Transfer Models Matthias Mrosek, Rolf Isermann ∗ ∗

Institute of Automatic Control, Technische Universit¨ at Darmstadt, Landgraf-Georg-Str.4, 64283 Darmstadt, Germany (e-mail: [email protected], [email protected])

Abstract: A semi-physical turbocharger power and heat transfer model built in a CR-Diesel engine is presented. The model is based on physical parameters described by Euler’s equation. For the dependency of the turbine power from the VGT-actuator a polynomial approach is chosen. A method for measurement and model parametrisation separating aerodynamic power from heat transfer is given. Further an analytical analysis for choosing ramp times of quasi stationary measurements is included. The model quality and the proportion of heat transfer in the measured temperatures is demonstrated with measured data from the testbed. Finally the parameters of the modelled GT1749MV turbocharger are included in this contribution. Keywords: Turbocharger modelling, heat transfer, quasi stationary measurement 1. INTRODUCTION A further increasing system complexity demands dynamic simulation models of turbocharged Diesel engines for design and testing of engine control structures, diagnosis and software or hardware simulation. The turbocharger model strongly influences the quality of the air and exhaust path model. Turbocharger performance maps supplied by the turbocharger manufacturer are mainly of limited ranges and have to be inter- and extrapolated. Other approaches parametrise efficiency maps with engine testbed measurements. Since the power measurement based on measured temperature differences is affected by several heat transfers in the turbocharger housing and the bearings (Rautenberg et al., 1983) the quality of the final maps is not satisfactorily. Further typical efficiency maps for the turbomachine’s power are not capable to separate aerodynamic work from the heat transfer. Thus, a parametric semiphysical model is presented which allows to overcome this problem. A general methodology to parametrise VGTturbocharger power models and the heat transfer models will be given in the following. The turbocharger modelled is the GT1749MV turbocharger of an Opel DTH-Z19 engine at the testbed of the institute of automatic control.

In which the fraction qc,1 represents the specific heat input at the compressor inlet, while qc,2 represents the specific heat added after the adiabatic compression process. The specific heat input in the compressor results from the heat input Q˙ c and the compressor mass flow rate m ˙ c . Assuming the heat addition taking place at constant pressure and the compression process to be irreversible adiabatic ∆hc,adi , the real and the modelled change of state are shown in a h-s-diagram in Fig. 1 a). Similar assumptions for the heat transfer Q˙ t on the turbine side of the turbocharger can be made resulting in the fraction of heat transfer before the turbine qt,3 , after the turbine qt,4 and the drop in total enthalpy by adiabatic expansion ∆ht,adi . This results in the modelled and real states shown in the h-s-diagram in Fig. 1 b). With the turbine mass flow rate m ˙ t the specific heat transfer qt of the turbine follows as qt = Q˙ t /m ˙ t = qt,3 + qt,4 .

(2)

2. TURBOCHARGER MODEL Zahn and Isermann (2008) present a parametric model based on the mean line theory for the aerodynamic turbocharger power and the heat transfer. Caused by the heat transfer from the turbine via the housing, the compression can not be assumed to be adiabatic. In order to avoid a complex modelling of the heat transfer over the flow path, the specific heat input qc is divided into two fractions qc = Q˙ c /m ˙ c = qc,1 + qc,2 . 978-3-902661-72-2/10/$20.00 © 2010 IFAC

(1) 210

Fig. 1. Schematic h-s diagram of non-adiabatic compression and expansion Fig. 2 shows a schematic thermal and aerodynamic model of a VGT-turbocharger with the sensor positions. Zahn (2010) shows that there are many heat flows in the turbocharger housing, see also Rautenberg et al. (1983). The aerodynamic compressor power Pc and the turbine 10.3182/20100712-3-DE-2013.00075

AAC 2010 Munich, Germany, July 12-14, 2010

power Pt result in enthalpy changes described by the nonmeasurable temperatures Ti′ . Directly measurable are the turbocharger speed ntc and the position of the VGT-vanes st . In the thermal system possible heat flows by conduction Q˙ i,d are illustrated. The heat transfers from the turbine into the housing by conduction are described by Q˙ 3,d and Q˙ 4,d . Q˙ 1,d and Q˙ 2,d depict the heat addition on the compressor side. Furthermore Q˙ in,d and Q˙ out,d denote the heat transfers between the engine block, the oil and the turbocharger housing. Heat transfer from the housing to the environment can be modelled by convection Q˙ tc,v and radiation Q˙ tc,r . Measured quantities are the temperatures Ti comprising all heat transfers, the pressures pi , the compressor and turbine mass flow rates m ˙ c and m ˙ t . For simplification the pressures pi are regarded as the pressures directly measured at the compressor and the turbine.

a)

impeller inlet with prewhirl without prewhirl

b)

impeller outlet ideal actual

Fig. 3. Compressor velocity triangles cc2,u = µcc2,u,th = µ (uc2 + cc2,m cot(βc2,b )) .

(6)

Using the continuity equation and neglecting blade blockage the meridional component cc2,m in (6) results in cc2,m =

Fig. 2. VGT-turbocharger with measured variables and a schematic illustration of the heat flows

The adiabatic irreversible change in enthalpy ∆hc,adi represents the energy transferred to the impeller ac . It can be derived from Euler’s turbo-machine equation (3)

Here uc1 and uc2 are the circumferential velocities at compressor inlet and outlet. cc1,u and cc2,u are the corresponding circumferential components of the flow velocities. The flow entering the impeller has normally no prewhirl (cc1,u = 0). Hence Eq. (3) simplifies to (4)

The circumferential velocity at the compressor impeller outlet results from uc2 = πdc2 ntc

µ=

cc2,u cslip =1− cc2,u,th cc2,u,th

(8)

and

2.1 Compressor Model

ac = uc2 cc2,u .

(7)

Here ρ′2 is the gas density at impeller exit and bc2 the width at impeller exit. The slip factor µ in (6) is a widely used method to describe the deviation between the flow angle βc2 and the blade angle βc2,b . It is modelled by a classical approach proposed by Stodola (1945). The slip factor can be derived from the velocity triangles by

In the following sections the turbocharger model of Zahn and Isermann (2008) is briefly summarised.

ac = uc2 cc2,u − uc1 cc1,u .

m ˙c . ′ ρ2 πdc2 bc2

(5)

were dc2 is the impeller’s outer diameter. From the compressor velocity triangles in Fig. 3 the peripheral component of the absolute velocity at the impeller outlet is given as 211

cc2,u,th = uc2 + cc2,m cot (βc2,b ) .

(9)

According to Stodola (1945) cslip can be expressed as cslip = kslip uc2 in which kslip is a constant depending on the impeller geometry. Hence applying (7), (8) and (9) the slip factor results in µ=1−

kslip uc2 uc2 +

m ˙ c cot(βc2,b ) ρ′2 πdc2 bc2

.

(10)

Consequently substituting (5), (6) and (7) into (4), the aerodynamic power of the compressor results in   m ˙ c ntc cot(βc2,b ) 2 Pc = m ˙ c a c = µm ˙ c (πdc2 ntc ) + . ρ′2 bc2 (11) Thus the parameters to identify are dc2 ,

cot(βc2,b ) bc2

and kslip .

2.2 Turbine Model The aerodynamic work at the turbine is also determined by Euler’s equation and results summarised in

AAC 2010 Munich, Germany, July 12-14, 2010

Pt = m ˙ t at = −

m ˙ 2t ntc cot(αt3 ) . ρ′3 bt3

(12)

ρ′3 is the gas density and bt3 the blade width at the rotor inlet. Here αt3 is the absolute flow angle at the rotor inlet. It is determined by the angle of the guiding vanes. cot(αt3 )/bt3 is modelled as a third degree polynomial of the VGT-actuator position st as 3

cot(αt3 ) X = wt,i sit . bt3 i=0

(13)

2.3 Shaft and Heat Transfer Models According to Newton’s second law of motion the shaft acceleration results with inertia of the shaft Itc as

n˙ tc = −



1 2π

2

Pt + Pc + Pf . ntc Itc

(14)

Pf is the friction power of the journal bearing modelled as viscous friction 2

Pf = (2πntc ) kf .

(15)

In conclusion the additional aerodynamic model parameters wt,i for the turbine power and kf for the friction have to be identified. According to Fig. 1 the drop in total enthalpy for the turbine results from ht,3 − ht,4 = qt,3 + ∆ht,adi + qt,4 .

(16)

Rearranging gives the turbine exit enthalpy including all heat flows ht,4 = ht,3 −

Q˙ t − Pt . m ˙t

(17)

Similarly the compressor’s gain in enthalpy is calculated. For an ideal gas on the compressor side c or turbine side t its enthalpy results from h·,i = u(Ti , λi ) + RspecificTi ,

(18)

in which u represents the internal energy, Rspecific the specific gas constant and Ti the temperature of the gas. The internal energy u depending on gas temperature Ti and air mass fraction λi is calculated by the polynomial approach based on the tables published by Justi (1938), Lutz and Wolf (1938). According to Pischinger et al. (2002) R can be treated as constant for temperatures below 1500 K. The turbine enthalpy is calculated depending on the temperature and the gas composition. The compressor enthalpy is calculated assuming constant cP of air. Stationary heat transfer by conduction between location a and location b with the temperatures Ta and Tb is modelled by Newton’s law of cooling Q˙ ·,a/b,d = k·,a/b A·,a/b (Ta − Tb ) ,

(19) 212

with k·,a/b being the heat transmission coefficient and A·,a/b the surface area. c denotes heat transfer into the compressor, while t stands for heat transfer from the turbine. Since the surface area is hard to determine, the product of both coefficients is identified. Q˙ c is the sum of all heat transfers into the compressor housing, while Q˙ t represents all heat transfers from the turbine housing. 3. TESTBED SETUP AND DATA ACQUISITION Mrosek et al. (2010) investigate detailed the appropriate sensor and testbed setup for the parametrisation of air path models. Additionally to the setup of the measurement chain all external heat sources like λ-sensors or NOx sensors have to be disabled. A dismounting of the exothermic oxidation catalytic converter influencing T4 can be an advantage for measuring the real turbine exit temperature. For a better separation of aerodynamic power and the effects of heat transfer two different measurements are made. As proposed in Zahn (2010) the first measurement is a cold measurement with no injection and a hot measurement with injection, both covering a wide engine operation range. The next section will show that quasi stationary measurements are not suitable for the identification of processes with a slow process dynamics. Thus the turbocharger has to be parametrised with stationary measurements. Dwell times of 5 . . . 7 min have shown good results for the settling of all heat exchanges. A longer dwell is necessary when the engine speed or the injection quantity changes. Here settling times of 10 . . . 15 min showed good results. Finally a sensor offset compensation according to Mrosek et al. (2010) can further reduce the effect of measurement errors. 3.1 Limitations of quasi stationary excitation signals Contrary to stationary measurement quasi stationary excitation signals require a continuously measuring of data points (Fehl, 1992; Schwarte et al., 2004). Using only a weak excitation of the process’s dynamic the static relationship between process input u and output y can be identified saving measurement time. Fig. 4 shows a Hammerstein model assumption with an non-linearity u ˜ = f (u). According to Schwarte et al. (2004) the process dynamics can be simplified as a first lag with τP being the Porder n summed n time constants i=1 τi of the process.

Fig. 4. Hammerstein model assumption for quasistationary identification This non-linearity can be identified by a ramp excitation. Fig. 5a) shows the ramp excitation signal u ˜(t) with the ramp-up time τR /τP in respect to the process time constant. The dynamic error e decreases by mirroring the down-ramp at τR /τP and averaging it with the upramp (Schwarte et al., 2004). However, still a dynamic error occurs. Parametrising physical or semi-physical models with data affected by measurement error can have severe consequences for the model quality (Mrosek et al., 2010). In Fehl

AAC 2010 Munich, Germany, July 12-14, 2010

With the mean value of the up-ramp and the mirrored down-ramp the dynamic error can be expressed as

a) uniformly distributed in u ˜ and y b) uniformly distributed u

u ˜ [-]

1 0.5

e(˜ u) = u˜ −

0 0 1 y [-]

e(˜ u) = 0.5

e [-]

0 0.2 0.1 0

−0.1 0

= 10

τR 2 τP

=5

τR 2 τP

= 10

0.5 u ˜ [-]

τR 2 τP τR 2 τP

1 0

= 10

=5 τR 2 τP

0.5 u ˜ [-]

= 10 1

Fig. 5. Comparision of the dynamic error e caused by excitation signals beeing uniformly distributed a) in output y and b) in input u (1992) an investigation of the dynamic error resulting from an unidirectional ramp excitation is presented. Analysis of the dynamic error caused by bidirectional ramp excitation are not known to the author. Hence this error is investigated analytically here. For the analysis the process gain and unit are added to the non-linearity. With the process input u ˜ in Laplace domain τ ˜ (s) = 1 τP − 2 τP e− τRP U s2 τR s2 τR

(20)

the process response follows as Y (s) =

τ 2τP /τR τP /τR − τR P . − e s2 (τP s + 1) s2 (τP s + 1)

(21)

A transformation of (21) to the time domain yields the time response y´(t) of the up-ramp   τP  τR − t :0
(22)

   t−τR /τP 2τP τR − τP y`(t) = − t− − τP 1 − e τR τP   τ 2τR τP  R − t
(23)

and y`(t) of the down-ramp

Mirroring the down-ramp response y`(t) at τR /τP and mapping the time t to the appropriate process input u ˜ results in the image  τ P  t τR t→u ˜: τP  2 − t τR

τR :0
(25)

which finally results in

τR /τP 2τR /τP 0 τR /τP 2τR /τP Ratio of time constants Ratio of time constants

τR 2 τP

y´(˜ u) + y`(˜ u) , 2

(24)

213

τP2 − e τR

(1−u)τ ˜ R τ2 P

−

τP2 − e 2τR

uτ ˜ R τ2 P

−

τP2 − e 2τR

(2−u)τ ˜ R τ2 P

. (26)

Obviously the highest magnitude of the dynamic error occurs when the ramp signal reaches its peak at u ˜ = 1. Here with τR > τP2 the error can be approximated by e(1) ≈ τP2 /τR , which gives an easy tuning rule to design excitation signals. Interestingly the dynamic error resulting from quasi stationary excitation signals rises quadratically with the time constant of the process. This limits quasi stationary identification to processes with a fast dynamic response. For slow dynamic processes like thermal processes, stationary measurement is favourable. Further a process with constant dynamics is considered here. Regarding processes with non-linear time constants the slowest dynamics has to be considered for the choice of the ramp time. Fig. 5 b) shows the ramp excitation u on the static non-linearity u ˜ = u2 . The data samples are uniformly distributed over time in the input u. Shaped by the non-linearity the input u˜ of the linear model changes its distribution. This causes an excitation of the linear model at different rates. Regarding the model output y and the dynamic error e, an error twice as large as in Fig. 5 a) occurs. Practically there are three ways to avoid this increase in the dynamic error. Obviously an increase in the ramp time reduces the error. Secondly an approximate shaping of the process input u to reach a uniformly distributed u ˜ can lessen the error. Finally adjusting the output y in a ramp leads to a uniformly distributed u˜. For the last solution a closed loop control of y is necessary. 4. MODEL PARAMETRISATION Fig. 6 shows the parametrisation process of the aerodynamic compressor power (11) and the heat transferred into the compressor housing. γP,c denotes the parameters for the aerodynamic model, while γQ,c are the parameters ˙ of the heat transfer model. Initially the heat transfer is neglected. Knowing the impeller diameter dc2 from measurements or the model number of the turbocharger, Eq. (11) can be expressed as a linear optimisation problem. A least-squares estimation delivers the initial parameters for γP,c . The parameters are estimated by applying a nonlinear constrained optimisation with the Matlab function fmincon. In the first step the aerodynamic temperatures T1′ and T2′ are estimated from the heat transfer model. Then the density at the compressor outlet ρ′2 as well as the aerodynamic enthalpy difference ∆hc,adi are calculated. Evaluation of the compressor power and taking the L2norm with ∆hc,adi m ˙ c delivers the quality criterion for the optimisation. The directions of heat transfer are constrained to physical meaningful values. All possible heat transfers from T3 , T4 , Teng and Toil to T1 and T2 are tested.

AAC 2010 Munich, Germany, July 12-14, 2010

The best results showed the model with the summed two heat transfers Q˙ c,eng/1 and Q˙ c,3/1 .

Neglecting the heat transfer a least-squares estimation delivers the initial compressor power parameters γP,t . The constrained non-linear optimisation is largely common to the compressor power optimisation. After estimating the heat transfer, the pure aerodynamic enthalpy ∆ht,adi and the density at turbine entry ρ′3 are estimated. For the turbine there are two quality criteria JEuler and Jenthalpy weighted by wEuler and wenthalpy . Jenthalpy is the criterion for the heat transfers. The aerodynamic power and the friction are described by JEuler. For simulation the matching of turbine power, friction and compressor power is directly affecting the turbocharger speed and is therefore crucial for the model quality. Jenthalpy only affects the exhaust gas temperature. Therefore the weighting is chosen here as wEuler = 2 and wenthalpy = 1. The directions of heat transfer are constrained to physical meaningful values.

Fig. 6. Optimisation algorithm to determine the compressor power model parameters γPc and heat transfer parameters γQc ˙ Fig. 7 shows the results for the compressor model and the heat transfer model Q˙ c . The measurements represented by ∆hc,adi m ˙ c are relieved from the heat transfer. Measurements 1-74 represent the cold measurement, data points 75-176 are measured with injection. Even in the zoomed section with low turbocharger power the model fits the data well. The averaged root mean squared error of a seven fold cross validation is 39.6 W . At operation points with low compressor load the modelled heat transfer into the compressor housing has a large effect with respect to the aerodynamic power. Especially for the data points 1 to 15 the contribution of enthalpy by the heat transfer is large. The heat transfer into the compressor housing ranges from 100 W to 300 W.

Fig. 8. Optimisation algorithm to determine the turbine power model parameter γP,t the heat transfer parameter γQ,t and the shaft friction coefficient kf . ˙ Tab.1 shows the averaged root mean square error RM SE of a seven-fold cross-validation. It determines the optimal model structure. Modelling the turbine power Pt by (12) and cot (αt3 ) /bt3 as polynomial of order 3 showed good validation results. A heat transfer between T3 and T1 could be estimated. Further modelling showed that there is a dependency to the exhaust gas density ρ′3 of fourth order and a term ρ′3 st described by 3

cot(αt3 ) X wt,i sit + wt,ρ′3 ,1 st ρ′3 + wρ3 ,4 ρ′4 = 3. bt3 i=0 Fig. 7. Model of the compressor power Pc with the corrected aerodynamic temperatures T1′ and T2′ and the modelled heat transfer into the compressor housing Q˙ c , zoom into data with low power Next Fig. 8 shows the optimisation of the turbine power parameters γP,t , the friction coefficient kf and the heat transfers γQ,t out of the turbine housing. First the aero˙ dynamic power of the turbine is approximated by (13). 214

(27)

For this reason the engine cold measurement helps to reach densities in the exhaust gas between 3 . . . 4 kg/m3 to ease the parametrisation of the of the ρ′3 dependency in the turbine power model. Considering an additional heat transfer between T3 and T4 results in the best model. However, the heat transfer Q˙ t,3/4 seems not to be physical. Eq. (19) shows that the heat transfer depends only on the differences of two temperatures. So T4 might be substitutional for a not measured temperature rising and falling in the same manner. Possibly it might be the not locally measured oil temperature cooling the

AAC 2010 Munich, Germany, July 12-14, 2010

turbocharger housing. A heat transfer by radiation could not be modelled in a physical meaningful direction. Table 1. Turbine model cross-validation Model structure Pt = f (st ); Q˙ t = f (T1 , T3 ) Pt = f (st ); Q˙ t = f (T1 , T3 , T4 ) Pt = f (st , ρ′3 ); Q˙ t = f (T1 , T3 ) Pt = f (st , ρ′ ); Q˙ t = f (T1 , T3 , T4 ) 3

RM SE Euler

RM SE enthalpy

171 W 168 W 99 W 84 W

263 W 235 W 232 W 158 W

Fig. 9 shows the results for the turbine power model Pt and the heat transfer model Q˙ t . The turbine power is compared to the sum of friction losses Pf and the compressor power Pc . Here the model quality for the aerodynamic turbine power is slightly worse that the aerodynamic compressor power model. Depending on the engine operation point the modelled heat transfer varies from 100 W to 2000 W.

Fig. 9. Model of the turbine power Pt with the corrected aerodynamic temperatures T3′ and the modelled heat transfer Q˙ t , zoom into data section with low power Finally Tab. 2 shows the identified model parameters for the turbocharger power and the heat transfer models. Whereas st is normalised between 0 and 1. Table 2. Model parameters Compressor parameters dc2 kslip

cot(βc2,b ) bc2

kc,3/1 Ac,3/1 kc,eng/1 Ac,eng/1

wt,0 wt,1

254.04 m−1 155.02 m−1

−270.07 m−1 0.3415 W/K 2.01 W/K

wt,2 wt,3 wt,ρ′ ,1 3 wρ′ ,4 3 kt,3/1 At,3/1

−66.18 m−1 426.08 m−1 115.39 m2 /kg 3.83 m11 /kg4 0.82 W/K

kt,3/4 At,3/4

13.38 W/K

Friction parameter kf

Turbine parameters

0.049 m 0.1686

3.23e-6 kgm2 /s

5. CONCLUSION Parametrising the turbocharger power model with this methodology results in excellent models. Validation results for the complete air path model in the time domain can be found in Mrosek et al. (2010) and Mrosek and Isermann (2009). It is analytically shown that quasi stationary excitation signals are limited to dynamical fast processes, 215

since the measurement error rises quadratically with the proportion of the process time constant. Mainly depending on dynamical slow temperatures, the turbocharger power models should be parametrised by stationary measurements. Here a cold and hot engine measurement in conjunction with the optimisation methodology allow to separate the heat transfer and aerodynamic power. It is shown that especially at low turbocharger power the heat transfer has a large impact on the measured enthalpy difference. Extending the turbine power model of Zahn (2010) by the model input ρ′3 and the heat transfer Q˙ t,3/4 significantly improves the model quality. REFERENCES Fehl, G. (1992). Entwicklung einer Meßstrategie zur kontinuierlichen Erfassung von Kennlinien auf Motorenpr¨ ufst¨ anden. Dissertation, Technische Hochschule Darmstadt. Justi, E.W. (1938). Spezifische W¨ arme – Enthalpie, Entropie und Dissoziation technischer Gase. Springer, Berlin. Lutz, O. and Wolf, F. (1938). IS-Tafel f¨ ur Luft und Verbrennungsgase. Springer, Berlin. Mrosek, M. and Isermann, R. (2009). System properties and control of turbocharged Diesel engines with highand low-pressure EGR. In 2009 IFAC Workshop on Engine and Powertrain Control, Simulation and Modeling. Rueil-Malmaison, France. Mrosek, M., Zahn, S., and Isermann, R. (2010). Parameter estimation for physical based air path models of turbocharged Diesel engines - An experience based guidance. SAE International Journal of Engines, 2(2), 570–583. Pischinger, R., Klell, M., and Sams, T. (2002). Thermodynamik der Verbrennungskraftmaschine (Der Fahrzeugantrieb). Springer, Wien, New York, 2. edition. Rautenberg, M., Mobarak, A., and Malobabic, M. (1983). Influence of heat transfer between turbine and compressor on the performance of small turbochargers. In International Gas Turbine Congress, JSME Paper 83TOKYO-IGTC-73. Tokyo. Schwarte, A., Hack, L., Isermann, R., Nitzke, H.G., Jeschke, J., and Piewek, J. (2004). Automated calibration of engine control units with continuous engine data recording. In Steuerung und Regelung von Fahrzeugen und Motoren - AUTOREG 2004, 651–663. Wiesloch. Stodola, A. (1945). Steam and Gas Turbines. McGrawHill, New York. Reprinted by Peter Smith. Zahn, S. (2010). Elektronisches Management motorischer Fahrzeugantriebe: Elektronik-Architektur, Modellbildung, Regelung und Diagnose von Verbrennungsmotoren, Getrieben und Elektroantrieben, chapter Mittelwert- und Arbeitstaktsynchrone Simulation von Dieselmotoren. Vieweg+Teubner. Zahn, S. and Isermann, R. (2008). Crank angle synchronous modelling and real-time simulation of diesel engines for ECU function development and testing. In 9th Int. Symp. On Advanced Vehicle Control. Kobe.