Estimation of Internal Heat Flux on Pulsating Heat Pipes using Kalman Filter: Numerical and Experimental Results

Pulsating Heat Pipes (PHPs) are two-phase passive thermal devices characterized by the presence of significant fluid oscillations inside the tubes that permit to fast and efficiently transfer heat from a hot region to a cold one. They are present in many engineering applications, e.g., electronics, sustainable energies, aerospace. They are very attractive due to their high heat removal capability, flexibility, and low manufacturing cost. Although they have been widely studied their working principles are still not fully understood. To better comprehend their thermal behaviour recent works presented different approaches to assess the internal heat flux in PHPs. However, all the adopted approaches are based on whole-domain techniques that require a higher computational cost compared with sequential estimations. Therefore, to allow heat flux estimation with lower computational effort, in this work a procedure based on the Kalman filter has been adopted. The heat flux has been estimated by solving the inverse heat conduction problem in the PHP’s wall adopting as input data the temperature measurement on the external surface of the pipe acquired by an infrared camera. Firstly, the procedure has been validated adopting synthetic data. Then, experimental data referring to actual operating conditions have been employed to estimate the internal heat flux. The Kalman filter has been adopted as technique to solve the classical instability intrinsically present in inverse problems. It is relatively simple to implement and requires only previous time information, requiring low computational costs. Moreover, the Kalman filter allows the real time estimation of the heat flux: assessing the heat flux on PHPs during operation could be a useful instrument to provide information about on-time working conditions and thus avoid any possible malfunctioning.


Introduction
Pulsating Heat Pipes (PHPs) are devices that transfer heat from a hot region to a cold one by fluid oscillations in a two-phase passive flow.Due to their efficient and high heat transfer rate, geometrical flexibility, and low manufacturing cost, these devices are promising for many engineering areas, such as sustainable energy, electronics, fuel cells and battery, and aerospace applications [1].Usually, PHPs thermal performance are evaluated by computing an overall device thermal resistance [2].However, such approach evaluates only globally the thermal performance and only in pseudo-steady state condition.Many studies try to obtain deeper evaluation and comprehension of the thermal behavior along the tubes to characterize the fluid working conditions [3].In this sense, recent works showed strong relation of the fluid motion inside the tube and the heat flux exchanged between the fluid and the pipe wall.In these works, the internal heat flux was non-intrusively assessed considering the temperature measurements on the external pipe wall and solving the inverse heat conduction problem (IHCP) in the PHP's wall.Cattani et al. [4] used the Tikhonov Regularization to estimate the heat flux on a PHP portion made of Sapphire and filled with ethanol, providing stunning results related to the fluid flow pattern inside the tube.Pagliarini et al. [5] proposed the use of Gaussian filter to evaluate the heat flux exchanged between the fluid and the pipe wall on a multi-turn PHP of aluminum under microgravity conditions.The work considered a transient one-dimensional heat conduction model and showed remarkable results regarding the fluid motion inside the tubes.Recently, Colaço et al. [6] employed the Gradient Conjugated Method with Adjoint Operator to solve the IHCP with the same experimental setup from [4] and [5], and it showed good agreement with the heat flux obtained when compared to the previous works.Recently, some works (see [7,8]) provided a deep insight of PHP mechanism using the Gaussian filter to estimate the heat flux in PHPs.These works presented a whole-domain approach to solve the inverse problem, nevertheless it is attractive to obtain the thermal performance of such devices during operation.Such information could predict the device malfunction and failure, and it is crucial to avoid undesired system crashing due to overheating.For online control purposes, the sequential estimation with Bayesian filters is a suitable choice that correlates information of the desired states and measurements by Bayes' theorem [9].Recursive filtering provides the estimation of the state sequentially, and its advantage is that it does not require storing all the complete dataset from the past.Furthermore, the use of Kalman Filter (KF) provides the optimal estimation when the problem under question is linear-Gaussian [10], that it is the case of the current study.Tuan et al. [11] first used the KF to estimate the heat flux with the Recursive Least-Square Estimator (RLSE).Daouas et al. [12] proposed the augmented vector technique to reduce the two-steps calculation from the previous approach that considers the temperature and heat flux on the target estimated state vector.Due to its simplicity and low computational effort, such approach gained attention and will be used in this work to assess the internal heat flux in PHPs.Therefore, the current work aims at evaluating the capability of the internal heat flux estimation in PHPs with KF for online monitoring purpose.The suggested technique for the internal heat flux estimation solves an IHCP considering available the temperature measurements on the external surface of the pipe wall.The obtained results showed good accordance when compared with the synthetic data and with the experimental data from previous references.The proposed technique is promising to real time performance monitoring since it is non-intrusive and needs lower computational effort compared to the whole-domain techniques since it performs the estimation sequentially.

Forward Problem
The physical problem is represented by a hollow cylinder (Figure1) indicating one of the pipes that constitute the PHP.The domain is a solid wall, defined by Ω, with a length of L and inner and outer radii rint and rext, respectively.The inner surface, Γint, is submitted to a heat flux q(z,t), that represents the heat exchanged between the pipe wall and the working fluid.The external surface, Γext, exchanges heat by natural convection with the environment with a convective heat transfer coefficient of h∞ and environmental temperature of T∞.The other surfaces Γ0 and ΓL are considered adiabatic.The boundary conditions are set according to the PHP operation conditions and were experimentally verified in [4].
The study presumes the thermal properties constant and neglect the temperature variation over the angular direction to model the physical problem, which has been verified experimentally from [4].Therefore, a 2D heat conduction model in radial coordinate rules the physical description as follows: (, ,  = 0) =  ∞ in Ω for  = 0 (6) This work adopts the Finite Difference Method to solve the two-dimensional model.We consider the explicit scheme for time discretization, and we restrict the discretization according to the equation (7) to guarantee the stability of the solution [13].

Inverse Problem
The purpose of solving the inverse problem is to estimate the internal heat flux in PHP devices with temperature measurement available on the external surface.The inverse problem under question is also known as nonstationary inverse problem [10] since the desired quantities in this work vary over the time.Its solution in this work is the state estimation approach that assess the states quantities correlating them with measurements.Firstly, it is defined the discrete-time stochastic evolution and observation models as: where the subscript k and k+1 represent the k th and (k+1) th instant, respectively, f(•) and h(•) are the evolution and observation functions, possibly nonlinear, which describe the evolution of the states and return the measurable states, respectively.The state vector x contains all the information about the states to be estimated, and in this work, it contains the temperature and heat flux nodes, while u is known as input vector.The vector z is the measurement vector, and v and w are mutually independent and represent the process and measurement noises, respectively.When the problem is linear and it considers the Gaussian error distribution, equation ( 8) and ( 9) provide the following system:  +1 =   +  +   (10)  +1 =  +1 +  +1 (11) where F is the evolution matrix; both F and u are constructed with the finite difference method coefficients.H is the observation matrix and returns the measurable states.The evolution and observation models do not vary over time and the subscripts k and k+1 are neglected for matrices F and H and input vector u.The sequential estimation with Kalman Filter consists in two basic steps: prediction and update.The first step consists of evolving the temperature and the heat flux in the nodes from instant k to k+1 with the evolution model, represented in equation (10).When the measurement for instant k+1 is available, represented by equation (11), the predicted vector is corrected.Then, with the corrected state, a new prediction of the state can be performed for the next instant and consequently be updated.Figure 2 presents the flowchart of the described sequential estimation.
+1 = ( −  +1 ) +1 − (16) where P represents the covariance matrix of the state vector, Qk and Rk+1 are the covariance matrices related to wk and vk+1, respectively.The matrix I is the identity matrix and Kk+1 is known as Kalman Gain matrix, which correct the predict states by correlating them with the measurement.The superscript " T " denotes the transpose of the matrix while "-" means the predicted state or covariance matrix.To estimate the temperature and the heat flux simultaneously, the evolution model considers the augmented state vector approach [12], where the estimation for both temperature and heat fluxes adopt them as a state.In this sense, the state vector is constructed by: (17) where T comprises all the temperature nodes originated from the finite difference method discretization and q contains all the heat flux nodes on the inner surface.Since the evolution model of the heat flux is unknown, it is common to describe unknown functions by a black-box random model [12].The equation for such approach follows as: where ω is a Gaussian random number and the parameter σrw is the standard deviation of the heat flux which controls the heat flux search window.High values of σrw admit higher range of the heat flux, but on the other hand it generates higher estimation instabilities.Therefore, this parameter should be analysed for each studied case.

Results and Discussions
To study the estimation accuracy of the here proposed procedure, the working conditions were considered equal to the ones of the experimental investigation performed in Ref. [4].The tube is made of Sapphire with density ρ of 3970 kg/m 3 , specific heat cp of 419 J/kg•K and thermal conductivity λ of 27 W/m•K.The tube length  is 40mm and inner and outer radius are  int = 1 mm and  ext = 2 mm respectively.The environmental temperature  ∞ is 21°C and the convective heat transfer coefficient ℎ ∞ is considered to be 5 W/m 2 •K.Temperature was acquired with a frequency of 50Hz, which means 0.02 s between two measurements, while on the spatial point of view values were taken every 1 mm.

Synthetic data
Firstly, the procedure was validated considering synthetic data.They comes from the solution of the forward problem with COMSOL ® Multiphysics software.The different method for forward and inverse problem alleviates the inverse crime [15].To emulate the heat exchanged between the fluid and the wall, the heat flux function varies over time and space, and it is defined as [4]: where p rules the heat flux frequency variation over time and A is the heat flux amplitude.Then, a random Gaussian noise sequence was added to exact data.Therefore, the definition of the synthetic temperature follows the equation below: where  is a random Gaussian sequence and the value of the standard deviation σmeas is 0.1K.The equation below evaluates the estimation error between the estimated and the exact heat flux: Figures 3 and 4 present the exact and estimated heat flux distributions for the two cases of A equal to 2000 and 5000 W/m 2 respectively.In both the cases the frequency parameter p was considered equal to of 1.2.Figures 3a and 4a represent the numerical distribution of q (Eq.19) that have been used as input data in the solution of the forward problem.Solving forward problem with COMSOL® Multiphysics software it was possible to obtain distributions of temperature that after being spoiled with random noise they have been used as input data for the inverse problem to obtain the distributions of q showed in Figs.3b and 4b.It is possible to notices that the KF presented a small delay when compared to the exact solution, which is intrinsic to the filter [10].In addition, it can be seen that the estimation was less accurate near to the boundaries.Table 1 presents the estimation errors for the Kalman Filter and Tikhonov Regularization for the two considered cases.The parameter for the evolution model of the heat flux from Eq. ( 18) σrw was 20 W/m 2 for A equal to 2000 W/m 2 and 50 W/m 2 for A equal to 5000 W/m 2 .Even if the estimation error shows better estimation with the Tikhonov Regularization technique than with KF, the Tikhonov Regularization considers all the time window of the experiment, and consequently it should provide better results than the sequential estimation with KF.Even so, the KF was able to recover the heat flux with good accuracy and the computation effort was lower.The computational code for the inverse problem with KF was implemented on FORTRAN language with LAPACK package [15].The CPU used for the simulations was an Apple Silicon M1 chip.The computational time, considering the smoothing procedure for the temperature measurement of each instant and the KF calculation, was about 0.12s.
Moreover, in both techniques the estimation error reduced when the amplitude parameter A increases.This can be explained since the difference of temperature between two measurement instants is higher when the value of A increases, and consequently the measurement errors influence less during the estimation.

Experimental data
The experiment was performed according to reference [4].The experiment consisted in a single loop PHP where the evaporator and condenser were made of copper with inner and outer diameters of

Conclusions
The present work evaluated the use of Kalman Filter to estimate the internal heat flux on Pulsating Heat Pipes by solving an inverse heat conduction problem.A high-speed and high-resolution infrared camera acquired the temperature the adiabatic section.A two-dimensional heat conduction model was used for the forward problem.First the study evaluated the heat flux estimation by considering synthetic data and, consequently, experimental data.The sequential estimation is very sensitive to errors and to improve the accuracy of the results the estimation procedure adopted the smoothed temperature with Cubic Smoothing Spline with Cross Validation following the Discrepancy Principle.The results showed good agreement with the exact solution when considering the synthetic data.The estimation errors were higher with the Kalman Filter when compared to the Tikhonov Regularization method due to the different time domain window considered, but on the other hand, lower computational cost was required due to the non-necessity of storing all time domain.Then, the results with Kalman Filter considering experimental data were in good correspondence with the Tikhonov Regularization, whose results were a reference to evaluate the accuracy.The proposed methodology based on Kalman Filter encourages its application to heat flux estimation during Pulsating Heat Pipes operation to control the correct working of such devices and avoid undesired malfunctioning.

Figure 1 .
Figure 1.PHP geometry (left) and the two-dimensional axisymmetric PHP wall domain (right).

Figure 2 .
Figure 2. Flowchart of the sequential estimation with the Kalman Filter (Adapted from [10]).Considering the described calculation procedure, the prediction and update steps comprehend therefore the equations below Prediction: 40th UIT International Heat Transfer Conference (UIT 2023) Journal of Physics: Conference Series 2685 (2024) 012071

Figure 5 .Figure 6 .
Figure 5. Results with experimental data (a) Estimated heat flux with Tikhonov and (b) Estimated heat flux with Kalman filter for the Oscillating flow.