© AFM, EDP Sciences 2017

1 Introduction

As engineering systems become increasingly sophisticated, the complexity of the problems faced by engineers also increases. The lack of efficient physical models is one consequence of this complexity. Empirical studies must be developed to understand the behavior of systems, to enhance models and allow for design optimization. This is particularly true for engine design, where typically, the engine is controlled by a large number of parameters in order to meet multiple performance objectives [1]. It is complicated to take all of them into account in an optimization process.

The lack of good models also has an impact on optimization processes because new designs need to be tested in order to define their performances. In the field of engines, engineers usually have two ways to validate the optimal design found: to run tests with a prototype or to analyze the results of CFD simulations. Both solutions are expensive so a third option is proposed in this paper based on a kriging method.

The kriging method is already used in various fields, including the field of engines, and for multi-objective optimization issues. In particular, kriging is fully integrated into the design space exploration phase using Response Surface Models. Multiple cylinder designs have already been improved using this method, such as Jeong et al. [2], Castric [3] or Brahmi et al. [1]. In the design optimization process presented in this paper, kriging is used in a different way: it is not part of the exploration of the design space, but it used after this phase. Indeed, after finding the optimal cylinder design, the behavior of gas flows inside the optimal chamber is interpolated from previous results thanks to the kriging technique and the gas flow distribution can be visualized. This visualization of the interpolated flow behavior is called the “Fast Kriging Simulation”.

To validate the method, it has been applied to the optimization of a 2-stroke Diesel engine with ports. The complexity of the engine cycle led us to focus on scavenging. Less studied than combustion, the scavenging process also plays a major role in the formation of pollutants. Indeed, it is during scavenging that fresh air enters the cylinder pushing out the burnt gases. As with all fluid flows, any change in the cylinder design has consequences on fluid dynamics. Thus, improving cylinder design cannot be done without studying scavenging behavior.

After defining the cylinder variables influencing the scavenging process [4], a new model was developed: it explicitly integrates these variables. Based on this model, the design space was explored converging to different design solutions depending on the form of the fitness function. For the designer, the last step in design optimization is to finally test the different designs which have been defined.

First, the methodology we developed is introduced. Then, we present the β-NTF and kriging methods. β-NTF is used to establish the data set needed for kriging interpolation. Finally, the numerical results are given followed by a discussion underlining the interest of the method.

2 Overall methodology and case of application

Kriging visualization is a part of the overall design optimization method we developed. It is used as the final step in our methodology. The method was applied to enhance the cylinder design of a 2-stroke Diesel engine. The particularity of our engine is the use of ports for both the intake and emission of gases. This creates a singular flow in the cylinder during the scavenging process which justifies us focusing on this part of the engine cycle.

The whole methodology is represented in Figure 1 and is briefly summarized in this section.

The first step is to define the design parameters influencing the flow and their range of values [4]. Then 27 cylinder configurations were simulated thanks to CFD models. The evolution of distribution of gas in the chamber during the scavenging process for each configuration was extracted [5]: these data will be used for kriging approximation as the known points. Based on CFD results, certain scavenging models were developed and used to determine the most suitable cylinder configuration.

After finding an optimal cylinder configuration, the designer needs to evaluate its new performances in order to validate the new design. Experimental tests on prototypes or numerical simulations are the most common ways to do so. However, they are both expensive: CFD is generally time-consuming; the experiments are expensive unique models and require special instruments (e.g. transparent cylinder). Before investing time and money in simulations or experiments, the designer should be certain to have the best design. In this approach, we have developed a transitional solution based on the kriging method. It enables in a very short time for the distribution of gases to be visualized during the scavenging process whatever the cylinder design. Before going further, the designer can do a qualitative analysis and evaluate the relevance of the optimal solution. This is the purpose of the last step of the methodology presented in this article.

Fig. 1 Method of model reduction and optimization.

3 NTF reduced model and Kriging

The kriging method implies the use of a data set. Twenty-seven engine configurations have been simulated, the simulations last 174 crankshaft degrees. Based on CFD simulations and results, several quantities (fields of pressure, temperature, speed, etc.) have been extracted at each crankshaft angle degree. This represents over 30 million data items. Not all of them can be used with the kriging method; it is too much data to deal with. A model reduction was carried out using the Non-negative Tensor Factorization (NTF) algorithm to prepare the data for kriging exploitation.

3.1 The β-NTF reduction

NTF was first proposed by Shashua et al. [7]. This algorithm is a generalization of Non-negative Matrix factorization (NMF) presented by Lee et al. [8]. NTF is based on a Parallel Factor Analysis (PARAFAC) analysis. Its particularity is that it imposes non-negative constraints on tensor and factor matrices. PARAFAC is a tensor (multiway array) factorization method which allows us to find hidden factors (component matrices) from multidimensional data.

We consider as a 3D tensor, . We have no negative value y_ijk ≥ 0   ∀ i, j, k. The NTF concept is illustrated in Figure 2. (1)

With the complete model, the reduced model and the error or the noise, the decomposition generates a set of three unknown matrices U⁽ⁿ⁾ such as , with n the number of modes and l = (I, J, K). The NTF algorithm is very attractive because of its ability to take into account spatial and temporal correlations between variables more accurately than 2D NMF. As is the case with NMF, NTF also provides greater stability and a unique solution, as well as meaningful latent (hidden) components or features with physical or physiological meaning and interpretation [6]. Finally, the NTF algorithm is very simple to implement in python and provides powerful implementation with multi-array data. The NTF method also usually provides sparse common factors or hidden (latent) components.

In this study, in addition to its efficiency and speed, the capability of the NTF algorithm to deal with spatial and temporal correlations between variables is highly useful. Indeed, the correlations between variables have to be identified depending on numerous parameters, such as spatial localization in the cylinder, temporal moment during the engine cycle and the cylinder's design configuration. Moreover, the spatial, temporal and parametric dimensions are totally independent of each other.

Fig. 2 NTF concept.

3.2 Application to the data set

The data are collected from the CFD results of the 27 cylinders tested (defined by an experiment design). The data extracted from the CFD results characterize gas distribution, temperature, pressure, etc. Only gas distribution will be presented here, but the method has been applied to all the fields extracted from the CFD simulations.

Because of the cylinder's piston motion during scavenging, some cells are added or deleted between two time steps. That is why all the results have been projected onto the same regular mesh (with 6561 square cells). Each cell is associated with its spatial location in the cylinder (from the regular mesh), the configuration of the cylinder (i.e. the value of each design parameter), the crankshaft angle at which the distribution is represented and the mass fraction of fresh gases (the output we want to forecast). The scavenging lasts 171 crankshaft angles: the results are collected from θ = 95 crank degrees to θ = 266 crank degrees after Top Dead Center (TDC). The piston motion in the cylinder is considered between 2 extreme positions: the TDC represents the highest position, whereas the Bottom Dead Center (BDC) is the lowest, as illustrated in Figure 3.

Due to the number of cylinder configurations, the duration of the scavenging process and the number of mesh cells, the number of known points is equal to 30 823 578. The data set therefore had to be reduced: this was done thanks to the β-NTF method.

To use the NTF algorithm, the data should be organized so as to obtain the model. The three dimensions were respectively associated with space, time and input parameter combinations (cf. Fig. 4)

The NTF decomposition was done with 160 modes, the best compromise between efficiency and accuracy. The average relative error is 2.5% considering all the data. The average error reaches 8% at 180 crankshaft degrees (when the ports are fully opened and the gas exchanges are greatest).

The maximum error of 19.6% is represented in Figure 5. In this figure we can observe the gas distribution in a cylinder before and after β-NTF reduction.

The data set now consists of 1 081 920 points, which represents a reduction of 96.5% compared to the CFD data thanks to the NTF method.

Now the NTF model is established. The 1 081 920 points included in the model are associated with the two spatial coordinates, the temporal moment during the engine cycle and the seven parameters of the cylinder design. All of this constitutes the data set of known points necessary to apply the kriging method. Without the reduction, the number of points is too high to apply the kriging method: the memory is full and cannot deal with all the data points, and solve the kriging algorithm at the same time.

Fig. 3 TDC and BDC piston positions.

Fig. 4 NTF model of engine quantity.

Fig. 5 Example of fresh gas distribution from (a) CFD results (b) β-NTF reduced model.

3.3 Kriging principles

Kriging is an optimal interpolation technique based on regression against observed Z values of surrounding data points, weighted according to spatial covariance values. The method was first developed by Krige [9] for the spatial interpolation of soil properties, but Matheron [10] was the first person to formalize the method. The kriging principle is to interpolate the value of an unknown function Z at a located point, given the values of function Z at other points. The interest of the method lies in its capability to take into account not only the distance from the evaluated point to the known points, but also the distance between two known points.

Many kriging algorithms have been developed since 1973 [10], but they are all based on the same basic linear regression estimator Z*(u) defined as [11]: (2)

With u the point where the Z function is estimated and u_α the neighborhood points; n(u) is the number of local points used to estimate Z*(u); m(u) and m(u) are the expected values of Z respectively for points u and u_α and λ_α(u) the kriging weights associated to u_α estimate the point u.

m(u) is the trend component associated with the random field Z(u). The aim of the kriging method is to evaluate and to lower the residual R(u) = Z(u) – m(u) as the weighted sum of the surrounding residuals calculated for the u_α points. Therefore, the key point of this technique is the way the kriging weights are evaluated: they are derived from covariance function or semi-variogram.

3.4 The semi-variogram

The semi-variogram is a function representing the spatial dependency of the known points and is obtained from the stationarity definition. Indeed, the stationarity hypothesis is an essential condition for using kriging method. So the variation of a data set is only dependent on the distance r between two locations where variables values are Z(u_i + h) and Z(u_i) with r = |h|, can be given by the following semi-variogram: (3)

With N(r) ={ (i, j)such as|u_i − u_j| = r } is the pair number of Z(u_i + h) and Z(u_i) and is the experimental semi-variogram.

The experimental semi-variogram presented in equation (3) estimates the theoretical semi-variogram only for a finite number of distances. It is assumed that two spatially close points should have similar values and the similarity decreases as the distance increases.

Then, the variogram is established as a model adjustment based on the semi-variogram points. Various forms of variogram model are available. If the semi-variogram represents the spatial auto-correlation between the known points, it does not indicate the possible directions or distances which justify the need for the model's adjustment. This is usually done by linear regression or by the least squares method.

Finally the weights λ_α in (2) are computed solving the following system: (4)

With (5) (6) (7) λ is the Lagrange multiplier.

The variance of the square of the standard error each point is obtained thanks to the relationship: (8) indicates the incertainity of the estimated value. The accuracy of the prediction value depends directly on the distance between the evaluated point and the known points. If we consider the estimation errors normally distributed around the true value, the probability that the true value will be in Z(u_i) ± s_i is 68%, while the probability that the true value will be in Z(u_i) ± 2s_i is 95% (Davis [12]).

4 Application to the cylinder design

4.1 The data set

The use of the kriging method is motivated by its capability to take into account all the information relating to the data points in order to determine the weights λ. Not only the spatial coordinates are compared but also the temporal parameter (crankshaft angle at which the quantity is evaluated) and the 7 variables characterizing the cylinder design [4]:

• the height of the exhaust port θ_{end_exh_port} which is the crankshaft angle (before the BDC) the exhaust port is fully uncovered;

• the advance of opening of intake and exhaust ports, θ_{in_advance} and θ_{exh_advance}, are respectively the crankshaft angles at which the piston begins to uncover the intake and exhaust ports;

• the boost pressure P_boost, pressure of intake gases;

• the difference between intake and exhaust pressures ΔP;

• the angles β_in and β_exh are respectively the angle between the piston plane at BDC and the intake or exhaust ports (Fig. 6).

Other parameters can be chosen or can be removed depending on the designer's framework.

The study is focused on the scavenging process. During the phase, the fresh gases enter the cylinder pushing out the burnt gases. During the gas exchanges, backflow and short-cutting should be avoided, which means no burnt gases go through the intake ports (this happens when the pressure in the chamber is higher than the intake pressure) and no fresh gases go through the exhaust ports.

Fig. 6 β angle.

4.2 Fast kriging simulation

From the database, the designer is able to evaluate any quantity in the cylinder for any configuration at any moment during scavenging: the kriging method defines the quantity depending on the spatial location in the cylinder and the value of the chamber variables at any moment. Indeed, the “distance” between the evaluated point and the known points is also linked to the cylinder design in this approach.

Coded in Matlab^©, the kriging method is put into practice with the optimal cylinder found in [13] thanks to genetic algorithms. Kriging is used to determine the distribution of gas in the optimal design to illustrate the method's results.

4.3 Results of kriging fast simulations

The optimal cylinder is defined by the following variable values [13]:

• β_in = 60°;

• β_exh = 5.88°;

• θ_{end_exh_port} = 9.81 crank degrees before BDC;

• θ_{in_advance} = 49.02 crank degrees before BDC;

• θ_{exh_advance} = 76.43 crank degrees before BDC;

• P_boost = 2.16 bar;

• ΔP = 0.29 bar.

Concerning the execution time, kriging interpolation requires around 20 min to map the gas distribution for a given crankshaft angle. Hence, 2.5 days are needed to map the whole scavenging process for one cylinder configuration. To compare, 10 h are needed to simulate the whole process of scavenging with a 2D numerical model.

Prima facie, kriging appears quite long to model the process. Because of the previous model reductions, we would expect less time. Indeed, it is just a way to globally observe the flow in the optimal geometry in order to validate it, so the execution time should be as short as possible. However, the interpolation done with the kriging method is a loop process which requires time. For each point, the kriging algorithm determines the weights comparing spatial coordinates and cylinder design of the data set points, and then approximates the expected value for the new design. Moreover, it is still faster than to get results from CFD simulations or any experiment.

In addition, there is no need to krige the whole process: we can focus on certain significant crankshaft angles (intake opening angle, BDC position, etc.) and, with the kriging method, there is no need to solve the previous time step. Therefore, it only takes 20 min to obtain the gas distribution in the cylinder after closing the ports, for example. With a numerical model, the whole process needs to be simulated to get the same information. Considering a global quantitative validation, only few crankshaft angles are needed to validate (or not) the design and to choose if the optimal design should be fully tested.

Concerning the error, Figure 7 shows the distribution of the fresh gas for the optimal configuration previously defined. Three crankshaft angles were selected:

• 150 crank angles after TDC, ten crankshaft angles after we see the first fresh gases entering the cylinder;

• 180 crank angles after TDC (the piston is at the TDC position);

• 256 crank angles after TDC, all ports are closed.

These 3 moments are the most significative piston positions: opening of the ports, ports fully opened and closing of the ports. To get the 3 pictures of the gas distribution, a few minutes are needed with the kriging method, only these moments are interpolated. With the CFD simulation, it is necessary to run the whole process in order to get the data about the 3 moments, around 10 h of calculation are needed.

To evaluate the error of the kriged interpolation, the CFD model for the optimal configuration was used. Figure 7 illustrates the difference between the kriged (up) and numerical results (down).

Some differences can be observed. One striking difference is the gap between the expected opening of intake ports and the kriged opening. With the fast simulation, the fresh gas enters the cylinder 139 crank degrees after TDC whereas it is supposed to open at 131 crank degrees after TDC. There is a difference of 8 crank degrees which cannot be justified by the blackflow of burnt gases in the intake duck blocking the entry of fresh gases. Indeed, all simulations show that the backflow lasts less than 5° whatever the design. More realistic, fast kriging is an interpolation of values which implies an error in terms of estimation. In addition, the data set used to interpolate the value comes from a reduced model which also implies some reduction errors as illustrated by Figure 5. The data set points should be carefully chosen and used.

The first observation is that the global path of fluids is well represented. In both simulations, the entering fresh gases follow the cylinder wall in a loop path: gases enter on the bottom left side of the cylinder, touch the cylinder head and go out on the bottom right side. This motion is due to the inclination of the intake duct. So, the kriging method well represents the general expected flow path of gases.

In addition, the evolution of the average mass fractions of fresh gases in the cylinder was drawn comparing the CFD and kriging results. Figure 8 shows that the two curves are close except at the beginning. As explained before, there is a gap between the opening of the intake ports and the entry of the first fresh gases with the kriging method. Nevertheless, it does not seem to influence the overall proportion of fresh gases (except for the beginning). At the end of the two simulations, the averaged mass fraction is almost the same with both methods.

The fast kriging simulation provides the flow trend, as we can see in Figure 7 (its path, the mixing of gases, etc.) with a reduced cost compared to CFD simulations or prototype testing. The fast simulation therefore allows a first validation of the optimal configuration found before going further or looking for another configuration.

Fig. 7 Comparison of fresh gases distribution in the cylinder between the kriged results (up) and the CFD results (down).

Fig. 8 Comparison of average mass fraction of fresh gases between CFD and kriging results.

5 Conclusion

The kriging method is a good alternative option to globally reproduce the flow trend. It offers initial results at a very low cost (in terms of time and resources) compared to CFD simulations or experiments. If the data are interpolated at several well selected moments, the fast kriging simulation provides a good observation of optimal design behavior. Comparing CFD and kriging results, the flow path and the average mass fraction of fresh gases in the cylinder are close enough to validate the method.

However, use of the method is mostly restricted by the data set needed for the interpolation. Indeed, the number of known points impacts the time and the precision of the interpolation. Thanks to CFD simulations, a lot of data items were extracted (more than 30 million) which motivates model reduction with the β-NTF algorithm. However, model reduction induces an error which is added to the kriging interpolation error. That is why the kriging technique is a very interesting alternative way to validate a design, but does not replace a CFD simulation or tests on real prototypes.

Nevertheless, the approach presented here is also a new way to apply the powerful kriging method in a design optimization process. Instead of using it to explore the design space and search for optimal solutions, it is used to interpolate the behavior of the system during its use and validate the optimal design early on.

References

All Figures

	Fig. 1 Method of model reduction and optimization.
In the text

	Fig. 2 NTF concept.
In the text

	Fig. 3 TDC and BDC piston positions.
In the text

	Fig. 4 NTF model of engine quantity.
In the text

	Fig. 5 Example of fresh gas distribution from (a) CFD results (b) β-NTF reduced model.
In the text

	Fig. 6 β angle.
In the text

	Fig. 7 Comparison of fresh gases distribution in the cylinder between the kriged results (up) and the CFD results (down).
In the text

	Fig. 8 Comparison of average mass fraction of fresh gases between CFD and kriging results.
In the text