Control Oriented Gas Exchange Models for CVCP Engines and their Transient Sensitivity

Résumé — Modèles d’échange de gaz contrôlé pour les moteurs CVCP et leur sensibilité — Le présent article analyse un ensemble de modèles de contrôle de la phase d’échange des gaz dans les moteurs à synchronisation des cames variable et continue (CVCP). Ces modèles décrivent le ﬂux massif de gaz nouveaux et de gaz résiduels pris dans le cylindre au cours de la phase d’échange des gaz. Des simulations à l’aide de CVCP provisoires sont également effectuées afin d’analyser la performance des modèles au cours des essais provisoires. Abstract — Control Oriented Gas Exchange Models for CVCP Engines and their Transient Sensitivity — The paper analyzes a set of control oriented models for the gas exchange phase in engines with continuously variable cam phasing ( CVCP ). These models describe the mass ﬂow of fresh gases and the residual gases caught in the cylinder during the gas exchange phase. Simulations with CVCP transients are also performed to analyze the models performance during transients.


INTRODUCTION
There are high demands on the engine control units (ecu) since they are central for achieving good performance such as stable combustion, good torque response, good fuel economy and low emissions.New technologies are introduced to meet the increasing demands from legislators and customers and these must also be handled and properly controlled by the ecu.
Continously variable cam phasing (cvcp) provides an interesting possibility to reduce the fuel consumption but it also influences the amount of air and residual gases in the engine.
Air charge is important for engine fuel-and torque control while residual mass fraction is a crucial factor that limits stable engine operation since it influences the combustion variability.Therefore it is essential for the control system to know the air charge and residual mass fraction but these are not measured and hence models are required.When cam timing is used to its full extent the existing control oriented models for residual mass fraction and cylinder air charge have the shortcoming of not covering the entire engine operating region.Another issue is that transient effects introduced by the cam phasing mechanism have to be handled correctly.
Three models [2,6,8] that all predict residual mass fraction are investigated in this paper.The model in [6] also predicts air charge.
These models have previously been studied in [12] under static conditions.In particular their sensitivity to input data and model parameters were analyzed.In the paper it is concluded that while the more advanced models in [8] and [6] perform well on data from simulations they are sensitive to input data and model parameters, and that while the model in [2] is less sensitive it has other shortcomings.For example it is concluded that Models B and C are not very sensitive to errors in residual gas temperature but very sensitive to errors in mass flow and fresh gas temperature.
There are two main additions made in this paper.The first addition is the investigation of how much transients in cam phasing affect the air charge as well as residual mass fraction and the models ability to capture these effects.Additions have also been made to the models of the intake-and exhaust ducts which are here 1-D models, that captures the ram effects in the manifolds, instead of the static pressure models used in the earlier publication.
There are other publications, e.g.[3][4][5], that propose control and or estimation algorithms for different types of vvt systems, but the focus has been on the air charge and not on the residual gas fraction.In [4] and [5] air charge for dual equal and intake only cvcp systems, that have moderate valve overlap, is studied.In [3] the focus is on fuel injection for the same type of engines as mentioned earlier.

The Models
Model A [2] is a simple generalized flow restriction model created with physical insight.Two model parameters are tuned to measurements or simulations of residual mass fractions.
Model B [8], and Model C [6], are both based on thermodynamic relations, e.g. an energy balance at ivc.They both require measurement of cylinder pressure.
The reference model that is used is a crank angle based multi-zone model.A reference model is necessitated by the fact that it is very difficult to directly measure the actual residual mass fraction in a production engine.

The Investigation
An important issue in engine control are transients and the focus of this investigation is on how sensitive these three models are to transients in cvcp position.Therefore simulations have been performed to obtain cylinder pressure and reference data for static operating points and for transients in cam phasing.Only effects due to the different effective valve-lifts that is the results of transients in cvcp position are studied and the intake and exhaust manifold pressures are assumed to be measured.

REFERENCE MODEL
As reference cylinder model a crank angle based, single cylinder, zero-dimensional, cf.[9], multi-zone model with 1-D flow models for intake and exhaust ducts.The incylinder model is tailored to track the flows of fresh and burned gases during the gas exchange period.It upon the method presented in [7] that solves the in-cylinder differential equations.
In the reference model the engine cycle is divided into three phases, gas exchange, compression and combustion.During the gas exchange and the combustion phase two zones are used.During the compression phase only one zone is used.
The two zones in the gas exchange phase are used to separate the residual gas and the fresh charge.When entering the compression phase the gases are instantaneously mixed into one zone.In the combustion phase one zone is used to track the unburned gases and the other is used to track the burned gas.Combustion is modeled as a mass transfer from the unburned to the burned zone.
Flow into and out of the cylinder is modeled by standard flow restrictions with constant discharge coefficients and the valve opening areas are modeled using the geometries and valve lifts from the engine.Intake and exhaust manifold dynamics are modeled using a pipe connected to the valve in one end and to a constant pressure source in the other.The pipe is divided into a configurable number of sections that all have their own pressure, temperature and mixture.
All boundaries are associated with a plug-flow which defines the mass transfer between the sections.The velocity of the plug is modeled using Newtonian physics.This way the ram effects in the ducts are captured but further improvement could be made by using a multi cylinder model and including wave reflection for example from the supercharger.
Heating/cooling by the exhaust/intake ports is not considered.The fuel that enters the model is regarded as well mixed and no consideration is taken to charge cooling by evaporation.
In-cylinder heat transfer is modeled using [11] with constant wall temperature.For burn ratio the standard function [10] is used.
Thermochemical data for the model is obtained by using polynomials fitted to calculations using the chemical equilibrium program package, CHEPP, [1].For the specific heat values, c p and c v , dissociation is considered but the specific gas constant, R, is constant.The fuel used is isooctane.
Geometries and valve timings for the model are those of a 2.0 dm 3 turbocharged SAAB 9 3 Sport Sedan engine (SAAB B207R) equipped with a cvcp mechanism.A complete list of parameters and constants used in the reference model can be found in Table 5 in the Appendix.

Model A
Model A, in [2], defines an overlap factor (OF) that measures the amount of overlap between intake and exhaust valves.The overlap factor is then used as effective area in a generalized flow restriction model.The model has two parameters that need tuning.These parameters are tuned to residual mass fractions from measurements or simulations.The following quantities are needed when using the model: -Intake manifold pressure, p im -Exhaust manifold pressure, p em -Engine speed, N -Fuel air ratio, Φ -Valve positions (which gives OF)

Computational Scheme
The key element in Model A is the overlap factor which is defined as: where D i and D e are the inner seat diameters of the intake and exhaust valves respectively and V d is the displacement volume of the engine.A i and A e are the areas under the valve-lift/crank-angle curves and are defined as: where L i and L e are the intake and exhaust valve lifts respectively.In this implementation L i and L e have been calculated using the assumption that they are proportional to the valves effective flow areas, i.e.L i ∝ A eff,i Since the valve lifts in the calculation of A i and A e are small this has little effect on the calculation of OF.
Finally the residual mass fraction is calculated using the following expression: where C 1x rg and C 1x rg are constants that have to be fitted to data from measurements or simulations.These constants have been tuned to data from the reference model and the results are shown in Figure 9 in the Appendix.

Extension for Turbocharging
When using Model A for a turbocharged engine the pressure difference in (3), p em − p im , becomes negative for some operating points.To avoid taking the square root of negative numbers the absolute value of the pressure difference is used.The sign is then moved outside of the square root to obtain: where sgn is the sign operator.
Since the first term of the equation may now become negative it is possible to obtain negative residual mass fractions.This would indicate flow through of fresh gas to the exhaust manifold.The residual mass fraction in the cylinder is zero under these conditions.Therefore the residual mass fraction is set to zero when negative values are obtained.

Model B
Model B is presented in [8] and is a simplification of Model C. The model was originally used for an engine with external egr and vvt.External egr is omitted in this implementation since the engine used lacks external egr capabilities.The key element in this method is an energy balance for the mixing of ideal gases at ivc: (5) where: m fc is mass of fresh charge m rg is mass of residual gas c v,fc is specific heat at constant volume for fresh charge c v,rg is specific heat at constant volume for residual gas T fc is temperature for fresh charge T rg is temperature for residual gases T ivc is mean gas temperature at ivc Using this equation the residual mass can be calculated as long as T ivc , T fc and T rg are known.Estimates of T fc and T ivc are calculated using measurements of the following quantities: -Intake manifold pressure, p im -In cylinder pressure, p cyl -Intake manifold temperature, T im -Air mass entering the cylinder, m fc -Engine speed, N A good estimate of T rg , however, is trickier and is therefore obtained using a correlation between T rg and m Tot • N.

Algorithm
The first step of the algorithm in [8] is to adjust the cylinder pressure to the correct level.The reason for this is that the pressure measurement system used have a slowly varying offset.This adjustment is done by considering the cylinder pressure during the compression phase as a polytropic process for which p • V γ is constant.Here p is the pressure, V the volume and γ the ratio of specific heats.A special variant of this technique is used in [8] but since the model is validated on simulated data no compensation is needed and therefore the description of this technique is omitted.
The temperatature at ivc is calculated using cylinder pressure at ivc, measured air charge, m fc , an estimate of the residual mass fraction, x rg , and the following equations: In the next step the temperature of the fresh charge is estimated.In [6] the estimation of the fresh charge temperature for this model is discussed in detail.It is stated that the heat flux from the walls to the fresh charge is partly compensated for the fuel evaporating process.Therefore, the fresh charge temperature can be estimated using a polytropic compression from manifold conditions to in-cylinder conditions: where γ fc is the ratio of specific heat for the fresh charge.
When it comes to the temperature for the residual gas other measures have to be taken.In [6] a correlation between mass flow and burned gas temperature at ivc is used.The model is: (10) where the three constants C 1T rg , C 2T rg and C 3T rg are estimated using simulations.Here these parameters have been tuned to data from the reference model and results from the tuning are shown in Figure 10 in the Appendix.
Finally the new residual mass fraction can be estimated using the energy balance at ivc in (5): where the specific heat values at constant volume c v,fc , c v,rg and c v,ivc are the same as in the reference model.Given an initial estimate of x rg and using ( 6)-( 11) a new x rg can thus be calculated.The problem is solved by fixed point iteration which gives the final estimate of x rg .Unfortunately the original fixed point iteration does not always converge.In [8] this has been addressed and a stabilizing scheme has been suggested.Instead of updating the residual estimate with the new value in each iteration the following weighting is performed: The iterative process is summarized by the flowchart in Figure 1.

Extension for Charge Heating
In [6] charge cooling due to evaporation is said to be partly compensated by charge heating by heat transfer from the cylinder wall.In the reference model charge cooling due to evaporation is not included.Therefore charge heating has to be added to compensate for this.A crude model, proposed in [12], is that the fresh charge estimate is heated about 10% by the cylinder walls.Equation ( 9) is therefore replaced by:

Model C
Model C is presented in [6] and its key elements are energy balance equations at two different points during the engine cycle.Using these equations both cylinder air charge and residual mass fraction can be calculated.In addition to well known thermodynamic relations two approximations are used to get estimates of the residual mass temperature at ivc and the relative heat loss to cylinder walls at x b = 50%.

Yes
No ivc temperature evaluation Ideal gas law (Equation ( 8)) x rg evaluation Energy balance at IVC (Equation ( 11)) Output: x RG

Total mass calculation
Equation (6) Initial estimation of x rg x rg update x rg,i+1 = x rg,i Steplength adjustment (Equation ( 12)) Figure 1 Flowchart of Model B. Using an initial estimate of x rg a new value can be obtained using ( 6), ( 8) and (11).Repeating this in an iterative manner a final value of x rg is obtained.
The model requires on-line measurements of -Intake manifold pressure, p im -In cylinder pressure, p cyl -Intake manifold temperature, T im -Engine speed, N

Original Algorithm
Model C uses much of the same principles as Model B. As mentioned before Model B is in fact a simplification of Model C. The difference is that the fresh charge, m fc , is considered known in Model B but is estimated in Model C. Equations ( 7)-(13) are thus the same for Model C as for Model B.
One extra equation is needed to be able to estimate the total mass in the cylinder, m Tot , energy balance at 50% mass fraction burned, x b = 50%, is therefore introduced.
The position for x b = 50% is estimated using a burn rate analysis as follows: 1.The end of combustion is found by using the following relation ca eoc = argmax(p • V 1.15 ) + 10, where ca eoc is crank angle in degrees at end of combustion.2. Weighting values are calculated so that the final x b will be 0% at start of combustion and 100% at end of combustion.The weighting values are: where soc , eoc , tdc is start of combustion, end of combustion and top dead center respectively and γ comp and γ exp are the polytropic indices for compression and expansion phases respectively.The values of the latter are set to 1.32 and 1.27, as suggested in [6].3. x b is then calculated as: where γ(ca) is equal to γ comp for ca ≤ ca tdc and γ exp for ca > ca tdc .4. Finally the crank angle for x b = 50% is obtained from the estimated x b .The energy balance equation at x b = 50% is stated and solved to give: where T50 is the mean gas temperature at x b = 50%, q cool is a factor describing the heat loss to the cylinder wall, x c is a compensating factor for incomplete combustion, AF s is the stoichiometric air fuel ratio, λ is the normalized air fuel ratio and q HV is the higher heating value of the fuel.
Since q cool is unknown it has to be estimated and this is done using the following correlation: where C 1q cool and C 2q cool are constants that have to be tuned to data from simulations.Values obtained for the reference model can be found in the Appendix where Figure 11 shows the result of the fit.

Algorithm Modification
Using the c v,50 calculated from the gas composition, i.e. calculated in the same way as for the reference model, will not give accurate results when using the energy balance in ( 16) above.This is because: for the typical temperature interval (16) has been modified in this implementation.Another issue is that ca 50 ca ivc pdV is negative up till tdc since the volume is decreasing.Since work is added to the gas during this period a sign shift is necessary for the equation to be correct.The final expression is therefore: Note, once again, that c v,50 and c v,ivc are functions of their respective temperature.
Finally, when T50 is known a new total mass estimate can be performed.In the model no dissociation is considered.Therefore the gas constant for the composition is straight forward to calculate and is equal to R ivc .Hence the total mass in the cylinder can be calculated using the relations:

Complete Algorithm
The residual mass fraction, x rg , can be calculated using an initial estimate of x rg and the total mass in the cylinder, m Tot , much in the same way as in Model B. A flowchart of the process is shown in Figure 2.

INVESTIGATION SETUP
The cvcp mechanism on the engine has a cam phase shifting speed of 100 deg/s measured in terms of crank angle revolution.This corresponds to a change of about 3 degrees in intake valve duration if a transient in cam position is triggered at intake valve opening when running at 800 rpm's.The purpose of this investigation is to capture the effect that this difference in duration has on the air charge and residual mass fraction and to evaluate the performance of Models A -C under these conditions.Simulations have been performed to obtain cylinder pressure and reference data for static operating points and for transients in cam phasing.The simulations are designed to capture the difference in air charge and mass flow between The transients are simulated by steps up and down in cam phasing for all static operating points.The steps have been large enough not to finish within the actual cycle and all combinations of positive and negative steps have been performed for both the intake and exhaust cam.Altogether this yields 8 different types of steps.For transients in intake cam position the transient is triggered at ivo and for transients in exhaust cam position the transient is triggered at evo.
The static operating points are shown in Table 1.However, operating points that give residual mass fractions above 50% are removed since the engine does not run in these points.Altogether 36 of the 2430 simulations are removed because of this.

DEVIATIONS FOR TRANSIENT CYCLES
In Figure 4 the differences in residual mass fraction and total mass in the cylinder between static and transient cycles are shown.The relative difference in mass flow can be as large as 80% and the difference in residual mass fraction can be as large as 20%.Because of the low number of operating points for each variable it is not advantageous to plot the differences in mass flow and residual mass fraction against any of them.It is however possible to see a strong dependence on engine speed which occurs because the relative speed of the cvcp mechanism is higher for lower engine speeds.

EVALUATION OF MODEL PERFORMANCE
The transient and the static cycles are used to evaluate the performance of Models A -C.To be able to quantify the performance of the evaluated models the term detrended errors is introduced.The detrended error is defined as follows:

Static cycle
Transient cycle Comparison between static cycles and transient cycles.As can be seen both in (a) and (b) the differences are larger for low engine speeds.For 1000 rpm the maximum difference in mass flow is about 80% and for 3000 rpm the maximum difference is about 40%.The same effect is visible for the residual mass fraction where the maximum difference at 1000 rpm is about 20% while the maximum difference at 3000 rpm is about 10%.The reason for the higher impact at low speeds this is that the relative speed of the cvcp mechanism is higher for lower engine speeds.
where k and m are the values that minimize Σ i |e d,i | 2 and x is the estimate of x.The standard deviation of the detrended error may now be used to quantify the spreading of data around a thought straight-line fit of a plotted data set from an evaluated model.

Model A
The main advantage with Model A is that it does not require measurement of in cylinder pressure.Another advantage is that it directly yields the residual mass fraction by the use of simple calculations.As long as it is possible to precalculate everything that is not directly measurable this is very efficient.This is the case for static operating points since the overlap factor, OF, only depends on the difference in intake and exhaust maximum opening position imop and emop.For transients however the overlap factor has to be calculated for the actual effective valve lifts which may be different for every single cycle during the transient.Another issue with Model A is that while it is feasible to obtain data for the calculation of the constants C 1x rg and C 2x rg in static operating points it is problematic to include transient cycles in the model tuning.This might be a problem if the residual mass fractions is not only dependent of overlap factor but also the shape of the valve lift curves.
Altogether it is interesting to see how Model A performs if the overlap factors of the static cycles are used in contrast to the performance when the overlap factors are recalculated with the effective valve lifts of the transient cycles.
In Figure 5 the performance when used on static cycles (a), transient cycles with recalculated overlap factor (b) and transient cycles with overlap factor from the static cycle (c) are shown.Since there are 8 times as many different operating points for the transient cycles than for the static cycles it is hard to spot any difference between Figure 5(a) and (b).Therefore the standard deviation of the detrended errors have been used to quantize the differences.The standard deviations of the detrended errors are shown in Table 2.The increase in standard deviation is about 17% between the results in Figure 5(a) and (b).The difference between 5(b) and (c) is another 33% relative to 5(a).The increase in error in Figure 5(c) looks smaller at first sight but remember that the number of different operating points are the same as in Figure 5(b).

Model B
For Model B it might look like there is no reason why the performance should be different between the transient cycles and the static cycles since the estimation is based on measurable quantities that is not affected by the cam phasing.There are however some caveats.As shown in Figure 4(a) the difference in mass flow might be as large as 80% between two consecutive cycles.Model B relies on the measured air flow and since there are delays and dynamics in air flows it is likely that Model B is fed the wrong air mass flow during transients.For model B this gives two possible problems.
Firstly, the residual gas temperature of the first transient cycle is that of the last static cycle.When calculating the residual gas temperature it is necessary to use the mass flow of the last cycle.However, the model that predicts the residual gas temperature is a rather crude model and in [12] it is shown that Model B is not very sensitive to this type of error.Secondly, in [12] it is also concluded that the residual mass fraction estimate of Model B (and C) is sensitive to errors in its input quantities coupled to the fresh gas and that small errors in mass flow estimates may result in large errors in residual mass fraction estimates.
In Figure 6 the residual fraction estimates are shown for both transient and static cycles.In Figure 6(a) and (b) it can be seen that residual mass fraction estimates are almost as good for transient cycles as for static cycles.This is also supported by the standard deviation of the detrended errors in Table 3 where it can be seen that the increase in standard deviation is about 13%.Note that the standard deviations for Model B are one order of magnitude lower than the standard deviations for Model A. Using the wrong mass esti- mate does however have a detrimental effect on the residual mass fraction as can be seen in Figure 6(c).The increase in standard deviation is about 287% which takes it is close to the standard deviation of Model A.

Model C
In Figure 7 it can, in the same way as with Model B, be seen that the performance of the residual mass estimation for Model C does not degrade for transient cycles.The increase in standard deviation is as small as 3% as can be seen in Table 4.An interesting point is that while the residual mass estimates are worse than for Model B with the correct mass flow they are better than for Model B using the wrong mass flow.This is because Model C has the benefit of using its own mass flow estimate in the estimation process.
The performance of the mass flow estimation is shown in Figure 8.The difference in error between the mass flows estimates is hardly visible and the standard deviations differs with about 3.5%.

CONCLUSION
A cvcp engine has been simulated both in static operating points and with transients in cvcp position.It has been shown that transients in cam positions have a large impact on air charge and residual mass fraction.Three models, that predict residual mass fraction of which one also predicts air charge, have been investigated and the models' sensitivity to transients have been studied.
Model A is very simple and it does not require measurement of the cylinder pressure but it is also the model with the largest errors and highest sensitivity to transients.The errors for the transient cycles are of the same magnitude as the errors for the static cycles even if the method is not adjusted to handle the transients.It is however, amongst the evaluated models, the model that gives the largest errors both for static cycles and for transient cycles.Models B and C handles the transient cycles well when used correctly.The errors for the transient cycles are of the same magnitude as the errors for the static cycles for both models.
For static cycles it is advantageous to use an external measurement of air mass, as in Model B, since it makes the model less sensitive to input data.During transients however, if the external measurement is delayed, it is better to use Model C. Burned gas temperature at IVC The correlation between mass flow and temperature of residual gas at ivc is weak for a cvcp engine.
The parameters can be found by using a least square fit or by studying the data in Figure 9 Heat loss vs. mean gas temperature for a cvcp engine.
Looking at Figure 10 there is little chance of getting a good fit between the data and the model used.This is not surprising since a cvcp engine has been used instead of a conventional engine as in [6].For a cvcp engine the residual gas temperature, T rg , dependence on cam timing is much stronger than the dependence on mass flow.
The model has been tuned to the data in Figure 10  Model for q cool q cool is introduced in (17), Section 2.3.The equation is repeated here for reference: For q cool the results are somewhat better than for the T rg model.This can be seen in Figure 11.The model has been tuned to the data with the parameters: Copyright © 2007 Institut français du pétrole Permission to make digital or hard copies of part or all of this work for personal or classroom use is granted without fee provided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice and the full citation on the first page.Copyrights for components of this work owned by others than IFP must be honored.Abstracting with credit is permitted.To copy otherwise, to republish, to post on servers, or to redistribute to lists, requires prior specific permission and/or a fee: Request permission from Documentation, Institut français du pétrole, fax.+33 1 47 52 70 78, or revueogst@ifp.fr.

Figure 2 Flowchart
Figure 2Flowchart of Model C. Using an initial estimate of x rg and m Tot , new values can be obtained using (8),(11), (18) and (20).Repeating this in an iterative manner final values of x rg and m Tot is obtained.

Figure 3 Definition
Figure 3 Definition of transient and static cycle.The last cycle before the transient is referred to as the static cycle while the first cycle of the transient is referred to as the transient cycle.

Figure 5 Figure 6
Figure 5Residual gas fraction estimates, x rg , from Model A. The spread in x rg for the transient cycles are larger than for the static cycles but remember that there are 8 times as many operating points in (b) compared to (a)).Using the overlap factor from the static cycles and comparing with the residual mass fractions from the transient cycles decreases the performance even more (c).

Figure 7 Residual
Figure 7 Residual gas fraction estimates for Model C. The errors for the transient cycles in (a) are of the same size as for the static cycles in (b).Note once again that there are 8 times as many different transient cycles than static cycles.

8
Figure 8 Total mass estimates for Model C. The errors for the transient cycles in (a) are of the same size as for the static cycles in (b).Even though there are 8 times as many transient cycles as static cycles the difference is hardly visible.

Figure 9
Figure 9The constants C 1x rg and C 2x rg can be taken directly from the figure or fitted using a least square fit.
. The parameters have been fitted to:C 1x rg = 0.1578 C 2x rg = 0.4712Model for T RGThe model for T rg in (10), Section 2.2, is used both in Model B and Model C. The equation is repeated here for reference:T rg = −(C 1T rg • (m Tot • N)) C 2T rg + C 3T rg with the following values:C 1T rg = 0.02456 [-] C 2T rg = −0.6960[-] C 3T rg = 1101.2[K]

TABLE 2
Standard deviations of detrended errors for Model A

TABLE 4
Standard deviations of detrended errors for Model C