- Email: [email protected]
Filtered reciprocity functional approach to estimate internal heat transfer coefficients in 2D cylindrical domains using infrared thermography
International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Contents lists available at ScienceDirect
International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt
Filtered reciprocity functional approach to estimate internal heat transfer coefficients in 2D cylindrical domains using infrared thermography Andrea Mocerino a, Marcelo J. Colaço b, Fabio Bozzoli a,⇑, Sara Rainieri a a b
Department of Industrial Engineering, University of Parma, 181/A Parco Area delle Scienze, I-43124 Parma, Italy Department of Mechanical Engineering, Federal University of Rio de Janeiro, Cid. Universitaria, Cx. Postal 68503, Rio de Janeiro, RJ 21941-972, Brazil
a r t i c l e
i n f o
Article history: Received 22 February 2018 Accepted 18 April 2018
Keywords: Inverse heat conduction problem Reciprocity functional Filtering technique Local heat transfer coefficient
a b s t r a c t Performance estimation of heat exchanger apparatuses requires a proper knowledge of the thermo-fluid dynamic interaction between fluid and device. One of the aspects to take into account is related to the estimation of the internal convective heat transfer coefficient. This requires the knowledge of the internal temperature and the internal wall heat flux, although these quantities are complicated to be measured. Therefore, the use of inverse problem techniques may be very useful for estimating these quantities. This approach deals with the estimation of the local internal properties, given some external temperature measurements, possibly by means of contactless experimental methodologies (i.e. infrared camera imaging). The solution strategy here presented, for a 2D model, is based on the Reciprocity Functional approach (RF), which requires the solution of two auxiliary problems that are solved, in this paper, by means of the Classical Integral Transform Technique (CITT). The RF approach coupled with CITT, presents some advantages over other techniques since its solution is completely analytical and therefore computationally fast, and it does not require any numerical simulation. Moreover, it deals with the noise content of the experimental data, filtering it and overcoming the documented limits of the classical RF approach. This original method, named Filtered Reciprocity Functional (FRF) approach, was applied to different test cases using simulated measurements, and the results were compared to both the exact and the reconstruct solutions obtained by the Truncated Singular Value Decomposition (TSVD) method, showing a good agreement. Finally, the FRF approach was applied to an experimental case in order to test its robustness when dealing with real data. Ó 2018 Elsevier Ltd. All rights reserved.
1. Introduction: Estimating the local heat flux on surfaces is crucial in multiple industrial applications. While measuring temperature is usually an easy task, direct measurement of heat flux is not always straightforward and may not be even possible. For instance, in commercial heat exchangers, for savings in materials and energy use, passive heat transfer enhancement techniques are usually adopted (e.g., treated surfaces, rough surfaces, displaced enhancement devices, swirl-flow devices, surface-tension devices, coiled tubes, or flow additives), [1]. These techniques originate a significant variation of the heat transfer coefficient at the fluid-wall interface along the surfaces and, in
⇑ Corresponding author. E-mail address: [email protected] (F. Bozzoli). https://doi.org/10.1016/j.ijheatmasstransfer.2018.04.089 0017-9310/Ó 2018 Elsevier Ltd. All rights reserved.
some applications in which uniform thermal processes are required, this irregular distribution behaviour impacts negatively on the efficiency of these devices, [2]. For this reason, for optimising these devices, it is fundamental to determine what happens on the internal surfaces in terms of local surface temperature and heat fluxes. The hint for succeeding in the estimation of the local behaviour comes from some research papers in which the space-variable heat transfer coefficient on the inner/outer pipe surface was determined, given measurements of temperature at some points of the pipe wall. Some of these authors [3,4] used a specific experimental setup, based on a thin metallic layer, applied to the tube wall, which was heated by Joule effect. Other authors [5–9] adopted the inverse heat conduction problem (IHCP) solution approach to estimate the local convective heat transfer coefficient starting from the temperature distribution without the need of a special heating system.
1182
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Nomenclature r
H C h qg X n K Y G,F S
w N
c m
k
radial coordinate (m) circumferential coordinate (rad) domain’s boundary convective heat transfer coefficient (W/m2 K) heat source per unit volume (W/m3) cross sectional domain unity normal vector pointing outward the domain thermal conductivity (W/m K) vector of measured temperatures available at the external surface (K) auxiliary functions numbers of external surface temperature measurements generic orthonormal basis function total number of basis trace of the auxiliary function at the internal boundary numbers of harmonics of the eigenfunction number of harmonic of the orthonormal basis
Although the IHCP approach is more promising because it can be applied to every kind of heating configuration and pipe geometry, it presents some complications due to the fact that IHCPs are ill-posed and, consequently, are very sensitive to small perturbations in the input data. Different methods have been developed for solving IHCPs using analytical and numerical techniques. Some of these methods are the conjugate gradient iterative method; Laplace transform method; sequential function specification method; regularization methods, such as Tikhonov regularization; mollification method; truncated singular value decomposition method; and filtering technique method [10,11]. Particularly promising for this kind of application is the Reciprocity Function approach (RF), which needs low computational resources because it avoids the use of iterative techniques. Two important works by Andrieux and Ben Abda [12,13] showed the concept and use of reciprocity functional. From these, other studies based on the reciprocity functional began to emerge in different areas. Delbary et al. [14] developed a qualitative method for breast cancer detection by combining the reciprocity functional method with the linear sampling method. Colaço and Alves [15] estimated spatial variation of the thermal contact conductance by using a reciprocity functional approach with the method of fundamental solutions and non-intrusive temperature measurements. Shifrin and Shushpannikov [16] developed a method for identifying small defects in an anisotropic elastic body based on the reciprocity functional. Other studies, regarding the estimation of the thermal contact conductance through reciprocity functional, using non-intrusive measures, can be found in [15,17]. In recent times, Colaço et al. [18] presented a methodology, based on the reciprocity functional approach, to estimate internal convective heat transfer coefficients in ducts, using only data available at an exterior boundary and the solutions of two auxiliary problems that depends only on the system geometry. By this approach, here referred as Numerically integrated Reciprocity Function (NRF) method, the two auxiliary problems are solved numerically and the unknown function can be estimated solving a linear system, where the solution vector is composed of integrals of the measured boundary data. The procedure was analysed and compared with traditional techniques showing very good results. However, some limits were found: in particular, it was not clear how to choose the regularization term (i.e., the number of
a,b X U; V; R
coefficient of the Fourier’s series expansion of the measurement sensitivity matrix matrices of the singular value decomposition
Subscripts int internal ext external b bulk env environmental avg average Acronyms RF Reciprocity Functional NRF Numerical Integrated Reciprocity Functional FRF Filtered Reciprocity Functional CITT Classical Integral Transform Technique TSVD Truncated Singular Value Decomposition
orthonormal functions used for the two auxiliary problems) by traditional criteria (e.g., discrepancy principle, L-curve, fixed point). The present paper starts from the findings of Colaço et al. [18] and overcomes the limit of the previous work solving the auxiliary problems not by numerical approaches, but instead by the Classical Integral Transform Technique (CITT), a well-known analytical method [19]. In this way, thanks to a fully analytical approach, it is possible to understand how the RF approach deals with the noise content in the experimental data and therefore the regularization factor can be tuned in order to improve the robustness of the approach. The new solution methodology, named Filtered Reciprocity Functional (FRF), for its ability to filter out the unwanted noisy content, is firstly tested using synthetic data, considering different profiles for the internal heat transfer coefficient, and then it is applied to an experimental case with real measurements. 2. Direct problem The physical problem, represented in Fig. 1, consists in a crosssection of a circular duct with internal radius rint and external radius rext. The thermal conductivity of the material, K, is assumed to be constant and uniform such as the environmental temperature Tenv and the bulk temperature Tb of the fluid that flows inside the tube. Both the external and internal surfaces are subjected to Robin’s boundary condition: there is a prescribed convective heat transfer coefficient distribution hint at the inner surface Cint of the pipe, while the exterior surface Cext is subjected to an overall convective heat transfer coefficient hext. The domain X is also subjected to an internal heat generation per unit of volume qg. The equations that describe the physical problem are listed below:
r2 T ¼
qG K
in X
ð1:aÞ
K
@T ¼ hint ðT T b Þ on Cint @nCint
ð1:bÞ
K
@T ¼ hext ðT T env Þ on Cext @nCext
ð1:cÞ
1183
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
on Cext
Gk ¼ wk
ð2:cÞ
where wk is a set of orthogonal functions, which, for a 2D case, can be written as a standard orthonormal Fourier’s basis (3.a)–(3.c):
1 wk ¼ pffiffiffiffiffiffiffipffiffiffiffiffiffiffi r ext 2p
for k ¼ 1
ð3:aÞ
1 k h for k ¼ 2; 4; 6; . . . N 1 wk ¼ pffiffiffiffiffiffiffipffiffiffiffi cos 2 r ext p
ð3:bÞ
1 k1 wk ¼ pffiffiffiffiffiffiffipffiffiffiffi sin h for k ¼ 3; 5; 7; . . . N 2 r ext p
ð3:cÞ
Defining the reciprocity functional and the integral of the internal heat source, respectively, as follows
@T @Gk Y wk dCext @n @n Cext Z henv ðT T env Þ @Gk Y ¼ wk dCext K @n Cext Z
RG;k ¼
Q G;k ¼
Z qg dX Gk K X
ð4Þ
ð5Þ
it is possible to write the following identity:
Fig. 1. Physical model.
@T dCint ¼ RG;k Q G;k Gk @n Cint
Z
3. Inverse problem The direct problem presented in Eqs. (1.a)–(1.c) concerns the determination of the temperature field T within the O domain given the distribution of the convective heat transfer coefficient on both the internal and external surfaces. On the other hand, the inverse problem here described is addressed to restore the internal heat transfer coefficient hint given some extra temperature measurements Y taken on the external surface Cext. The methodology used to solve this problem is based on the reciprocity functional approach [12,13]; this technique, under the assumption of linear inverse problem, is able to restore the internal heat transfer coefficient avoiding the use of iterative procedures. This approach, as it will be shown, estimates the unknown function by solving a linear system described by a diagonal matrix, which speeds up the estimation procedure, that involves the integrals of the external measurements. The determination of the unknown internal heat transfer coefficient is performed by solving two auxiliary problems [18]: the first one concerns the determination of the normal derivative of the temperature @T/@n at the internal surface Cint, which is related to the heat flux per unit of surface; while the second one estimates the wall temperature at the same surface. According to Eq. (1.b), and assuming that Tb and K are known values, the local convective heat transfer coefficient distribution hint can be easily determined on Cint. 3.1. Determining the normal derivative of the temperature at the internal boundary @T/@n As suggested by Colaço et al. [18], the following well-posed auxiliary problem, given by Eqs. (2.a)–(2.c), can be used to obtain @T/@n|Cint:
r2 Gk ¼ 0 in X
ð2:aÞ
@Gk ¼ 0 on Cint @n Cint
ð2:bÞ
ð6Þ
More details about the equation above can be found in [18]. Using the orthonormal basis defined in Eqs. (3.a)–(3.c), and defining ck as the trace of the solution Gk on Cint, the last expression can be manipulated in order to get:
Z
RG;k Q G;k ¼
Z @T @T Gk ck dCint ¼ dCint @n @n Cint Cint
ð7Þ
The normal derivative of the temperature on the internal surface can be approximated by using the induced basis c1, . . ., cN as follows:
@T ¼ a1 c1 þ þ aN cN @nCint
ð8Þ
Therefore, truncating the expansion with N terms, the system of Eq. (7) is reduced to the following linear system N Z X k¼1
Cint
ck aj cj dCint ¼
N Z X k¼1
Cint
ck cj aj dCint
¼ RG Q G for j ¼ 1; 2; 3; . . . ; N
ð9Þ
Eq. (9) can be written in matrix form as:
Ma ¼ RG Q G
ð10Þ
The system is invertible since c1,. . .,cN are linearly independent and, therefore, the aj coefficients can be found as:
a ¼ M1 ðRG Q G Þ
ð11Þ
After that, the estimation of @T/@n|Cint is obtained by simply using the expansion presented in Eq. (8). The auxiliary problem presented in Eqs. (2.a)–(2.c) is then solved analytically by using the Classical Integral Transform Technique [19]. Initially, Eq. (2.a) is expressed in the cylindrical coordinate system:
@ 2 Gk 1 @Gk 1 @ 2 Gk þ ¼0 þ 2 r @r r @h2 @r 2 Then, defining the integral transform as
ð12Þ
1184
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Gk;m Gk ðr; mÞ ¼
Z 2p h¼0
Gk ðr; h0 Þcos½mðh h0 Þdh0
ð13Þ
and the inversion formula as
Gk Gk ðr; hÞ ¼
1 X 1 1 Gk ðr; m ¼ 0Þ þ G ðr; mÞ 2p p k m¼1
ð14Þ
Eq. (12), multiplied by cos[m(h-h’)] and integrated, becomes: 2
d Gk;m 1 dGk;m 1 2 þ þ 2 m Gk;m ¼ 0 r dr r dr2
ð15Þ
With the same procedure, Eq. (2.b), which expresses the boundary conditions imposed on the internal boundary Cint, becomes:
dGk;m ¼0 dr
ð16Þ
and the external boundary condition can be re-written as follows:
Gk;m ¼ wk;m
ð17Þ
where the transformed potential, present in the first auxiliary problem given by Eq. (2.c), is given as:
wk;m ¼
Z 2p h0 ¼0
wk cos ½mðh h0 Þdh0
ð18Þ
Using the above equations, the transformed first auxiliary problem is finally obtained as: 2
d Gk;m 1 dGk;m 1 2 þ þ 2 m Gk;m ¼ 0 in X r dr r dr2 dGk;m ¼ 0 on Cint dr
ð19:aÞ
ð19:bÞ
Cint
on Cext
Gk;m ¼ wk;m
ð19:cÞ
whose solution can be found as:
Gk;m ¼ wk;m cosh
n h
m ln
r io r int int sech m ln r r ext
ð20Þ
In order to obtain the solution of the first auxiliary problem, the above equation has to be inverted using the definition given by Eq. (14):
Gk ¼
1 w 2p k;m¼0
1 n h r io X 1 r int int þ wk;m cosh m ln sech m ln p r r ext m¼1
2
1
6 6 6 6 6 6 6 6 6 rint 6 6 M¼ r ext 6 6 6 6 6 6 6 6 6 4
2
h
sech 1 ln
r int r ext
ð21Þ
In this paper, the measurements Y are decomposed by the Fourier series expansion, in order to take advantage of the harmonic decomposition of the Classical Integral Transform Technique and efficiently approach the inverse problem:
YðhÞ ¼
N a0 X ½am cos ðmhÞ þ bm sin ðmhÞ þ 2 m¼1
where ai and bi are real numbers and m is an integer number. Thanks to this, it is possible to distinguish between the signal and the noise. In almost every practical situation [7], the signal content is mainly present at low frequencies, while noise is present over the whole spectrum. This approach, whose details are reported in Appendix A, enables to conveniently express the terms of Eq. (10) in the following way:
3
2
am¼1 6a 7 6 m¼2 7
7 a¼6 6 .. 7 5 4
ð23Þ
.
am¼m pffiffiffi h 2 enK v ða0 T env Þ
2
3
6 6 n h i o 6 henv 1 int 6 tanh ln rrext rext K 6 6 6 n h i o 6 henv 1 int 6 tanh ln rrext rext K 6 6 6n h i o 6 6 henv tanh 2 ln rint 2 6 K r ext rext pffiffiffiffiffiffiffiffiffiffi6 RG ¼ pr ext 6 6n h i o 6 henv r int 2 6 K tanh 2 ln rext rext 6 6 6 6 .. 6 . 6 6 6n h i o 6 henv rint m 6 K tanh m ln rext r ext 6 6 4n h i o henv m int tanh m ln rrext K r ext 2 q qffiffiffiffiffiffiffi g
6 6 6 QG ¼ 6 6 4
K
p
2r ext
r 2ext r 2int 0 .. .
7 7 7 a1 7 7 7 7 7 b1 7 7 7 7 7 a2 7 7 7 7 7 7 b2 7 7 7 7 7 7 7 7 7 7 am 7 7 7 5 bm
ð24Þ
3 7 7 7 7 7 5
ð25Þ
0
3
0
i 2
h
sech 1 ln
rint r ext
i h i 2 int sech 2 ln rrext
h i 2 int sech 2 ln rrext ..
. sech
0
ð22Þ
2
h
m ln
r int r ext
i sech
2
h
m ln
r int r ext
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 i 5
ð26Þ
1185
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
For this reason, these terms work as operators that amplifies the high frequency components of RG more than the low ones of a given input data set Y. In Fig. 2, the values of M1components are plotted for a representative ratio rint/rext = 0.9 as a function of the harmonic order. For cases where both the signal and the noise have frequency components equally distributed over the spectrum of interest, the M1 matrix does not deteriorate the signal content in comparison to the original data. Unfortunately, in almost all practical situations, the signal has frequency components concentrated in the low frequency range of the spectrum, while the spectral components of noise is expected to be uniformly distributed over the entire frequency domain, [7]. Therefore the M1 operator makes very difficult the extraction of some useful information from the RG vector. For this reason, some regularization technique is needed. In the present work this problem is tackled by conveniently truncating the M1 matrix, considering only the first m harmonics as suggested by Bozzoli et al. in a analogous inverse problem [22]. The number of harmonics to keep is chosen by the classical discrepancy principle, originally formulated by Morozov [21]. More efficient approaches in terms of computational costs are available, but this particular aspect goes beyond the scope of the present paper and it will be investigated in upcoming works.
Fig. 2. Diagonal elements of the M1 matrix.
2
3
p1ffiffi 2
h i 7 6 6 cos ð1hÞsech 1 ln rint 7 7 6 r ext 6 h i 7 7 6 r int 7 6 sin ð1hÞsech 1 ln r 7 6 ext 6 h i 7 7 6 rint 1 6 cos ð2hÞsech 2 ln rext 7 h i 7 c ¼ pffiffiffiffiffiffiffiffiffiffi 6 7 6 rint 7 prext 6 6 sin ð2hÞsech 2 ln rext 7 7 6 7 6 .. 7 6 .h 7 6 i 7 6 6 cos ðmhÞsech m ln rint 7 7 6 r ext 4 h i 5 int sin ðmhÞsech m ln rrext
3.2. Determining the internal wall temperature T|Cint
ð27Þ
To estimate the internal wall temperature, the reciprocity functional approach is used. As suggested by Colaço et al. [18], the following auxiliary problem, given by Eqs. (29.a)–(29.c), can be used to obtain T|Cint:
The values of the aj coefficients, which are the goal of the first auxiliary problem, are easily obtained from a mathematical point of view by Eq. (11). Due to the ill-posed nature of the problem and the presence of the noise in the measurement vector Y {Y1, Y2, . . ., Yn}, the inversion of the matrix M has to be approached by inverse problem solution techniques in order to find a realistic estimation of the unknown heat flux. Since the matrix M is diagonal, its inverse can be easily determined as follows:
2
1
M
6 6 6 6 6 6 6 6 6 6 rext 6 6 ¼ r int 6 6 6 6 6 6 6 6 6 6 4
r2 F k ¼ 0 in X
ð29:aÞ
F k ¼ 0 on Cint
ð29:bÞ
F k ¼ wk
on Cext
ð29:cÞ
This problem is well posed and, for simplicity, we used the same basis defined in the first auxiliary problem, Eqs. (3.a)–(3.c). Following the same steps used in the first auxiliary problem, the second reciprocity functional and the integral of the heat generation terms can be defined as:
1
3
0
1
sech2 1 ln
r int r ext
1
sech2 1 ln
r int rext
1
sech2 2 ln
r int rext
1
sech2 2 ln
r int r ext
..
. sech2
0
The analytical expression achieved clearly highlights the destructive effect of the noise on the reconstruction of the aj coefficient: the values of M1 matrix components, following 1/sech2[m ln(rint/rext)], increase with the number of the harmonic.
1
m ln
r int r ext
sech2
1
m ln
r
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5
ð28Þ
int r ext
@T @F k dCext wk Y @n @n Cext Z henv ðT T env Þ @F k Y ¼ wk dCext K @n Cext Z
RF;k ¼
ð30Þ
1186
Q F;k
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Z qg ¼ Fk dX K X
ð31Þ
Therefore, it is possible to write the following expression:
@F k T dCint ¼ RF;k þ Q F;k @n Cint
Z
ð32Þ
r2 F k;m ¼ 0 in X
ð44:aÞ
F k;m ¼ 0 on Cint
ð44:bÞ
on Cext
F k;m ¼ wk;m
ð44:cÞ
whose solutions are:
r io r int int csch m ln r r ext
Using the orthonormal basis defined in Eqs. (3.a)–(3.c) and defining n as the trace of the solution Fk on Cint, the last expression can be manipulated in order to get:
F k;m ¼ wk;m sinh
RF;k þ Q F;k ¼
In order to obtain the solution of the second auxiliary problem, the above equations need to be inverted using the definition given in Eq. (40):
Z
Z @F k T ðnk T ÞdCint dCint ¼ @n Cint Cint
ð33Þ
The temperature on the internal surface can be approximated by using the induced basis n1, . . ., nN as follows:
TjCint ¼ b1 n1 þ þ bN nN
ð34Þ
Therefore, truncating the expansion with N terms, Eq. (33) is reduced to the following linear system (35): N Z X k¼1
Cint
N X
nk bj nj dCint ¼ k¼1
Z
Cint
nk nj bj dCint
¼ RF;k Q F;k for j ¼ 1; 2; 3; . . . ; N
r ext
þ
h¼0
F k ðr; h0 Þcos½mðh h0 Þdh0
ð38Þ
ð39Þ
and the inversion formula as (40):
F k F k ðr; hÞ ¼
1 X 1 1 F k ðr; m ¼ 0Þ þ F k ðr; mÞ 2p p m¼1
ð40Þ
it is possible to transform Eq. (38) as follows: 2
d F k;m 1 dF k;m 1 2 þ þ 2 m F k;m ¼ 0 r dr r dr 2
ð41Þ
By an analogous elaboration, the transformation of the boundary conditions imposed on the internal boundary Cint, Eq. (29.b) results in:
F k;m ¼ 0
ð42Þ
In the same way, the external boundary condition can be rewritten as:
F k;m ¼ wk;m
p
n h
m ln
r io r int int csch m ln r r ext
ð46Þ
YðhÞ ¼
Defining the integral transform as:
Z 2p
m¼1
wk;m sinh
ð36Þ
After that, the estimate of T|Cint is obtained simply by using the expansion presented in Eq. (34). Analogously to the approach adopted in the previous problem, the second auxiliary problem presented in Eqs. (29.a)–(29.c) is solved analytically by using the Classical Integral Transform Technique. For 2D problems, it is possible to write the Laplacian defined in Eq. (29.a) as:
F k;m F k ðr; mÞ ¼
1 X 1
ð35Þ
ð37Þ
@ 2 F k 1 @F k 1 @ 2 F k þ ¼0 þ r @r r 2 @h2 @r2
ð45Þ
Analogously to the previous problem, in order to take advantage of the harmonic decomposition of the Classical Integral Transform Technique and to obtain a fully analytical expression, it is useful to represent the measurements Y as a Fourier series expansion:
The system is invertible since n1, . . ., nN are linearly independent. Eq. (35) is easy to solve having the bj coefficient as unknown:
b ¼ N1 ðRF Q F Þ
m ln
ln rint 1 r Fk ¼ w 2p k;m¼0 ln rint
which can be written in the following matrix form:
Nb ¼ RF Q F
n h
ð43Þ
where wk;m is the transformation of the potential imposed in the second auxiliary problem, Eq. (29.c), according to Eq. (18). Using the above equations, the transformed second auxiliary problem is obtained as:
N a0 X ½am cos ðmhÞ þ bm sin ðmhÞ þ 2 m¼1
ð47Þ
This approach, whose details are reported in Appendix B, enables to conveniently express the terms of Eq. (36) in the following way:
3 bm¼1 7 6 7 6 6 bm¼2 7 7 6 7 b¼6 6 . 7 6 .. 7 7 6 5 4 bm¼m 2
ð48Þ
2 pffiffiffi henv K 2 1 T a 0 en v r K 6 henv rext ln r int 6 ext 6 6 n h i o 6 henv 1 int 6 coth ln rrext K r ext 6 6 6 n h i o 6 henv 1 int 6 coth ln rrext r ext K 6 6 6 n h i o 6 rint henv 2 coth 2 ln pffiffiffiffiffiffiffiffiffiffi6 r ext K r ext 6 RF ¼ pr ext 6 6 n h i o 6 henv 2 int 6 coth 2 ln rrext K r ext 6 6 6 6 .. 6 . 6 6 6 n h i o 6 henv m int coth m ln rrext 6 K r ext 6 6 4 n h i o henv m int coth m ln rrext K r ext 2 qg
6K 6 6 6 6 QF ¼ 6 6 6 6 4
qffiffiffiffiffiffiffi p
2r ext
r 2ext þ 12 0 .. . 0
ðr2ext r2int Þ ln
r int rext
3 7 7 7 7 7 a1 7 7 7 7 7 b1 7 7 7 7 7 a2 7 7 7 7 7 b2 7 7 7 7 7 7 7 7 7 7 am 7 7 7 5 bm
ð49Þ
3 7 7 7 7 7 7 7 7 7 5
ð50Þ
1187
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
2 6 6 6 6 6 6 6 6 6 6 1 6 6 N¼ rint r ext 6 6 6 6 6 6 6 6 6 6 4
ln2
1r
3
0
int r ext
h i 2 int csch 1ln rrext
h i 2 int csch 1ln rrext
h i 2 int 22 csch 2 ln rrext
h i 2 int 22 csch 2 ln rrext ..
.
h
m2 csch2 m ln
r int rext
i h
m2 csch2 m ln
0
h
r int r ext
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 ii 5 ð51Þ
2
3
1 p1ffiffi 2 ln rint r int r
ext h i 7 7 7 int cos ð1hÞcsch 1 ln rrext 7 h i 7 7 r int sin ð1hÞcsch 1 ln rext 7 7 h i 7 7 r int cos ð2hÞcsch 2 ln rext 7 h i 7 7 7 int sin ð2hÞcsch 2 ln rrext 7 7 7 .. 7 . h i 7 7 int 7 cos ðmhÞcsch m ln rrext 7 h i 5 m sin ðmhÞcsch m ln rint r int rext
6 6 6 1 6 rint 6 6 1 6r 6 int 6 6 2 1 6 rint n ¼ pffiffiffiffiffiffiffiffiffiffi 6 prext 6 6 2 6 rint 6 6 6 6 6 m 6 6 rint 4
ð52Þ
The values of bj coefficients, which are the goal of the second auxiliary problem, are obtained from a mathematical point of view by Eq. (37). Due to the ill-posed nature of the problem and the presence of the noise in the measurement vector Y, the inversion of the matrix N has to be approached by inverse problem solution techniques in order to find a realistic estimation of the unknown heat flux. The diagonal behaviour of the matrix of the inner product N allows us to easily calculating its inverse as follow:
2
1
N
6 6 6 6 6 6 6 6 6 6 6 6 ¼ r int r ext 6 6 6 6 6 6 6 6 6 6 6 6 4
2
ln
r int r ext
Fig. 3. Diagonal elements of the N1 matrix.
3
0 2
csch
1 1 ln
r int rext
1
csch2 1 ln
r int rext
1 r
22 csch2 2 ln
int rext
1 r
22 csch2 2 ln
int rext
..
. 1 r
m2 csch2 m ln
0
int rext
1
m2 csch2 m ln
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 r 5 int r ext
ð53Þ
1188
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
This operator has many points in common with the operator defined in Eq. (28) for the previous auxiliary problem. It means that, also in this case, the presence of noise has a negative impact on the estimation of the b coefficients because the values of N1 matrix components, following 1/m2 csch2[m ln(rint/rext)], increase with the number of the harmonic. In Fig. 3, the values of N1 components are plotted for a representative ratio rint/rext = 0.9 as a function of the harmonic order. Also in this case, the values of the diagonal terms of the N1 matrix increase as the harmonic order increases; for this reason the matrix is regularized by keeping the same numer of harmonics employed in the solution of the first auxiliary problem. 4. Results The solution methodology above presented is firstly tested using synthetic data, considering different profiles for the internal heat transfer coefficient, and then it is applied to an experimental case with real measurements. The result given by the Filtered Reciprocity Functional (FRF) approach are compared to the result given by the Numerically integrated Reciprocity Functional (NRF) [18] and then to the well-known TSVD [20] method in order to qualify the goodness of the new proposed approach. Since the direct problem described by Eq. (1) is linear with respect to the internal heat flux q, in the discrete domain, the following linear system can be written:
T ¼ Xq þ Tq¼0
ð54Þ
where T is the vector of the temperature distribution at the external tube surface Cext given by the internal wall heat flux q at the fluidwall interface, Tq=0 is a constant term and X is the sensitivity matrix. More details about this approach can be found in [8]. Defining a relative temperature as follows:
Trel ¼ T Tq¼0
ð55Þ
from Eq. (54), the internal wall heat flux q can be estimated:
q ¼ X1 Yrel ¼ X1 Y Tq¼0
ð56Þ
The TSVD approach here used, requires the decomposition of the X matrix into three different matrices U, R and V:
X ¼ URVT ¼
n X ui ri v Ti
ð57Þ
i¼1
where U = [u1, . . ., un] and V = [v1, . . ., vn] are matrices with orthogonal columns, and the matrix R has diagonal form:
R ¼ diag ðr1 ; r2 ; . . . ; rn Þ
ð58Þ
The matrix R has non-negative diagonal elements, called the singular values of X, whose magnitude decrease as their order increase. The purpose of the TSVD is to force the solution of the
inverse problem, given by Eq. (56), to have a small norm in order to damp the contributions from the errors present in the measurement contained in the Trel term [20]. This is obtained by the truncation of the solution components that correspond to the smallest singular values, because they bring the largest contribution to the noise amplification in the inversion of the matrix. The truncation parameter tp can be effectively chosen by the well-known Morozov’s discrepancy principle [21]. Therefore, the TSVD of X becomes:
Xt ¼ U Rt V with Rt ¼ diag
r1 ; r2 ; . . . rtp 1 ; 0; . . . ; 0
ð59Þ
where tp is the truncation parameter adopting the Morozov’s discrepancy principle [21]. The internal heat flux is given by the following expression:
1 T T q ¼ VR1 t U Y rel ¼ VRt U Y Tq¼0
ð61Þ
which can be seen also like the solution of the following leastsquares problem:
min q Xt q þ Tq¼0 Y
ð62Þ
Once the heat-flux distribution at the internal wall surface has been determined, the local convective heat-transfer coefficient can be easily determined, as follows:
hint ¼
q TjCint T b
ð62Þ
where T|Cint is the internal wall temperature determined by solving the direct problem once the internal wall heat flux is estimated and Tb is the bulk-fluid temperature on the test section. 4.1. Synthetic data The physical problem considered in the present work, according to the conditions reported in the available literature [5–9], consist in a cross-section of a circular duct with internal radius rint = 0.007 m and external radius rext = 0.008 m. The thermal conductivity of the pipe K = 15 W/mK was assumed to be constant such as the environmental temperature Tenv = 298.15 K and the bulk temperature Tb = 292.15 K of the fluid that flows inside the tube. The exterior surface Cext of the tube was subjected to an overall convective heat transfer coefficient henv = 5W/m2K with the environment. The domain X was assumed to be homogeneous and isotropic, and subjected to an internal heat generation per unit of volume qg = 106 W/m3. The Filtered Reciprocity Functional approach (FRF) was first compared to Numerically integrated Reciprocity Functional (NRF) [18] and then to the well-known TSVD approach [20] within the MatlabÒ environment by adopting synthetic data. The finite element method implemented in Comsol MultiphysicÒ was used in order to generate the synthetic measurements (sampled at the
Table 1 hint profiles used as test functions. Case h1
h2
hint
h
i
W m2 K
p < h < 23 p
0
i
W m2 K
hmin
h
i
W m2 K
hmax
0
600
4000
0
6000
p < h < p
400
0
600
p < h < p
4000
0
6000
23 p < h < 23 p
0
2 3
0
p < h < 23 p
p
23 p < h < 23 p
2 3
0
h4
h
400
600
6000
h3
have
h2
600 p2 þ 600 h2
6000 p2 þ 6000
p
h
i
W m2 K
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
1189
Fig. 4. Reconstruction of the internal heat transfer coefficient for (a) test case h1, (b) test case h2, (c) test case h3 and (d) test case h4 (r = 0.1 K).
external surface of the pipe) by imposing a known distribution of h at the internal wall surface of the tube. In the present analysis, to compare straightforwardly the three considered approaches, the criterion provided by the discrepancy principle [21] was adopted for all of them. Since it is well known that, the effectiveness of the estimated solution of the inverse problem depends on the shape of the function being estimated, different test cases were implemented in order to perform a robust comparison. The shapes of the convective heat transfer coefficients that were used as test functions are reported in Table 1. Test cases h1 and h2 are the classical step function with different amplitudes, chosen in order to test the limits of the estimation approaches to estimate the discontinuities, while cases h3 and h4 are instead examples of extremely smooth functions. Moreover, test cases h1-h3 and h2-h4 were assumed to have the same maximum hmax, minimum hmin and average have in order to make the comparison of different test functions more objective, minimizing the effect of the signal’s amplitude. In order to perform a trustworthy comparison among the techniques, the synthetic measurement temperatures Y were deliberately spoiled by random noise and then used as input data for both the mathematical models. In particular, a Gaussian white noise characterized by a standard deviation ranging from r = 0.01 K to r = 5K was introduced according to:
Y ¼ T exact jCext þ r
ð63Þ
where e is a random Gaussian variable with zero mean and unit variance. As shown in Figs. 4 and 5, which compare the estimated functions to the exact ones for some representative cases, the three approaches show different degrees of accuracy depending on the noise level and the test case considered. In order to perform a quantitative comparison among the three considered approaches, the global relative estimation error between the exact and the estimated internal heat transfer coefficient was calculated as:
E¼
khreconstructed hexact k2 khexact k2
ð64Þ
The added noise depends intrinsically on the random sequence generated by MatlabÒ. Therefore the estimation procedure was repeated 50 times for different random noise sequences and an average value Eavg was calculated for each noise level. As highlighted by Fig. 6, FRF estimates far better than the classical NRF, which shows always the worst behaviour. It has to be highlighted that the stopping criterion used affects the quality of the reconstruction as already shown by Bazàn et al. [24] that compared the fixed-point criteria with the Morozov’s discrepancy principle; moreover, Colaço et al. [15] already showed that the reconstruction due the NFR could be better if a different stopping criterion is used. It means that the harmonic decomposition of the Classical Integral Transform Technique enables to efficiently remove the unwanted information from the measured data. Unex-
1190
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Fig. 5. Reconstruction of the internal heat transfer coefficient for (a) test case h1, (b) test case h2, (c) test case h3 and, (d) test case h4 (r = 1K).
pectedly, FRF works also better than TSVD in almost all conditions, giving this approach interest for future applications. As observed in all inverse problems [23], the presence of the noise in the measurements corrupts the reconstruction of the unknown heat flux distribution. In almost every practical condition, the noise is uniformly distributed on the whole spectrum while the signal is almost always present in the lower harmonic orders; under this statement, the magnification effect of the M1 and N1 matrices, which increases as the harmonics order increase, requires the truncation of their diagonal elements. This kind of truncation can be performed according to the Morozov’s discrepancy principle in the same fashion as it was already adopted for the TSVD approach. This behaviour is clarified by Fig. 7, which shows the reciprocity functional values for both the first and the second auxiliary problem as a function of the harmonic order, for the two representative test case h1 and h3. In order to improve the readability of these plots, the reciprocity functionals are expressed in a more compact polar form, following the classical notation:
RG ¼
rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 R2G þ R2G
ð65:aÞ
RF ¼
rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 R2F k þ R2F k þ1
ð65:bÞ
ðk Þ
ð Þ
ðk
ð
þ1Þ
Þ
where: k⁄ = 2, 4, 6, . . ., N 1.
In all the graphs present in Fig. 7, the reciprocity functional of the noiseless signal is compared to the reciprocity functional of a representative set of white noise, having standard deviation of r = 0.1 K, that were used to spoil the synthetic data. They confirm that, in the reciprocity functionals, the signal prevails for the lower harmonics, while the noise overcomes the signal for the higher harmonics. Lastly, Table 2 reports t that is the number of elements kept in the truncation of M1 and N1 matrix for FRF methodology, for the four test functions and for different noise levels. This number decreases as far as the noise increases; less harmonics are kept for the continuous functions (i.e., case h3 and h4) than for the discontinuous ones (i.e., test case h1 and h2). This behaviour is similar to what was observed for TSVD and NRF in [10] and [18], respectively.
4.2. Experimental data An experimental data set [8] was used to validate the proposed technique. A sketch of the experimental setup is reported in Fig. 8: the solid domain is internally heated by Joule effect and it is thermally insulated, except for a small portion at the end of the tube; in this section the thermographic images of the external wall surface were taken. The experimental problem consists in a stainless steel coiled tube with circular cross-section having an internal radius rint = 0.007 m and an external radius rext = 0.008 m. The thermal conductivity of the pipe K = 15 W/mK is assumed constant as well
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
1191
Fig. 6. Comparison between the NRF, FRF and the TSVD methods for different test case as a function of the dimensionless noise level: (a) h1, (b) h2, (c) h3, and (d) h4.
as the environmental temperature Tenv = 296.8 K and the bulk temperature Tb = 294.7 K of the fluid that flows inside the tube. The exterior surface Cext of the tube was subjected to an overall convective heat transfer coefficient henv = 5 W/m2K with the environment while the domain X, that is assumed to be homogeneous and isotropic, was subjected to an internal heat generation per unit of volume qg = 4.8 106 W/m3. The result given by the FRF approach, here developed and presented in Fig. 9, was compared to the result given by the TSVD and the NRF where the regularization parameter (i.e., the number of bases considered) was determined using Morozov’s Discrepancy Principle. The results, reported in Fig. 9, show that the FRF, NFS and TSVD give the same results confirming the robustness of the new approach here presented. 5. Conclusions The approach here developed deals with the estimation of the local internal heat transfer coefficient in a convective heat transfer problem, given the external temperature measurements, by solving an Inverse Heat Conduction Problem. The solution strategy presented was developed for the 2D model, and was based on the Reciprocity Functional approach (RF), solved by mean of the Classical Integral Transform Technique (CITT). This approach, named Filter Reciprocity Functional (FRF) approach, presents some practical advantages over the more clas-
sical inverse problem solution techniques since it is completely analytical and, for this reason, it is computationally inexpensive. The application of FRF to numerical and experimental benchmarks highlighted its effectiveness and robustness, suggesting the application to other challenging inverse problems. Moreover, the obtained formulation enables, from the theoretical point of view, to highlight the ill-posed nature of the inverse problem and, from the practical point of view to filter efficiently the noisy content thanks to the harmonic nature of CITT and the Fourier series expansion. Conflict of interest The authors declare that there is no conflict of interest. Appendix A. Derivation of the solution of the first auxiliary problem According to the definition given in Eqs. (3.a)–(3.c), the solution has to be divided in at least in three parts for different values of the harmonic frequencies k of the orthonormal basis function: (a) k = 1 According to Eqs. (3.a) and (18), for k = 1, we have wm–0;k¼1 ¼ 0, qffiffiffiffiffi and wm¼0;k¼1 ¼ r2extp . Substituting these results in Eq. (21) the auxil-
1192
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Fig. 7. RG and RF for two representative test case: (a) RG for test case h1, (b) RG for test case h3, (c) RF for test case h1 and (d) RF for test case h3.
Table 2 Number of harmonic components t used in FRF method as a function of the standard deviation of the noise. Case
r = 0.01[K]
r = 0.05[K]
r = 0.1[K]
r = 0.5[K]
r = 1[K]
r = 5[K]
h1 h2 h3 h4
35 34 24 22
30 28 18 17
23 21 15 13
16 15 9 8
10 8 6 5
8 5 4 3
1 ffi. Therefore, its normal derivaiary function becomes Gk¼1 ¼ pffiffiffiffiffiffiffiffiffi 2pr ext k¼1 ¼ 0. tive, evaluated on the external boundary, becomes @G@n
Cext
The trace of the auxiliary function at the internal boundary Cint 1 ffi. Thus, the reciprocity functional is then obtained as ck¼1 ¼ pffiffiffiffiffiffiffiffiffi 2pr ext
defined in Eq. (4) becomes:
1 henv RG;k¼1 ¼ r ext pffiffiffiffiffiffiffiffiffiffiffiffiffi 2pr ext K
Fig. 8. Experimental set-up.
Z 2p h0 ¼0
Ydh0
Z 2p h0 ¼0
T env dh0
ðA1Þ
R 2p According to Eq. (22), we can write h0 ¼0 Ydh0 ¼ 2pa0 being a0 the first coefficient of the Fourier expansion. The environmental temperature, Tenv, is supposed to be constant, thus its integral R 2p could be re-written as h0 ¼0 T env dh0 ¼ 2pT env . According to these results, the reciprocity functional becomes:
1193
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
The trace of the auxiliary function at the internal boundary Cint is then:
ck¼even ¼
sffiffiffiffiffiffiffiffiffiffiffi 1 1 X
r int cos ðmhÞsech m ln prext m¼1 rext
ðA7Þ
Finally, the reciprocity functional, defined in Eq. (4), becomes:
1 henv rint r ext pffiffiffiffiffiffiffipffiffiffiffi sech m ln K r p r ext ext m¼1
Z 2p rint m Ycosðmh0 Þdh0 : sinh m ln rext r ext h0 ¼0
RG;k¼even ¼
Fig. 9. Convective heat transfer coefficient estimated with the TSVD, NRF and FRF approaches.
pffiffiffiffiffiffiffiffiffiffiffiffiffi henv ða0 T env Þ RG;k¼1 ¼ 2pr ext K
ðA2Þ
In agreement with Eq. (5), the integral of the internal heat generation is given as:
Q G;k¼1
qg ¼ K
p r2ext r2int pffiffiffiffiffiffiffiffiffiffiffiffiffi 2prext
hc1 ; cj¼1 i ¼
r int rext
hc1 ; cj–1 i ¼ 0
RG;k¼even ¼
D
m¼1
am
pffiffiffiffiffiffiffiffiffiffiffi henv rint m prext tanh m ln K rext r ext
ðA9Þ
E
ck¼even ; cj¼k ¼
r int 2 sech r ext
m ln
rint rext
E
ck¼even ; cj–k ¼ 0
ðA10:1Þ ðA10:2Þ
ðA4:1Þ
(d) k = 3, 5, 7, . . ., N and k = 2m + 1 For these values of k, according to Eqs. (3.c) and (18), we can qffiffiffiffiffiffi p sin ðmhÞ, and w obtain wk;m–0 ¼ rext k;m¼0 ¼ 0.
ðA4:2Þ
Substituting these results in Eq. (21) the auxiliary function becomes:
(b) k = 2, 4, 6, . . ., N 1 and k – 2m; k = 3, 5, 7, . . ., N and k – 2m +1 It is easy to show that, according to Eqs. (3.b) and (18), both wk;m–0 and wk;m¼0 are null as well as the auxiliary function Gk. Therefore, these frequencies of the orthonormal basis function do not participate in the solution. (c) k = 2, 4, 6, . . ., N 1 and k = 2m For these values of k, according to Eqs. (3.b) and (18), we can qffiffiffiffiffiffi p cos ðmhÞ, and w obtain wk;m–0 ¼ rext k;m¼0 ¼ 0.
Gk¼odd
Therefore, its normal derivative, evaluated at the external boundary, becomes:
sffiffiffiffiffiffiffiffiffiffiffi
1 @Gk¼odd 1 X r ¼ sin ð m h Þsech m ln int @n Cext prext m¼1 r ext
r int m sinh m ln r ext rext
Gk¼even
ck¼odd
ðA5Þ Therefore, its normal derivative, evaluated at the external boundary, becomes:
ðA6Þ
ðA12Þ
The trace of the auxiliary function at the internal boundary Cint is:
sffiffiffiffiffiffiffiffiffiffiffi
1 n h r io 1 X r int int sech m ln ¼ cos ðmhÞ cosh m ln prext m¼1 r r ext
sffiffiffiffiffiffiffiffiffiffiffi
1 n h r io 1 X r int int ¼ sin ½mh cosh m ln sech m ln prext m¼1 r r ext ðA11Þ
Substituting these results in Eq. (19), the auxiliary function becomes:
sffiffiffiffiffiffiffiffiffiffiffi
1 @Gk¼even 1 X r ¼ cos ð m h Þsech m ln int @n Cext prext m¼1 r ext
r int m sinh m ln r ext rext
1 X
agreement with Eq. (5), the integral of the internal heat generation is null, while according to the definition given by Eq. (9) the elements of the matrix even rows that contain the inner product between the trace of the solutions become:
ðA3Þ
while the elements of the first row of the matrix given by Eq. (9), which contains the inner product between the trace of the solutions, becomes:
ðA8Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement at the external boundary Cext can be R 2p written as: h0 ¼0 Ycosðmh0 Þdh0 ¼ pam being am the coefficient of the Fourier series expansion referring to the m-th harmonic. According to this result, the reciprocity functional becomes:
D
1 X
sffiffiffiffiffiffiffiffiffiffiffi
1 1 X r int ¼ sin ðmhÞsech m ln prext m¼1 r ext
ðA13Þ
The reciprocity functional defined in Eq. (4) becomes:
1 henv rint r ext pffiffiffiffiffiffiffipffiffiffiffi sech m ln K r p r ext ext m¼1
Z 2p rint m Ysinðmh0 Þdh0 sinh m ln r ext r ext h0 ¼0
RG;k¼odd ¼
1 X
ðA14Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement at the external boundary Cext can be
1194
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
R 2p written as: h0 ¼0 Ysinðmh0 Þdh0 ¼ pbm being bm the coefficient of the Fourier series expansion referring to the m-th harmonic. According to this result, the reciprocity functional becomes:
RG;k¼odd
1 X pffiffiffiffiffiffiffiffiffiffiffi henv rint m ¼ bm pr ext tanh m ln K rext r ext m¼1
ðA15Þ
In agreement with Eq. (5), the integral of the internal heat generation is null, while according to the definition given by Eq. (9), the elements of the matrix odd rows that contain the inner product between the trace of the solutions becomes:
n1 ; nj–1 ¼ 0
ðB4:2Þ
(b) k = 2, 4, 6, . . ., N 1 and k – 2m; k = 3, 5, 7, . . ., N and k – 2m +1 It is easy to show that, according to Eqs. (3.b), (3.c) and (18), both wk;m–0 and wk;m¼0 are null as well as the auxiliary function Fk. Therefore, these frequencies of the orthonormal basis function do not participate in the solution. (c) k = 2, 4, 6, . . ., N 1 and k = 2m
D E r r ck¼odd ; cj¼k ¼ int sech2 m ln int rext r ext
ðA16:1Þ
For these values of k, according to Eqs. (3.b) and (18), we can qffiffiffiffiffiffi p cos ðmhÞw 1ffipffiffiffi cos ½mh, and wk;m¼0 ¼ 0. obtain wk;m–0 ¼ rext k;m–0 ¼ pffiffiffiffiffi rext :p
hck¼odd ; cj–k i ¼ 0
ðA16:2Þ
Substituting these results in Eq. (46), the auxiliary function becomes:
Appendix B. Derivation of the solution of the second auxiliary problem According to the definition given in Eqs. (3.a)–(3.c), the solution has to be divided in at least three steps for different values of harmonic frequencies k of the orthonormal basis function: (a) k = 1 According to Eqs. (3.a) and (18), for k = 1, we have wm–0;k¼1 ¼ 0, qffiffiffiffiffi and wm¼0;k¼1 ¼ r2extp . Substituting these results in Eq. (46) the auxil1 ffi iary function becomes F k¼1 ¼ pffiffiffiffiffiffiffiffiffi
ln ð
2pr ext ln
Þ . Therefore, its normal r int r r int r ext
derivative, evaluated on the external boundary, becomes 1 ffi r 1 . The trace of the auxiliary function at ¼ pffiffiffiffiffiffiffiffiffi int
@F k¼1 @n
2pr ext ln
Cext
the internal 1 ffi nk¼1 ¼ pffiffiffiffiffiffiffiffiffi
2pr ext ln
rext
r ext
Cint is then obtained as boundary r 1 . Thus, the reciprocity functional defined int r ext
r int
F k¼even ¼
RF;k¼1
ðB2Þ
r ext
while Eq. (31) can be written as:
Q F;k¼a
2 3 rffiffiffiffiffiffiffiffiffiffi qg p 4 2 1 r2ext r 2int 5 ¼ r þ K 2rext ext 2 ln rint r ext
ðB3Þ
n1 ; nj¼1 ¼
2
ln
1
rint rext
rint rext
m ln
r io rint int csch m ln r rext
ðB4:1Þ
ðB5Þ
ðB6Þ
The trace of the auxiliary function at the internal boundary Cint is:
nk¼even
sffiffiffiffiffiffiffiffiffiffiffi
1 1 X r int m ¼ cos ðmhÞcsch m ln prext m¼1 r ext r int
ðB7Þ
Finally, the reciprocity functional defined in Eq. (30) becomes:
sffiffiffiffiffiffiffiffiffiffiffi
1 henv r int csch m ln p r K r ext ext m¼1
Z 2p r int m Ycosðmh0 Þdh0 cosh m ln r ext rext h0 ¼0
RF;k¼even ¼
1 X
r ext
ðB8Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement on the external boundary Cext can R 2p be written as h0 ¼0 Ycosðmh0 Þdh0 ¼ pam being am the coefficient of the Fourier’s series expansion referring to the m-th harmonic. According to this result, the reciprocity functional becomes:
RF ¼
1 X
m¼1
am
pffiffiffiffiffiffiffiffiffiffiffi henv r int m prext coth m ln K r ext rext
ðB9Þ
In agreement with Eq. (31), the integral of the internal heat generation is null, while according to the definition given by Eq. (35), the elements of the matrix even rows that contain the inner product between the trace of the solutions becomes:
hnk¼even ; nj¼k i ¼
m2 r ext r int
2
csch
m ln
rint r ext
hnnk¼even ; nj–k i ¼ 0
According to Eq. (35), the first raw of the matrix of the inner product becomes:
n h
sffiffiffiffiffiffiffiffiffiffiffi
1 @F k¼even 1 X r ¼ cos ð m h Þcsch m ln int @n Cext prext m¼1 r ext
r int m cosh m ln r ext rext
ðB1Þ
20 1 3 pffiffiffiffiffiffiffiffiffiffiffiffiffi henv K 4@1 Aa0 T env 5 ¼ 2pr ext K ln rint r ext henv
p rext
Therefore, its normal derivative, evaluated at the external boundary, becomes:
r ext
Once again, in order to make the above expression fully analytical, we propose, in this paper, to represent the measurements, Y, as a Fourier series expansion. Therefore, can be write R 2p Ydh0 ¼ 2pa0 being a0 the first coefficient of the Fourier’s series h0 ¼0 expansion. As previously observed, the environmental temperature is supposed to be constant, and its integral can be written R 2p as h0 ¼0 T env dh0 ¼ 2pT env . According to these results, the reciprocity functional becomes:
m¼1
cos ðmhÞ sinh
in Eq. (30) becomes:
82 9 3 Z 2p Z = rext <4henv 1 henv 2p 0 0 5 RF;k¼1 ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffi Ydh T env dh r 0 0 : ; K K int 2prext h ¼0 ln r ext h ¼0
rffiffiffiffiffiffiffi 1 X 1 p
ðB10:1Þ ðB10:2Þ
(d) k = 3, 5, 7, . . ., N and k = 2m + 1 For these values of k, according to Eqs. (3.c) and (18), we can qffiffiffiffiffiffi p sin ðmhÞ, and w obtain wk;m–0 ¼ rext k;m¼0 ¼ 0 .Substituting those results in Eq. (46) the auxiliary function becomes:
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
F k¼odd
sffiffiffiffiffiffiffiffiffiffiffi
1 n h r io 1 X r int int ¼ sin ðmhÞ sinh m ln csch m ln prext m¼1 r r ext ðB11Þ
Therefore, its normal derivative, evaluated at the external boundary, becomes: sffiffiffiffiffiffiffiffiffiffiffi
1 @F k¼odd 1 X r int r int m ¼ sin ðmhÞcsch m ln cosh m ln @n Cext prext m¼1 r ext r ext r ext ðB12Þ
The trace of the auxiliary function at the internal boundary Cint is:
nk¼odd
sffiffiffiffiffiffiffiffiffiffiffi
1 1 X rint m ¼ sin ðmhÞcsch m ln prext m¼1 rext r int
ðB13Þ
The reciprocity functional defined in Eq. (30) becomes:
RF;k¼odd
sffiffiffiffiffiffiffiffiffiffiffi
1 henv r int ¼ rext csch m ln prext K r ext m¼1
Z 2p r int m cosh m ln Ysinðmh0 Þdh0 r ext r ext h0 ¼0 1 X
ðB14Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement at the external boundary Cext in the R 2p last expression can be written as h0 ¼0 Y sin ðmh0 Þdh0 ¼ pbm being bm the coefficient of the Fourier’s series expansion referring to the mth harmonic. According to this result, the reciprocity functional becomes:
RF;k¼odd ¼
1 X
m¼1
bm
pffiffiffiffiffiffiffiffiffiffiffi henv rint m prext coth m ln K rext r ext
ðB15Þ
In agreement with Eq. (31), the integral of the internal heat generation is null, while according to the definition given by Eq. (35), the elements of the matrix odd rows that contain the inner product between the trace of the solutions becomes:
hnk¼odd ; nj¼k i ¼
m r int r ext
hnk¼odd ; nj–k i ¼ 0
csch
2
r m ln int r ext
ðB16:1Þ ðB16:2Þ
References [1] R. Webb, N.H. Kim, Principles of Enhanced Heat Transfer, second ed., Taylor & Francis, New York, NY, 2005.
1195
[2] A.E. Bergles, Techniques to enhance heat transfer, in: Handbook of Heat Transfer, McGraw-Hill, New-York, NY, 1998. [3] W.H. Giedt, Investigation of variation of point unit heat-transfer coefficient around a cylinder normal to an air stream, ASME Trans. 71 (1949) 375–381. [4] S. Aiba, Y. Yamazaki, An experimental investigation of heat transfer around a tube in a bank, J. Heat Transf. 98 (3) (1976) 503–508. [5] J. Taler, Nonlinear steady-state inverse heat conduction problem with spacevariable boundary conditions, J. Heat Transf. 114 (4) (1992) 1048–1051. [6] J. Taler, Determination of local heat transfer coefficient from the solution of the inverse heat conduction problem, Forsch. Ingenieurwes. 71 (2) (2007) 69–78. [7] F. Bozzoli, G. Pagliarini, S. Rainieri, Experimental validation of the filtering technique approach applied to the restoration of the heat source field, Exp. Therm. Fluid Sci. 44 (2013) 858–867. [8] F. Bozzoli, L. Cattani, S. Rainieri, F.S.V. Bazán, L.S. Borges, Estimation of the local heat-transfer coefficient in the laminar flow regime in coiled tubes by the Tikhonov regularisation method, Int. J. Heat Mass Transf. 72 (2014) 352–361. [9] G.M. Carlomagno, G. Cardone, Infrared thermography for convective heat transfer measurements, Exp. Fluids 49 (2010) 1187–1218. [10] J.V. Beck, B. Backwell, C.R.S. Clair, Inverse Heat Conduction – Ill-posed Problems, Wiley-Interscience, New-York, NY, 1985. [11] H.R.B. Orlande, F. Olivier, D. Maillet, R.M. Cotta, Thermal Measurements and Inverse Techniques, Taylor and Francis, New York, NY, 2011. [12] S. Andrieux, A.B. Abda, The reciprocity gap: a general concept for flaws identification problems, Mech. Res. Commun. 20 (1993) 415–420. [13] S. Andrieux, A.B. Abda, Identification of planar cracks by complete overdetermined data: Inversion formulae, Inverse Prob. 12 (1996) 553–563. [14] F. Delbary, R. Aramini, G. Bozza, M. Brignone, M. Piana, On the use of the reciprocity gap functional in inverse scattering with near-field data: an application to mammography, J. Phys. Conf. Ser. 135 (2008). [15] M.J. Colaço, C.J.S. Alves, A fast non-intrusive method for estimating spatial thermal contact conductance by means of the reciprocity functional approach and the method of fundamental solutions, Int. J. Heat Mass Transf. 60 (2013) 653–663. [16] E.I. Shifrin, P.S. Shushpannikov, Identification of a spheroidal defect in an elastic solid using a reciprocity gap functional, Inverse Prob. 26 (5) (2010) 055001. [17] M.J. Colaço, C.J. Alves, H.R. Orlande, Transient non-intrusive method for estimating spatial thermal contact conductance by means of the reciprocity functional approach and the method of fundamental solutions, Inverse Prob. Sci. Eng. 23 (4) (2015) 688–717. [18] M.J. Colaço, C.J. Alves, F. Bozzoli, The reciprocity function approach applied to the non-intrusive estimation of spatially varying internal heat transfer coefficients in ducts: numerical and experimental results, Int. J. Heat Mass Transf. 90 (2015) 1221–1231. [19] M. Necati Ozisik, Heat Conduction, second ed., John Wiley & Sons, Inc., 1993. [20] P.C. Hansen, The truncated SVD as a method for regularization, BIT Numerical Math. 27 (4) (1987) 534–553. [21] A.N. Tikhonov, V.Y. Arsenin, Solution of Ill-Posed Problems, Winston & Sons, Washington DC, USA, 1977. [22] F. Bozzoli, S. Rainieri, G. Pagliarini, Experimental validation of the filtering technique approach applied to the restoration of the heat source field, Exp. Therm. Fluid Sci. 44 (2013) 858–867. [23] James V. Beck, Kenneth J. Arnold, Parameter Estimation in Engineering and Science, John Wiley & Sons, 1997. [24] Fermín S.V. Bazán, J.B. Francisco, Koung Hee Leem, G. Pelekanos, A maximum product criterion as a Tikhonov parameter choice rule for Kirsch’s factorization method, J. Comput. Appl. Math. 236 (2012) 4264–4275.
Contents lists available at ScienceDirect
International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt
Filtered reciprocity functional approach to estimate internal heat transfer coefficients in 2D cylindrical domains using infrared thermography Andrea Mocerino a, Marcelo J. Colaço b, Fabio Bozzoli a,⇑, Sara Rainieri a a b
Department of Industrial Engineering, University of Parma, 181/A Parco Area delle Scienze, I-43124 Parma, Italy Department of Mechanical Engineering, Federal University of Rio de Janeiro, Cid. Universitaria, Cx. Postal 68503, Rio de Janeiro, RJ 21941-972, Brazil
a r t i c l e
i n f o
Article history: Received 22 February 2018 Accepted 18 April 2018
Keywords: Inverse heat conduction problem Reciprocity functional Filtering technique Local heat transfer coefficient
a b s t r a c t Performance estimation of heat exchanger apparatuses requires a proper knowledge of the thermo-fluid dynamic interaction between fluid and device. One of the aspects to take into account is related to the estimation of the internal convective heat transfer coefficient. This requires the knowledge of the internal temperature and the internal wall heat flux, although these quantities are complicated to be measured. Therefore, the use of inverse problem techniques may be very useful for estimating these quantities. This approach deals with the estimation of the local internal properties, given some external temperature measurements, possibly by means of contactless experimental methodologies (i.e. infrared camera imaging). The solution strategy here presented, for a 2D model, is based on the Reciprocity Functional approach (RF), which requires the solution of two auxiliary problems that are solved, in this paper, by means of the Classical Integral Transform Technique (CITT). The RF approach coupled with CITT, presents some advantages over other techniques since its solution is completely analytical and therefore computationally fast, and it does not require any numerical simulation. Moreover, it deals with the noise content of the experimental data, filtering it and overcoming the documented limits of the classical RF approach. This original method, named Filtered Reciprocity Functional (FRF) approach, was applied to different test cases using simulated measurements, and the results were compared to both the exact and the reconstruct solutions obtained by the Truncated Singular Value Decomposition (TSVD) method, showing a good agreement. Finally, the FRF approach was applied to an experimental case in order to test its robustness when dealing with real data. Ó 2018 Elsevier Ltd. All rights reserved.
1. Introduction: Estimating the local heat flux on surfaces is crucial in multiple industrial applications. While measuring temperature is usually an easy task, direct measurement of heat flux is not always straightforward and may not be even possible. For instance, in commercial heat exchangers, for savings in materials and energy use, passive heat transfer enhancement techniques are usually adopted (e.g., treated surfaces, rough surfaces, displaced enhancement devices, swirl-flow devices, surface-tension devices, coiled tubes, or flow additives), [1]. These techniques originate a significant variation of the heat transfer coefficient at the fluid-wall interface along the surfaces and, in
⇑ Corresponding author. E-mail address: [email protected] (F. Bozzoli). https://doi.org/10.1016/j.ijheatmasstransfer.2018.04.089 0017-9310/Ó 2018 Elsevier Ltd. All rights reserved.
some applications in which uniform thermal processes are required, this irregular distribution behaviour impacts negatively on the efficiency of these devices, [2]. For this reason, for optimising these devices, it is fundamental to determine what happens on the internal surfaces in terms of local surface temperature and heat fluxes. The hint for succeeding in the estimation of the local behaviour comes from some research papers in which the space-variable heat transfer coefficient on the inner/outer pipe surface was determined, given measurements of temperature at some points of the pipe wall. Some of these authors [3,4] used a specific experimental setup, based on a thin metallic layer, applied to the tube wall, which was heated by Joule effect. Other authors [5–9] adopted the inverse heat conduction problem (IHCP) solution approach to estimate the local convective heat transfer coefficient starting from the temperature distribution without the need of a special heating system.
1182
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Nomenclature r
H C h qg X n K Y G,F S
w N
c m
k
radial coordinate (m) circumferential coordinate (rad) domain’s boundary convective heat transfer coefficient (W/m2 K) heat source per unit volume (W/m3) cross sectional domain unity normal vector pointing outward the domain thermal conductivity (W/m K) vector of measured temperatures available at the external surface (K) auxiliary functions numbers of external surface temperature measurements generic orthonormal basis function total number of basis trace of the auxiliary function at the internal boundary numbers of harmonics of the eigenfunction number of harmonic of the orthonormal basis
Although the IHCP approach is more promising because it can be applied to every kind of heating configuration and pipe geometry, it presents some complications due to the fact that IHCPs are ill-posed and, consequently, are very sensitive to small perturbations in the input data. Different methods have been developed for solving IHCPs using analytical and numerical techniques. Some of these methods are the conjugate gradient iterative method; Laplace transform method; sequential function specification method; regularization methods, such as Tikhonov regularization; mollification method; truncated singular value decomposition method; and filtering technique method [10,11]. Particularly promising for this kind of application is the Reciprocity Function approach (RF), which needs low computational resources because it avoids the use of iterative techniques. Two important works by Andrieux and Ben Abda [12,13] showed the concept and use of reciprocity functional. From these, other studies based on the reciprocity functional began to emerge in different areas. Delbary et al. [14] developed a qualitative method for breast cancer detection by combining the reciprocity functional method with the linear sampling method. Colaço and Alves [15] estimated spatial variation of the thermal contact conductance by using a reciprocity functional approach with the method of fundamental solutions and non-intrusive temperature measurements. Shifrin and Shushpannikov [16] developed a method for identifying small defects in an anisotropic elastic body based on the reciprocity functional. Other studies, regarding the estimation of the thermal contact conductance through reciprocity functional, using non-intrusive measures, can be found in [15,17]. In recent times, Colaço et al. [18] presented a methodology, based on the reciprocity functional approach, to estimate internal convective heat transfer coefficients in ducts, using only data available at an exterior boundary and the solutions of two auxiliary problems that depends only on the system geometry. By this approach, here referred as Numerically integrated Reciprocity Function (NRF) method, the two auxiliary problems are solved numerically and the unknown function can be estimated solving a linear system, where the solution vector is composed of integrals of the measured boundary data. The procedure was analysed and compared with traditional techniques showing very good results. However, some limits were found: in particular, it was not clear how to choose the regularization term (i.e., the number of
a,b X U; V; R
coefficient of the Fourier’s series expansion of the measurement sensitivity matrix matrices of the singular value decomposition
Subscripts int internal ext external b bulk env environmental avg average Acronyms RF Reciprocity Functional NRF Numerical Integrated Reciprocity Functional FRF Filtered Reciprocity Functional CITT Classical Integral Transform Technique TSVD Truncated Singular Value Decomposition
orthonormal functions used for the two auxiliary problems) by traditional criteria (e.g., discrepancy principle, L-curve, fixed point). The present paper starts from the findings of Colaço et al. [18] and overcomes the limit of the previous work solving the auxiliary problems not by numerical approaches, but instead by the Classical Integral Transform Technique (CITT), a well-known analytical method [19]. In this way, thanks to a fully analytical approach, it is possible to understand how the RF approach deals with the noise content in the experimental data and therefore the regularization factor can be tuned in order to improve the robustness of the approach. The new solution methodology, named Filtered Reciprocity Functional (FRF), for its ability to filter out the unwanted noisy content, is firstly tested using synthetic data, considering different profiles for the internal heat transfer coefficient, and then it is applied to an experimental case with real measurements. 2. Direct problem The physical problem, represented in Fig. 1, consists in a crosssection of a circular duct with internal radius rint and external radius rext. The thermal conductivity of the material, K, is assumed to be constant and uniform such as the environmental temperature Tenv and the bulk temperature Tb of the fluid that flows inside the tube. Both the external and internal surfaces are subjected to Robin’s boundary condition: there is a prescribed convective heat transfer coefficient distribution hint at the inner surface Cint of the pipe, while the exterior surface Cext is subjected to an overall convective heat transfer coefficient hext. The domain X is also subjected to an internal heat generation per unit of volume qg. The equations that describe the physical problem are listed below:
r2 T ¼
qG K
in X
ð1:aÞ
K
@T ¼ hint ðT T b Þ on Cint @nCint
ð1:bÞ
K
@T ¼ hext ðT T env Þ on Cext @nCext
ð1:cÞ
1183
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
on Cext
Gk ¼ wk
ð2:cÞ
where wk is a set of orthogonal functions, which, for a 2D case, can be written as a standard orthonormal Fourier’s basis (3.a)–(3.c):
1 wk ¼ pffiffiffiffiffiffiffipffiffiffiffiffiffiffi r ext 2p
for k ¼ 1
ð3:aÞ
1 k h for k ¼ 2; 4; 6; . . . N 1 wk ¼ pffiffiffiffiffiffiffipffiffiffiffi cos 2 r ext p
ð3:bÞ
1 k1 wk ¼ pffiffiffiffiffiffiffipffiffiffiffi sin h for k ¼ 3; 5; 7; . . . N 2 r ext p
ð3:cÞ
Defining the reciprocity functional and the integral of the internal heat source, respectively, as follows
@T @Gk Y wk dCext @n @n Cext Z henv ðT T env Þ @Gk Y ¼ wk dCext K @n Cext Z
RG;k ¼
Q G;k ¼
Z qg dX Gk K X
ð4Þ
ð5Þ
it is possible to write the following identity:
Fig. 1. Physical model.
@T dCint ¼ RG;k Q G;k Gk @n Cint
Z
3. Inverse problem The direct problem presented in Eqs. (1.a)–(1.c) concerns the determination of the temperature field T within the O domain given the distribution of the convective heat transfer coefficient on both the internal and external surfaces. On the other hand, the inverse problem here described is addressed to restore the internal heat transfer coefficient hint given some extra temperature measurements Y taken on the external surface Cext. The methodology used to solve this problem is based on the reciprocity functional approach [12,13]; this technique, under the assumption of linear inverse problem, is able to restore the internal heat transfer coefficient avoiding the use of iterative procedures. This approach, as it will be shown, estimates the unknown function by solving a linear system described by a diagonal matrix, which speeds up the estimation procedure, that involves the integrals of the external measurements. The determination of the unknown internal heat transfer coefficient is performed by solving two auxiliary problems [18]: the first one concerns the determination of the normal derivative of the temperature @T/@n at the internal surface Cint, which is related to the heat flux per unit of surface; while the second one estimates the wall temperature at the same surface. According to Eq. (1.b), and assuming that Tb and K are known values, the local convective heat transfer coefficient distribution hint can be easily determined on Cint. 3.1. Determining the normal derivative of the temperature at the internal boundary @T/@n As suggested by Colaço et al. [18], the following well-posed auxiliary problem, given by Eqs. (2.a)–(2.c), can be used to obtain @T/@n|Cint:
r2 Gk ¼ 0 in X
ð2:aÞ
@Gk ¼ 0 on Cint @n Cint
ð2:bÞ
ð6Þ
More details about the equation above can be found in [18]. Using the orthonormal basis defined in Eqs. (3.a)–(3.c), and defining ck as the trace of the solution Gk on Cint, the last expression can be manipulated in order to get:
Z
RG;k Q G;k ¼
Z @T @T Gk ck dCint ¼ dCint @n @n Cint Cint
ð7Þ
The normal derivative of the temperature on the internal surface can be approximated by using the induced basis c1, . . ., cN as follows:
@T ¼ a1 c1 þ þ aN cN @nCint
ð8Þ
Therefore, truncating the expansion with N terms, the system of Eq. (7) is reduced to the following linear system N Z X k¼1
Cint
ck aj cj dCint ¼
N Z X k¼1
Cint
ck cj aj dCint
¼ RG Q G for j ¼ 1; 2; 3; . . . ; N
ð9Þ
Eq. (9) can be written in matrix form as:
Ma ¼ RG Q G
ð10Þ
The system is invertible since c1,. . .,cN are linearly independent and, therefore, the aj coefficients can be found as:
a ¼ M1 ðRG Q G Þ
ð11Þ
After that, the estimation of @T/@n|Cint is obtained by simply using the expansion presented in Eq. (8). The auxiliary problem presented in Eqs. (2.a)–(2.c) is then solved analytically by using the Classical Integral Transform Technique [19]. Initially, Eq. (2.a) is expressed in the cylindrical coordinate system:
@ 2 Gk 1 @Gk 1 @ 2 Gk þ ¼0 þ 2 r @r r @h2 @r 2 Then, defining the integral transform as
ð12Þ
1184
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Gk;m Gk ðr; mÞ ¼
Z 2p h¼0
Gk ðr; h0 Þcos½mðh h0 Þdh0
ð13Þ
and the inversion formula as
Gk Gk ðr; hÞ ¼
1 X 1 1 Gk ðr; m ¼ 0Þ þ G ðr; mÞ 2p p k m¼1
ð14Þ
Eq. (12), multiplied by cos[m(h-h’)] and integrated, becomes: 2
d Gk;m 1 dGk;m 1 2 þ þ 2 m Gk;m ¼ 0 r dr r dr2
ð15Þ
With the same procedure, Eq. (2.b), which expresses the boundary conditions imposed on the internal boundary Cint, becomes:
dGk;m ¼0 dr
ð16Þ
and the external boundary condition can be re-written as follows:
Gk;m ¼ wk;m
ð17Þ
where the transformed potential, present in the first auxiliary problem given by Eq. (2.c), is given as:
wk;m ¼
Z 2p h0 ¼0
wk cos ½mðh h0 Þdh0
ð18Þ
Using the above equations, the transformed first auxiliary problem is finally obtained as: 2
d Gk;m 1 dGk;m 1 2 þ þ 2 m Gk;m ¼ 0 in X r dr r dr2 dGk;m ¼ 0 on Cint dr
ð19:aÞ
ð19:bÞ
Cint
on Cext
Gk;m ¼ wk;m
ð19:cÞ
whose solution can be found as:
Gk;m ¼ wk;m cosh
n h
m ln
r io r int int sech m ln r r ext
ð20Þ
In order to obtain the solution of the first auxiliary problem, the above equation has to be inverted using the definition given by Eq. (14):
Gk ¼
1 w 2p k;m¼0
1 n h r io X 1 r int int þ wk;m cosh m ln sech m ln p r r ext m¼1
2
1
6 6 6 6 6 6 6 6 6 rint 6 6 M¼ r ext 6 6 6 6 6 6 6 6 6 4
2
h
sech 1 ln
r int r ext
ð21Þ
In this paper, the measurements Y are decomposed by the Fourier series expansion, in order to take advantage of the harmonic decomposition of the Classical Integral Transform Technique and efficiently approach the inverse problem:
YðhÞ ¼
N a0 X ½am cos ðmhÞ þ bm sin ðmhÞ þ 2 m¼1
where ai and bi are real numbers and m is an integer number. Thanks to this, it is possible to distinguish between the signal and the noise. In almost every practical situation [7], the signal content is mainly present at low frequencies, while noise is present over the whole spectrum. This approach, whose details are reported in Appendix A, enables to conveniently express the terms of Eq. (10) in the following way:
3
2
am¼1 6a 7 6 m¼2 7
7 a¼6 6 .. 7 5 4
ð23Þ
.
am¼m pffiffiffi h 2 enK v ða0 T env Þ
2
3
6 6 n h i o 6 henv 1 int 6 tanh ln rrext rext K 6 6 6 n h i o 6 henv 1 int 6 tanh ln rrext rext K 6 6 6n h i o 6 6 henv tanh 2 ln rint 2 6 K r ext rext pffiffiffiffiffiffiffiffiffiffi6 RG ¼ pr ext 6 6n h i o 6 henv r int 2 6 K tanh 2 ln rext rext 6 6 6 6 .. 6 . 6 6 6n h i o 6 henv rint m 6 K tanh m ln rext r ext 6 6 4n h i o henv m int tanh m ln rrext K r ext 2 q qffiffiffiffiffiffiffi g
6 6 6 QG ¼ 6 6 4
K
p
2r ext
r 2ext r 2int 0 .. .
7 7 7 a1 7 7 7 7 7 b1 7 7 7 7 7 a2 7 7 7 7 7 7 b2 7 7 7 7 7 7 7 7 7 7 am 7 7 7 5 bm
ð24Þ
3 7 7 7 7 7 5
ð25Þ
0
3
0
i 2
h
sech 1 ln
rint r ext
i h i 2 int sech 2 ln rrext
h i 2 int sech 2 ln rrext ..
. sech
0
ð22Þ
2
h
m ln
r int r ext
i sech
2
h
m ln
r int r ext
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 i 5
ð26Þ
1185
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
For this reason, these terms work as operators that amplifies the high frequency components of RG more than the low ones of a given input data set Y. In Fig. 2, the values of M1components are plotted for a representative ratio rint/rext = 0.9 as a function of the harmonic order. For cases where both the signal and the noise have frequency components equally distributed over the spectrum of interest, the M1 matrix does not deteriorate the signal content in comparison to the original data. Unfortunately, in almost all practical situations, the signal has frequency components concentrated in the low frequency range of the spectrum, while the spectral components of noise is expected to be uniformly distributed over the entire frequency domain, [7]. Therefore the M1 operator makes very difficult the extraction of some useful information from the RG vector. For this reason, some regularization technique is needed. In the present work this problem is tackled by conveniently truncating the M1 matrix, considering only the first m harmonics as suggested by Bozzoli et al. in a analogous inverse problem [22]. The number of harmonics to keep is chosen by the classical discrepancy principle, originally formulated by Morozov [21]. More efficient approaches in terms of computational costs are available, but this particular aspect goes beyond the scope of the present paper and it will be investigated in upcoming works.
Fig. 2. Diagonal elements of the M1 matrix.
2
3
p1ffiffi 2
h i 7 6 6 cos ð1hÞsech 1 ln rint 7 7 6 r ext 6 h i 7 7 6 r int 7 6 sin ð1hÞsech 1 ln r 7 6 ext 6 h i 7 7 6 rint 1 6 cos ð2hÞsech 2 ln rext 7 h i 7 c ¼ pffiffiffiffiffiffiffiffiffiffi 6 7 6 rint 7 prext 6 6 sin ð2hÞsech 2 ln rext 7 7 6 7 6 .. 7 6 .h 7 6 i 7 6 6 cos ðmhÞsech m ln rint 7 7 6 r ext 4 h i 5 int sin ðmhÞsech m ln rrext
3.2. Determining the internal wall temperature T|Cint
ð27Þ
To estimate the internal wall temperature, the reciprocity functional approach is used. As suggested by Colaço et al. [18], the following auxiliary problem, given by Eqs. (29.a)–(29.c), can be used to obtain T|Cint:
The values of the aj coefficients, which are the goal of the first auxiliary problem, are easily obtained from a mathematical point of view by Eq. (11). Due to the ill-posed nature of the problem and the presence of the noise in the measurement vector Y {Y1, Y2, . . ., Yn}, the inversion of the matrix M has to be approached by inverse problem solution techniques in order to find a realistic estimation of the unknown heat flux. Since the matrix M is diagonal, its inverse can be easily determined as follows:
2
1
M
6 6 6 6 6 6 6 6 6 6 rext 6 6 ¼ r int 6 6 6 6 6 6 6 6 6 6 4
r2 F k ¼ 0 in X
ð29:aÞ
F k ¼ 0 on Cint
ð29:bÞ
F k ¼ wk
on Cext
ð29:cÞ
This problem is well posed and, for simplicity, we used the same basis defined in the first auxiliary problem, Eqs. (3.a)–(3.c). Following the same steps used in the first auxiliary problem, the second reciprocity functional and the integral of the heat generation terms can be defined as:
1
3
0
1
sech2 1 ln
r int r ext
1
sech2 1 ln
r int rext
1
sech2 2 ln
r int rext
1
sech2 2 ln
r int r ext
..
. sech2
0
The analytical expression achieved clearly highlights the destructive effect of the noise on the reconstruction of the aj coefficient: the values of M1 matrix components, following 1/sech2[m ln(rint/rext)], increase with the number of the harmonic.
1
m ln
r int r ext
sech2
1
m ln
r
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5
ð28Þ
int r ext
@T @F k dCext wk Y @n @n Cext Z henv ðT T env Þ @F k Y ¼ wk dCext K @n Cext Z
RF;k ¼
ð30Þ
1186
Q F;k
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Z qg ¼ Fk dX K X
ð31Þ
Therefore, it is possible to write the following expression:
@F k T dCint ¼ RF;k þ Q F;k @n Cint
Z
ð32Þ
r2 F k;m ¼ 0 in X
ð44:aÞ
F k;m ¼ 0 on Cint
ð44:bÞ
on Cext
F k;m ¼ wk;m
ð44:cÞ
whose solutions are:
r io r int int csch m ln r r ext
Using the orthonormal basis defined in Eqs. (3.a)–(3.c) and defining n as the trace of the solution Fk on Cint, the last expression can be manipulated in order to get:
F k;m ¼ wk;m sinh
RF;k þ Q F;k ¼
In order to obtain the solution of the second auxiliary problem, the above equations need to be inverted using the definition given in Eq. (40):
Z
Z @F k T ðnk T ÞdCint dCint ¼ @n Cint Cint
ð33Þ
The temperature on the internal surface can be approximated by using the induced basis n1, . . ., nN as follows:
TjCint ¼ b1 n1 þ þ bN nN
ð34Þ
Therefore, truncating the expansion with N terms, Eq. (33) is reduced to the following linear system (35): N Z X k¼1
Cint
N X
nk bj nj dCint ¼ k¼1
Z
Cint
nk nj bj dCint
¼ RF;k Q F;k for j ¼ 1; 2; 3; . . . ; N
r ext
þ
h¼0
F k ðr; h0 Þcos½mðh h0 Þdh0
ð38Þ
ð39Þ
and the inversion formula as (40):
F k F k ðr; hÞ ¼
1 X 1 1 F k ðr; m ¼ 0Þ þ F k ðr; mÞ 2p p m¼1
ð40Þ
it is possible to transform Eq. (38) as follows: 2
d F k;m 1 dF k;m 1 2 þ þ 2 m F k;m ¼ 0 r dr r dr 2
ð41Þ
By an analogous elaboration, the transformation of the boundary conditions imposed on the internal boundary Cint, Eq. (29.b) results in:
F k;m ¼ 0
ð42Þ
In the same way, the external boundary condition can be rewritten as:
F k;m ¼ wk;m
p
n h
m ln
r io r int int csch m ln r r ext
ð46Þ
YðhÞ ¼
Defining the integral transform as:
Z 2p
m¼1
wk;m sinh
ð36Þ
After that, the estimate of T|Cint is obtained simply by using the expansion presented in Eq. (34). Analogously to the approach adopted in the previous problem, the second auxiliary problem presented in Eqs. (29.a)–(29.c) is solved analytically by using the Classical Integral Transform Technique. For 2D problems, it is possible to write the Laplacian defined in Eq. (29.a) as:
F k;m F k ðr; mÞ ¼
1 X 1
ð35Þ
ð37Þ
@ 2 F k 1 @F k 1 @ 2 F k þ ¼0 þ r @r r 2 @h2 @r2
ð45Þ
Analogously to the previous problem, in order to take advantage of the harmonic decomposition of the Classical Integral Transform Technique and to obtain a fully analytical expression, it is useful to represent the measurements Y as a Fourier series expansion:
The system is invertible since n1, . . ., nN are linearly independent. Eq. (35) is easy to solve having the bj coefficient as unknown:
b ¼ N1 ðRF Q F Þ
m ln
ln rint 1 r Fk ¼ w 2p k;m¼0 ln rint
which can be written in the following matrix form:
Nb ¼ RF Q F
n h
ð43Þ
where wk;m is the transformation of the potential imposed in the second auxiliary problem, Eq. (29.c), according to Eq. (18). Using the above equations, the transformed second auxiliary problem is obtained as:
N a0 X ½am cos ðmhÞ þ bm sin ðmhÞ þ 2 m¼1
ð47Þ
This approach, whose details are reported in Appendix B, enables to conveniently express the terms of Eq. (36) in the following way:
3 bm¼1 7 6 7 6 6 bm¼2 7 7 6 7 b¼6 6 . 7 6 .. 7 7 6 5 4 bm¼m 2
ð48Þ
2 pffiffiffi henv K 2 1 T a 0 en v r K 6 henv rext ln r int 6 ext 6 6 n h i o 6 henv 1 int 6 coth ln rrext K r ext 6 6 6 n h i o 6 henv 1 int 6 coth ln rrext r ext K 6 6 6 n h i o 6 rint henv 2 coth 2 ln pffiffiffiffiffiffiffiffiffiffi6 r ext K r ext 6 RF ¼ pr ext 6 6 n h i o 6 henv 2 int 6 coth 2 ln rrext K r ext 6 6 6 6 .. 6 . 6 6 6 n h i o 6 henv m int coth m ln rrext 6 K r ext 6 6 4 n h i o henv m int coth m ln rrext K r ext 2 qg
6K 6 6 6 6 QF ¼ 6 6 6 6 4
qffiffiffiffiffiffiffi p
2r ext
r 2ext þ 12 0 .. . 0
ðr2ext r2int Þ ln
r int rext
3 7 7 7 7 7 a1 7 7 7 7 7 b1 7 7 7 7 7 a2 7 7 7 7 7 b2 7 7 7 7 7 7 7 7 7 7 am 7 7 7 5 bm
ð49Þ
3 7 7 7 7 7 7 7 7 7 5
ð50Þ
1187
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
2 6 6 6 6 6 6 6 6 6 6 1 6 6 N¼ rint r ext 6 6 6 6 6 6 6 6 6 6 4
ln2
1r
3
0
int r ext
h i 2 int csch 1ln rrext
h i 2 int csch 1ln rrext
h i 2 int 22 csch 2 ln rrext
h i 2 int 22 csch 2 ln rrext ..
.
h
m2 csch2 m ln
r int rext
i h
m2 csch2 m ln
0
h
r int r ext
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 ii 5 ð51Þ
2
3
1 p1ffiffi 2 ln rint r int r
ext h i 7 7 7 int cos ð1hÞcsch 1 ln rrext 7 h i 7 7 r int sin ð1hÞcsch 1 ln rext 7 7 h i 7 7 r int cos ð2hÞcsch 2 ln rext 7 h i 7 7 7 int sin ð2hÞcsch 2 ln rrext 7 7 7 .. 7 . h i 7 7 int 7 cos ðmhÞcsch m ln rrext 7 h i 5 m sin ðmhÞcsch m ln rint r int rext
6 6 6 1 6 rint 6 6 1 6r 6 int 6 6 2 1 6 rint n ¼ pffiffiffiffiffiffiffiffiffiffi 6 prext 6 6 2 6 rint 6 6 6 6 6 m 6 6 rint 4
ð52Þ
The values of bj coefficients, which are the goal of the second auxiliary problem, are obtained from a mathematical point of view by Eq. (37). Due to the ill-posed nature of the problem and the presence of the noise in the measurement vector Y, the inversion of the matrix N has to be approached by inverse problem solution techniques in order to find a realistic estimation of the unknown heat flux. The diagonal behaviour of the matrix of the inner product N allows us to easily calculating its inverse as follow:
2
1
N
6 6 6 6 6 6 6 6 6 6 6 6 ¼ r int r ext 6 6 6 6 6 6 6 6 6 6 6 6 4
2
ln
r int r ext
Fig. 3. Diagonal elements of the N1 matrix.
3
0 2
csch
1 1 ln
r int rext
1
csch2 1 ln
r int rext
1 r
22 csch2 2 ln
int rext
1 r
22 csch2 2 ln
int rext
..
. 1 r
m2 csch2 m ln
0
int rext
1
m2 csch2 m ln
7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 r 5 int r ext
ð53Þ
1188
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
This operator has many points in common with the operator defined in Eq. (28) for the previous auxiliary problem. It means that, also in this case, the presence of noise has a negative impact on the estimation of the b coefficients because the values of N1 matrix components, following 1/m2 csch2[m ln(rint/rext)], increase with the number of the harmonic. In Fig. 3, the values of N1 components are plotted for a representative ratio rint/rext = 0.9 as a function of the harmonic order. Also in this case, the values of the diagonal terms of the N1 matrix increase as the harmonic order increases; for this reason the matrix is regularized by keeping the same numer of harmonics employed in the solution of the first auxiliary problem. 4. Results The solution methodology above presented is firstly tested using synthetic data, considering different profiles for the internal heat transfer coefficient, and then it is applied to an experimental case with real measurements. The result given by the Filtered Reciprocity Functional (FRF) approach are compared to the result given by the Numerically integrated Reciprocity Functional (NRF) [18] and then to the well-known TSVD [20] method in order to qualify the goodness of the new proposed approach. Since the direct problem described by Eq. (1) is linear with respect to the internal heat flux q, in the discrete domain, the following linear system can be written:
T ¼ Xq þ Tq¼0
ð54Þ
where T is the vector of the temperature distribution at the external tube surface Cext given by the internal wall heat flux q at the fluidwall interface, Tq=0 is a constant term and X is the sensitivity matrix. More details about this approach can be found in [8]. Defining a relative temperature as follows:
Trel ¼ T Tq¼0
ð55Þ
from Eq. (54), the internal wall heat flux q can be estimated:
q ¼ X1 Yrel ¼ X1 Y Tq¼0
ð56Þ
The TSVD approach here used, requires the decomposition of the X matrix into three different matrices U, R and V:
X ¼ URVT ¼
n X ui ri v Ti
ð57Þ
i¼1
where U = [u1, . . ., un] and V = [v1, . . ., vn] are matrices with orthogonal columns, and the matrix R has diagonal form:
R ¼ diag ðr1 ; r2 ; . . . ; rn Þ
ð58Þ
The matrix R has non-negative diagonal elements, called the singular values of X, whose magnitude decrease as their order increase. The purpose of the TSVD is to force the solution of the
inverse problem, given by Eq. (56), to have a small norm in order to damp the contributions from the errors present in the measurement contained in the Trel term [20]. This is obtained by the truncation of the solution components that correspond to the smallest singular values, because they bring the largest contribution to the noise amplification in the inversion of the matrix. The truncation parameter tp can be effectively chosen by the well-known Morozov’s discrepancy principle [21]. Therefore, the TSVD of X becomes:
Xt ¼ U Rt V with Rt ¼ diag
r1 ; r2 ; . . . rtp 1 ; 0; . . . ; 0
ð59Þ
where tp is the truncation parameter adopting the Morozov’s discrepancy principle [21]. The internal heat flux is given by the following expression:
1 T T q ¼ VR1 t U Y rel ¼ VRt U Y Tq¼0
ð61Þ
which can be seen also like the solution of the following leastsquares problem:
min q Xt q þ Tq¼0 Y
ð62Þ
Once the heat-flux distribution at the internal wall surface has been determined, the local convective heat-transfer coefficient can be easily determined, as follows:
hint ¼
q TjCint T b
ð62Þ
where T|Cint is the internal wall temperature determined by solving the direct problem once the internal wall heat flux is estimated and Tb is the bulk-fluid temperature on the test section. 4.1. Synthetic data The physical problem considered in the present work, according to the conditions reported in the available literature [5–9], consist in a cross-section of a circular duct with internal radius rint = 0.007 m and external radius rext = 0.008 m. The thermal conductivity of the pipe K = 15 W/mK was assumed to be constant such as the environmental temperature Tenv = 298.15 K and the bulk temperature Tb = 292.15 K of the fluid that flows inside the tube. The exterior surface Cext of the tube was subjected to an overall convective heat transfer coefficient henv = 5W/m2K with the environment. The domain X was assumed to be homogeneous and isotropic, and subjected to an internal heat generation per unit of volume qg = 106 W/m3. The Filtered Reciprocity Functional approach (FRF) was first compared to Numerically integrated Reciprocity Functional (NRF) [18] and then to the well-known TSVD approach [20] within the MatlabÒ environment by adopting synthetic data. The finite element method implemented in Comsol MultiphysicÒ was used in order to generate the synthetic measurements (sampled at the
Table 1 hint profiles used as test functions. Case h1
h2
hint
h
i
W m2 K
p < h < 23 p
0
i
W m2 K
hmin
h
i
W m2 K
hmax
0
600
4000
0
6000
p < h < p
400
0
600
p < h < p
4000
0
6000
23 p < h < 23 p
0
2 3
0
p < h < 23 p
p
23 p < h < 23 p
2 3
0
h4
h
400
600
6000
h3
have
h2
600 p2 þ 600 h2
6000 p2 þ 6000
p
h
i
W m2 K
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
1189
Fig. 4. Reconstruction of the internal heat transfer coefficient for (a) test case h1, (b) test case h2, (c) test case h3 and (d) test case h4 (r = 0.1 K).
external surface of the pipe) by imposing a known distribution of h at the internal wall surface of the tube. In the present analysis, to compare straightforwardly the three considered approaches, the criterion provided by the discrepancy principle [21] was adopted for all of them. Since it is well known that, the effectiveness of the estimated solution of the inverse problem depends on the shape of the function being estimated, different test cases were implemented in order to perform a robust comparison. The shapes of the convective heat transfer coefficients that were used as test functions are reported in Table 1. Test cases h1 and h2 are the classical step function with different amplitudes, chosen in order to test the limits of the estimation approaches to estimate the discontinuities, while cases h3 and h4 are instead examples of extremely smooth functions. Moreover, test cases h1-h3 and h2-h4 were assumed to have the same maximum hmax, minimum hmin and average have in order to make the comparison of different test functions more objective, minimizing the effect of the signal’s amplitude. In order to perform a trustworthy comparison among the techniques, the synthetic measurement temperatures Y were deliberately spoiled by random noise and then used as input data for both the mathematical models. In particular, a Gaussian white noise characterized by a standard deviation ranging from r = 0.01 K to r = 5K was introduced according to:
Y ¼ T exact jCext þ r
ð63Þ
where e is a random Gaussian variable with zero mean and unit variance. As shown in Figs. 4 and 5, which compare the estimated functions to the exact ones for some representative cases, the three approaches show different degrees of accuracy depending on the noise level and the test case considered. In order to perform a quantitative comparison among the three considered approaches, the global relative estimation error between the exact and the estimated internal heat transfer coefficient was calculated as:
E¼
khreconstructed hexact k2 khexact k2
ð64Þ
The added noise depends intrinsically on the random sequence generated by MatlabÒ. Therefore the estimation procedure was repeated 50 times for different random noise sequences and an average value Eavg was calculated for each noise level. As highlighted by Fig. 6, FRF estimates far better than the classical NRF, which shows always the worst behaviour. It has to be highlighted that the stopping criterion used affects the quality of the reconstruction as already shown by Bazàn et al. [24] that compared the fixed-point criteria with the Morozov’s discrepancy principle; moreover, Colaço et al. [15] already showed that the reconstruction due the NFR could be better if a different stopping criterion is used. It means that the harmonic decomposition of the Classical Integral Transform Technique enables to efficiently remove the unwanted information from the measured data. Unex-
1190
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Fig. 5. Reconstruction of the internal heat transfer coefficient for (a) test case h1, (b) test case h2, (c) test case h3 and, (d) test case h4 (r = 1K).
pectedly, FRF works also better than TSVD in almost all conditions, giving this approach interest for future applications. As observed in all inverse problems [23], the presence of the noise in the measurements corrupts the reconstruction of the unknown heat flux distribution. In almost every practical condition, the noise is uniformly distributed on the whole spectrum while the signal is almost always present in the lower harmonic orders; under this statement, the magnification effect of the M1 and N1 matrices, which increases as the harmonics order increase, requires the truncation of their diagonal elements. This kind of truncation can be performed according to the Morozov’s discrepancy principle in the same fashion as it was already adopted for the TSVD approach. This behaviour is clarified by Fig. 7, which shows the reciprocity functional values for both the first and the second auxiliary problem as a function of the harmonic order, for the two representative test case h1 and h3. In order to improve the readability of these plots, the reciprocity functionals are expressed in a more compact polar form, following the classical notation:
RG ¼
rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 R2G þ R2G
ð65:aÞ
RF ¼
rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 R2F k þ R2F k þ1
ð65:bÞ
ðk Þ
ð Þ
ðk
ð
þ1Þ
Þ
where: k⁄ = 2, 4, 6, . . ., N 1.
In all the graphs present in Fig. 7, the reciprocity functional of the noiseless signal is compared to the reciprocity functional of a representative set of white noise, having standard deviation of r = 0.1 K, that were used to spoil the synthetic data. They confirm that, in the reciprocity functionals, the signal prevails for the lower harmonics, while the noise overcomes the signal for the higher harmonics. Lastly, Table 2 reports t that is the number of elements kept in the truncation of M1 and N1 matrix for FRF methodology, for the four test functions and for different noise levels. This number decreases as far as the noise increases; less harmonics are kept for the continuous functions (i.e., case h3 and h4) than for the discontinuous ones (i.e., test case h1 and h2). This behaviour is similar to what was observed for TSVD and NRF in [10] and [18], respectively.
4.2. Experimental data An experimental data set [8] was used to validate the proposed technique. A sketch of the experimental setup is reported in Fig. 8: the solid domain is internally heated by Joule effect and it is thermally insulated, except for a small portion at the end of the tube; in this section the thermographic images of the external wall surface were taken. The experimental problem consists in a stainless steel coiled tube with circular cross-section having an internal radius rint = 0.007 m and an external radius rext = 0.008 m. The thermal conductivity of the pipe K = 15 W/mK is assumed constant as well
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
1191
Fig. 6. Comparison between the NRF, FRF and the TSVD methods for different test case as a function of the dimensionless noise level: (a) h1, (b) h2, (c) h3, and (d) h4.
as the environmental temperature Tenv = 296.8 K and the bulk temperature Tb = 294.7 K of the fluid that flows inside the tube. The exterior surface Cext of the tube was subjected to an overall convective heat transfer coefficient henv = 5 W/m2K with the environment while the domain X, that is assumed to be homogeneous and isotropic, was subjected to an internal heat generation per unit of volume qg = 4.8 106 W/m3. The result given by the FRF approach, here developed and presented in Fig. 9, was compared to the result given by the TSVD and the NRF where the regularization parameter (i.e., the number of bases considered) was determined using Morozov’s Discrepancy Principle. The results, reported in Fig. 9, show that the FRF, NFS and TSVD give the same results confirming the robustness of the new approach here presented. 5. Conclusions The approach here developed deals with the estimation of the local internal heat transfer coefficient in a convective heat transfer problem, given the external temperature measurements, by solving an Inverse Heat Conduction Problem. The solution strategy presented was developed for the 2D model, and was based on the Reciprocity Functional approach (RF), solved by mean of the Classical Integral Transform Technique (CITT). This approach, named Filter Reciprocity Functional (FRF) approach, presents some practical advantages over the more clas-
sical inverse problem solution techniques since it is completely analytical and, for this reason, it is computationally inexpensive. The application of FRF to numerical and experimental benchmarks highlighted its effectiveness and robustness, suggesting the application to other challenging inverse problems. Moreover, the obtained formulation enables, from the theoretical point of view, to highlight the ill-posed nature of the inverse problem and, from the practical point of view to filter efficiently the noisy content thanks to the harmonic nature of CITT and the Fourier series expansion. Conflict of interest The authors declare that there is no conflict of interest. Appendix A. Derivation of the solution of the first auxiliary problem According to the definition given in Eqs. (3.a)–(3.c), the solution has to be divided in at least in three parts for different values of the harmonic frequencies k of the orthonormal basis function: (a) k = 1 According to Eqs. (3.a) and (18), for k = 1, we have wm–0;k¼1 ¼ 0, qffiffiffiffiffi and wm¼0;k¼1 ¼ r2extp . Substituting these results in Eq. (21) the auxil-
1192
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
Fig. 7. RG and RF for two representative test case: (a) RG for test case h1, (b) RG for test case h3, (c) RF for test case h1 and (d) RF for test case h3.
Table 2 Number of harmonic components t used in FRF method as a function of the standard deviation of the noise. Case
r = 0.01[K]
r = 0.05[K]
r = 0.1[K]
r = 0.5[K]
r = 1[K]
r = 5[K]
h1 h2 h3 h4
35 34 24 22
30 28 18 17
23 21 15 13
16 15 9 8
10 8 6 5
8 5 4 3
1 ffi. Therefore, its normal derivaiary function becomes Gk¼1 ¼ pffiffiffiffiffiffiffiffiffi 2pr ext k¼1 ¼ 0. tive, evaluated on the external boundary, becomes @G@n
Cext
The trace of the auxiliary function at the internal boundary Cint 1 ffi. Thus, the reciprocity functional is then obtained as ck¼1 ¼ pffiffiffiffiffiffiffiffiffi 2pr ext
defined in Eq. (4) becomes:
1 henv RG;k¼1 ¼ r ext pffiffiffiffiffiffiffiffiffiffiffiffiffi 2pr ext K
Fig. 8. Experimental set-up.
Z 2p h0 ¼0
Ydh0
Z 2p h0 ¼0
T env dh0
ðA1Þ
R 2p According to Eq. (22), we can write h0 ¼0 Ydh0 ¼ 2pa0 being a0 the first coefficient of the Fourier expansion. The environmental temperature, Tenv, is supposed to be constant, thus its integral R 2p could be re-written as h0 ¼0 T env dh0 ¼ 2pT env . According to these results, the reciprocity functional becomes:
1193
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
The trace of the auxiliary function at the internal boundary Cint is then:
ck¼even ¼
sffiffiffiffiffiffiffiffiffiffiffi 1 1 X
r int cos ðmhÞsech m ln prext m¼1 rext
ðA7Þ
Finally, the reciprocity functional, defined in Eq. (4), becomes:
1 henv rint r ext pffiffiffiffiffiffiffipffiffiffiffi sech m ln K r p r ext ext m¼1
Z 2p rint m Ycosðmh0 Þdh0 : sinh m ln rext r ext h0 ¼0
RG;k¼even ¼
Fig. 9. Convective heat transfer coefficient estimated with the TSVD, NRF and FRF approaches.
pffiffiffiffiffiffiffiffiffiffiffiffiffi henv ða0 T env Þ RG;k¼1 ¼ 2pr ext K
ðA2Þ
In agreement with Eq. (5), the integral of the internal heat generation is given as:
Q G;k¼1
qg ¼ K
p r2ext r2int pffiffiffiffiffiffiffiffiffiffiffiffiffi 2prext
hc1 ; cj¼1 i ¼
r int rext
hc1 ; cj–1 i ¼ 0
RG;k¼even ¼
D
m¼1
am
pffiffiffiffiffiffiffiffiffiffiffi henv rint m prext tanh m ln K rext r ext
ðA9Þ
E
ck¼even ; cj¼k ¼
r int 2 sech r ext
m ln
rint rext
E
ck¼even ; cj–k ¼ 0
ðA10:1Þ ðA10:2Þ
ðA4:1Þ
(d) k = 3, 5, 7, . . ., N and k = 2m + 1 For these values of k, according to Eqs. (3.c) and (18), we can qffiffiffiffiffiffi p sin ðmhÞ, and w obtain wk;m–0 ¼ rext k;m¼0 ¼ 0.
ðA4:2Þ
Substituting these results in Eq. (21) the auxiliary function becomes:
(b) k = 2, 4, 6, . . ., N 1 and k – 2m; k = 3, 5, 7, . . ., N and k – 2m +1 It is easy to show that, according to Eqs. (3.b) and (18), both wk;m–0 and wk;m¼0 are null as well as the auxiliary function Gk. Therefore, these frequencies of the orthonormal basis function do not participate in the solution. (c) k = 2, 4, 6, . . ., N 1 and k = 2m For these values of k, according to Eqs. (3.b) and (18), we can qffiffiffiffiffiffi p cos ðmhÞ, and w obtain wk;m–0 ¼ rext k;m¼0 ¼ 0.
Gk¼odd
Therefore, its normal derivative, evaluated at the external boundary, becomes:
sffiffiffiffiffiffiffiffiffiffiffi
1 @Gk¼odd 1 X r ¼ sin ð m h Þsech m ln int @n Cext prext m¼1 r ext
r int m sinh m ln r ext rext
Gk¼even
ck¼odd
ðA5Þ Therefore, its normal derivative, evaluated at the external boundary, becomes:
ðA6Þ
ðA12Þ
The trace of the auxiliary function at the internal boundary Cint is:
sffiffiffiffiffiffiffiffiffiffiffi
1 n h r io 1 X r int int sech m ln ¼ cos ðmhÞ cosh m ln prext m¼1 r r ext
sffiffiffiffiffiffiffiffiffiffiffi
1 n h r io 1 X r int int ¼ sin ½mh cosh m ln sech m ln prext m¼1 r r ext ðA11Þ
Substituting these results in Eq. (19), the auxiliary function becomes:
sffiffiffiffiffiffiffiffiffiffiffi
1 @Gk¼even 1 X r ¼ cos ð m h Þsech m ln int @n Cext prext m¼1 r ext
r int m sinh m ln r ext rext
1 X
agreement with Eq. (5), the integral of the internal heat generation is null, while according to the definition given by Eq. (9) the elements of the matrix even rows that contain the inner product between the trace of the solutions become:
ðA3Þ
while the elements of the first row of the matrix given by Eq. (9), which contains the inner product between the trace of the solutions, becomes:
ðA8Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement at the external boundary Cext can be R 2p written as: h0 ¼0 Ycosðmh0 Þdh0 ¼ pam being am the coefficient of the Fourier series expansion referring to the m-th harmonic. According to this result, the reciprocity functional becomes:
D
1 X
sffiffiffiffiffiffiffiffiffiffiffi
1 1 X r int ¼ sin ðmhÞsech m ln prext m¼1 r ext
ðA13Þ
The reciprocity functional defined in Eq. (4) becomes:
1 henv rint r ext pffiffiffiffiffiffiffipffiffiffiffi sech m ln K r p r ext ext m¼1
Z 2p rint m Ysinðmh0 Þdh0 sinh m ln r ext r ext h0 ¼0
RG;k¼odd ¼
1 X
ðA14Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement at the external boundary Cext can be
1194
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
R 2p written as: h0 ¼0 Ysinðmh0 Þdh0 ¼ pbm being bm the coefficient of the Fourier series expansion referring to the m-th harmonic. According to this result, the reciprocity functional becomes:
RG;k¼odd
1 X pffiffiffiffiffiffiffiffiffiffiffi henv rint m ¼ bm pr ext tanh m ln K rext r ext m¼1
ðA15Þ
In agreement with Eq. (5), the integral of the internal heat generation is null, while according to the definition given by Eq. (9), the elements of the matrix odd rows that contain the inner product between the trace of the solutions becomes:
n1 ; nj–1 ¼ 0
ðB4:2Þ
(b) k = 2, 4, 6, . . ., N 1 and k – 2m; k = 3, 5, 7, . . ., N and k – 2m +1 It is easy to show that, according to Eqs. (3.b), (3.c) and (18), both wk;m–0 and wk;m¼0 are null as well as the auxiliary function Fk. Therefore, these frequencies of the orthonormal basis function do not participate in the solution. (c) k = 2, 4, 6, . . ., N 1 and k = 2m
D E r r ck¼odd ; cj¼k ¼ int sech2 m ln int rext r ext
ðA16:1Þ
For these values of k, according to Eqs. (3.b) and (18), we can qffiffiffiffiffiffi p cos ðmhÞw 1ffipffiffiffi cos ½mh, and wk;m¼0 ¼ 0. obtain wk;m–0 ¼ rext k;m–0 ¼ pffiffiffiffiffi rext :p
hck¼odd ; cj–k i ¼ 0
ðA16:2Þ
Substituting these results in Eq. (46), the auxiliary function becomes:
Appendix B. Derivation of the solution of the second auxiliary problem According to the definition given in Eqs. (3.a)–(3.c), the solution has to be divided in at least three steps for different values of harmonic frequencies k of the orthonormal basis function: (a) k = 1 According to Eqs. (3.a) and (18), for k = 1, we have wm–0;k¼1 ¼ 0, qffiffiffiffiffi and wm¼0;k¼1 ¼ r2extp . Substituting these results in Eq. (46) the auxil1 ffi iary function becomes F k¼1 ¼ pffiffiffiffiffiffiffiffiffi
ln ð
2pr ext ln
Þ . Therefore, its normal r int r r int r ext
derivative, evaluated on the external boundary, becomes 1 ffi r 1 . The trace of the auxiliary function at ¼ pffiffiffiffiffiffiffiffiffi int
@F k¼1 @n
2pr ext ln
Cext
the internal 1 ffi nk¼1 ¼ pffiffiffiffiffiffiffiffiffi
2pr ext ln
rext
r ext
Cint is then obtained as boundary r 1 . Thus, the reciprocity functional defined int r ext
r int
F k¼even ¼
RF;k¼1
ðB2Þ
r ext
while Eq. (31) can be written as:
Q F;k¼a
2 3 rffiffiffiffiffiffiffiffiffiffi qg p 4 2 1 r2ext r 2int 5 ¼ r þ K 2rext ext 2 ln rint r ext
ðB3Þ
n1 ; nj¼1 ¼
2
ln
1
rint rext
rint rext
m ln
r io rint int csch m ln r rext
ðB4:1Þ
ðB5Þ
ðB6Þ
The trace of the auxiliary function at the internal boundary Cint is:
nk¼even
sffiffiffiffiffiffiffiffiffiffiffi
1 1 X r int m ¼ cos ðmhÞcsch m ln prext m¼1 r ext r int
ðB7Þ
Finally, the reciprocity functional defined in Eq. (30) becomes:
sffiffiffiffiffiffiffiffiffiffiffi
1 henv r int csch m ln p r K r ext ext m¼1
Z 2p r int m Ycosðmh0 Þdh0 cosh m ln r ext rext h0 ¼0
RF;k¼even ¼
1 X
r ext
ðB8Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement on the external boundary Cext can R 2p be written as h0 ¼0 Ycosðmh0 Þdh0 ¼ pam being am the coefficient of the Fourier’s series expansion referring to the m-th harmonic. According to this result, the reciprocity functional becomes:
RF ¼
1 X
m¼1
am
pffiffiffiffiffiffiffiffiffiffiffi henv r int m prext coth m ln K r ext rext
ðB9Þ
In agreement with Eq. (31), the integral of the internal heat generation is null, while according to the definition given by Eq. (35), the elements of the matrix even rows that contain the inner product between the trace of the solutions becomes:
hnk¼even ; nj¼k i ¼
m2 r ext r int
2
csch
m ln
rint r ext
hnnk¼even ; nj–k i ¼ 0
According to Eq. (35), the first raw of the matrix of the inner product becomes:
n h
sffiffiffiffiffiffiffiffiffiffiffi
1 @F k¼even 1 X r ¼ cos ð m h Þcsch m ln int @n Cext prext m¼1 r ext
r int m cosh m ln r ext rext
ðB1Þ
20 1 3 pffiffiffiffiffiffiffiffiffiffiffiffiffi henv K 4@1 Aa0 T env 5 ¼ 2pr ext K ln rint r ext henv
p rext
Therefore, its normal derivative, evaluated at the external boundary, becomes:
r ext
Once again, in order to make the above expression fully analytical, we propose, in this paper, to represent the measurements, Y, as a Fourier series expansion. Therefore, can be write R 2p Ydh0 ¼ 2pa0 being a0 the first coefficient of the Fourier’s series h0 ¼0 expansion. As previously observed, the environmental temperature is supposed to be constant, and its integral can be written R 2p as h0 ¼0 T env dh0 ¼ 2pT env . According to these results, the reciprocity functional becomes:
m¼1
cos ðmhÞ sinh
in Eq. (30) becomes:
82 9 3 Z 2p Z = rext <4henv 1 henv 2p 0 0 5 RF;k¼1 ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffi Ydh T env dh r 0 0 : ; K K int 2prext h ¼0 ln r ext h ¼0
rffiffiffiffiffiffiffi 1 X 1 p
ðB10:1Þ ðB10:2Þ
(d) k = 3, 5, 7, . . ., N and k = 2m + 1 For these values of k, according to Eqs. (3.c) and (18), we can qffiffiffiffiffiffi p sin ðmhÞ, and w obtain wk;m–0 ¼ rext k;m¼0 ¼ 0 .Substituting those results in Eq. (46) the auxiliary function becomes:
A. Mocerino et al. / International Journal of Heat and Mass Transfer 125 (2018) 1181–1195
F k¼odd
sffiffiffiffiffiffiffiffiffiffiffi
1 n h r io 1 X r int int ¼ sin ðmhÞ sinh m ln csch m ln prext m¼1 r r ext ðB11Þ
Therefore, its normal derivative, evaluated at the external boundary, becomes: sffiffiffiffiffiffiffiffiffiffiffi
1 @F k¼odd 1 X r int r int m ¼ sin ðmhÞcsch m ln cosh m ln @n Cext prext m¼1 r ext r ext r ext ðB12Þ
The trace of the auxiliary function at the internal boundary Cint is:
nk¼odd
sffiffiffiffiffiffiffiffiffiffiffi
1 1 X rint m ¼ sin ðmhÞcsch m ln prext m¼1 rext r int
ðB13Þ
The reciprocity functional defined in Eq. (30) becomes:
RF;k¼odd
sffiffiffiffiffiffiffiffiffiffiffi
1 henv r int ¼ rext csch m ln prext K r ext m¼1
Z 2p r int m cosh m ln Ysinðmh0 Þdh0 r ext r ext h0 ¼0 1 X
ðB14Þ
Following the same procedure as before, we propose to represent the measurements Y as a Fourier expansion. Therefore, the integral of the measurement at the external boundary Cext in the R 2p last expression can be written as h0 ¼0 Y sin ðmh0 Þdh0 ¼ pbm being bm the coefficient of the Fourier’s series expansion referring to the mth harmonic. According to this result, the reciprocity functional becomes:
RF;k¼odd ¼
1 X
m¼1
bm
pffiffiffiffiffiffiffiffiffiffiffi henv rint m prext coth m ln K rext r ext
ðB15Þ
In agreement with Eq. (31), the integral of the internal heat generation is null, while according to the definition given by Eq. (35), the elements of the matrix odd rows that contain the inner product between the trace of the solutions becomes:
hnk¼odd ; nj¼k i ¼
m r int r ext
hnk¼odd ; nj–k i ¼ 0
csch
2
r m ln int r ext
ðB16:1Þ ðB16:2Þ
References [1] R. Webb, N.H. Kim, Principles of Enhanced Heat Transfer, second ed., Taylor & Francis, New York, NY, 2005.
1195
[2] A.E. Bergles, Techniques to enhance heat transfer, in: Handbook of Heat Transfer, McGraw-Hill, New-York, NY, 1998. [3] W.H. Giedt, Investigation of variation of point unit heat-transfer coefficient around a cylinder normal to an air stream, ASME Trans. 71 (1949) 375–381. [4] S. Aiba, Y. Yamazaki, An experimental investigation of heat transfer around a tube in a bank, J. Heat Transf. 98 (3) (1976) 503–508. [5] J. Taler, Nonlinear steady-state inverse heat conduction problem with spacevariable boundary conditions, J. Heat Transf. 114 (4) (1992) 1048–1051. [6] J. Taler, Determination of local heat transfer coefficient from the solution of the inverse heat conduction problem, Forsch. Ingenieurwes. 71 (2) (2007) 69–78. [7] F. Bozzoli, G. Pagliarini, S. Rainieri, Experimental validation of the filtering technique approach applied to the restoration of the heat source field, Exp. Therm. Fluid Sci. 44 (2013) 858–867. [8] F. Bozzoli, L. Cattani, S. Rainieri, F.S.V. Bazán, L.S. Borges, Estimation of the local heat-transfer coefficient in the laminar flow regime in coiled tubes by the Tikhonov regularisation method, Int. J. Heat Mass Transf. 72 (2014) 352–361. [9] G.M. Carlomagno, G. Cardone, Infrared thermography for convective heat transfer measurements, Exp. Fluids 49 (2010) 1187–1218. [10] J.V. Beck, B. Backwell, C.R.S. Clair, Inverse Heat Conduction – Ill-posed Problems, Wiley-Interscience, New-York, NY, 1985. [11] H.R.B. Orlande, F. Olivier, D. Maillet, R.M. Cotta, Thermal Measurements and Inverse Techniques, Taylor and Francis, New York, NY, 2011. [12] S. Andrieux, A.B. Abda, The reciprocity gap: a general concept for flaws identification problems, Mech. Res. Commun. 20 (1993) 415–420. [13] S. Andrieux, A.B. Abda, Identification of planar cracks by complete overdetermined data: Inversion formulae, Inverse Prob. 12 (1996) 553–563. [14] F. Delbary, R. Aramini, G. Bozza, M. Brignone, M. Piana, On the use of the reciprocity gap functional in inverse scattering with near-field data: an application to mammography, J. Phys. Conf. Ser. 135 (2008). [15] M.J. Colaço, C.J.S. Alves, A fast non-intrusive method for estimating spatial thermal contact conductance by means of the reciprocity functional approach and the method of fundamental solutions, Int. J. Heat Mass Transf. 60 (2013) 653–663. [16] E.I. Shifrin, P.S. Shushpannikov, Identification of a spheroidal defect in an elastic solid using a reciprocity gap functional, Inverse Prob. 26 (5) (2010) 055001. [17] M.J. Colaço, C.J. Alves, H.R. Orlande, Transient non-intrusive method for estimating spatial thermal contact conductance by means of the reciprocity functional approach and the method of fundamental solutions, Inverse Prob. Sci. Eng. 23 (4) (2015) 688–717. [18] M.J. Colaço, C.J. Alves, F. Bozzoli, The reciprocity function approach applied to the non-intrusive estimation of spatially varying internal heat transfer coefficients in ducts: numerical and experimental results, Int. J. Heat Mass Transf. 90 (2015) 1221–1231. [19] M. Necati Ozisik, Heat Conduction, second ed., John Wiley & Sons, Inc., 1993. [20] P.C. Hansen, The truncated SVD as a method for regularization, BIT Numerical Math. 27 (4) (1987) 534–553. [21] A.N. Tikhonov, V.Y. Arsenin, Solution of Ill-Posed Problems, Winston & Sons, Washington DC, USA, 1977. [22] F. Bozzoli, S. Rainieri, G. Pagliarini, Experimental validation of the filtering technique approach applied to the restoration of the heat source field, Exp. Therm. Fluid Sci. 44 (2013) 858–867. [23] James V. Beck, Kenneth J. Arnold, Parameter Estimation in Engineering and Science, John Wiley & Sons, 1997. [24] Fermín S.V. Bazán, J.B. Francisco, Koung Hee Leem, G. Pelekanos, A maximum product criterion as a Tikhonov parameter choice rule for Kirsch’s factorization method, J. Comput. Appl. Math. 236 (2012) 4264–4275.