|
|
Article October 2014 Vol. 57 No. 10: 1–8
doi: 10.1007/s11433-014-5469-2
Interval analysis of transient temperature field with uncertain-but-bounded parameters
Institute of Solid Mechanics, Beihang University, Beijing 100191, China
Received October 22, 2013; accepted March 13, 2014
Based on the traditional finite volume method, a new numerical technique is presented for the transient temperature field
prediction with interval uncertainties in both the physical parameters and initial/boundary conditions. New stability theory
applicable to interval discrete schemes is developed. Interval ranges of the uncertain temperature field can be approximately
yielded by two kinds of parameter perturbation methods. Different order Neumann series are adopted to approximate the interval
matrix inverse. By comparing the results with traditional Monte Carlo simulation, a numerical example is given to demonstrate
the feasibility and effectiveness of the proposed model and methods.
finite volume method, temperature field prediction, interval uncertainties, stability theory, parameter perturbation method
PACS number(s):
Citation:
1 Introduction
Thermal analysis has undergone a rapid development in engineering, especially in the field of aeronautics and astronautics,
where the coupled interaction of structure and heat conduction is playing a more and more significant role. Traditional thermal
analysis has been conducted under the assumption that the physical properties and initial/boundary conditions are deterministic.
But in the actual engineering, due to the model inaccuracies, physical imperfections and system complexities, uncertainties
in material properties, geometric dimensions and initial/boundary conditions are unavoidable, which will lead to the uncertainty
of the structural temperature field [1–3]. Popular approaches to these uncertain problems are probabilistic methods, where
the probability density functions are defined unambiguously. Hien and Kleiber suggested the stochastic variational principle
and stochastic finite element method for transient heat problem with random parameters [4]. By combining boundary element
method with orthogonal expansion theory, a new numerical method was proposed for solving stochastic heat transfer problems
in the literature [5]. On the basis of generalized polynomial chaos, Xiu presented a random spectral decomposition method
for the solution of transient heat conduction subjected to random inputs [6].
The fuzzy set theory, introduced by Zadeh [7] in 1965, is another efficient category to represent the uncertain parameter
whose subjective probability based on the expert opinions is available [8]. Unfortunately, for many engineering problems,
it’s often too difficult or costly to collect enough information about the uncertainty to construct the probability distribution
function and membership function, which leads to a certain extent limit applicability of the stochastic and fuzzy methods.
In overcoming the shortcomings, some non-probabilistic approaches such as convex model and interval analysis are worked out
to deal with the uncertain problems without sufficient information [9,10]. Neumaier investigated the hypercube approximation
for the united solution set of the interval equations by the Gaussian elimination scheme [11], which will be extremely conservative
due to a large number of elimination operations [12]. The exact response intervals can be achieved by the vertex method using
all possible combinations of the interval parameters [13]. However, the computational effort by this method increases exponentially
with the increase in the number of uncertain parameters. Compared with above interval approaches, the interval perturbation
method, proposed by Qiu et al. [14,15], has been widely applied in the structural response analysis out of its simplicity
and efficiency [16,17]. Compared with other interval analysis approaches, the computational cost of the interval parameter
perturbation method is smaller, and the condition of convergence which is related to the ranges of interval parameters is
more easily guaranteed. Nevertheless, the unpredictable effect of neglecting the higher order term of Neumann series is the
main disadvantage of traditional perturbation methods.
Hence, establishing effective numerical analyzing methods for the transient heat conduction problem with interval parameters
has great significance in both theory and engineering. Based on the thermal finite volume method [18] and interval definition,
the interval equilibrium equation for the transient temperature field prediction is established in sect. 2. Subsequently,
new stability theory applicable to interval discrete schemes is developed. The first interval analysis method named IPFVM
is presented in sect. 4, where the interval matrix and vector are expanded by the first-order Taylor series and the interval
matrix inverse is approximated by the first-order Neumann series. Considering higher order terms in the Neumann expansion,
sect. 5 proposes a modified parameter perturbation method to approximately predict the interval transient temperature field.
Numerical results on a 2D thermal plate are given to verify the accuracy of the proposed methods in sect. 6, and we draw the
conclusions in the last section.
2 Transient heat conduction problem
The governing equation of one-dimensional transient heat conduction problem with a heat source can be written as:
where u(x, t) denotes the temperature field; r, c are the density, specific heat capacity of the material respectively; k(x) is the thermal conductivity related to the position x; f(x) stands for the intensity of the heat source.
For an interior domain bounded by G, an initial condition and three sets of boundary conditions exist as shown in Figure 1.
Initial condition: the distribution law j (x) of initial temperature field is known, i.e.
Dirichlet boundary condition: the temperature u on D is given.
Neumann boundary condition: the heat flux q on N is determined.
where n is the normal vector of the boundary.
Robin boundary condition: the ambient temperature uf and heat transfer coefficient h are known
In the standard finite volume method (FVM), space step D x and time step D t are adopted to discrete the spatial domain
and time domain, and the temperature at grid node x j and time t k is denoted as u j k .
Figure 1 Three kinds of boundary conditions.
A transition variable W ( x ) = k ( x ) & # x2202; u ( x , t ) & # x2202; x & # x2202; u ( x , t ) & # x2202; x = W ( x ) k ( x ) , [ x i & # x2212; 1 , x i ]
In terms of the common integral rectangular formula, we can derive
Suppose that the temperature and its derivative in time step t are constant and equal to the values at the forward time t=kt. Eq. (6) can be simplified to be
which is called the forward discrete scheme and can be rewritten as:
Similarly, if the values at the backward time t=(k+1)t are adopted in the time step t, we can obtain the so-called backward discrete scheme:
By extracting the arithmetic mean of above two equations, the central scheme is easily derived
Based on integral interpolation, the original differential equation (1) can be transformed into the conservation format as
shown in eq. (6), which drops the differential terms from the second order to the first order and weakens the smoothness requirement
of the term k(x). Besides, the discrete format derived from the integral conservation equation is easier to be extended to any type of meshes
and boundary conditions, which is the unique advantage of FVM.
By means of matrices and vectors, the commonly used two-layer finite volume formats can be expressed as:
where Uk stands for the temperature field at time layer t=(k+1)t; A, B are the coefficient matrices; F is the equivalent nodal heat flow vector.
For the actual engineering problem, due to the limitations of production and measurement error, uncertainties in material
properties, external loads and boundary conditions are unavoidable. In this paper, the uncertainties whose lower and upper
bounds can be determined by the limited information are quantitatively described as the interval parameters and assumed to
belong to an uncertain-but-bounded parameter vector as follows:
where & # x03b1; _ i & # x03b1; & # x00af; i & # x03b1; i c = ( & # x03b1; & # x00af; i + & # x03b1; _ i ) / 2 & # x0394; & # x03b1; i = ( & # x03b1; & # x00af; i & # x2212; & # x03b1; _ i ) / 2 & # x03b4; I = [ & # x2212; 1 , 1 ]
Therefore, the coefficient matrices A, B and load vector F become interval matrices and vector with respect to the interval parameter vector aI. The recurrence equation (13) can be rewritten as:
Theoretical solution set of eq. (15) is defined as:
where & # x03a9; k + 1 ( U k + 1 ) I = [ U _ k + 1 , U & # x00af; k + 1 ] , & # x03a9; k + 1 .
Similarly, by extending the coefficient matrices and vectors, the same discrete equation also can be derived to high-dimensional
interval transient heat conduction problems.
3 Interval stability analysis
For a discrete equation with an initial condition, if the error introduced in one time layer cannot be constantly amplified
in the later time, then the discrete format is claimed to be stable. Besides the interval parameter vector a I , A ( a I
), B ( a I ) and F ( a I ) in eq. (17) are also dependent on the time step D t and space step D x . Suppose that for
any & # x03b1; & # x2208; & # x03b1; I ,
Accordingly, given the initial condition U 0 = & # x03a8; ,
In this paper, the error e 0 in the initial condition is taken into account. Denoting the corresponding solution by U & # x02dc;
By using ek to denote the error in the (k+1)th time layer and subtracting eq. (20) from eq. (19), we obtain
If there exists a constant K >0, for any U 0 & # x2208; R N , & # x03b1; & # x2208; & # x03b1; I ,
and then the discrete format is considered to be stable.
In other words, the necessary and sufficient condition to ensure the stability is
which means the interval matrix family
is uniformly bounded.
Based on the traditional von Neumann qualification, the spectral radius of the coefficient matrix C ( a , D t , D x ) is
required to be less than 1 for any & # x03b1; & # x2208; & # x03b1; I ,
which is called the interval von Neumann condition suitable for stability analysis of interval finite volume schemes.
For instance, in order to ensure the stability of eq. (10)tx
When the interval uncertainties in system parameters k, r, c are considered, the corresponding interval stability condition turns to be
4 Interval perturbation finite volume method (IPFVM)
As we know, the first-order Taylor expansion for the linear function is accurate, whose accuracy for the non-linear function
is also acceptable if the ranges of the interval parameters are narrow enough [19]. In this paper, we just discuss the problem
of the temperature field prediction with narrow uncertain-bout-bounded parameters. Based on the first-order Taylor series
expansion, the coefficient matrices and load vector in eq. (17) can be developed as:
By substituting eq. (28) into eq. (17), we can obtain
where (Uk) and D(Uk)I stand for the midpoint and deviation interval of the (k+1)th layer temperature field, respectively.
By pre-multiplying both sides of Eq. (29) with (Ac+ AI), one gets
where
If the spectral radius of (Ac)A is less than 1, (Ac+ AI) can be expanded by the Neumann series
By substituting eq. (32) into eq. (30), one gets
Without considering the higher order terms, we can derive the iterative formulas
Furthermore, by substituting Eqs. (28) and (31) into Eq. (34) and neglecting the cross small terms, one can specify above
formulas to be
It is obvious that D ( U k +1 ) I is a monotonic function with respect to & # x03b4; i I & # x03b4; u k I
where | & # x2022; |
Hence, the upper and lower bounds of the transient temperature field can be expressed as:
5 Modified interval perturbation finite volume method (MIPFVM)
The unpredictable effect of neglecting the higher order terms in the interval matrix inverse is the inherent disadvantage
of the traditional interval perturbation method. Based on the work of Impollonia [10], in this section, a modified Neumann
expansion is proposed to calculate the inversion of the non-singular interval matrix.
By substituting eq. (28) into eq. (32), one can get
where A i = & # x2202; A ( & # x03b1; c ) & # x2202; & # x03b1; i & # x0394; & # x03b1; i ( A c ) & # x2212; 1
where the expression & # x00ac; i = j
According to eq. (40), eq. (39) can be converted into
If the norm & # x2016; & # x0394; & # x03b1; i A i & # x2016; < 1 & # x2211; r = 1 & # x221e; ( & # x2212; & # x03b4; i I A i ) r
where E i I = & # x2212; & # x03b4; i I A i I + & # x03b4; i I A i ;
Based on the algorithms of interval numbers, it is easy to get the calculating formulas of median and radius
By substituting eqs. (42) and (43) into eq. (30), we can obtain
Furthermore, by substituting eqs. (28) and (31) into eq. (44) and neglecting the cross small terms, one gets
If & # x03b4; e i I & # x03b4; u k I & # x03b4; i I & # x0394; ( U k + 1 ) I & # x03b4; e i I & # x03b4; u k I & # x03b4; i I ( U k + 1 ) I
The upper and lower bounds of the transient temperature field can be calculated by
It is obvious that the interval analyses of IPFVM and the inversion of the non-singular interval matrix in IPFVM is approximated
by the first-order Neumann expansion, while the inversion in MIPFVM is calculated by retaining part of higher order terms.
Therefore, the accuracy of the temperature bounds obtained by MIPFVM is higher than that calculated by IPFVM. Meanwhile, additional
numerical burden suffering from matrix inverse will be introduced in MIPFVM.
6 Numerical example
In order to verify the performance of the proposed methods, consider a complex thermal plate with the thickness of 0.1 m as
depicted in Figure 2. The circular part and the rectangular-sectioned part are discretized by 188 triangular elements and
100 quadrilateral elements, respectively. A volumetric heat changing from to is generated in the shaded portion of the plate,
while a heat flux changing from TDirichlet boundary at the left side is also regarded as an interval parameter with the range from ambient temperature range
from initial temperature is between heat conductivity, specific heat capacity, heat transfer coefficient and density are considered
interval parameters kchr.
In the conventional interval analysis, different parameters are supposed to be independent of each other [20]. Time step t=5 s is selected to ensure the stability of the interval discrete scheme in this numerical example. On this basis, simulations
of IPFVM and MIPFVM for the transient heat conduction problem are carried out by Matlab R2012 on a
Figure 2 Finite volume mesh model of the complex thermal plate.
3.00 GHz Petium(R) 4 CPU computer. The lower bound (LB) and upper bound (UB) of the interval temperature responses at two
feature points are listed in Tables 1 and 2. Assuming that the uncertain parameters satisfy the uniform distribution in the
given intervals, the response intervals obtained by Monte Carlo method (MCM) with 10 samples are used as referenced results
for validating the accuracies and efficiencies of the proposed perturbation methods. It is obvious that the accuracy of MIPFVM
in which the higher order terms are retained is much higher than that of IPFVM approximated by the first-order Neumann expansion.
From the computational execution time listed in the tables, we can see that the execution time of MIPFVM is a little longer,
and the exceeded time is implemented to calculate the higher order terms.
Time evolutions of the interval temperature responses at the feature points are plotted in Figure 3 and 4. Compared with the
mean values, it is clear that the uncertainties in structural parameters and initial/boundary conditions have a profound impact
on the transient temperature field. Besides, we can easily see that the temperature intervals obtained by MIPFVM are much
smaller than those obtained by IPFVM, and match the referenced ranges yielded by
Table 1 Interval temperature response at feature point 1
|
Bounds
|
|
MCM (102.4s)
(C)
|
|
IPFVM (5.52s)
|
|
MIPFVM (6.73s)
|
|
Value (C)
|
|
Error (%)
|
|
Value (C)
|
|
Error (%)
|
|
100
|
|
LB
|
|
27.27
|
|
25.19
|
|
7.63
|
|
26.62
|
|
2.38
|
|
UB
|
|
31.31
|
|
33.84
|
|
8.08
|
|
31.82
|
|
1.63
|
|
200
|
|
LB
|
|
32.73
|
|
30.54
|
|
6.69
|
|
32.02
|
|
2.17
|
|
UB
|
|
36.78
|
|
39.11
|
|
6.33
|
|
37.33
|
|
1.50
|
|
500
|
|
LB
|
|
39.45
|
|
37.55
|
|
4.82
|
|
38.73
|
|
1.83
|
|
UB
|
|
43.39
|
|
45.53
|
|
4.93
|
|
44.00
|
|
1.41
|
|
1000
|
|
LB
|
|
44.25
|
|
42.39
|
|
4.20
|
|
43.43
|
|
1.85
|
|
UB
|
|
48.44
|
|
50.97
|
|
5.22
|
|
49.22
|
|
1.61
|
|
2000
|
|
LB
|
|
49.85
|
|
47.59
|
|
4.53
|
|
48.75
|
|
2.21
|
|
UB
|
|
54.66
|
|
57.61
|
|
5.40
|
|
55.69
|
|
1.88
|
|
3000
|
|
LB
|
|
53.69
|
|
51.20
|
|
4.64
|
|
52.44
|
|
2.33
|
|
UB
|
|
58.68
|
|
61.34
|
|
4.53
|
|
59.80
|
|
1.91
|
Table 2 Interval temperature response at feature point 2
|
Bounds
|
|
MCM (102.4s)
(C)
|
|
IPFVM (5.52s)
|
|
MIPFVM (6.73s)
|
|
Value (C)
|
|
Error (%)
|
|
Value (C)
|
|
Error (%)
|
|
100
|
|
LB
|
|
19.60
|
|
19.46
|
|
0.71
|
|
19.53
|
|
0.36
|
|
UB
|
|
21.81
|
|
22.21
|
|
1.83
|
|
21.90
|
|
0.41
|
|
200
|
|
LB
|
|
20.45
|
|
20.00
|
|
2.20
|
|
20.24
|
|
1.03
|
|
UB
|
|
23.19
|
|
24.58
|
|
5.99
|
|
23.45
|
|
1.12
|
|
500
|
|
LB
|
|
24.94
|
|
22.98
|
|
7.86
|
|
24.39
|
|
2.21
|
|
UB
|
|
29.44
|
|
31.80
|
|
8.02
|
|
30.23
|
|
2.68
|
|
1000
|
|
LB
|
|
32.76
|
|
30.29
|
|
7.54
|
|
31.98
|
|
2.38
|
|
UB
|
|
38.95
|
|
41.51
|
|
6.57
|
|
40.02
|
|
2.75
|
|
2000
|
|
LB
|
|
44.58
|
|
41.53
|
|
6.84
|
|
43.30
|
|
2.87
|
|
UB
|
|
52.38
|
|
55.98
|
|
6.87
|
|
53.82
|
|
2.75
|
|
3000
|
|
LB
|
|
52.95
|
|
49.35
|
|
6.80
|
|
51.52
|
|
2.70
|
|
UB
|
|
61.15
|
|
65.13
|
|
6.51
|
|
62.69
|
|
2.52
|
Figure 3 Time evolution of temperature response at feature point 1.
Figure 4 Time evolution of temperature response at feature point 2.
MCM perfectly. The numerical results obtained by MIPFVM may be conservative compared with MCM. However, the differences are
so small that the proposed MIPFVM is completely reliable and greatly reduces the huge computational cost of MCM.
7 Conclusions
By combining the interval analysis theory with the finite volume method, this paper proposes a new technique named interval
finite volume method to predict the transient temperature field with uncertain-but-bounded parameters. The interval uncertainties
in structural parameters and initial/boundary conditions are all fully considered, which makes the calculation more objective.
A new stability analysis theory suitable for transient interval discrete schemes is developed. On the basis of parameter perturbation
theory, two interval analysis methods, named IPFVM and MIPFVM, are presented to predict the upper and lower bounds of the
temperature field effectively. Compared with IPFVM, MIPFVM obtains more accurate results by retaining higher order terms and
greatly reduces the huge computational cost of Monte Carlo simulation. Numerical results on a 2D thermal plate fully validate
the feasibility and superiority of the suggested model and method to deal with heat conduction problems with many kinds of
interval uncertainties. It should be noted that the additional computational burden to calculate matrix inverses is the disadvantage
of MIPFVM, especially for the case with a large number of uncertain parameters. Besides, the proposed perturbation methods
are more suitable for solving the problems with narrow parameter intervals due to the unpredictable effect of neglecting the
higher order terms of Taylor series.
This work was supported by the National Special Fund for Major Research Instrument Development (Grant No. 2011YQ140145), 111
Project (Grant No. B07009), National Natural Science Foundation of China (Grant No. 11002013) and Defense Industrial Technology
Development Program (Grant Nos. A2120110001 and B2120110011).
Efe M O, Ozbay H. Multi input dynamical modeling of heat flow with uncertain diffusivity parameter. Math Comp Model Dyn, 2003,
9: 437–450
Nicolai B M, Egea J A, Scheerlinck N. Fuzzy finite element analysis of heat conduction problems with uncertain parameters.
Int J Food Eng, 2011, 103: 38–46
Lu Q. Some results on the controllability of forward stochastic heat equations with control on the drift. J Funct Anal, 2011,
260: 832–851
Hien T D, Kleiber M. Stochastic finite element modeling in linear transient heat transfer. Comput Methods Appl Mech Eng, 1997,
144: 111–124
Emery A F. Solving stochastic heat transfer problems. Eng Anal Bound Elem, 2004, 28: 279–291
Xiu D B. A new stochastic approach to transient heat conduction modeling with uncertainty. Int J Heat Mass Tran, 2003, 46:
4681–4693
Zadeh L A. Fuzzy sets. Inf Control, 1965, 8: 338–353
Hanss M, Turrin S. A fuzzy-based approach to comprehensive modeling and analysis of systems with epistemic uncertainties.
Struct Saf, 2010, 32: 433–441
Ben-Haim Y, Elishakoff I. Convex Models of Uncertainty in Applied Mechanics. Amsterdam: Elsevier Science Publishers, 1990
Impollonia N, Muscolino G. Interval analysis of structures with uncertain-but-bounded axial stiffness. Comput Methods Appl
Mech Eng, 2011, 220: 1945–1962
Neumaier A. Interval Methods for Systems of Equations. Cambridge: Cambridge University Press, 1990
Berkel L, Rao S S. Analysis of uncertain structure systems using interval analysis. AIAA J, 1997, 35: 727–773
Qiu Z P, Xia Y Y, Yang J. The static displacement and the stress analysis of structures with bounded uncertainties using the
vertex solution theorem. Comput Methods Appl Mech Eng, 2007, 196: 4965–4984
Qiu Z P, Chen S H, Elishakoff I. Bounds of eigenvalues for structures with an interval description of uncertain-but-non-random
parameters. Chaos Soliton Fract, 1996, 7(3): 425–434
Qiu Z P, Elishakoff I. Anti-optimization of structures with large uncertain-but-non-random parameters via interval analysis.
Comput Methods Appl Mech Eng, 1998, 152: 361–372
Sliva G, Brezillon A, Cadou J M, et al. A study of the eigenvalue sensitivity by homotopy and perturbation methods. J Comput
Appl Math, 2010, 234: 2297–2302
Xia B Z, Yu D J. Modified sub-interval perturbation finite element method for 2D acoustic field prediction with large uncertain-but-bounded
parameters. J Sound Vib, 2012, 331: 3774–3790
Versteeg H K, Malalsekera W. An Introduction to Computational Fluid Dynamics: The Finite Volume Method. Essex: Longman Scientific
& Technical, 1995
Qiu Z P, Wang X J. Parameter perturbation method for dynamic response of structures with uncertain-but-bounded parameters
based on interval analysis. Int J Solids Struct, 2005, 42(18-19): 4958–4970
Degrauwe D, Lombaert G, De Roeck G. Improving interval analysis in finite element calculations by means of affine arithmetic.
Comput Struct, 2010, 88: 247–254