Mathematical Problems in Engineering

Volume 2015, Article ID 680823, 15 pages

http://dx.doi.org/10.1155/2015/680823

## Inverse Transient Radiative Analysis in Two-Dimensional Turbid Media by Particle Swarm Optimizations

School of Energy Science and Engineering, Harbin Institute of Technology, Harbin 150001, China

Received 30 September 2014; Revised 4 December 2014; Accepted 9 December 2014

Academic Editor: Subhash Chandra Mishra

Copyright © 2015 Yatao Ren et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

Three intelligent optimization algorithms, namely, the standard Particle Swarm Optimization (PSO), the Stochastic Particle Swarm Optimization (SPSO), and the hybrid Differential Evolution-Particle Swarm Optimization (DE-PSO), were applied to solve the inverse transient radiation problem in two-dimensional (2D) turbid media irradiated by the short pulse laser. The time-resolved radiative intensity signals simulated by finite volume method (FVM) were served as input for the inverse analysis. The sensitivities of the time-resolved radiation signals to the geometric parameters of the circular inclusions were also investigated. To illustrate the performance of these PSO algorithms, the optical properties, the size, and location of the circular inclusion were retrieved, respectively. The results showed that all these radiative parameters could be estimated accurately, even with noisy data. Compared with the PSO algorithm with inertia weights, the SPSO and DE-PSO algorithm were demonstrated to be more effective and robust, which had the potential to be implemented in 2D transient radiative transfer inverse problems.

#### 1. Introduction

The problem of transient radiative transfer in participating media has attracted increasing interest over the last two decades due to the emergence of ultrashort pulse laser and its application to the picosecond level or higher accuracy level of time-resolved techniques [1]. Lasers with pulse durations of pico- to femtoseconds have been widely utilized in the areas of laser light-scattering flame diagnose, the prediction of the temperature distribution in combustion chambers, atmospheric science, remote sensing, nondestructive testing, laser machining, and medical diagnostics, to name a few [2–4]. In recent years, the medical applications of short pulse laser, such as the laser tissue welding [5, 6], the laser tissue ablation [7, 8], the optical tomography in medical imaging [9], and the photodynamic therapy [10], have drawn more and more attentions all over the world. In particular, the short pulse laser of near-infrared spectrum is extensively used in noninvasive techniques to retrieve the characteristics of the semitransparent media such as the biological tissues and turbid media [11–14]. Theoretically speaking, the research emphases of noninvasive reconstruction technology are focused on using the model-based iterative imaging reconstruction techniques, which employ a precise and efficient numerical forward model based on solving complete transient radiative transfer equation (TRTE) and then establish a proper inverse algorithm to retrieve the optical property parameters or internal geometry of the medium by using the measured time-resolved transmittance and reflectance signals. Compared with the traditional measurement methods, its biggest advantages are noninvasive, in situ measurement, containing a lot of time-resolved information about internal composition and characteristic of the participating media.

For solving the inverse transient radiation problems in the noninvasive reconstruction, the first step is to solve the direct TRTE model accurately and efficiently. Nowadays, several numerical strategies have been developed to solve the forward model, including the discrete ordinate method (DOM) [15–18], the finite volume method (FVM) [18–21], the integral equation (IE) method [22], the finite-element method FEM [23, 24], and the Monte Carlo method [25]. Each of these methods has its own relative advantages and disadvantages and none of them is superior to others in all aspects. However, many numerical complications may arise from incompatibilities between the radiative transport model and the discrete formulation employed for the fluid dynamics and other heat transfer modes. The FVM for RTE avoids many of these complications and has the advantage of greater compatibility with existing finite volume based heat transfer software. Besides, the FVM also could simplify the problems of complex geometries using the unstructured meshes. Through decades of development, FVM has become a sophisticated technology in mechanics and thermal analysis. There have been several investigations focused on the interaction of ultrafast radiation with participating media using FVM [19, 20]. In these studies, both the homogeneous and inhomogeneous media have been considered. However, most of the studies of inverse radiation problems concentrate on one-dimensional (1D) cases and there are very few researches about multidimensional problems, especially the inverse transient radiative problems.

To solve the inverse problems, a wide variety of inverse techniques have been successfully employed in the inverse radiation analyses which can be grouped into two categories, that is, the traditional algorithm based on gradient and the intelligent optimizations [26]. To date, many gradient-based techniques have been employed in the inverse radiation analysis such as the constrained least-squares method [27], the conjugate gradient (CG) method [28], and Levenberg-Marquardt method [29]. However, all these traditional methods depend on the initial value or the derivatives and gradients which are difficult to be solved accurately by numerical simulation in some cases. Furthermore, if the initial value is not chosen properly, the solution may be infeasible in the original domain, or in a worse case, not be convergent. In addition, these conventional methods using derivatives are generally restricted to the conditions of continuity, sensitivity, convexity, monotonicity, and nonlinearity of both objective functions and constraints [30]. Compared with traditional gradient-based methods, the intelligent optimization algorithms have a couple of outstanding characteristics; that is, the ill-posed inverse problem could also be solved, the derivative of the objective function is not necessary, and the a priori information is not needed. A typical characteristic of the intelligent optimization methods is that they can solve the global optimal problem reliably and obtain high quality global solutions with enough computational time, especially for higher problem dimensions [31]. This is primarily due to the fact that the random search technique is used to modify locally examined regions, which gives the algorithm a global optimization capability. The Particle Swarm Optimization (PSO), as an intelligent swarm optimization, first introduced in 1995 by Eberhart and Kennedy [32], has been studied extensively by many researchers in recent years. It is a potential heuristic bionic evolutionary algorithm originally inspired by the locking behavior of birds. Generally speaking, PSO is characterized to be simple in concept, easy to implement, and computationally efficient. Early studies have found the implementation of PSO to be effective and robust in solving problems featuring nonlinearly, nondifferentiability, and high dimensionality. Many modifications have been made to improve the convergence rate of the original PSO algorithm. More recently, our research group has done several studies on the PSO algorithm to solve the one-dimensional inverse radiation problems [33–37]. However, to the best of our knowledge, few researches are presented to solve the multidimensional inverse transient radiation problems by the intelligent optimization methods.

The objective of present work is to apply the PSO-based algorithms to solve the inverse transient radiation problem in two-dimensional (2D) participating medium. The optical properties, the size, and location of the circular inclusion were retrieved, respectively. The remainder of the paper is organized as follows. The detailed mathematical formulation and computational steps of the direct model are described in Section 2. The theoretical models of the standard PSO, the Stochastic PSO, and the hybrid Differential Evolution and Particle Swarm Optimization (DE-PSO) algorithms are described in Section 3. The solving model for direct problem and the inverse transient radiation analysis by the PSO-based algorithms are examined in Section 4. In Section 5, the main conclusions and perspectives are provided.

#### 2. The Direct Model

In the present work, the thermal effect caused by the incident laser was ignored and only the transient radiative transfer was considered in the numerical model. Thus, the transient radiative transfer equation can be expressed as [38]where is the length along the direction , radiative intensity is the function of position , direction , and time , , , and denotes extinction, absorption, and scattering coefficients, respectively, and . represents the scattering phase function between incoming direction and scattering direction . is the solid angle in the direction . In the 2D transient radiative transfer problem caused by collimated irradiation, the radiative intensity within the medium is separated into two parts [35]: (a) , the remnant of the collimated beam after partial extinction, by absorption and scattering, along its path; (b) , the diffuse part, which is the result of emission from the boundaries, emission from the medium inside, and the radiation scattered away from the collimated irradiation. Thus, the radiative intensity in the media can be written as

For the square pulse laser, the intensity within the participating media can be expressed as [35]where is Heaviside function [35]; denotes Dirac function [35]; is the heat flux of incident laser at boundary ; is the incident direction; and represents the geometric distance along the direction .

For the 2D participating media, the discrete equation could be obtained using the FVM model:where

In (4), denotes the total source term. The source terms and result from the collimated radiation and the diffused radiation , respectively. For the diffuse reflection boundary, the boundary condition is given as [35]where and are the outgoing and the incoming intensities at the boundary, respectively; is the outward unit normal vector at the boundary; is the emissivity of the boundary. Substituting the length-dimensional time term into (4), we can get

The total source term can be depicted aswhere denotes the direction of incident radiation. Using the fully implicit scheme, the diffusion intensity for the control volume in the direction at time can be written aswhere represents the velocity of light; is the time step; represents the direction weight; and are set as and , respectively. and denote the -direction radiative intensities of node at the upstream boundaries of control volume which are along axis and axis , respectively. The time-domain thermal signals, that is, the transmittance signal and reflectance signal , can be expressed aswhere is the direction cosine and denotes the direct transmission flux of the collimated light. The FVM was chosen to solve the equation of TRT. For the sake of simplicity, the details of FVM are available in [20] and are not repeated here.

#### 3. The Inverse Model

Commonly, there are three ways to obtain the internal information of media, that is, the continuous wave method, the time-domain method, and the frequency-domain method [39]. All these methods have the similar procedure in the experimental studies and the numerical simulations. Taking the time-domain method as an example, the reconstruction scheme consists of three major parts [40]: (1) a forward model that predicts the detector signals based on the solution of the transient radiative transfer equation; (2) an objective function that provides criterion of the differences between the detected and the predicted data; (3) the reconstruction technique which can minimize the objective function to get new guesses of the estimated parameters such as the optical properties or the geometrical parameters. Based on the new guesses of the optical properties, a new forward calculation is performed to get the corresponding detector predictions. The reconstruction process is completed when the value of the objective function is less than a preset value or the number of iterations exceeds the maximum number of iterations. Foremost, the forward model must be solved precisely enough so that the measurement data obtained by detectors could be simulated correctly. Consequently, the forward model based on the complete transport-theory TRTE should be utilized [40]. The details of the inverse optimization algorithms are shown as follows.

##### 3.1. Theoretical Model of Basic PSO Algorithm

The PSO algorithm was first introduced by Eberhart and Kennedy [32] inspired by social behavior simulations of flocking birds. The flocking population is called a swarm, and the individuals are called the particles. The basic idea of PSO can be described simply as follows. Each particle in a swarm represents a potential solution which changes its position and updates its velocity based on its own and its neighbors’ experience. The velocity and position of the th particle at generation could be expressed as [32]where and represent two positive constants called acceleration coefficients; and are uniformly distributed random numbers in the interval ; is the inertia weight factor which is used to control the impact of the previous velocities on the current velocity. In the present work, is defined as , where represents the maximum generation number; denotes the present location of particle , which represents a potential solution; is the present velocity of particle , which is based on its own and its neighbors’ flying experience; is the “local best” in the th generation of the swarm. The global best position of the swarm is expressed as .

##### 3.2. Theoretical Model of DE-PSO Algorithm

In order to avoid premature convergence of the standard PSO and ensure the overall convergence, a considerable number of modified PSO algorithms were proposed in recent years, like the Stochastic PSO (SPSO), the Multiphase PSO (MPPSO), the Quantum-Behaved PSO (QPSO), the hybrid Ant Colony Optimization and PSO (ACO-PSO), and so forth [33, 36, 41]. All these improved PSO-based algorithms are proved to be adaptive and robust parameter-searching techniques. In the present work, the SPSO and the hybrid Differential Evolution and PSO (DE-PSO) are utilized in the inverse procedure [42]. The DE-PSO is the combination of two stochastic algorithms and is expected to accelerate the convergence process by introducing the DEA which has strong global searching ability to the PSO. The detailed model description of SPSO is available in [41] and will not be repeated here. The theoretical model of DE-PSO is described as follows.

The differential evolution algorithm (DEA) is a novel parallel, direct, and stochastic searching method which was proposed to solve the continuous domain problems in 1995 [43]. The basic idea of DEA is to use the differential of two individuals in the current generation, namely, vectors and , to combine with a third vector to generate a new particle, which is called the trial vector. Then, the objective function was calculated to decide whether the trial vector is to replace vector or not. The trial vector is generated according to [43]where is a real and constant factor which controls the amplification of the differential variation , and in this paper it was set as 0.5. The integers , , and are chosen randomly from the interval [] and are different from running index . Other parameters are defined as same as those of the PSO. Figure 1 illustrates the vector-generation process defined by (13).