Temperature Distribution Measurement Using the Gaussian Process Regression Method

The temperature distribution in real-world industrial environments is often in a three-dimensional space, and developing a reliable method to predict such volumetric information is beneficial for the combustion diagnosis, the understandings of the complicated physical and chemical mechanisms behind the combustion process, the increase of the system efficiency, and the reduction of the pollutant emission. In accordance with the machine learning theory, in this paper, a new methodology is proposed to predict three-dimensional temperature distribution from the limited number of the scattered measurement data. The proposed prediction method includes two key phases. In the first phase, traditional technologies are employed tomeasure the scattered temperature data in a large-scale three-dimensional area. In the second phase, the Gaussian process regression method, with obvious superiorities, including satisfactory generalization ability, high robustness, and low computational complexity, is developed to predict threedimensional temperature distributions. Numerical simulations and experimental results from a real-world three-dimensional combustion process indicate that the proposed prediction method is effective and robust, holds a good adaptability to cope with complicated, nonlinear, and high-dimensional problems, and can accurately predict three-dimensional temperature distributions under a relatively low sampling ratio. As a result, a practicable and effectivemethod is introduced for three-dimensional temperature distribution.


Introduction
Three-dimensional (3D) temperature distribution plays an important role in combustion diagnosis tasks.The acquisition of the rich and accurate combustion temperature distribution details is of a great significance to the combustion adjustment and control.Currently, two kinds of approaches are available for achieving such task, for example, the numerical simulation technique and the experiment method.Owing to the challenges, such as high computational cost and complexity, the inaccurate properties of initial conditions, boundary conditions, geometrical conditions, and physical property parameters, it is hard for the former to achieve combustion diagnosis tasks in real-world applications.With the development of modern measurement technologies, the latter has attracted more and more attentions.
Conventional flame monitoring systems are only used to judge whether there is a flame in a combustion space, and it is difficult to realize the quantitative measurement of the combustion process parameters.A lot of technologies based on the radiation image and radiation energy signal processing had been developed for 3D temperature distribution measurements [1][2][3][4][5][6].
Optical CCD methods are based on the radiation transfer characteristics of the flue gas in a combustion system, which focus on the mathematical relationship between the radiation accumulation image and 3D radiation energy distributions, and thus the temperature distribution profiles are reconstructed via an appropriate algorithm [7][8][9][10][11].Though this kind of method has been accepted for industrial measurements, the accuracy is relatively low, the measurement systems are complex and costly, and the final result is sensitive to exterior conditions.Particularly, the increase of the spatial scale of the measurement area will lead to the rapid degeneration of the precision.
The acoustic tomography (AT) method employs an appropriate algorithm to realize the reconstruction of the temperature distribution profiles from the given acoustic wave time-of-flight (TOF) data and has gained extensive acceptances in a variety of areas [12][13][14].The effectiveness of the inversion method has a great influence on practical applications of the AT technology.Unfortunately, owing to the illposed nature in the AT inverse problem, seeking for an efficient inversion method with a good numerical stabilizing to ensure a high-precision reconstruction remains a formidable challenge.Particularly, with the augment of the spatial scale of the measurement domain, the above challenges become more serious.
Are there other measurement strategies for 3D temperature distributions?The development of modern measurement instruments makes it possible to measure scattered temperature values.Is it possible to predict 3D temperature distribution from the finite number of the temperature measurement data?If the answer is positive, a natural problem will appear, that is, how to achieve it?With such insight in mind, in this work, we put forward a new methodology to predict 3D temperature distributions from the finite local temperature measurement values, and the main highlights of the study can be summarized as follows.
(1) A new methodology is proposed to predict 3D temperature distribution from the limited number of the scattered temperature measurement data.The proposed method includes two key phases.In the first phase, traditional measurement technologies are employed to measure the local temperature values in a large-scale 3D measurement area.In the second phase, the GPR method is developed to predict 3D temperature distribution details.Numerical simulation results indicate that the execution the GPR method is easy and can accurately predict 3D temperature distributions under a relatively low SA.Furthermore, the GPR method holds satisfactory robustness and good adaptability of coping with complicated, nonlinear, and high-dimensional problems.
(2) A 3D combustion experiment system is constructed, and a series of combustion experiments are conducted.The effectiveness and robustness of the GPR method are further validated by experimental results.As a result, a practicable and efficient method is introduced for 3D temperature distribution measurements.
In accordance with the research target, the rest of this paper is arranged as follows.In Section 2, we propose a new methodology for 3D temperature distribution predictions, and the appealing superiorities are summarized.Section 3 introduces the principle of the GPR method.In Section 4, we specify the prediction procedure, and attractive superiorities are summarized.Section 5 performs numerical simulations and discusses numerical results.In Section 6, a 3D experiment system is constructed; a series of the combustion experiments are carried out to validate the practicability and the effectiveness of the proposed prediction method.Finally, we draw conclusions of the study in Section 7.

Methodology
The development of modern measurement instruments makes it possible to measure scattered temperature values.In this study, the 3D temperature distribution measurement problem will predict unknown temperature distributions in terms of the provided finite observations.
For convenience, 3D temperature distribution predictions from the limited number of the measurement data can be formulized to be a tensor completion (TC) problem.Mathematically speaking, a 3D measurement area can be represented as a tensor X ∈   1 × 2 × 3 by means of the discretization process, in which some entries are missing.The TC task tries to recovery the tensor X with missing data in terms of a given sampling set X ,, , (, , ) ∈ Ω, where Ω is a subset of the complete set.For the convenience of calculation, the solution of the TC problem can be written as the following optimization problem with the introduction of the low rank constraint [15]: where the elements of X Ω and T Ω in the set Ω are given while the remnant elements are missing.
It is noticeable that the TC problem has attracted the increasing attention over the past several years and has been successfully applied to signal and image processing area.A variety of numerical methods have been developed for solving the problem, and the interested readers are referred to [16,17] for more details.But, it is essential to stress that seeking an efficient method to solve the TC problem is a formidable task.Besides, numerous applications indicate that the TC method requires a high SA to ensure a satisfactory solution.Obviously, in practical applications, it is hard to satisfy such requirement due to the limitation of measurement conditions and costs.
We find that if the mapping of the measurement point positions and the corresponding temperature values are abstracted, the temperature distribution at other positions in a 3D measurement area can be predicted.In a mathematical notation, the mapping between the measurement point positions and the corresponding temperature values can be formulized by where   means the spatial coordinates of the th measurement positions;   is the temperature value corresponding to the th measurement position; (⋅) is a mapping describing the correlation between the spatial coordinates and the temperature values.In (2), if the mapping (⋅) is obtained, we can estimate the temperature distribution at other positions in a 3D measurement area.
According to the analysis presented above, in this work, we put forward an alternative methodology to overcome the drawbacks of the TC technique to predict 3D temperature distributions, which can be divided into the following five steps.
Step 1.In accordance with practical demands, a 3D measurement area is appropriately determined.
Step 2. Acquire finite temperature measurement data using one of conventional measurement technologies.
Step 3. The raw measurement data is refined by an appropriate method to ameliorate the data quality.
Step 4. Abstract the mapping between the temperature value and the measurement position according to the finite measured temperature values.
Step 5. Use the mapping abstracted in Step 4 to predict 3D temperature distributions at other locations in a measurement area.

Gaussian Process Regression Method
The above discussions result in a fact that, in order to successfully predict 3D temperature distributions, we must seek an effective method to abstract the mapping between the temperature information and the position information in terms of the finite temperature measurement values.Naturally, in the rest of the paper, we will answer the above problems to achieve the goal of predicting 3D temperature distributions.
A variety of methods are available for the estimation of the mapping function, (⋅), including the multivariate linear regression (MLR) method, the multivariate nonlinear regression (MNR) method, the robust estimation (RE) method [18], the regularized MLR (RMLR) method, the artificial neural network (ANN) technique [19,20], and the GPR method [21][22][23][24][25][26][27][28][29][30].The least squares based MLR (LSMLR) method is the simplest and holds a closed solution, but the estimation results are far from satisfactory.Considering the inaccurate attributes of the measurement data, the RE method is introduced to improve the robustness of the estimation.The Mestimator is one kind of the RE methods, which has been successfully applied to different fields.Popular M-estimation functions include the Huber function, the Cauchy function, and the Fair function.In practical applications, a key task for the M-estimator is to determine the form of the M-estimation function.The RMLR technique is developed to improve the numerical stability of the LSMLR method.For the complicated regression tasks, the results estimated by the MNR often exceed the LSMLR method, but the determination of the form of the regression function is formidable in real-world applications.Owing to the superior fitness capacity, the ANN technique has found wide applications in different areas.In real-world applications, how to select a suitable network structure is full of challenges.Due to the prominent superiorities, for example, satisfactory generalization ability, high robustness, and low computational complexity, the GPR method is successfully introduced to a variety of areas, including image and signal processing, and wide speed prediction, to cope with hard prediction tasks.Inspired by successes of the GPR method, in this study, the method is introduced to abstract the mapping from the position information to the temperature values.
The GPR method is a kernel based learning machine and holds a good adaptability to deal with complicated, nonlinear, and high-dimensional problems.The execution of the GPR method includes two key phases, the training process and the prediction process, and more details can be found in [21][22][23][24][25][26][27][28][29][30].
3.1.Prediction Process.Given a training set  = {(  ,   ) |  = 1, . . ., }, where   ∈   means the spatial coordinate of the measurement positions in a measurement area;   ∈  denotes the temperature corresponding to the measurement positions.In the finite set of the given , as a set of random variables, ( (1) ), ( (2) ), . . ., ( () ), obeys the joint Gaussian distribution.The statistical characters of the GP can be clearly represented via the mean function () and the covariance function (,   ): where  and   mean random variables.
If the observation value, , is perturbed by noises, the model of the GPR problem can be expressed as where  stands for the independent random variables, which obeys the Gaussian distribution: where  2  stands for the variance of the noise.Similarly, the prior distribution of  can be specified as Eventually, we can express the joint prior distribution of  and the prediction  * as follows: where (, ) = (  ) = (  ,   ) represents a symmetric and positive definite covariance matrix with the dimensionality  × ; (,  * ) = ( * , )  stands for a covariance matrix between the training set and the test point  * , and its dimensionality  × 1; ( * ,  * ) defines the covariance of  * .Meanwhile, we can write the joint posterior distribution of the prediction as follows: where  * and cov( * ) define the prediction mean and covariance corresponding to  * , respectively.It is necessary to mention that the covariance matrix is positive definite for the limited data set, which is consistent with the property of the kernel in Mercer's theorem, and thus the covariance function and kernel function are equivalent.Consequentially, we can rewrite (9) as follows: where   = ((, ) +  2    ) −1 .

Training Process.
The selection of various covariance functions is crucial for practical applications of the GPR method.Currently, a variety of covariance functions are available.The squared exponential function has been successfully applied to different areas, and it has following form: where  = {,  2  ,  2  } means the set of the hyperparameters and  2  defines the variance of the kernel function.The optimal hyperparameter  can be solved by the maximum likelihood estimate.The negative log marginal likelihood function () can be expressed as In order to obtain the optimal solution of (13), the partial derivatives with respect to  can be specified as follows: where  =   +  2    and  = (  +  2    ) −1  =  −1 .Finally, we can compute the prediction mean and variance corresponding to a new input  * provided that the optimal hyperparameter is determined.
The execution details of the GPR technique can be divided into the following three steps.
Step 1.According to the temperature measurement values and the measurement position coordinates, a set of training samples, (, ), in which  represents the coordinates of the measurement positions in a 3D measurement area and  stands for the temperature values, is determined.
Step 2. Solve (13) using a suitable method to get the optimal hyperparameters.
Step 3. The mean and variance of other positions can be computed via solving ( 9) and (10).

Prediction Procedure of 3D Temperature Distributions
Owing to the conspicuous superiorities, for example, satisfactory generalization ability, easy numerical execution, low computational complexity and cost, and high robustness, in this study, the GPR method is employed to predict 3D temperature distributions.The prediction details can be divided into the following five steps.
Step 1.A 3D measurement area is determined according to the requirements of practical measurements.
Step 2. Measure the finite temperature data via traditional measurement technologies according to practical measurement conditions.
Step 3. Refine the raw measurement data to improve the data quality.
Step 4. The GPR model is used to abstract the mapping between the temperature values and the position information.
Step 5. Use the trained GPR model to predict 3D temperature distributions of other positions at a predetermined measurement area.

Numerical Simulations
In this section, we use the temperature distribution function to generate the temperature distribution of the whole measurement area firstly.The known temperature data is obtained by random sampling with different sampling ratio, and other temperature data is considered unknown.Then, we perform numerical simulations to make a fair assessment for the GPR method, and the results are compared with the MNR method, the generalized regression neural network (GRNN) method, and the TC methods [17].For computational convenience, all algorithms are implemented by the MATLAB software.The mean relative error (MRE) is used to make a quantitative evaluation for the prediction accuracy of the GPR method, which is defined by where x  and x  denote the true values and the estimated values, respectively.In order to simulate a real-world measurement environment, the simulation data is perturbed with the normal distribution random number with different noise variances (NVs), which can be formulated as where r =  1 ⋅rand ;  1 represents the standard deviation and rand  means a normal distribution random number with the mean of 0 and the standard deviation of 1; the original and noisy measurement data are defined as y  and y  , respectively.We use the SA to make a quantitative evaluation for the competing methods, which is shown in (17), that is, where the number of the samples and the number of the unknown variables are defined as SN and VN, respectively.

Case 1.
In order to fairly assess the effectiveness and robustness, the following temperature distribution is reconstructed: (, , ) = 1200 (0.02 ( 2 +  2 +  2 ) + 1) , where −5 ≤  ≤ 5, −5 ≤  ≤ 5, and −5 ≤  ≤ 5.The temperature distribution is shown in Figure 1.For a fair comparison, the measurement data is randomly sampled.In the GPR method, the kernel function is the squared exponential kernel, and the quasi-Newton algorithm is employed to solve the hyperparameters.In the MNR method, the polynomial basis function is used, the highest order is two, and the coefficients are estimated by means of the least squares method.When the SA and the NV are 9.80% and 25, the temperature distributions predicted by the GPR method, the MNR method, the GRNN algorithm, and the TC algorithm are shown in Figures 2-5.With the SA of 7.84%, Figure 6 illustrates the MREs for the compared prediction techniques under different NVs, for example, 1, 25, 100, 225, and 400.Similarly, when the NV of the measurement data is 25, in Figure 7, we present the MREs for the compared methods under different SAs, that is, 1.83%, 3.93%, 5.88%, 7.84%, and 9.80%.
Figures 2-5 show the results predicted by the GPR method, the MNR technique, the GRNN algorithm, and the TC method, respectively.Numerical simulation results  indicate that the numerical execution of the GPR method is easy and computational complexity and cost are low.Furthermore, unlike inverse heat transfer problems and conventional numerical simulation approaches, the GPR method does not solve complicated governing equations, without initial conditions, boundary conditions, geometrical conditions, physical property parameters, and so forth.Moreover, the GPR method is able to cope with complicated, nonlinear, and high-dimensional problems, thereby conducing to practical measurements.In particular, it can be found from Figures 2-5 that the GPR method can get the higher prediction accuracy, and the temperature distribution forecasted is in a better agreement with the original temperature distribution shown  in Figure 1 in comparison with the other competitors such as the MNR technique, the GRNN algorithm, and the TC method, and the details of the temperature distribution can be clearly revealed.
Under a specific SA level, that is, 7.84%, the MREs for the compared algorithms under varying NVs, that is, 1, 25, 100, 225, and 400, are shown in Figure 6.Evidently, the GPR method is robust, and can cope with the noises contained in the measurement data, and the prediction accuracies under different NVs are superior to other competitors.In Figure 6, under all NVs, the MREs of the GPR method are the smallest among all competitors.In particular, the MRE of  the GPR method is only 0.015%, and the prediction result is encouraging.
The SA is a good criterion for assessing the effectiveness of a prediction method.Under a specific NV level, that is, 25, Figure 7 shows the MREs for the compared prediction techniques under varying SAs, that is, 1.83%, 3.93%, 5.88%, 7.84%, and 9.80%.We can find that the MREs for all compared methods generally decrease with the increase of SAs, and the MREs of the TC method are the largest among all competitors.Furthermore, in Figure 3, when the SA is low, the MREs of the GPR method are the smallest among all compared methods.For example, when the SA is 1.96%, the MRE is only 0.092%.Such precision level can meet the accuracy requirements of practical measurements, which indicates that the GPR method is appropriate for 3D temperature measurement tasks with a low SA.The attribute is highly conducive for practical engineering measurements since the number of the measurement data is always small due to the restrictions from the measurement cost and condition.

Case 2.
In this section, we use another temperature distribution, which is described in (19), to further assess the effectiveness and robustness of the GPR method.
where −5 ≤  ≤ 5, −5 ≤  ≤ 5, and −5 ≤  ≤ 5.The temperature distribution is visualized in Figure 8.The numerical simulation conditions and algorithmic parameters are identical to Section 5.1.When the SA and the NV are 9.80% and 25, the temperature distributions predicted by the compared methods are shown in Figures 9-12.With the SA of 9.8%, Figure 13 illustrates the MREs under different NVs, that is, 1, 25, 100, 225, and 400, for the quantitative assessment.Similarly, when the NV of the measurement data is 25, in Figure 14, we present the MREs under different SAs, that is, 1.96%, 3.92%, 5.88%, 7.84%, and 9.80%, for the quantitative appraisal.
When the SA and the NV are 9.80% and 25, Figures 9-12 show the temperature distributions predicted by the  GPR method, the MNR technique, the GRNN algorithm, and the TC method, respectively.Obviously, the prediction results of the GPR method are better than those of the other competitors and are in a good agreement with the original temperature distribution shown in Figure 8, and the details of the temperature distribution can be well reconstructed.With a specific SA, that is, 9.8%, the MREs for all compared prediction techniques under varying NVs, for example, 1, 25, 100, 225, and 400, are shown in Figure 13.We can observe that the MREs of the GPR method under varying NVs are the smallest among all prediction methods, and the prediction accuracies are encouraging.The property is highly  conducive for practical engineering measurements since the measurement data always contains noises.
When the NV level of the measurement data is 25, Figure 14 shows the MREs for the compared methods under different SAs, that is, 1.96%, 3.92%, 5.88%, 7.84%, and 9.80%.For all methods, we can find that the MREs generally decrease with the increase of the SAs.Furthermore, we can observe in Figure 14 that the MER for the GPR method is still very small even when the SA is very low.For example, under the SA of 1.96%, the MRE for the GPR method is only 0.32%.We thus can draw a conclusion from these results, the GPR method is competent in predicting 3D temperature distributions under a low SA.However, in practical applications, it is difficult to  directly determine the lowest sampling ratio that can ensure a reliable reconstruction result.In fact, the reconstruction accuracy is influenced by many factors, such as the complexity of the temperature distribution and the sampling ratio.In real-world applications, a high sampling ratio will be helpful for ensuring a reliable reconstruction result.When the temperature distribution is simple and uniform, the required sampling ratio is low.However, when the temperature is complicated a high sampling ratio is needed.In the future study, we will further investigate this issue.

Experiment Evaluation
In order to verify the practicability and effectiveness of the proposed prediction method, that is, the GPR method, an experiment study was carried out in this section.The experimental procedure is summarized as the following three steps.
Step 1.The scattered temperature values at different positions are measured by traditional thermocouples.
Step 2. The measurement data is used to train the GPR model to abstract the mapping between the temperature values and the measurement positions.
Step 3. The trained GPR model is used to predict 3D temperature distributions.

Experiment
System.An experiment system is designed and constructed, which contains three main parts, that is, a combustion system, a temperature measurement system, and a data acquisition module.The combustion system consists of a diesel supply system, burners, a gas/air system, a water cooling system, and a control system.The temperature measurement system automatically achieves the temperature measurement.The data acquisition module achieves the data storage and transmission.The layout of the experiment system is shown in Figure 15.

Combustion System.
The test chamber is designed as an adiabatic combustion system.Geometric dimensions of the chamber, whose walls are made of light mullite bricks, are 900 mm () × 900 mm () × 1500 mm ().Nonpremixed burners are used.Diesel is mixed with air, and the ratio of the air/fuel is 15.The flow rate of the fuel (under the operating pressure of 7 bar) is about 0.0011 kg/s; the air velocity is about 4 m/s.The diesel supply system consists of a VSC63A5-2 electromagnetic pump, 1.25 Danfoss atomizers (spray angle is 60 degrees), diesel filters, a high pressure ignition bags, ignition needles, and stable flame discs.The combustion air is provided by a CX-75 5.5 medium-pressure fan, and the air is distributed into the header to support the combustion.Different conditions can be accomplished via changing the number of the burners.6.1.3.Data Acquisition.Four DAM-3038 acquisition modules are employed to collect the measurement data.Each module can simultaneously collect 8 channel electric signals, and the acquisition frequency is 10 Hz, and the electric signal is transferred into the corresponding temperature signal in the module.Through the automatic temperature acquisition system, the temperature under different positions and time instances can be obtained.

Results and Discussions.
In this section, we use the experiment data to evaluate the practicability and effectiveness of the GPR method.The number of the known temperature measurement data is 33.In the GPR method, the squared exponential kernel is used, and the quasi-Newton algorithm is employed to solve the hyperparameters.
Figure 18 shows 3D temperature distribution predicted by the GPR method.For a fair quantitative comparison, Figure 20 presents the temperature values from the thermocouples and the results predicted by the GPR method; Figure 19 is the positions of the measurement points for the quantitative comparison; and Figure 21 illustrates the relative errors.
Figure 18 shows the 3D temperature distribution predicted by the GPR method.We find that the predicted temperature distribution is consistent with experimental observations.As the diesel-liquid is vertically injected into the combustion area, the high temperature region of the combustion is close to the walls.Therefore, when number 1 and number 2 burners are under operation, the temperature distribution predicted by the GPR method is in accordance with the numerical calculation and the laws of the combustion and heat transfer [31][32][33], and the temperature values at the bottom of the test chamber are higher than at the top.
Figure 20 shows the temperature distributions from the GPR method and the thermocouple measurement, it can be seen that the forecasting result is encouraging, the temperature distribution predicted by the GPR method is in a good agreement with the measured temperature values, and the temperature difference is very small.Furthermore, we can observe from Figure 21 that the minimum and maximum relative errors are 5.64 × 10 −4 and 0.021, respectively.We can draw a conclusion from the above results; for example, the GPR method is successful in predicting temperature distributions in a large-scale 3D area.

Conclusions
High-precision temperature distribution information is important for various industrial processes.In this study, a new methodology is put forth to predict 3D temperature  practicable and efficient method is introduced for temperature distribution measurements in a 3D area.

Figure 2 :
Figure 2: Temperature distribution predicted by the GPR method when the SA and the NV are 9.80% and 25.

Figure 3 :
Figure 3: Temperature distribution predicted by the MNR method when the SA and the NV are 9.80% and 25.

Figure 4 :
Figure 4: Temperature distribution predicted by the GRNN method when the SA and the NV are 9.80% and 25.

Figure 5 :
Figure 5: Temperature distribution predicted by the TC method when the SA and the NV are 9.80% and 25.

Figure 7 :
Figure 7: Mean relative errors under different SAs.

Figure 9 :
Figure 9: Temperature distribution predicted by the GPR method when the SA and the NV are 9.80% and 25.

Figure 10 :
Figure 10: Temperature distribution predicted by the MNR method when the SA and the NV are 9.80% and 25.

Figure 11 :
Figure 11: Temperature distribution predicted by the GRNN method when the SA and the NV are 9.80% and 25.

Figure 12 :
Figure 12: Temperature distribution predicted by the TC method when the SA and the NV are 9.80% and 25.

Figure 14 :
Figure 14: Mean relative errors under different SAs.

6. 1 . 2 .
Temperature Measurement System.HT-116 K thermocouples are above the burners, which are arranged around the front and back walls.The temperature can be measured via changing measurement locations.The arrangement of the measurement point is presented in Figures16 and 17 .