Heat Mass Transfer (2017) 53:1657–1666 DOI 10.1007/s00231-016-1923-1

ORIGINAL

Identification of the heat transfer coefficient in the two‑dimensional model of binary alloy solidification Edyta Hetmaniok1 · Jordan Hristov2 · Damian Słota1 · Adam Zielonka1 

Received: 7 November 2015 / Accepted: 4 October 2016 / Published online: 14 October 2016 © Springer-Verlag Berlin Heidelberg 2016

Abstract  The paper presents the procedure for solving the inverse problem for the binary alloy solidification in a twodimensional space. This is a continuation of some previous works of the authors investigating a similar problem but in the one-dimensional domain. Goal of the problem consists in identification of the heat transfer coefficient on boundary of the region and in reconstruction of the temperature distribution inside the considered region in case when the temperature measurements in selected points of the alloy are known. Mathematical model of the problem is based on the heat conduction equation with the substitute thermal capacity and with the liquidus and solidus temperatures varying in dependance on the concentration of the alloy component. For describing this concentration the Scheil model is used. Investigated procedure involves also the parallelized Ant Colony Optimization algorithm applied for minimizing a functional expressing the error of approximate solution.

* Edyta Hetmaniok [email protected] Jordan Hristov [email protected] Damian Słota [email protected] Adam Zielonka [email protected] 1

Institute of Mathematics, Silesian University of Technology, Kaszubska 23, 44‑100 Gliwice, Poland

2

Department of Chemical Engineering, University of Chemical Technology and Metallurgy, 8 Kliment Ohridsky blvd, Sofia 1756, Bulgaria



1 Introduction In technical modeling many identification, control or design problems take the form of the inverse problems characterized by the incomplete conditions defining them [1–3]. From engineering point of view the inverse problems are of great importance since they give possibilities to control various processes. In this paper we concentrate on the inverse problem connected with the process of binary alloy solidification [4–8]. The process of solidification can proceed in the constant temperature or in the temperature interval. If the solidification runs in the constant temperature, then we consider the so called Stefan problem or the problem of solidification with zero interval of solidification temperatures. In the Stefan problem the liquid phase is sharply separated from the solidified phase. Both phases are in contact by creating the solidification surface. Solidification of alloys proceeds in the so called solidification temperature intervals. In this case there is no sharp separation between the liquid phase and solid phase. Both phases are separated by the twophase zone, called the mushy zone, in which the liquid and solidified phases coexist. In a single component solidification the mushy zone can also appear if the non-equilibrium conditions occur at the interface. Width of the two-phase zone depends on the chemical composition of the solidified alloy and on the solidification velocity. Existence of the two-phase zone influences certainly the microstructure of the solidified metal. Solidification of the alloys is also influenced by the process of the alloy components segregation revealed for example in the dependance between the concentration of the alloy component and the values of liquidus and solidus temperatures. There are developed the following mathematical models, the most often used for describing the concentration: the

13

1658

model based on the lever arm rule and the Scheil model. In the first model the immediate equalization of chemical composition of the alloy in the liquid phase and solid phase is assumed [9–11], whereas the second one is based on the fact that diffusion coefficient in the solid phase is significantly lower than in the liquid phase. This leads to the assumption that the diffusion in the solid phase does not occur [9, 12– 15]. In this paper for describing the concentration we apply the Scheil model, whereas in paper [16] the similar problem was solved with the aid of the lever arm model. Thus, the current paper is a discussion with the results obtained in paper [16] and a trial of answering the question whether any of these two models can be judged as a better one for solving the problem of considered kind. Mathematical model used for describing the solidification in temperature interval is based on the heat conduction equation with the source element added, in which the latent heat of fusion and the volume contribution of solid phase are included. By assuming the function describing the solid mass fraction in some special form, we transform the heat conduction equation into the equation with the so called substitute thermal capacity and with the liquidus and solidus temperatures varying in dependance on the concentration of alloy component [17–19]. Such transformed differential equation defines the heat conduction in the entire domain (in solid phase, in two-phase zone (mushy zone) and in liquid phase). Aim of the paper is to solve the inverse problem of binary alloy solidification consisted in identification of the heat transfer coefficient occurring in the boundary condition defined on one of the boundaries and in reconstruction of the temperature distribution in the entire two-dimensional domain. Additional information, needed to solve such task, in our case is given by the measurements of temperature read in selected points of the domain. The current paper is a continuation of research undertaken in papers [20–22] where the similar problem was solved, but in one dimensional domain. Transformation of the discussed inverse problem to the two-dimensional space significantly increased the computational complexity of the procedure which resulted in a great extension of the computing time. To execute the calculations in a reasonable time, we decided to parallelize the solving procedure. Similarly as in the previous papers made by the authors (see for example [20–24]), the essential part of the procedure is based on the artificial intelligence optimization algorithm used for minimizing a functional expressing the differences between the calculated and measured values of temperature. Therefore we needed to choose an algorithm giving the possibility to parallelize the computations and we selected the Ant Colony Optimization algorithm, whereas, for example in [20– 22], the Artificial Bee Colony algorithm was used. The Ant Colony Optimisation (ACO) algorithm is the non-classical

13

Heat Mass Transfer (2017) 53:1657–1666

optimization algorithm belonging to the group of swarm intelligence algorithms inspired by the observations of creatures living in nature. The ACO algorithm imitates the behavior of ants looking for the food around the ant-hill and communicating between the swarm members with the aid of chemical substance left on the traversed paths [25– 28]. ACO algorithm has found many applications in solving a number of optimisation tasks because of its simplicity and relatively short time of working. The second essential method, the used approach is based on, is the generalized alternating phase truncation method dedicated for determining the approximate values of temperature. This method serves for solving the multiphase direct heat conduction problems and its idea lies in the proper transformation of multi-phase domain such that each step of the procedure requires few solutions of onephase direct problems - so many as many phases we consider [29–32].

2 Two‑dimensional problem formulation Let us consider the solidifying material occupying region Ω = [0, b] × [0, d] and varying in time. Simplified version of this region for the selected moment of time ¯t is shown in Fig. 1. In this figure two subregions taken by the liquid and solid phase, respectively, are separated by the intermediate two-phase zone, that is the mushy zone. According to the figure, boundary of region Ω = Ω × [0, t ∗ ] is divided into the following parts

Γ0 = {(x, y, 0); x ∈ [0, b], y ∈ [0, d]},   Γ1 = (0, y, t); y ∈ [0, d], t ∈ [0, t ∗ ] ,   Γ2 = (x, 0, t); x ∈ [0, b], t ∈ [0, t ∗ ] ,   Γ3 = (b, y, t); y ∈ [0, d], t ∈ [0, t ∗ ] ,   Γ4 = (x, d, t); x ∈ [0, b], t ∈ [0, t ∗ ] .

(1)

In region Ω the distribution of temperature, described by means of function T, is modeled with the aid of heat conduction equation [17–19, 29]:

c̺

∂fs (x, y, t) ∂T (x, y, t) =  ∇ 2 T (x, y, t), +L ̺ ∂t ∂t

(2)

where c, ̺ and  are the specific heat, mass density and thermal conductivity coefficient, respectively, L denotes the latent heat of solidification, fs describes the volumetric solid state fraction, T is the temperature, t means the time and finally x and y refer to the spatial coordinates. The volumetric solid state fraction depends on the temperature, thus we present it in the form

∂fs ∂T ∂fs . = ∂t ∂T ∂t

(3)

Heat Mass Transfer (2017) 53:1657–1666

1659

y

Finally, on boundaries Γ3 and Γ4 condition of the third kind has to be fulfilled

4

d

−

mushy zone

3

solidus

liquidus liquid phase

2

b

x

Fig. 1  Domain of the problem for the selected moment of time ¯t (Γ¯i = Γi ∩ {¯t })

Substituting relation (3) to equation (2), after some transformations we get   ∂T (x, y, t) ∂fs =  ∇ 2 T (x, y, t). ̺ c−L (4) ∂T ∂t Defining the so called substitute thermal capacity

C =c−L

∂fs , ∂T

(5)

equation (4) can be written in the form

C̺

∂ T (x, y, t) =  ∇ 2 T (x, y, t). ∂t

(6)

The above equation must be completed by the initial and boundary conditions defined on the respective parts of boundary described by sets (1). Thus, on boundary Γ0 the following initial condition must be satisfied

T (x, y, 0) = T0 ,

(7)

where T0 denotes the initial temperature. On boundaries Γ1 and Γ2 the homogeneous condition of the second kind is defined

−

∂ T (x, y, t) = 0. ∂n

(9)

where α denotes the heat transfer coefficient and T∞ expresses the ambient temperature. The substitute thermal capacity is variable in dependence on temperature and is defined as  T > TL (ZL ), c   l L T ∈ [TS (ZL ), TL (ZL )], C = cmz +  TL (ZL ) − TS (ZL )  cs T < TS (ZL ), (10)

solid phase

1

∂ T (x, y, t) = α(x, y, t) (T (x, y, t) − T∞ ), ∂n

(8)

where cl, cmz and cs denote, respectively, the specific heat of liquid phase, mushy zone and solid phase(value cmz is calculated as the arithmetic average of values cl and cs), L describes the latent heat of fusion and finally TL (ZL ) and TS (ZL ) refer to the liquidus and solidus temperatures varying while the concentration of the alloy component ZL is changing. Definition (10) has been obtained by assuming the linear dependence of the solid state fraction on the temperature in the mushy zone [17, 18, 29]. In consequence, the following form of function fs can be given

fs (T ) =

TL − T TL − T S

for T ∈ [TS (ZL ), TL (ZL )]

(11)

which leads to (10). It should be noticed that the above assumption about the dependance of the solid state fraction only on temperature is valid for low concentration of the dissolved species. Values of density and the thermal conductivity coefficient in equation (6) vary as well in dependence on temperature   ̺l T > TL (ZL ), ̺ = ̺mz T ∈ [TS (ZL ), TL (ZL )], (12)  ̺s T < TS (ZL ),

  l T > TL (ZL ),  = mz T ∈ [TS (ZL ), TL (ZL )],  s T < TS (ZL ).

(13)

In the above equations the values of density and the thermal conductivity coefficient for the mushy zone are calculated as the arithmetic averages of the values of these parameters for the liquid and solid phases. While creating the model of considered process we take into account one more phenomenon affecting the solidification, that is the macrosegregation referring to the variations in composition of the cast. For describing this phenomenon we decided to use the Scheil model [9, 12–15]. In this approach, with respect to the fact that the diffusion coefficient in the solid phase is significantly lower than in

13

1660

Heat Mass Transfer (2017) 53:1657–1666

the liquid phase, we assume that Ds = 0. It means that the diffusion in the solid phase does not occur. From the other side, the convection in the liquid phase causes the equalization of the alloy fraction concentration in the liquid phase, thus, in this connection, we assume that Dl → ∞. Let us discretize interval [0, t ∗ ] with nodes ti, i = 0, 1, . . . , p∗ , and assume that the values of concentration in moments ti, for i = 1, 2, . . . , p∗, are known. By applying the mass balance of the alloy component in region of the cast for the moment of time tp+1, we receive equation

m0 Z0 = mL (tp+1 ) ZL (tp+1 ) +

p+1 

mS (ti ) ZS (ti ),

(14)

i=1

where m0 means the mass of alloy, Z0 denotes the initial concentration of the alloy component, ZL (tp+1 ) and ZS (tp+1 ) are the concentrations of the alloy in liquid phase and solid phase in moment tp+1, respectively, whereas mL (tp+1 ) and mS (tp+1 ) describe mass of the alloy in liquid state and solid state in moment tp+1. By introducing the partition coefficient k = ZZLS (t) (t), we get

ZL (tp+1 ) =

p m0 Z0 − i=1 mS (ti ) ZS (ti ) . k mS (tp+1 ) + mL (tp+1 )

(15)

Investigated problem is the inverse problem, therefore for solving it we need some additional information. In the discussed case this information is given in the form of temperature measurements in selected points of the domain, that is in the form of values

T (xi , yi , tj ) = Uij , i = 1, . . . , N1 , j = 1, . . . , N2 ,

mS,j (t) =Vj ̺s fj (t),

(16)

mL,j (t) =Vj ̺l (1 − fj (t)).

(17)

By using the above relations equation (15) is transformed to the form ZL (tp+1 )      p n b d ̺l Z0 − ̺s i=1 ZS (ti ) j=0 ∆rj fj (ti ) − fj (ti−1 )     ,     = k ̺s nj=0 ∆rj fj (tp+1 ) − fj (tp ) + ̺l nj=0 ∆rj 1 − fj (tp+1 )

(18)

where b and d refer to the length and width of twodimensional region. Values ZS (ti ) of concentration of the alloy component in the solid phase in moments ti , for i = 1, . . . , p, appearing in the above equation, are calculated from the relation k = ZZLS (t(tii )) for the known value of k and already calculated ZL (ti ).

(19)

where N1 denotes the number of sensors and N2 means the number of measurements taken from each sensor. On the basis of these data we intend to identify the values of heat transfer coefficient α on boundaries Γ3 and Γ4 as well as the distribution of temperature in region Ω. If we assume the fixed form of function α, then problem (6)–(9) can be considered as the direct solidification problem. By solving the direct problem we may determine the values of temperature Tij = T (xi , yi , tj ) corresponding with the assumed form of function α. By using the calculated values of temperature Tij and the given values Uij we formulate the functional

J(α) =

N1  N2   i=1 j=1

Considered region of the alloy is divided into the control volumes Vj based on the area ∆rj, j = 0, . . . , n, resulting from the mesh assumed in the two-dimensional region of the problem. If the solid mass fraction in volume Vj in moment t is denoted by fj (t), then the masses of metal in the solid state and in the liquid state contained in volume Vj in moment t are represented by relations

13

3 Solving procedure

Tij − Uij

2

(20)

expressing the error of approximate solution. While solving the inverse tasks, the problem with stability of the solution often occurs. For avoiding this problem the regularization procedures are frequently used, like for example the Tikhonov regularization or the regularization by means of the sensitivity coefficients. However, in functional (20) we do not use any regularization term. It is our intention since, according to our experiences, in the problems of considered kind the addition of any regularization techniques makes the computations longer without any significant improvement of the results. Instead, in the investigated case the ill-posed problem is regularized by reducing the problem to the determination of a small number of parameters [33]. Perhaps, if we increase the number of determined parameters, it may appear necessary to apply some regularization techniques. In future we plan to execute the suitable analysis in order to investigate this issue more carefully. Presented method of solution is based on two techniques. The first one is needed for solving the direct problem obtained for the assumed form of the sought heat transfer coefficient. Solution of the direct problem is determined on the way of solving equation (6) accompanied with initial condition (7) and boundary conditions (8)–(9). Next, by using formula (18) the concentration of alloy component in the liquid phase is calculated which enables to determine

Heat Mass Transfer (2017) 53:1657–1666

the new values of the liquidus and solidus temperatures TL and TS. For solving the direct problem we use the finitedifference method combined with the generalized alternating phase truncation method. Idea of connecting these two methods lies in introducing the enthalpy in place of the temperature. Next, the calculations are performed in three stages associated with three states—at first the entire domain is reduced to the liquid phase, next to the mushy zone and at the end to the solid phase. In this way, each step of calculations consists of three stages, but each one is reduced to the solution of one-phase direct heat conduction problem. More detailed description of the basic and generalized alternating phase truncation method can be found in papers [29–32]. The second procedure is needed for minimizing functional (20) and is realized with the aid of Ant Colony Optimization (ACO) algorithm—the algorithm of artificial intelligence imitating the behavior of natural ants. Detailed description of ACO algorithm can be found in papers [25, 26] and the mathematical foundations of this algorithm, together with some results concerning its convergence, have been presented in works [27, 28]. The Ant Colony Optimization is an heuristic optimization technique based on the collective behavior of simple individuals composing the self-organized system. The ants are almost or totally blind. Therefore the only chances for communication between the swarm members lie in some chemicals, called pheromones, left on the ground by moving ants to mark the traversed trail. The bigger number of ants use a trail, the more intense is the pheromone smell of this trail. Thus, by using this mechanism the swarm of ants is able to find the best way leading from the ant-hill to the source of food, that is the way which is the shortest and omits the possible obstacles. Below we present the general description of the procedure imitating the behavior of ants and serving for minimizing any function J : D → R, where D ⊂ Rn . Let the vectors x = (x1 , . . . , xn ) ∈ D, dispersed in domain D, be considered as the ants and let the value of minimized function, denoted by J(x), be interpreted as the distance between ant x and the desired location of the food source (for which the value of minimized function is minimal). To initialize the procedure we need to fix the number m of individuals in one population, the number I of iterations and the value of narrowing parameter β (which serves for simulating the evaporation of pheromone trace). Next, we randomly select the initial population of ants xk = (x1k , . . . , xnk ), where xk ∈ D, k = 1, 2, . . . , m, and we determine the best located ant xbest, that is the ant for which the value of minimized function is the lowest. The main part of the ACO algorithm is composed of the following steps.

1661

1. Updating of the ant locations: • random selection of vector x¯ such that

−βi ≤ x¯ j ≤ βi ,

j = 1, . . . , n,

where i denotes the number of current iteration; • creation of the new ant population:

xk = xbest + x¯ ,

k = 1, 2, . . . , m.

2. Determination of the best located ant xbest in the current ant population. 3. Points 1 and 2 are repeated I 2 times. 4. Narrowing of the ant dislocations range—evaporation of the pheromone trail—according to relation

βi+1 = 0.1βi . 5. Points 1–4 are repeated I times. Presented list of steps is designed such that the ants are forced to gather more and more densely around the best solution. The process is terminated after the assumed number of iterations and the obtained solution is the best, but the best obtained in the given run of the algorithm. Every other execution of the procedure can give slightly different result which characterises the algorithms of heuristic nature. Therefore, to avoid the mistake and to ensure the best solution, the procedure based on the heuristic algorithm suppose to be repeated some number of times. Another thing, which should be mentioned here, is that in our case each execution of the ACO algorithm requires to solve for many times the direct heat conduction problem (which is necessary to calculate the value of functional [20]). Therefore, thanks to the way of generating the new population of ants such that each new ant is created independently, it is possible to parallelize the procedure. Since for each found solution we need to solve the associated direct problem, we may execute the calculations independently for each tested solution which significantly decreases the working time of the procedure.

4 Numerical verification The numerical calculations verifying the discussed procedure will be carried out for the following values of parameters: b = d = 0.08 [m], t ∗ = 250 [s], l = 104 [W/(m K)], s = 162 [W/(m K)], cl = 1275 [J/(kg K)], cs = 1077 [J/(kg K)], ̺l = 2498 [kg/m3], ̺s = 2824 [kg/m3], L = 390000 [J/kg], T∞ = 300 [K] and T0 = 940 [K],

13

1662

initial concentration of the alloy component Z0 = 0.02 and partition coefficient k = 0.125, describing the solidification process defined by Eqs. (6)–(9). Dependance of the liquidus and solidus temperatures TL and TS on the concentration of alloy component in the liquid phase ZL is taken in the following form

TL (ZL ) =933.37 − 259.32 · ZL , TS (ZL ) =933.37 − 18076.32 · ZL , derived on the basis of [34]. The unknown element is the value of heat transfer coefficient appearing in the formulation of boundary condition (9) defined on boundaries Γ3 and Γ4. Thus, goal of the task is to reconstruct function α describing this coefficient in the following form    α1 for t ≤ t1 ∧ x ∈ [0, b] ∧ y = d,     α3 for t ∈ (t1 , t2 ] ∧ x ∈ [0, b] ∧ y = d,  α5 for t > t2 ∧ x ∈ [0, b] ∧ y = d, α(x, y, t) = α2 for t ≤ t1 ∧ y ∈ [0, d] ∧ x = b,     α  4 for t ∈ (t1 , t2 ] ∧ y ∈ [0, d] ∧ x = b,   α6 for t > t2 ∧ y ∈ [0, d] ∧ x = b,

where t1 = 70, t2 = 140 [s], on the basis of temperature measurements taken from four thermocouples (N1 = 4 ) located 4mm away from boundary of the region, in points of coordinates (0.0076, 0.0068),  (0.0076, 0.0072),  (0.0068, 0.0076) and (0.0072, 0.0076) [m]. From each sensor we have 250 and 50 values of temperature, since the temperatures were measured in two rounds, at every 1 and 5 s, respectively. The elaborated procedure will be investigated with respect to the precision of the obtain results and their stability. Therefore we know the exact values of the sought αi [W/(m2 K)]:

α1 = 800, α3 = 600, α5 = 250, α2 = 600, α4 = 400, α6 = 200, and the computations were executed for exact values of temperature as well as for the values burdened by pseudorandom 1 and 2 % error of normal distribution. Functional (20) was minimized with the aid of ACO algorithm performed for 18 ants in one population and 6 external iterations (within each iteration there is executed 36 internal iterations—see the algorithm description). The initial values of the sought coefficients αi , i = 1, . . . , 6, were selected from the intervals adequate for the expected values, that is from the intervals: [500, 2000],  [300, 2000],  [300, 1500],  [100, 1000],  [100, 1000] and [100, 600],  respectively. The initial value of the narrowing parameter β was also adjusted to the expected value of the retrieved coefficient αi. Thus, for the respective αi the value of β0 was equal to 800, 900, 600, 500, 500, 300. Next, during the run of the procedure the value of β decreased with coefficient

13

Heat Mass Transfer (2017) 53:1657–1666

0.1 at the beginning (see Step 4 of the ACO algorithm) and with coefficient 0.2 starting from the second iteration. We decided to this change for better concentration of the solutions around the best one. Thanks to this modification and thanks to the concept of parallelizing the ACO algorithm we were able to achieve the quite reasonable time of calculations. Each execution of the procedure took about 2 weeks, since for computations we used the system of ordinary PC computers with processor Intel Core i7. Additionally, because of the heuristic nature of the ACO algorithm, in each considered case of initial data perturbation and the number of measurements the procedure was run for several times to avoid uncertainty of the obtained results. For solving the direct problem, with the aid of finite difference method combined with the generalized alternating phase truncation method, we used the discretization grid of steps ∆t = 0.001, ∆x = b/100 and ∆y = d/100 (by the way we have observed that a reasonable change of the grid density does not influence significantly the obtained results). For verifying the proposed approach we performed a numerical experiment, in which the initial data were not the real measurements, only the values simulating them. For generating the values treated as the measurements we solved the direct problem for the known exact values of the sought αi, but for the grid of different density. In this way we avoided the inverse crime [35, 36]. The certain way for verifying correctness of elaborated algorithm is to investigate its behavior in case of the exact input data, that is for the data simulating the values of temperature in the measurement points obtained directly from the solution of the proper direct problem, not perturbed with any random error. Obviously, the errors of results should be close to zero in that case. In Fig. 2 there is presented the distribution of objective functional (20) values changing in dependence on the number of iterations executed in the ACO algorithm. As expected, the values visibly converge to zero (a little bit faster than in case of approach using the lever arm model presented in paper [16]). Consequence of this is the decrease of the mean and maximal relative errors in reconstruction of values αi calculated according to the formulas 6

1 δi , δ max = max δi , i∈{1,...,6} 6 i=1  approx exact  −α α  where δi =  i α exact i , which is the goal of our task and

δ mean =

i

can be observed in the same figure. Presented results are obtained for measurements taken at every 1 s. Similar collection of results, for the same frequency of measurements but burdened with 2 % error, is displayed in Fig. 3. Values of minimized functional, as well as the reconstruction errors, significantly decrease, however not to zero which

Heat Mass Transfer (2017) 53:1657–1666

1663 28 500

objective function

objective function

250 200 150 100 50

28 450 28 400 28 350 28 300 28 250 28 200 28 150 150

0

150

160

170

180

190

200

160

170

210

180

iterations

190

200

210

160

180

200

iterations 100

80

δ max

80

60

error

error

δ max 40 δ

mean

40

0

100

120

140

160

180

80

100

120

140

iterations

200

iterations Fig. 2  Values of the objective function (20) of the ACO algorithm (upper figure) together with the mean and maximal relative errors in reconstructing all values αi, i = 1, . . . , 6, (lower figure) in dependance on the iteration number obtained for measurements taken at every 1 s and exact input data

cannot be expected because of the input data perturbation. The decrease of objective function value is here a little bit slower than in case considered in [16], but the variation of errors presented in the current paper is smaller. Comparison of all mean errors of the heat transfer coefficient reconstruction, for all cases of perturbations on both sets of measurements, are shown in Fig. 4. It may be surprising that for more frequently taken measurements the reconstruction errors are higher in some cases, but this is exactly the specificity of heuristic algorithm. Some runs of this algorithm can give sometimes slightly better or slightly worse results. The most important is that all the mean errors are comparable with input errors, or even lower in some cases. By comparing the obtained results with the same collection of results shown in paper [16] one may observe that for 1 % input data error the results presented in the current paper are better, however for 2 % input data error the current results are worse. But still, the differences are not much significant.

Fig. 3  Values of the objective function (20) of the ACO algorithm (upper figure) together with the mean and maximal relative errors in reconstructing all values αi, i = 1, . . . , 6, (lower figure) in dependance on the iteration number obtained for measurements taken at every 1 s and 2 % perturbation of input data

5 4

1 noise 0 2 noise 1 3 noise 2

every 1 s 3

error

80

δ mean

20

20 0

60

every 5 s

3

3 2

2 1 0

1

2 1

Fig. 4  Mean relative errors in reconstructing all values αi, i = 1, . . . , 6, for various noises of input data and various frequencies of measurements

The last statement is confirmed by the results collected in Table 1 presenting the reconstructed values of the heat transfer coefficient, relative percentage errors of these reconstructions and standard deviations of results obtained in all runnings of the procedure. Most reconstruction errors are

13

1664

Heat Mass Transfer (2017) 53:1657–1666

Table 1  Results of the calculations for measurements taken at every 1 and 5 s and for various noises of input data (αi —reconstructed values of the heat transfer coefficient, δαi—relative percentage error, sαi —standard deviation) Measurements

Noise

i

αi

δαi [%]

sαi

Every 1 s

0 %

1 2 3 4 5 6 1 2 3 4 5 6 1 2 3 4 5 6 1 2 3 4 5 6 1 2 3 4 5 6 1 2 3 4 5

800.70 602.86 249.12 599.48 397.48 201.09 802.56 598.83 250.34 597.36 404.51 198.56 788.14 608.87 235.91 611.35 392.83 215.23 800.38 599.62 249.90 600.12 399.98 200.06 804.07 621.79 248.46 590.74 386.70 205.52 789.66 625.17 250.87 602.20 378.56

0.087 0.477 0.353 0.087 0.631 0.545 0.319 0.195 0.134 0.440 1.126 0.722 1.483 1.478 5.636 1.891 1.792 7.614 0.048 0.063 0.039 0.020 0.006 0.028 0.509 3.632 0.615 1.543 3.325 2.761 1.292 4.196 0.346 0.367 5.360

0.792 0.775 0.850 0.591 0.771 0.755 0.767 1.051 0.778 0.617 0.885 0.727 1.815 5.245 1.447 1.948 5.403 1.633 0.344 1.581 1.198 0.343 1.499 1.276 2.244 8.210 3.311 2.324 7.900 2.967 0.096 0.519 0.366 0.795 0.183

6

207.33

3.667

0.354

1 %

2 %

Every 5 s

0 %

1 %

2 %

comparable with the input data errors. The worst obtained value is about 7.6 % and was received for 2 % input data error and 250 measurements (the worst reconstruction error received by using the approach from paper [16] was higher—it was about 13.4 % obtained for 2 % input data error and 50 measurements). It is quite high value, however not enough high to disqualify the procedure. Especially since the small values of standard deviations of the results are small which indicates stability of the procedure.

13

Table 2  Relative and absolute errors of the temperature reconstruction in all four points of the sensor locations, for measurements taken at every 1 s and for various noises of input data Noise

Sensor

δ max (%)

δ mean (%)

max [K]

mean [K]

0 %

1 2 3 4 1 2 3 4 1 2 3

0.0183 0.0128 0.0069 0.0127 0.037 0.0441 0.0567 0.0578 0.0943 0.0737 0.0883

0.0059 0.0042 0.0025 0.0050 0.0160 0.0158 0.0183 0.0193 0.0381 0.0217 0.0389

0.1404 0.0977 0.0506 0.0967 0.2863 0.3313 0.4262 0.4410 0.8256 0.6478 0.6423

0.0471 0.0328 0.0193 0.0390 0.1249 0.1218 0.1428 0.1535 0.3057 0.1743 0.2997

4

0.1364

0.0574

0.9988

0.4476

1 %

2 %

Table 3  Relative and absolute errors of the temperature reconstruction in all four points of the sensor locations, for measurements taken at every 5 s and for various noises of input data Noise

Sensor

δ max (%)

δ mean (%)

max [K]

mean [K]

0 %

1 2 3 4 1 2 3 4 1 2 3

0.0071 0.0073 0.0066 0.0057 0.1584 0.1274 0.1262 0.1363 0.1406 0.1417 0.1562

0.0015 0.0015 0.0014 0.0013 0.0656 0.0639 0.0615 0.0693 0.0663 0.0603 0.0703

0.0577 0.0581 0.0526 0.0467 1.1949 0.9517 0.9179 0.9978 1.1393 1.0674 1.1353

0.0122 0.0122 0.0113 0.0108 0.5055 0.4889 0.4785 0.5478 0.5245 0.4722 0.5459

4

0.1662

0.0760

1.2163

0.5935

1 %

2 %

Relative and absolute errors of the temperature reconstruction, for measurements taken at every 1 s and for various noises of input data received for all four points of the sensors locations, are compiled in Table 2. All the presented errors, minimal and maximal, relative and absolute as well, are very small which shows that the reconstruction of temperature is very good. Next, Table 3 presents the similar set of results but for measurements read at every 5 s and in this case the reconstruction of temperature is also very good. Moreover, the distributions of relative error of the temperature reconstruction in all four control points for the case of measurements read at every 1 s and 2 % perturbation of input data, displayed in Fig. 5, confirm the above statements. Similarly, very good reconstruction of temperature was obtained by applying the approach described in paper [16].

Heat Mass Transfer (2017) 53:1657–1666

1665

Fig. 5  Errors of the temperature reconstruction in all four points of ▸ the sensors locations, for the case of measurements read at every 1 s and 2 % perturbation of input data

error

0.04 0.02 0.00

0

50

100

150

200

250

150

200

250

150

200

250

150

200

250

t s 0.07 0.06

error

0.05 0.04 0.03 0.02 0.01 0.00

0

50

100

t s 0.08 0.06

error

Scope of the paper was the procedure elaborated for solving the inverse solidification problem for binary alloy in the two-dimensional space. Aim of the task was to identify the values of heat transfer coefficient, representing the boundary condition on one of the boundaries, together with the temperature distribution in the entire domain. The task was possible to solve thanks to the measurements of temperature taken in two periods of time from four thermocouples located 4 mm away from boundary of the region. The values of sought coefficient were selected so that the differences between the measured and the reconstructed values of temperature are as small as possible. For obtaining this we minimized the appropriate functional representing the error of approximate solution. From mathematical point of view, for describing the process we used the heat conduction equation with the substitute thermal capacity completed by the appropriate initial and boundary conditions, for describing the macrosegregation we applied the Sheil model, for solving the direct problem we involved the finite difference method combined with the generalized phase truncation method and for solving the optimisation task we choose the Ant Colony Optimisation algorithm. Numerical verification of the approach showed that the elaborated procedure can be found as the effective tool for solving the inverse problem of considered kind. In most of investigated cases of input data the reconstruction errors were comparable with the input errors, the worst reconstruction error were obtained at the level of about 7.6 % which can be still accepted, especially in view of the satisfactorily small errors of temperature reconstruction. Also the standard deviations of the results, received in several executions of the algorithm, were low which indicates the desired stability of the procedure. Number of iterations and number of ants in the ACO algorithm, needed for obtaining the satisfactory results, were not big, 6 and 18, respectively. The calculations for the respective individuals were performed parallel, even so the computations within one execution of the procedure took about 2 weeks. The long time of working is the only disadvantage of the elaborated procedure. In future we plan to shorten this time by testing the compilation of the heuristic algorithm with some known analytical algorithm. The analytical algorithms work fast, but they need a good starting point, which can be determined exactly with the aid of an heuristic algorithm.

0.06

0.04 0.02 0.00

0

50

100

t s

0.12 0.10 0.08

error

5 Conclusions

0.08

0.06 0.04 0.02 0.00

0

50

100

t s

13

1666

And finally, the comparison of results delivered in the current paper by using the Scheil model with the results obtained by applying the approach based on the lever arm model led us to the conclusion that both models are comparatively good for approximating the macrosegregation process taking place in the solidifying alloys. With both approaches we obtained similarly good results (maybe with the tiny advantage of the Scheil model) in solving the problem of considered kind. However, to improve the mathematical model of the binary alloy solidification, we plan for the future to use the diffusion equation to determine the concentration of the alloy component. Acknowledgments This project has been financed from the funds of the National Science Centre granted on the basis of decision DEC-2011/03/B/ST8/06004. The Authors would like to thank the Reviewers for their valuable remarks and comments.

References 1. Tikhonov AN (1977) Solutions of Ill-posed problems. Halsted Press, Washington 2. Murio DA (1993) The mollification method and the numerical solution of ill-posed problems. Wiley, New York 3. Beck JV, Blackwell B, St. Clair CR (1985) Inverse heat conduction, Ill-posed problems. Wiley- Interscience, New York 4. Zabaras N, Ruan Y, Richmond O (1992) On the design of two-dimensional Stefan processes with desired freezing front motions. Numer Heat Transf B 21:307 5. Ang DD, Pham Ngoc Dinh A, Thanh DN (1996) An inverse Stefan problem: idenification of boundary value. J Comput Appl Math 66:75–84 6. Zabaras N, Kang S (1994) On the solution of an ill-posed inverse design solidification problem using minimization techniques in finite and infinite dimensional spaces. Int J Numer Methods Eng 36:3973 7. Kang S, Zabaras N (1995) Control of the freezing interface motion in two-dimensional solidification processes using the adjoint method. Int J Numer Methods Eng 38:63 8. Ren H-S (2007) Application of the heat-balance integral to an inverse Stefan problem. Int J Therm Sci 46:118–127 9. Mochnacki B, Majchrzak E, Szopa R (2008) Simulation of heat and mass transfer in domain of casting made from binary alloy. Arch Foundry Eng 8(4):121–126 10. Słota D (2011) Solution of the inverse solidification problems by using the genetic algorithms. Silesian University of Technology Press, Gliwice (in Polish) 11. Hetmaniok E, Słota D (2012) Numerical procedure of solving some inverse problem in solidification of the binary alloy. Comput Assisted Meth Eng Sci 19:393–402 12. Mochnacki B, Suchy JS, Praz˙mowski M (2000) Modelling of segregation in the process of Al–Si alloy solidification. Soldification Metals Alloys 2(44):229–234 13. Mochnacki B, Suchy JS (2006) Simplified models of macrosegregation. J Theor Appl Mech 44:367–379 14. Słota D (2011) Reconstruction of the boundary condition in the problem of the binary alloy solidification. Arch Metall Mater 56(2):279–285

13

Heat Mass Transfer (2017) 53:1657–1666 15. Hetmaniok E, Słota D (2012) Experimental verification of the procedure of reconstructing the boundary condition in the problem of binary alloy solidification. Steel Res Int, special edition Metal Forming:1043–1046 16. Hetmaniok E (2015) Solution of the two-dimensional inverse problem of the binary alloy solidification by applying the Ant Colony Optimization algorithm. Int Commun Heat Mass Transf 67:39–45 17. Majchrzak E, Mochnacki B (1995) Application of the BEM in the thermal theory of foundry. Eng Anal Bound Elem 16:99–121 18. Mochnacki B (2011) Numerical modeling of solidification process. In: Zhu J (ed) Computational Simulations and Applications. InTech, Rijeka, pp 513–542 19. Santos CA, Quaresma JMV, Garcia A (2001) Determination of transient interfacial heat transfer coefficients in chill mold castings. J Alloys Compd 319:174–186 20. Hetmaniok E (2016) Artificial bee colony algorithm in the solution of selected inverse problem of the binary alloy solidification. Thermal Sci. doi:10.2298/TSCI140715136H 21. Hetmaniok E (2015) Numerical procedure for the heat transfer coefficient identification in solidification of the binary alloy and its experimental verification. Numer Heat Transf B 68:93–114 22. Hetmaniok E (2016) Inverse problem for the solidification of binary alloy in the casting mould solved by using the bee optimization algorithm. Heat Mass Transf 52:1369–1379 23. Grzymkowski R, Hetmaniok E, Słota D, Zielonka A (2012) Application of the ant colony optimization algorithm in solving the inverse Stefan problem. Steel Res Int (special edition Metal Forming):1287–1290 24. Hetmaniok E, Słota D, Zielonka A (2015) Using the swarm intelligence algorithms in solution of the two-dimensional inverse Stefan problem. Comput Math Appl 69:347–361 25. Dorigo M, Stützle T (2004) Ant colony optimization. MIT Press, Cambridge 26. Toksari MD (2006) Ant colony optimization for finding the global minimum. Appl Math Comput 176:308–316 27. Dorigo M, Blumb Ch (2005) Ant colony optimization theory: a survey. Theor Comput Sci 344:243–278 28. Duan H (2011) Ant colony optimization: principle, convergence and application. Adapt Learn Optim 8:373–388 29. Mochnacki B, Suchy JS (1995) Numerical methods in computations of foundry processes. PFTA, Cracow 30. Rogers JCW, Berger AE, Ciment M (1979) The alternating phase truncation method for numerical solution of a Stefan problem. SIAM J Numer Anal 16:563–587 31. Słota D (2008) Solving the inverse Stefan design problem using genetic algorithms. Inverse Probl Sci Eng 16:829–846 32. Słota D (2011) Restoring boundary conditions in the solidification of pure metals. Comput Struct 89:48–54 33. Mera NS, Elliott L, Ingham DB (2004) Numerical solution of a boundary detection problem using genetic algorithms. Eng Anal Bound Elem 28:405–411 34. Xu D, Li Q (1991) Numerical method for solution of strongly coupled binary alloy solidification problems. Numer Heat Transf B 20:181–201 35. Fadale T, Nenarokomov A, Emery A (1995) Uncertainties in parameter estimation: the inverse problem. Int J Heat Mass Transfer 38:511–518 36. Ryfa A, Białecki R, Facchini B, Tarchi L (2009) Application of the inverse analysis for boundary condition retrieval. Inverse Probl Sci Eng 17:829–853