Inverse problem for porosity estimation during solidification of TNT

In the present study, the porosity formed during the solidification process is estimated by an inverse problem technique based on particle swarm optimization. The effective heat capacity method is adopted to model the heat transfer problem. The transient-diffusive heat transfer equation is solved numerically by the finite volume method with an explicit scheme, employing the central difference interpolation function. The solution of the direct problem is compared to reference solutions. The model is applied to trinitrotoluene (TNT) solidification process. The results show that the proposed procedure was able to estimate the porosity for different Stefan numbers. The analysis of the heat flux in the mold is indicated to predict the porosity formation during the casting process.


Introduction
The solidification phenomenon is present in nature and in industrial processes.In nature, it can be observed during the formation of ice as well as volcanic lava, while in industry, it is found in freezing, casting and welding (Tan, Li, & Li, 2011, Verma & Gómez-Arias, 2013, Yajun, Pengfei, Yongjun, Zhijun, & Yintao, 2013, Humphreys et al., 2013).
The Stefan problem is the classical mathematical model for the solidification process.It is based on the heat equation in a solid region where the phase change interface is regarded as a moving boundary (Stefan, 1891).
The method of effective heat capacity has been proposed to avoid the difficulties of the Stefan approach in complex geometries.In this method, the solid and liquid regions are solved in a single domain and an effective heat capacity takes into account the phase change phenomenon, considering the existence of a mushy region, where liquid and solid are present (Hu & Argyropoulos, 1996, Poirier & Salcudean, 1998, Mosaffa, Ferreira, Talati, & Rosen, 2013).
The porosity formed during the solidification process should be avoided since it is usually considered a manufacturing defect.Porosity changes the mechanical and thermophysical properties of the material (Zhang, Kim, & Lu, 2009).
The solidification process of casting explosives in artillery shells is used to produce high-explosive ammunitions (HE ammunitions).These ammunitions are filled with explosives that detonate to produce the terminal effect.The presence of pores in the explosive inside these devices is a dangerous problem, since this defect could lead to an undesirable initiation of the explosive in addition to changing the terminal effect, therefore reducing the efficiency of the ammunition (Zhang, Hu, & Liang, 2013, Kumar & Rao, 2014).
The Particle Swarm Optimization (PSO) is based on social behavior of various species and aims at balancing individuality and sociability, in order to optimize a group of parameters.This stochastic and evolutionary method has been employed in inverse problems for parameters estimation (Colaço, Orlande, & Dulikravich, 2006).
The present study deals with the estimation of porosity formed during the solidification process.A mathematical model with a numerical solution is provided for the direct problem.PSO is used to solve the inverse problem.Trinitrotoluene (TNT) solidification is studied because TNT casting process is often employed in HE ammunitions.Moreover, in such devices, the presence of pores is dangerous and must be avoided.

Material and methods
In this study, a one-dimensional transient heat solidification problem in rectangular coordinates is analyzed.Initially, the entire domain is at uniform temperature and only the liquid phase is present.Suddenly, a temperature lower than the solidification temperature is applied at the origin of the system of coordinates (y = 0) and thus a prescribed temperature is attained at this boundary.Then, the solidification process starts at y = 0 and the solidification front moves away from such boundary.The physical properties of liquid and solid phases are assumed constant at each phase.The effective heat capacity method is adopted (Hu & Argyropoulos, 1996).The mathematical model represented in Equation 1 describes heat transfer in the solid phase, in the liquid phase and in the mushy region.
In this model, C ef is the effective heat capacity, k ef is the effective thermal conductivity, T is the temperature, t is the time, y is the position, L is the length of the domain, T 0 is the initial temperature and T w is the wall temperature of the mold.
By applying the Finite Volume Method (FVM) (Maliska, 2004), the effective heat capacity is evaluated at each finite volume as, according Equation 2: (2) The usual nomenclature adopted in FVM is employed in Equation 2. The subscript n is related to the north face of the finite volume, while the subscript s is related to the south face of the finite volume.The apparent heat capacity, C ap , is Equation 3: In Equation 3, ∆H is the latent heat, T m is the melting temperature and ∆T is half of phase change temperature range, which is a parameter of the effective heat capacity method (Hu & Argyropoulos, 1996).This parameter defines a theoretical mushy region for the evaluation of the effective heat capacity.Furthermore, the subscripts solid and liquid indicate, respectively, the solid phase and the liquid phase.In addition, ρ is the density, c p is the specific heat and ε is the porosity.
The proposed model considers that the values of thermal conductivity, k, and those of specific heat, c p , of the TNT are constant and they do not depend on the phase (Sun & Garimella, 2007).Moreover, the specific heat is not a function of the porosity (Gibson & Ashby, 1997).The porosity is the ratio of the volume of the pores to the total volume of the porous medium.
It is common knowledge that porosity changes the thermophysical properties of the solidified region.Therefore, by considering a homogenous distribution for the pores, the effective thermophysical properties must be defined for the porous medium.The pores are also considered spherical in materials with porosity less than 70% (Zhang et al., 2009).In addition, assuming that thermal conductivity is greater in the solid medium than in the pores, the effective thermal conductivity is defined (Bauer, 1993, Wang, Carson, North, & Cleland, 2006, Zhang et al., 2009) as Equation 4: where: β=1 is used to spherical pores (Zhang et al., 2009).
The buoyancy effects can be neglected in TNT casting (Kumar & Rao, 2014), since TNT viscosity (Sun & Garimella, 2007)  the viscosity of water at 300 K.
In order to solve the problem described by Equation 1, the FVM is employed by using an explicit scheme and the central differences interpolation function (Maliska, 2004).The stability criterion is satisfied and described by Equation 5: Simulations were performed in a computer with Quad Core i7-860 (2.80 GHz) processor and 8 Gb RAM, using Δy = 0.0021 m and Δt = 10 s, which results in 0.14x10 -6 ≤ Γ ≤ 0.34.
The parameters for the model and the physical properties of TNT are presented, respectively, in Table 1 and 2. They are based on real ammunition.The proposed inverse problem estimates the porosity formed during the solidification process of TNT.The inverse problem is an optimization problem, which is solved by using PSO.The optimum solution determines the estimated porosity.
In PSO, each particle i of a given population P in a N dimensional space has a velocity j i v and a position j i x , which are updated in accordance with Equations 6 and 7: ( ) ( ) where: i=1, 2,…, P; j is the iteration of the optimization process; r and i 2 r are random numbers between 0 and 1 with uniform distribution; p i is the best parameters set of the particle i x ; p g is the best parameters set of the population; ω is the inertial constant within the range [0, 1]; c 1 and c 2 are constants within the range [1,2].
The PSO is used to estimate the porosity parameter, considering the objective function, according Equation 8: where: S is the position of the solidification front.The subscript exp refers to the experimental data and the subscript num refers to the proposed FVM solution.nt is the total number of evaluated times and q w is the heat flux at the wall (y = 0), which is expressed as Equation 9: Consequently, Equation 10: Thus, measurements of the heat flux at the wall (y = 0) (or of the position of the solidification front) could be introduced in the objective function to estimate porosity.However, in the present study, since experimental data were not available, an analytical solution (Zhang et al., 2009) has mimicked the experimental data (Iacono, Sluys, & Van Mier, 2011).
The analytical solution (Zhang et al., 2009) is Equations 11 at 14: where: ξ is the dimensionless position of the solidification front and λ is a theoretical parameter determined from Equation 12, τ is the Fourier number and Ste is the Stefan number.
It is important to note that usually experimental data or synthetic experimental data that could be obtained from a mathematical solution and corrupted by random noise, with a prescribed noise level, are used in inverse problems.Such approach is employed to avoid the inverse crime.However, in this paper, the analytical solution (Zhang et al., 2009) solves the solidification process in a semi-infinite domain, considering a Stefan model, while the numerical solution solves the solidification process in a finite domain, considering the effective heat capacity approach.In this sense, the inverse crime is not present in this study, because the analytical solution is used to simulate experimental data and the numerical solution solves the inverse procedure (Chávez, Alonzo-Atienza, & Álvarez, 2013).However, different levels of random noise were also introduced in the analytical solution (Zhang et al., 2009) by Equation 15: where: S analytical is the solidification front computed with the analytical solution (Zhang et al., 2009); ϖ is a random variable with normal distribution, zero mean and unitary standard deviation; E is an arbitrary noise level and; L is the maximum value of S analytical (Machado & Orlande, 1998, Ozisik & Orlande, 2000).
The stopping criterion of the optimization procedure considers that Equation 16 must be satisfied during 50 consecutive iterations.

Results and discussion
Figure 1 and 2 are mainly presented to verify the quality of the proposed numerical solution.Discrepancies between numerical and analytical results for the position of the solidification front are shown when ξ→1.They are caused by the finite domain used in the numerical solution and by the semi-infinite domain used in the analytical solution.Such discrepancies are reduced when the length of the finite domain is increased.The influence of the length of the finite domain on the position of the solidification front can be observed in Figure 1, where the results for the solidification of TNT (Ste = 0.58) without pores are presented.By increasing the length of the domain, the numerical solution approaches the analytical solution.Hence, the effective heat capacity method and proposed FVM solution represents correctly the evolution of the solidification front of TNT.found in the estimated porosity is also increased.For 10% of random noise level, 15% of relative error is obtained for the estimated porosity.Nevertheless, by increasing the noise, the CPU time and the number of iterations to reach the convergence showed non-monotonic behaviors.This can be explained by the adopted stopping criterion.The noise hinders the convergence and establishes a limit for the minimization process, as seen in Figure 6.By increasing the noise, the value of this limit is also increased, thus reducing the number of required iterations for convergence.As expected, the noiseless case shows the best results for the parameter estimated, spending less time, since the solidification model was built for this case.

Conclusion
In this study, a mathematical model based on an effective heat capacity approach was proposed.The Finite Volume Method was employed to solve this model by using the explicit scheme and the central differences interpolation function.
The numerical solution was verified and reference results for solidification processes with and without pores were reproduced.The disagreements found between those solutions were due to differences in the models, since the model proposed in this study considers a finite domain and the reference results were obtained from models for semi-infinite domains.
The results showed that porosity reduces the solidification velocity, because it diminishes the effective diffusivity of the solid.Therefore, by increasing the porosity, total solidification time is also increased.Thus, measurements of the total time of solidification and the heat flux at the boundary wall can be used to predict the porosity formation.
The particle swarm optimization was employed in the inverse problem to estimate the porosity.The effects of the population size on PSO performance were analyzed, indicating that 10 particles were enough to execute the computations.The objective function was built, considering synthetic experimental data of the position of the solidification front, which were related to the heat flux at the boundary wall.
The results of the inverse problem were verified, estimating correctly the porosity parameter for noiseless synthetic experimental data.However, when noise was increased, the porosity estimation became less accurate.
This study focused on the prediction of the porosity formed during the solidification process of TNT.The provided information can be used for planning the process of fabrication of high-explosive ammunitions in order to avoid pores and voids inside the solidified explosive.Furthermore, this study has shown that the level of porosity can be estimated by using an inverse problem technique.

Figure 2
Figure2shows the numerical results provided by the present analysis and the ones byZhang et al. (2009) for the evolution of the solidification front of a metal (Ste = 0.06), considering the cases with pores (ε = 0.7) and without pores (ε = 0.0).The proposed numerical results match the ones by

Table 3 .
Effect of the number of particles on the optimal solution (Ste = 0.58, ε = 0.7).