Band Heater Heat Flux Characterization Using Inverse Heat Conduction Problem Models

Peruchi Pacheco da Silva, Ramon; Woodbury, Keith; Samadi, Forooza; Carpenter, Joseph

doi:10.1115/1.4064731

Abstract

An experimental apparatus was constructed to correlate water flowrate and temperature rise under an external band heater. Due to the physical characteristics of the band heater, its transient heating behavior is unknown. This paper investigates the application of inverse heat conduction problem (IHCP) methods to characterize the heat flux from the band heater. Three experiments with different heating times (5, 10, and 20 s) and no flowrate were conducted to measure the transient temperature under the 400 W band heater. Type-T thermocouples measure surface temperature at the centerline of the band heater. The experimental results are computed with five different heat conduction models. The models are chosen to identify how the heat flux response varies from a simplified to a realistic model. Additionally, the results of the experimental heat flux are compared to the manufacturer band heater data (58.9 kW/m²) for each model. The minimum time needed for the heater to fully energize the system is from 10 to 12 s. The residuals for each model are analyzed and used to evaluate the appropriateness of the five different models. The results show that the use of simpler models can achieve results similar to those of complex models, with less time and computational cost.

View large Download slide

1 Introduction

Flow meters stand as pivotal instruments within various industries. Nonintrusive flow meters, offering benefits such as minimal maintenance costs, absence of pressure drop or flow disturbance, not stopping production during installation, and ease of manipulation have garnered significance. Nonetheless, nonintrusive flow meters, such as ultrasound flow meters, often incur high costs. Considering the potential issues with the simplicity and accuracy of common methods, the development of a novel economical nonintrusive flowmeter is undertaken by the authors of this paper. This flowmeter utilizes partial heating of a pipe and measuring the temperature response on its outer surface [1]. Partial heating is one of the methods used in the literature for analyzing heat conduction problems and finding temperature distribution. It is used by Beck et al. [2] for developing the analytical solution for a rectangular geometry at zero temperature at all locations except one boundary that was subject to a step change in the temperature. Beck [3] also uses partial heating on a surface of a flat plate to develop a solution to a three-dimensional heat conduction problem. Solutions are also developed by Woodbury et al. [4] and McMasters et al. [5], using partial heating for two-dimensional transient heat conduction problems on a plate with different boundary conditions. Samadi et al. [1] extended this method to cylindrical geometries by relating the temperature response to the pulse of partial heating on the pipe to the heat transfer coefficient inside the pipe. All these solutions use Green's functions available for one-dimensional (1D) heat conduction problems [6] and adapt them to their objective by adjusting the boundary conditions using the superposition principle. Similar to Laplace [7] and Fourier [8] transforms, Green's functions are powerful tools for solving linear differential equations, such as heat conduction equations, analytically. They build up the solution to any load using the system response to a pulse load [9].

Current step of developing the flowmeter involves a thorough characterization of the band heater to ensure precision in subsequent measurements. Inverse heat conduction problem (IHCP) solution methods [10] are used to estimate transient surface heat flux and can be used as an approach for validation of analytical solutions. These methods use the temperature response as the input and provide the unknown boundary condition, such as heat flux, as an output. IHCPs are ill-posed, meaning that a small measurement noise may result in a big error in the output. This necessitates regularization techniques as an inseparable piece of IHCP solution methods. Therefore, this approach minimizes an objective function defined using the difference between the temperatures estimated using the unknown function and the measured temperatures. Tikhonov regularization technique [11], which is used in this paper, does this job by amplifying the objective function by a penalty term and minimizing it afterward [12].

This paper investigates the effectiveness of different inverse heat conduction models for heat flux estimation during partial heating of a pipe, which is needed for characterizing the band heater. The cases investigated here include: X22B10T0, X23B10T0, R32B01T0, R0C12B1T00, and R32B0(x5)X21B00T0. The first four models are simplified and can be found on the EXACT website [13]. These models are listed here using the solution numbering system [13] and are described in more detail in Sec. 3.

An experimental arrangement has been conducted, comprising a pipe equipped with a band heater that applies heat flux to its outer surface. Temperature sensors are affixed to the pipe beneath the band heater. The temperatures recorded beneath the band heater, in conjunction with the pertinent direct solution, serve as essential inputs for solving the IHCP and, consequently, for estimating the heat flux.

2 Materials and Methods

The test section of the apparatus consists of a steel pipe 10 ft in length. The pipe is a nominal 1.5-in diameter schedule 40, and the assumed dimensions and thermal properties of the pipe are listed in Table 1. The system is filled with water, and Table 1 presents the water properties in this study, which are considered constant at 20 °C. An electric band heater is attached to the outside of the pipe. The manufacturer's data for the heater are listed in Table 2.

Table 1

Dimensions and properties of steel pipe and water

Parameter	Value	Units
d_inside	40.9	mm
d_outside	48.3	mm
Δr, wall thickness	3.70	mm
k, thermal conductivity of the steel pipe	60.5	W/m K
α, thermal diffusivity of the steel pipe	17.7 × 10⁻⁶	m²/s
k, thermal conductivity of the water	0.598	W/m K
C, volumetric heat capacity of the water	4.2	MJ/m³ K

Parameter	Value	Units
d_inside	40.9	mm
d_outside	48.3	mm
Δr, wall thickness	3.70	mm
k, thermal conductivity of the steel pipe	60.5	W/m K
α, thermal diffusivity of the steel pipe	17.7 × 10⁻⁶	m²/s
k, thermal conductivity of the water	0.598	W/m K
C, volumetric heat capacity of the water	4.2	MJ/m³ K

Table 2

Band heater characteristics

Parameter	Value	Units
d_inside	50.8	mm
Width	50.8	mm
Max. temperature	482.2	°C
Power	400	W
Power density	58.9	kW/m²

Three Type-T surface mount thermocouples (uncertainty of 1.0 °C or 0.75% of reading) are installed on the pipe: TC0 and TC1 are 45 deg off the centerline of the band heater as can be seen in Fig. 1, and TC2 is approximately 10 mm away from the centerline (see Fig. 2). The thermocouple data are acquired using MAX31856 Quad Digital Thermocouple Shield. Arduino UNO Rev3 controls the on-and-off cycle for the band heater during each test. The system is insulated as illustrated in Fig. 2 with a Fiberglass pipe insulation 25.4 mm thick (K = 0.311 W/mK).

Fig. 1

View large Download slide

The location of the thermocouples TC0 and TC1 in the pipe heated over (small) partial exterior surface

Fig. 2

View large Download slide

Position of the thermocouples and band heater on the pipe

2.1 Data Reduction.

The band heater was energized for 5, 10, and 20 s in three subsequent trials. During each trial, the steel pipe remained full of water, but the pump in the system remained off. In each trial, temperature data were gathered for 10 s before heating began and continued for a variable period following the heating. The average data obtained from the Arduino for TC0 and TC1 (at the band heater centerline) are displayed in Fig. 3. The data time interval are nominally 0.5 s, but, due to the data acquisition system, the timing is not precise and varies from about 0.57 to 0.68 with an average of 0.62 s. Because the data intervals for the measurements are not uniform, the spline() function of matlab is used to find the measured temperatures at uniform intervals using the mean() time-step.

Fig. 3

View large Download slide

Average temperature data from thermocouples “TC0” and “TC1” at the band heater centerline for three different heating times

2.2 Data Analysis.

The analysis relies on the superposition of solutions for linear problems as in Duhamel's theorem or application of Green's functions. A series of exact solutions will be used to represent the physics and geometry starting from very simple and progressing to highly realistic models. The IHCP solution analysis follows the “whole domain” Tikhonov regularization method [10].

2.2.1 Forward Problem.

The solution to the heat conduction problem for a time-varying q(t) can be constructed by dividing the overall time into small time intervals and using the superposition of many solutions based on a constant heat flux over each time interval. This approximate solution to a direct heat conduction problem with a time-varying heat flux at the surface can be represented as

T = Xq

(1)

where X is the sensitivity matrix defining the dependence of the time vector of temperatures, T, to the corresponding time vector of heat fluxes, q. Note that X depends on the geometry (rectangular, cylindrical), etc. and the boundary conditions of the problem.

The entries of the X matrix follow from superposition concepts and can be interpreted as a discrete form of Duhamel's theorem [10]. The columns of the matrix are all zero above the diagonal and populated with

X_{i, j} = Δ ϕ_{j - i} = {\begin{matrix} 0 & (j - i) < 0 \\ ϕ_{j - i + 1} - ϕ_{j - i} & (j - i) \geq 0 \end{matrix}

(2)

here ϕ_i = ϕ(t_i) are elements of the response function which gives the response of the assumed model at the measurement location for a step change in the unknown surface disturbance (the heat flux).

2.2.2 Measurement Errors.

In an experiment, measured temperatures, Y, are related to the true temperature, T, by

Y = T + ε

(3)

where ε is the vector of measurement errors. Classically, and herein, ε is assumed to be randomly distributed with zero mean and a known constant standard deviation σ_Y.

2.2.3 Inverse Problem.

When temperatures are measured, and heat flux is desired, it is tempting to invert the matrix X in Eq. to compute estimates for q as

\hat{q} = X^{- 1} Y = X^{- 1} (T + ε)

(4)

However, the matrix X generally has small entries and a large condition number, so the entries of X⁻¹ are relatively large and will amplify the measurement errors in Y. This instability in the solution results in “noisy” estimates for the heat flux (see discussion in Sec. 4.3.2 of Ref. [10]).

To address this instability, some form of regularization is needed. In Tikhonov regularization [10,11,14,15], a penalizing term is added to a sum of squares error expression as

S = {(Y - Xq)}^{T} (Y - Xq) + α q^{T} q

(5)

Zeroth-order regularization, α, directly penalizes the function value, and when α is set to a large value, the regularization term dominates the objective function. As a result, minimizing the function will cause the solution for

\hat{q}

to be biased toward zero, resulting in a smoother or more constant function. Other schemes are possible for the penalty term but the one in Eq. is simplest and often used. By minimizing the function S in Eq. with respect to the desired heat flux vector, q, the estimator for

\hat{q}

is found

\hat{q} = {(X^{T} X + α I)}^{- 1} X^{T} Y

(6)

The regularization parameter α in Eq. (6) is relatively small but must be chosen in relation to the amount of noise in the temperature data (see Sec. 4.8 of Ref. [10] for additional information about Tikhonov regularization).

2.2.4 Optimization Regularization.

The regularization parameter, α, must be selected to balance the bias introduced by the extra terms in Eq. against the sensitivity to noise associated with the ill-posed problem. The Morozov discrepancy principle [11,16] is useful and easy to apply with only experimental data at hand. The optimal value of the regularization parameter is that which causes the standard error between the estimated temperatures and the measured temperatures to equal the amount of noise in the measured temperature data

α_{opt} = \arg \min [{(\sqrt{\frac{\sum_{i = 1}^{n} {(\hat{T} (t_{i}; α) - Y_{i})}^{2}}{n}} - σ_{Y})}^{2}]

(7)

here $\hat{T} (t_{i}, α)$ are the temperatures at the measurement location computed from the model (Eq. (1)) using the estimated heat fluxes from Eq. (6).

The amount of noise in the data, σ_Y, can be determined in this experiment by computing the standard deviation of the temperature measurements in the data record during the time before the heater is energized.

Equation (7) is used for each case to find the optimal value for the regularization parameter. Consequently, the RMS error for each of the models applied to the data will be the same and equal to σ_Y.

3 Conduction Models

The entries of the X matrix are based on the temperature response of the assumed model at the measurement location for a step change in the surface heat flux. These five different models are investigated in increasing order of complexity to determine the simplest effective model. Five different models were implemented in matlab and were used to analyze the studied case: X22B10T0, X23B10T0, R32B01T0, R0C12B1T00, and R32B0(x5)X21B00T0. Figure 4 illustrates the boundary conditions for each heat conduction model.

Fig. 4

View large Download slide

Boundary conditions for investigated models

The models are presented from the most simplified to the most complex model. The first two models are the simplest where the problem is assumed to be 1D Cartesian. Although the third and fourth models are cylindrical coordinates, they do not have the most detailed representation of the system. Finally, the last model better represents the physics of the problem, representing a pipe with a heat flux on part of it. For the first four models, the EXACT website documents these solutions [13]. The two-dimensional cylindrical solution is based on Ref. [1]. The numbering notation is detailed and explained in Ref. [6], Chap. 6, and in Ref. [10], Appendix A. Each numbering notation demonstrates the boundary conditions and initial conditions. Sections 3.1–3.5 describe the governing equations, boundary conditions, initial conditions, and temperature solutions in the dimensionless form for each of the models. In all investigated models following dimensionless parameters are used:

\tilde{x} = \frac{x}{L}, \tilde{L} = \frac{L}{R_{1}}, \tilde{a} = \frac{a}{L}, \tilde{r} = \frac{r}{R_{1}}, \tilde{R} = \frac{R_{2}}{R_{1}}, \tilde{t} = \frac{α t}{R_{1}^{2}}, \tilde{T} (\tilde{x}, \tilde{r}, \tilde{t}) = \frac{T}{q_{0} R_{1} / k}

(8)

3.1 X22B10T0 Model.

X22B10T0 is a transient 1D Cartesian model with the following governing equation:

\frac{\partial^{2} \tilde{T}}{\partial {\tilde{x}}^{2}} = \frac{\partial \tilde{T}}{\partial \tilde{t}}

(9)

The boundary conditions on both sides are Neumann, which at the first boundary condition is constant, and the second is adiabatic (homogenous)

- \frac{\partial \tilde{T} (0, \tilde{t})}{\partial \tilde{x}} = 1

(10)

\frac{\partial \tilde{T} (1, \tilde{t})}{\partial \tilde{x}} = 0

(11)

Additionally, the initial condition is zero

\tilde{T} (\tilde{x}, 0) = 0

(12)

The temperature solution for

0 \leq \tilde{t} \leq {\tilde{t}}_{d}^{(1)}

⁠, where

{\tilde{t}}_{d}^{(1)} = {(2 - \tilde{x})}^{2} / (10 A)

is

\tilde{T} (\tilde{x}, \tilde{t}) = \sqrt{4 \tilde{t}} [ierfc (\frac{\tilde{x}}{\sqrt{4 \tilde{t}}})]

(13)

where ierfc(z) is the complementary error function integral. For

{\tilde{t}}_{d}^{(1)} \leq \tilde{t} \leq {\tilde{t}}_{d}^{(2)}

⁠, where

{\tilde{t}}_{d}^{(2)} = {(2 + \tilde{x})}^{2} / (10 A)

the temperature response is

\tilde{T} (\tilde{x}, \tilde{t}) = 2 \sqrt{\tilde{t}} ierfc (\frac{{\tilde{x}}^{2}}{2 \sqrt{\tilde{t}}}) + 2 \sqrt{\tilde{t}} ierfc (\frac{2 - \tilde{x}}{2 \sqrt{\tilde{t}}})

(14)

when

\tilde{t} > {\tilde{t}}_{d}^{(2)}

\begin{matrix} \tilde{T} (\tilde{x}, \tilde{t}) = (\tilde{t} + \frac{{\tilde{x}}^{2}}{2} - \tilde{x} + \frac{1}{3}) \\ - 2 \sum_{m = 1}^{m_{\max}} \frac{cos (β_{m} \tilde{x})}{β_{m}^{2}} exp (- β_{m}^{2} \tilde{t}) \end{matrix}

(15)

where the eigenvalues are β_m = mπ. Considering an accuracy of 10⁻^A, the maximum number of terms is defined by the equation below

m_{\max} = ceil [\frac{1}{π} {(\frac{A \ln 10}{\tilde{t}})}^{1 / 2}]

(16)

For the model X22B10T0, Eq. (14) evaluated at $\tilde{x} = 0$ provides the ϕ(t_i) elements of the response function in Eq. (2).

3.2 X23B10T0 Model.

Similar to the previous model, X23B10T0 is a one-dimensional model with the same governing equation (Eq. ). Additionally, the first boundary condition and the initial conditional are equal to the X22B10T0 model. However, the second boundary condition is a Robin condition

- \frac{\partial \tilde{T} (1, t)}{\partial \tilde{x}} = h_{L} \tilde{T} (1, \tilde{t})

(17)

The temperature response depends on the time. For

0 \leq \tilde{t} \leq {0.12 (2 + \tilde{x})}^{2} / A

⁠, the temperature solution is

\begin{matrix} \tilde{T} (\tilde{x}, \tilde{t}) \approx 2 \sqrt{\tilde{t}} ierfc (\frac{{\tilde{x}}^{2}}{2 \sqrt{\tilde{t}}}) + 2 \sqrt{\tilde{t}} ierfc (\frac{2 - \tilde{x}}{2 \sqrt{\tilde{t}}}) + \\ \frac{2}{B_{2}} [erfc (\frac{2 - \tilde{x}}{2 \sqrt{\tilde{t}}}) - U (2 - \tilde{x}, \tilde{t}, B_{2})] \end{matrix}

(18)

where

\begin{matrix} U (2 - \tilde{x}, \tilde{t}, B_{2}) = exp (- \frac{{(2 - \tilde{x})}^{2}}{4 \tilde{t}} erfcx (\frac{2 - \tilde{x}}{2 \sqrt{\tilde{t}}} + {Bi}_{L} \sqrt{\tilde{t}})) \end{matrix}

(19)

However, if

\tilde{t} > {0.12 (2 + \tilde{x})}^{2} / A

\begin{matrix} \tilde{T} (\tilde{x}, \tilde{t}) = (\frac{1}{{Bi}_{L}} + 1 - \tilde{x}) \\ - 2 \sum_{m = 1}^{m_{\max}} \frac{(β_{m}^{2} + {Bi}_{L}^{2}) cos (β_{m} \tilde{x})}{β_{m}^{2} (β_{m}^{2} + {Bi}_{L}^{2} + {Bi}_{L})} \end{matrix}

(20)

where

tan (β_{m}) = \frac{{Bi}_{L}}{β_{m}}, m = 1, 2, \dots

(21)

For accuracy 10⁻^A, the number of terms needed is

m_{\max} = ceil [\frac{1}{π} {(\frac{A \ln 10}{\tilde{t}})}^{1 / 2}]

(22)

The ϕ(t_i) elements of the response function in Eq. (2) are provided by Eq. (20) evaluated at $\tilde{x} = 0$ for the model X23B10T0.

3.3 R32B01T0 Model.

The governing equation is defined below. Note that this model is in a one-dimensional cylindrical coordinate

\frac{1}{\tilde{r}} \frac{\partial}{\partial \tilde{r}} (\tilde{r} \frac{\partial \tilde{T}}{\partial \tilde{r}}) = \frac{\partial \tilde{T}}{\partial \tilde{t}}, 1 < \tilde{r} < \tilde{R}, \tilde{t} > 0

(23)

The boundary condition for the inner and outer radius are shown in Eqs. and

k \frac{\partial \tilde{T} (\tilde{R}, \tilde{t})}{\partial \tilde{r}} = q_{0}

(24)

- \frac{1}{Bi} \frac{\partial \tilde{T} (1, \tilde{t})}{\partial \tilde{r}} = - \tilde{T} (1, \tilde{t})

(25)

The initial condition is zero for all the radius locations

\tilde{T} (\tilde{r}, 0) = 0

(26)

The dimensionless solutions are given by the following equation:

\begin{matrix} \tilde{T} (\tilde{r}, \tilde{t}) = - \tilde{R} [\ln (\tilde{r}) + {Bi}^{- 1}] \\ - \frac{π^{2}}{2} \tilde{R} \sum_{m = 1}^{m_{\max}} exp (- β_{m}^{2} \tilde{t}) \frac{U_{0}^{2} N_{m} (\tilde{r}) N_{m} (\tilde{R})}{U_{0}^{2} - (β_{m}^{2} + {Bi}^{2}) J_{0}^{2} (β_{m} \tilde{R})} \end{matrix}

(27)

where

N_{m} (\tilde{r}) = J_{0} (β_{m} \tilde{r}) Y_{1} (β_{m} \tilde{R}) - J_{1} (β_{m} \tilde{R}) Y_{0} (β_{m} \tilde{r})

(28)

The eigenvalues are given as

U_{0} Y_{0} (β_{m} \tilde{R}) - W_{0} J_{0} (β_{m} \tilde{R}) = 0

(29)

where

β_{m} = \frac{π (m - \frac{1}{2})}{x}

(30)

where

x = \tilde{R} - 1

⁠. The number of terms, m_max is

m_{\max} = ceil [\frac{1}{2} + \frac{\tilde{R} - 1}{π} {(\frac{A \ln 10}{\tilde{t}})}^{1 / 2}]

(31)

For the model R32B01T0, Eq. (26) evaluated at $\tilde{r} = \tilde{R}$ provides the ϕ(t_i) elements of the response function in Eq. (2).

3.4 R0C12B1T00 Model.

The R0C12B1T00 models 2 layers following the governing equations that are defined below:

\frac{\partial^{2} {\tilde{T}}_{1}}{\partial {\tilde{r}}^{2}} + \frac{1}{\tilde{r}} \frac{\partial {\tilde{T}}_{1}}{\partial \tilde{r}} = \frac{1}{\tilde{α}} \frac{\partial {\tilde{T}}_{1}}{\partial \tilde{t}}, 0 < \tilde{r} < {\tilde{r}}_{1}

(32)

\frac{\partial^{2} {\tilde{T}}_{2}}{\partial {\tilde{r}}^{2}} + \frac{1}{\tilde{r}} \frac{\partial {\tilde{T}}_{2}}{\partial \tilde{r}} = \frac{1}{\tilde{α}} \frac{\partial {\tilde{T}}_{2}}{\partial \tilde{t}}, {\tilde{r}}_{1} < \tilde{r} < 1

(33)

The boundary conditions are

{\tilde{T}}_{1} (0, \tilde{t}) = Finite

(34)

{\tilde{T}}_{1} ({\tilde{r}}_{1}, \tilde{t}) = {\tilde{T}}_{2} ({\tilde{r}}_{1}, \tilde{t})

(35)

\frac{\partial {\tilde{T}}_{1} ({\tilde{r}}_{1}, \tilde{t})}{\partial {\tilde{r}}_{1}} = {\tilde{k}}_{2} \frac{\partial {\tilde{T}}_{2} ({\tilde{r}}_{1}, \tilde{t})}{\partial {\tilde{r}}_{1}}

(36)

\frac{\partial {\tilde{T}}_{2} (1, \tilde{t})}{\partial \tilde{r}} = - \tilde{q}

(37)

The initial conditions are

{\tilde{T}}_{1} (\tilde{r}, 0) = 0

(38)

{\tilde{T}}_{2} (\tilde{r}, 0) = 0

(39)

The temperature solutions for each layer are

{\tilde{T}}_{1} (\tilde{r}, \tilde{t}) = \frac{P_{0} \tilde{t}}{N_{0}} + \sum_{n = 1}^{n_{\max}} \frac{P_{j} R_{1 j} (\tilde{r}) [1 - exp (- λ_{j}^{2} \tilde{t})]}{N_{j}}

(40)

{\tilde{T}}_{2} (\tilde{r}, \tilde{t}) = \frac{P_{0} \tilde{t}}{N_{0}} + \sum_{n = 1}^{n_{\max}} \frac{P_{j} R_{2 j} (\tilde{r}) [1 - exp (- λ_{j}^{2} \tilde{t})]}{N_{j}}

(41)

where the dimensionless functions are described as

R_{1 j} (\tilde{r}) = J_{0} (γ_{j} \tilde{r}), \tilde{r} \leq {\tilde{r}}_{1}

(42)

\begin{matrix} R_{2 j} (\tilde{r}) = A_{2 j} J_{0} (β_{j} \tilde{r}) + B_{2 j} Y_{0} (β_{j} \tilde{r}), {\tilde{r}}_{1} < \tilde{r} \leq {\tilde{r}}_{1} = 1 \end{matrix}

(43)

The coefficients A₂_j and B₂_j

A_{2 j} (\tilde{r}) = \frac{γ_{j} k_{1} Y_{0} (β_{j} {\tilde{r}}_{1}) J_{1} (γ_{j} {\tilde{r}}_{1}) / (β_{j} k_{2}) - J_{1} (β_{j} {\tilde{r}}_{1}) J_{0} (γ_{j} {\tilde{r}}_{1})}{J_{1} (β_{j} {\tilde{r}}_{1}) Y_{0} (β_{j} {\tilde{r}}_{1}) - J_{0} (β_{j} {\tilde{r}}_{1}) Y_{1} (β_{j} {\tilde{r}}_{1})}

(44)

B_{2 j} (\tilde{r}) = - \frac{γ_{j} k_{1} j_{0} (β_{j} {\tilde{r}}_{1}) J_{1} (γ_{j} {\tilde{r}}_{1}) / (β_{j} k_{2}) - J_{1} (β_{j} {\tilde{r}}_{1}) J_{0} (γ_{j} {\tilde{r}}_{1})}{J_{1} (β_{j} {\tilde{r}}_{1}) Y_{0} (β_{j} {\tilde{r}}_{1}) - J_{0} (β_{j} {\tilde{r}}_{1}) Y_{1} (β_{j} {\tilde{r}}_{1})}

(45)

The dimensionless parameters N_j and N₀ are

N_{0} = {\tilde{C}}_{1} π ({\tilde{r}}_{1}^{2} - {\tilde{r}}_{0}^{2}) - {\tilde{C}}_{2} π (1 - {\tilde{r}}_{1}^{2})

(46)

\begin{matrix} N_{j} = {\tilde{C}}_{1} π {{[{\tilde{r}}_{1} J_{0} (γ_{j} {\tilde{r}}_{1})]}^{2} + {[{\tilde{r}}_{1} J_{1} (γ_{j} {\tilde{r}}_{1})]}^{2}} \\ + π {[A_{2 j} J_{0} (β_{j}) + B_{2 j} Y_{0} (β_{j})]}^{2} \\ - π {[A_{2 j} {\tilde{r}}_{1} J_{0} (β_{j} {\tilde{r}}_{1}) + B_{2 j} {\tilde{r}}_{1} Y_{0} (β_{j} {\tilde{r}}_{1})]}^{2} \\ - π {[A_{2 j} {\tilde{r}}_{1} J_{1} (β_{j} {\tilde{r}}_{1}) + B_{2 j} {\tilde{r}}_{1} Y_{1} (β_{j} {\tilde{r}}_{1})]}^{2} \end{matrix}

(47)

where

\tilde{C_{i}}

s are the thermal capacitances in regions i = 1, 2. The dimensionless parameters P_j and P₀ are

P_{0} = - 2 π \tilde{q}

(48)

P_{j} = - 2 π \tilde{q} [A_{2 j} J_{0} (β_{j}) + B_{2 j} Y_{0} (β_{j})]

(49)

The discussion about the eigenvalues can be found in Ref. [13]. Equation (41) evaluated at $\tilde{r} = 1$ provides the ϕ(t_i) elements of the response function in Eq. (2) for the model R0C12B1T00.

3.5 R32B0(x5)X21B00T0 Model.

Samadi et al. [1] find an analytical transient solution for the temperature in a pipe with pulse heating on the outside and convective boundary condition from the inside. This problem can be described as R32B0(x5t5)X21B00T0. The model best represents the physics of the system due to its two-dimensional (r, x) cylindrical coordinate system and includes constant heating on a portion of the outer surface of the pipe for a finite time. The solution to this problem is found by solving the heat conduction problems in radial and axial directions and multiplying the solutions, which is a convenient feature provided by Green's function-based solutions.

Different from the previous models, R32B0(x5)X21B00T0 is a two-dimensional cylindrical coordinate with the following governing equation:

\frac{\partial^{2} \tilde{T}}{\partial {\tilde{x}}^{2}} + \frac{1}{\tilde{r}} \frac{\partial}{\partial \tilde{r}} (\tilde{r} \frac{\partial \tilde{T}}{\partial \tilde{r}}) = \frac{\partial \tilde{T}}{\partial \tilde{t}}, 1 < \tilde{r} < \tilde{R}, \tilde{t} > 0

(50)

The boundary conditions for the horizontal direction are the Neumann and Dirichlet conditions

\frac{\partial \tilde{T} (0, \tilde{r}, \tilde{t})}{\partial \tilde{x}} = 0

(51)

\tilde{T} (1, \tilde{r}, \tilde{t}) = 0

(52)

Additionally, the boundary conditions for the radius direction are

- \frac{1}{Bi} \frac{\partial \tilde{T} (1, \tilde{r}, \tilde{t})}{\partial \tilde{r}} = \tilde{T} (1, \tilde{r}, \tilde{t})

(53)

- k \frac{\partial \tilde{T} (\tilde{x}, 1, \tilde{t})}{\partial \tilde{r}} = {\begin{matrix} 1, 0 \leq \tilde{x} \leq a \\ 0, otherwise \end{matrix}

(54)

where a is half of the heated width, see Fig. 1. The initial condition for every horizontal and radius directions is null

\tilde{T} (\tilde{x}, \tilde{r}, 0) = 0

(55)

According to Ref. [1], the temperature distribution during the heating time (⁠

\tilde{t} < t_{1}

⁠) is

\begin{matrix} \tilde{T} (\tilde{x}, \tilde{r}, \tilde{t}) = π \tilde{R} \sum_{n = 1}^{n_{\max}} \sum_{m = 1}^{m_{\max}} \frac{1}{β_{m} (\frac{β_{m}^{2}}{{\tilde{L}}^{2}} + γ_{n}^{2})} sin (β_{m} \tilde{a}) cos (β_{m} \tilde{x}) \\ \times \frac{γ_{n}^{2} U_{0}^{2}}{U_{0}^{2} - (γ_{n}^{2} + {Bi}^{2}) J_{1}^{2} (γ_{n} \tilde{R})} \\ \times ϑ (\tilde{r}) ϑ (\tilde{R}) [\begin{matrix} exp (- (\frac{β_{m}^{2}}{{\tilde{L}}^{2}} + γ_{n}^{2}) (\tilde{t} + {\tilde{t}}_{1})) \\ - exp (- (\frac{β_{m}^{2}}{{\tilde{L}}^{2}} + γ_{n}^{2}) \tilde{t}) \end{matrix}] \end{matrix}

(56)

when t > t₁, the temperature solution is given by the following equation:

\begin{matrix} \tilde{T} (\tilde{x}, \tilde{r}, \tilde{t}) = π^{2} \tilde{R} {\tilde{L}}^{2} \sum_{n = 1}^{n_{\max}} \sum_{m = 1}^{m_{\max}} \frac{1}{β_{m} (β_{m}^{2} + γ_{n}^{2} {\tilde{L}}^{2})} sin (β_{m} \tilde{a}) cos (β_{m} \tilde{x}) \\ \times \frac{γ_{n}^{2} U_{0}^{2}}{U_{0}^{2} - (γ_{n}^{2} + {Bi}^{2}) J_{1}^{2} (γ_{n} \tilde{R})} \\ \times ϑ (\tilde{r}) ϑ (\tilde{R}) [1 - exp (- (\frac{β_{m}^{2}}{{\tilde{L}}^{2}} + γ_{n}^{2}) \tilde{t})] \end{matrix}

(57)

where the eigenvalues

β_{m} = (m - (1 / 2)) π

⁠, the eigenvalues for the radial coordinate are

U_{0} Y_{1} (γ_{n} \tilde{R}) - W_{0} J_{1} (γ_{n} \tilde{R}) = 0

(58)

where

U_{0} = - γ_{n} J_{1} (γ_{n}) - Bi J_{0} (γ_{n})

(59)

W_{0} = - γ_{n} Y_{1} (γ_{n}) - Bi Y_{0} (γ_{n})

(60)

For the model R32B0(x5)X21B00T0, Eq. (56) evaluated at $\tilde{r} = \tilde{R}$ and $\tilde{x} = 0$ provides the ϕ(t_i) elements of the response function in Eq. (2).

4 Results and Discussion

This section presents the heat flux results corresponding to each of the described models. For ease of reference, the designations X22B10T0, X23B10T0, R32B01T0, R0C12B1T00, and R32B0(x5)X21B00T0 are, respectively, denoted as X22, X23, R32, R0C12, and X21R32.

Figure 5 illustrates the results from the three heating durations, revealing several observations. First, while the nominal power density is anticipated to be approximately 58.9 kW/m², Fig. 5 indicates a maximum of about 20 kW/m² for the Δt_heat = 10 and 20 cases. Second, following the peak, the heat flux gradually diminishes until the power is turned off, suggesting the need for a more sustained, approximately constant value. Third, the heat flux exhibits negative values for a duration after the power cessation. Although it is plausible that, momentarily after power-off, the active surface cools predominantly due to heat loss to the water, it is essential to note that the X22 model disregards this aspect, considering the second boundary as insulated.

Fig. 5

View large Download slide

Results using the X22 model for sensitivity

To improve the modeling, the X23 model is used to compute the sensitivities. In this case, a heat transfer coefficient at the remote boundary dissipates heat to the environment at zero temperature (homogeneous Type 3 condition). The heat transfer coefficient is unknown, and a “reasonable” value should be selected.

Results are shown in Fig. 6 for h = 50 W/m² K. The heat transfer coefficient h is organically selected such that the heat flux after heating stops is non-negative. The peak values are similar to those for the insulated case, but the values remain high while the heater is energized. Also, the negative heat flux after heating stops is eliminated, owing to the selected value of h.

Fig. 6

Results using the X23 model for sensitivity with h = 50 W/m2 K

View large Download slide

Results using the X23 model for sensitivity with h = 50 W/m² K

Figure 7 shows the results for the R32 model to compute the sensitivity components for the matrix X. Note that the value for the heat transfer coefficient is 90 W/m² K, therefore the heat flux at the end of time is approximately zero and not negative. The curve behavior is similar to the X23 and X22, whereas the highest heat flux is lower than the previous models.

Fig. 7

Results using the R32 model for sensitivity with h = 90 W/m2 K

View large Download slide

Results using the R32 model for sensitivity with h = 90 W/m² K

R0C12 model was implemented to find the sensitivity components, and Fig. 8 presents the results for each test. The heat flux peak value is higher than the previous models, reaching approximately 25 kW/m². Additionally, the heat flux response after the band heater is turned off it decreases slower than in the previous models.

Fig. 8

View large Download slide

Results using the R0C12 model for sensitivity

Results of the X21R32 model with h = 50 W/m² K are shown in Fig. 9. The heat flux near the peak increases sharply for every case, and until t_heat = 20 test the heat flux does not appear to have achieved the maximum values, as in the previous models. In addition, the decay is similar to the R0C12 model.

Fig. 9

Results using the X21R32 model for sensitivity with h = 50 W/m2 K

View large Download slide

Results using the X21R32 model for sensitivity with h = 50 W/m² K

The peak heat fluxes for all the models fall below the maximum value of 58.9 kW/m², which assumes all the power of the heater is dissipated into the pipe. However, in reality, some energy is lost to the environment from the back and sides of the heater. The first four models indicate the time required to reach the peak heating rate is about 10–12 s, but the last model suggests that additional heating time is needed to attain a maximum value.

Figure 10 exhibits the residuals for each model. Note that the models have similar residuals. Moreover, the magnitude of the residuals when the band heater is turned on and off increases significantly.

Fig. 10

View large Download slide

Results for t_heat = 20 s for sensitivity

The similarity observed in the residuals is a consequence of the optimization of the regularization parameter. The RMS error for the residuals is identical for all cases, as the Morozov discrepancy principle (as expressed in Eq. (7)) is employed to ascertain the optimal regularization parameter for each model.

Ideally, the residuals should have zero mean value, thus the mean and median values of the residuals can be examined to quantitatively compare these models. As a two-dimensional cylindrical coordinate system better represents the problem, the radial cases are shown in Table 3 showing the mean and median values of the residuals. In most cases, the mean and median values of the residuals for the R32 cases are smaller than the R0C12 cases, and the X21R32 cases have the highest values. Table 3 also provides the regularization parameter. The RMS error, σ_Y, for each case was used to determine the optimal regularization parameter (α_opt) is approximately 0.04 °C.

Table 3

Evaluation parameters for each cylindrical solution

Heating time	Parameter	R32	R0C12	X21R32
5 (s)	α_opt	1.73 × 10⁻²	1.10 × 10⁻²	1.79 × 10⁻²
	μ (°C)	1.10 × 10⁻⁴	5.21 × 10⁻⁴	8.43 × 10⁻⁴
	M (°C)	4.68 × 10⁻³	3.68 × 10⁻³	5.35 × 10⁻³
10 (s)	α_opt	2.33 × 10⁻²	1.30 × 10⁻²	2.48 × 10⁻²
	μ (°C)	2.67 × 10⁻⁴	1.01 × 10⁻³	2.43 × 10⁻³
	M (°C)	2.59 × 10⁻⁴	−8.92 × 10⁻⁵	3.90 × 10⁻³
20 (s)	α_opt	2.93 × 10⁻²	1.54 × 10⁻²	2.94 × 10⁻²
	μ (°C)	1.33 × 10⁻³	2.38 × 10⁻³	5.36 × 10⁻³
	M (°C)	1.64 × 10⁻³	6.48 × 10⁻³	7.63 × 10⁻³

Heating time	Parameter	R32	R0C12	X21R32
5 (s)	α_opt	1.73 × 10⁻²	1.10 × 10⁻²	1.79 × 10⁻²
	μ (°C)	1.10 × 10⁻⁴	5.21 × 10⁻⁴	8.43 × 10⁻⁴
	M (°C)	4.68 × 10⁻³	3.68 × 10⁻³	5.35 × 10⁻³
10 (s)	α_opt	2.33 × 10⁻²	1.30 × 10⁻²	2.48 × 10⁻²
	μ (°C)	2.67 × 10⁻⁴	1.01 × 10⁻³	2.43 × 10⁻³
	M (°C)	2.59 × 10⁻⁴	−8.92 × 10⁻⁵	3.90 × 10⁻³
20 (s)	α_opt	2.93 × 10⁻²	1.54 × 10⁻²	2.94 × 10⁻²
	μ (°C)	1.33 × 10⁻³	2.38 × 10⁻³	5.36 × 10⁻³
	M (°C)	1.64 × 10⁻³	6.48 × 10⁻³	7.63 × 10⁻³

5 Future Work

Each of the models offers advantages and none can be concluded as providing the answer to the heat flux into the pipe. The R32 case has the lowest residuals in temperature but also has the lowest peak heat flux values. Both the R32 and the X21R32 require the specification of internal heat transfer coefficient (h), which is not known. The R0C12 solution provides a reasonable estimate for the heat flux without the need to specify an unknown internal condition. However, all the models provide insight into the general nature of the heat flux (magnitude and time variation) and are in general agreement with each other.

Because of this, additional work is ongoing to characterize the heat delivered from the band heater by incorporating additional temperature sensors on both sides of the band heater surface. Work to develop an external water flow measurement device continues in parallel.

6 Conclusion

Different heat conduction models are used to characterize the heating of a simple band heater. Five increasingly physically realistic models for conduction in the pipe wall are considered and utilized in the solution of the inverse heat conduction problem using Tikhonov regularization to determine transient heating action.

The experimental results are compared to the manufacturer's data, in which the maximum heat flux should be near 58.9 kW/m². The tests show that the peak values are approximately 25 kW/m², less than half of the expected. Another characteristic evaluated is the minimum time required to fully energize the band heater, which ranges from 10 to 12 s. The models are evaluated by analyzing the residuals. Due to the optimization of the regularization parameter, the residuals have similar results. The mean and median values of the residuals are used to compare the models. For the 20 s heating interval R32 model has the lowest mean and median values, followed by R0C12 and then X21R32. The results show that although the X21R32 is the model that better physically represents the system, the other model can achieve similar results.

Acknowledgment

Financial support from U.S. Department of Energy (DOE) under the grant DE-EE0009715 through the Industrial Assessment Center program is gratefully acknowledged.

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

Nomenclature

a =: half of the heated length
A =: desired accuracy in the calculations
Bi =: Biot number
C =: volumetric heat capacity
d =: diameter
h =: heat transfer coefficient
I =: identity matrix
k =: thermal conductivity
L =: length of the pipe
m =: eigenvalue index
M =: median value
n =: eigenvalue index
q =: heat flux
r =: radial coordinate
R =: radius of the pipe
t =: time
T =: temperature
x =: Cartesian coordinate
X =: sensitivity matrix
Y =: measurements of temperature
α =: thermal diffusivity
β =: eigenvalue, horizontal direction
γ =: eigenvalue, radial direction
ε =: measurement error
μ =: mean value
σ_Y =: standard deviation
ϕ =: elements of the response function∼ = dimensionless parameter

References

1.

Samadi

,

F.

,

Woodbury

,

K. A.

, and

Beck

,

J. V.

,

2020

, “

Analytical Solution for Partial Heating on the Exterior of the Pipe With Application to Measuring Fluid Flow Rate

,”

Int. J. Heat Mass Transfer

,

159

, p.

120146

.10.1016/j.ijheatmasstransfer.2020.120146

Google Scholar

Crossref

2.

Beck

,

J. V.

,

Haji-Sheikh

,

A.

,

Amos

,

D. E.

, and

Yen

,

D.

,

2004

, “

Verification Solution for Partial Heating of Rectangular Solids

,”

Int. J. Heat Mass Transfer

,

47

(

19–20

), pp.

4243

–

4255

.10.1016/j.ijheatmasstransfer.2004.04.021

3.

Beck

,

J. V.

,

2011

, “

Transient Three-Dimensional Heat Conduction Problems With Partial Heating

,”

Int. J. Heat Mass Transfer

,

54

(

11–12

), pp.

2479

–

2489

.10.1016/j.ijheatmasstransfer.2011.02.014

4.

Woodbury

,

K. A.

,

Najafi

,

H.

, and

Beck

,

J. V.

,

2017

, “

Exact Analytical Solution for 2D Transient Heat Conduction in a Rectangle With Partial Heating on One Edge

,”

Int. J. Therm. Sci.

,

112

, pp.

252

–

262

.10.1016/j.ijthermalsci.2016.10.014

Google Scholar

Crossref

5.

McMasters

,

R. L.

,

de Monte

,

F.

,

Beck

,

J. V.

, and

Amos

,

D. E.

,

2018

, “

Transient Two-Dimensional Heat Conduction Problem With Partial Heating Near Corners

,”

ASME J. Heat Mass Transfer-Trans. ASME

,

140

(

2

), p. 021301.10.1115/1.4037542

6.

Cole

,

K. D.

,

Beck

,

J. V.

,

Haji-Sheikh

,

A.

, and

Litkouhi

,

B.

,

2010

,

Heat Conduction Using Green's Functions

,

CRC Press

, Boca Raton, FL.

Google Scholar

Crossref

7.

Nezhad

,

Y. R.

,

Asemi

,

K.

, and

Akhlaghi

,

M.

,

2011

, “

Transient Solution of Temperature Field in Functionally Graded Hollow Cylinder With Finite Length Using Multi Layered Approach

,”

Int. J. Mech. Mater. Des.

,

7

(

1

), pp.

71

–

82

.10.1007/s10999-011-9151-9

Google Scholar

Crossref

8.

Hörmander

,

L.

,

2015

,

The Analysis of Linear Partial Differential Operators I: Distribution Theory and Fourier Analysis

,

Springer

, Berlin.

9.

Hartmann

,

F.

,

2012

,

Green's Functions and Finite Elements

,

Springer Science and Business Media

, Berlin.

10.

Woodbury

,

K. A.

,

Najafi

,

H.

,

De Monte

,

F.

, and

Beck

,

J. V.

,

2023

,

Inverse Heat Conduction: Ill-Posed Problems

,

Wiley

, Hoboken, NJ.

Google Scholar

Crossref

11.

Tikhonov

,

A. N.

, and

Arsenin

,

V. I.

,

1977

,

Solutions of Ill-Posed Problems

, V. H. Winston & Sons, Washington, DC; John Wiley & Sons, New York

.

12.

Samadi

,

F.

,

Woodbury

,

K.

, and

Kowsary

,

F.

,

2021

, “

Optimal Combinations of Tikhonov Regularization Orders for IHCPs

,”

Int. J. Therm. Sci.

,

161

, p.

106697

.10.1016/j.ijthermalsci.2020.106697

Google Scholar

Crossref

13.

Cole

,

K. D.

,

2014

, “

ExACT–Exact Analytical Conduction Toolbox

,” University of Nebraska-Lincoln (UNL), Lincoln, NE.

14.

Tikhonov

,

A. N.

,

1943

, “

The Stability of Inverse Problems

,”

Dokl. Akad. Nauk SSSR

,

39

(

5

), pp.

176

–

179

.

15.

Tikhonov

,

A. N.

,

1963

, “

On the Regularization of Ill-Posed Problems

,”

Dokl. Akad. Nauk SSSR

,

153

(

1

), pp.

49

–

52

.https://www.mathnet.ru/eng/dan/v153/i1/p49

16.

Tikhonov

,

A. N.

,

1963

, “

On the Solution of Ill-Posed Problems and the Method of Regularization

,”

Dokl. Akad. Nauk SSSR

,

151

(

3

), pp.

501

–

504

.https://www.mathnet.ru/eng/dan/v151/i3/p501

Band Heater Heat Flux Characterization Using Inverse Heat Conduction Problem Models

Abstract

1 Introduction

2 Materials and Methods

2.1 Data Reduction.

2.2 Data Analysis.

2.2.1 Forward Problem.

2.2.2 Measurement Errors.

2.2.3 Inverse Problem.

2.2.4 Optimization Regularization.

3 Conduction Models

3.1 X22B10T0 Model.

3.2 X23B10T0 Model.

3.3 R32B01T0 Model.

3.4 R0C12B1T00 Model.

3.5 R32B0(x5)X21B00T0 Model.

4 Results and Discussion

5 Future Work

6 Conclusion

Acknowledgment

Data Availability Statement

Nomenclature

References

Get Email Alerts

Cited By

ASME Journals

ASME Conference Proceedings

ASME eBooks

Resources

Opportunities

Band Heater Heat Flux Characterization Using Inverse Heat Conduction Problem Models

Abstract

1 Introduction

2 Materials and Methods

2.1 Data Reduction.

2.2 Data Analysis.

2.2.1 Forward Problem.

2.2.2 Measurement Errors.

2.2.3 Inverse Problem.

2.2.4 Optimization Regularization.

3 Conduction Models

3.1 X22B10T0 Model.

3.2 X23B10T0 Model.

3.3 R32B01T0 Model.

3.4 R0C12B1T00 Model.

3.5 R32B0(x5)X21B00T0 Model.

4 Results and Discussion

5 Future Work

6 Conclusion

Acknowledgment

Data Availability Statement

Nomenclature

References

Get Email Alerts

Cited By

Related Articles

Related Proceedings Papers

Related Chapters

ASME Journals

ASME Conference Proceedings

ASME eBooks

Resources

Opportunities

This Feature Is Available To Subscribers Only