Retrieval of Optical Constant and Particle Size Distribution of Particulate Media Using the PSO-Based Neural Network Algorithm
Participating medium is widely used in the petroleum chemical industry, bio-medicine, infrared detection, remote sensing, and many other fields [1, 2]. The particle system is one kind of participation media, such as liquids or gases, with small particles in microns or even nanometer level. Researches on the microphysical properties of particulate media, such as the optical properties and particle size distributions (PSDs), are very important for many engineering applications, such as power plant, coal-fired furnaces, combustion chambers, or other utilization systems, which will contribute to improve the combustion efficiency [3–5] and environment control [6, 7]. Especially, thermal radiation is the dominant mode of heat transfer in power plant [8], and the optical properties and PSD of fly ash particles are of significant importance in radiative heat transfer analysis [9]. Knowledge of the optical constants of the particulate fly ash is necessary for monitoring of temperature and condition, and prediction of heat transfer through fly ash formed in the pulverized fuel-fired furnace. However, determination of optical constants requires a priori information on the particle number density and size distribution.
Based on this fact, uncertainty in the PSD is the principal source of error in determining single particle optical properties. In many cases of practical interest, however, these quantities are not known. In spite of the widely use of high temperature fly ash particles, their refractive indices, particularly in the infrared range, are known insufficiently. This is especially true of the optical properties of infrared particles being considered in the high-efficiency combustors. Accurate values of the optical constants of particles, which are significant functions of wavelength, are of particular interest in radiative transfer analysis in power plant because of their significant influence on heat transfer in such applications. Light scattering analyses of the structure and refractive indices of fly ash particles have been proposed as a means of studying the interactions between fly ash particles, flue gases, and the surrounding atmosphere [10]. Moreover, fly ash is a very promising combustion by-product of the coal-fired power station. For example, it can be added into the cement pastes to improve their fluidity where the PSD is an important control parameter [11]. Generally, optical properties and PSD of participating particulate medium are of great importance in power plant systems and have evoked the wide interest of many researchers.
The optical constant of participating medium particle system cannot be directly measured by experiments [12]. Its direct calculation through the measurable physical parameters is also difficult. So its value is usually obtained by inversion algorithm combining with experimental measurements. Optical constant contains real part n and imaginary part k. These two parts, respectively, represent refraction behavior and absorption behavior of participating media. The optical constants play an indispensable role in determining the radiation properties of a medium. Accurate knowledge of the refractive index over a wide range of wavelength is indispensable for many engineering applications. For instance, combustion diagnostic techniques employing laser light scattering require accurate refractive indices in the visible wavelengths. Various environmental, atmospheric, astrophysical, and military applications rely on data for such particles in both the visible and the infrared wavelength. However, even for the same kind of participating particle system, its optical constants may also be variable because of some factors, such as temperature range, geometry size, and wavelength range. Because of its importance for figuring out the radiative transfer mechanism in the particle system, it has attracted significant attention of scholars all over the world [13–19]. In fact, the accurate determination of the refractive index of dispersed media should be considered today as an unsolved problem, open to research.
To date, the commonly used particle size measurement technique can be divided by its principle into the following categories: microscope method, sieving method, sedimentation, and ultrasonic method. For the practical problem of interest, the sieving method is rarely used recently because it has low precision and it is extremely easy to damage the samples in the measurement process [20]. With the development of the microscope, the accuracy of the microscope method improves. However, it is very time consuming. Ultrasonic method is more suitable for higher particle concentration due to its better ability of transmission in the participating medium compared to light. However, for some unstable particles systems, the ultrasonic method is not applicable. Therefore, spectral extinction method is more widely applied because of its wide measurement range, relatively simple measuring principle and device [21, 22]. Literature survey reveals that there is no reported work about estimating of optical constants and PSDs using artificial neural network (ANN). An important advantage of this method in comparison with other inverse heat transfer modeling approaches, such as the Gauss–Newton method [23], conjugate gradient method [24], genetic algorithm [25], and ant colony optimization [26, 27], to name a few, is that detailed knowledge of the geometrical and thermal properties of the system (such as wall conductivity and emissivity) is not necessary [28, 29]. In many cases, the measurement of such physical properties is extremely difficult or even impossible. Moreover, to use all conventional inverse methods, the direct problem must be solved firstly.
This may lead to significant computing errors and excessively time-consuming calculations. ANN solutions are based on experimental data rather than mathematical governing equations [30]. The direct model needs not be solved repeatedly, which will save considerable computation time. Furthermore, genetic algorithm and other evolutionary algorithms require large populations and involve long computing time for convergence, which make them being inappropriate for many inverse heat transfer problems. In ANN, the only iterative part is in training, which is containing simple mathematical relations. Therefore, the proposed method is much faster than conventional iterative methods. In this chapter, a particle swarm optimization (PSO) based back propagation (BP) neural network algorithm was developed to estimate the PSD and optical constant of a particle system. As one kind of intelligence optimization algorithms, the BP multi-layer feed-forward ANN algorithm was proposed by D.E. Rumelhart firstly [31]. In complicated systems with several effective input parameters, ANN can be used to predict output data. The BP neural network algorithm is a stochastic heuristic intelligent optimization technique, which was investigated thoroughly [32, 33]. However, there remain some unsolved problems. For example, as the initial weight coefficients are generated randomly, it may lead to the instability of the algorithm. Furthermore, the algorithm may fall into local optimum. Therefore, in the present work, to improve the performance of the BP neural network algorithm, the PSO algorithm was applied.The remainder of this chapter is organized as follows. First, the principle of the inverse method was introduced briefly. Afterwards, the direct model was presented. Subsequently, the PSDs and optical constant were retrieved separately using the PSO-based BP neural network algorithm. Finally, a technique was proposed to obtain the PSDs and optical constant simultaneously.
Inverse method
BP NEURAL NETWORK
An ANN model includes five parts, i.e. input values (x1, x2, … , xn), thresholds of network, weight (W[w1,1, w1,2,…, w1,N]), transfer function, and output values. The size of weight stands for the importance of input data for forward transmission by indicating the connection degree between the input information and related neurons. The basic unit of neural network is the neuron. In each neuron, the sum of input values is weighted and added with a parameter called bias, and the sum (see Equation 1) is passed through a function, which is called transfer function or activation function. The transfer function calculates the output of a neuron from its input.
Some ANN models include several layers. Each layer includes several neurons and performs a simple process on data. The actions of neurons do not directly depend on the summation but a certain threshold. Only when the summation netj is bigger than the threshold, neurons can be activated by activation function. Otherwise, the neuron will not send signals.
The working process of the BP neural network includes training and prediction. Training is to prepare the network to work for certain system, which is of highly importance. In the training process, the weights are adjusted to obtain better performance of the network. Prediction is to solve the inverse problem by using the prepared network [34]. The detailed principle of the BP neural network is described in Ref. [32], which will not be repeated here.
THE PSO ALGORITHM
The PSO, which has rooted in artificial life and social psychology as well as engineering and computer science, is a random searching technique inspired by bird foraging. It differs from evolutionary computation methods in that the population members, called particles, are flown through the problem hyperspace. The PSO finds the optimum value more quickly than traditional evolutionary algorithms due to the fact that PSO uses a combination of local and global searches with the sharing evolutionary information among the individual particles. Many modifications have been made to improve convergence speed of the basic PSO algorithm and to increase the diversity of particles [35]. The velocity and position of each particle in the basic PSO can be obtained as follows [36]:
where c1 and c2 represent two positive constants called acceleration coefficients; r1 and r2 are uniformly distributed random numbers in the interval [0, 1]. w is the inertia weight factor, which is used to control the impact of the previous velocities on the current velocity. In the present work, w is defined as w = 0.9 – 0.4 × t/tmax, where t and tmax represent the current and maximum generation number, respectively. Xi(t) = [Xi1(t), Xi2(t), …, Xin(t)]T denotes the present location of particle i, which represents a potential solution. Vi(t) = [Vi1(t), Vi2(t), …, Vin(t)]T denotes the present velocity of particle i, which is based on its own and its neighbors’ flying experience. Pi(t) = [Pi1(t), Pi2(t), …, Pin(t)]T is the ‘local best’ in the tth generation of the swarm. The global best position of the swarm is expressed as Pg.
In the present work, an improved multi-strategy PSO (MSPSO) was applied to BP neural network training [37]. Three improvements were made based on the standard PSO. First, the opposition learning [38] strategy is introduced into every iteration of the PSO algorithm. The idea of opposition learning is that when a new individual particle X = (X1, X2,…, Xn) is generated during the evolution process, its opposite position X¯¯¯=(X1¯¯¯¯,X2¯¯¯¯,…,Xn¯¯¯¯¯)X¯=X1¯X2¯…Xn¯ should be considered to test if it is better than X. If xi is in the interval of [ai, bi], Xi¯¯¯¯Xi¯ can be calculated by Xi¯¯¯¯=ai+bi−XiXi¯=ai+bi-Xi. Second, the acceleration coefficient c2 is set to vary with iterations instead of being constant. Generally speaking, c2 represents social ability of the swarm. In the later stage of the iteration, the society guide ability gradually strengthens, which makes the diversity of the swarm decline. It should properly weaken in the later stage of the algorithm. Therefore, c2 can be set as follows [37]:
where c2max and c2min are the maximum and minimum of c2, which are set as 2.05 and 0.4, respectively, in the present work [37].
Finally, it is known that the standard PSO may be trapped into local optimum as the decreasing of the population diversity. Therefore, a multi-start strategy is introduced to the MSPSO, which means that, when the optimal solution is not updated for Smax iterations, it will force the particle swarm to initialize, but the previous optimal solution is recorded, where Smax is the maximum number of stagnation. Note that during the particle initialization, it still reserves the individual best position and global best position. The detailed model flow chart of MSPSO is shown in Figure 1.
MSPSO-BP ALGORITHM
In the present work, the MSPSO is applied to optimize the threshold and weight of the BP neural network. The training process is required before prediction for BP neural network, which is introduced in Section 2.1. Before the optimization of BP neural network, the topological structure of the network should be decided. For a three-layered BP neural network, i.e. to choose the node numbers in input layer, hidden layer, and output layer. The structure of input layer and output layer is mainly determined by the number of input and output data. The node number in hidden layer is of significant importance, which highly affects the training speed and training accuracy. Too many nodes may lead to slow convergence rate and over learning. However, if the node number is not large enough, the network may have small training precision. The values of the node number can be determined using many different ways [39]. After the node number of the hidden layer is determined, the next step is to decide the weight and threshold values, which are optimized by the MSPSO in the present work. Each particle stands for a possible solution of the optimal weight and threshold.
Direct model
A one-dimensional (1-D) absorbing, scattering, and non-emitting particle system was under consideration in the present work (see Figure 2). The left side of the system was exposed to continuous wave laser beam of different wavelengths. For the short interaction time and low laser intensity, the media were assumed to be cold without considering the thermal effect caused by the interaction of the laser and media. The radiative transfer equation in a 1-D particle system can be expressed as follows [40]:
|
∂I(x,θ)∂x=−(κλ+σsλ)I(x,θ)+σsλ2∫0πI(x,θ′)Φλ(θ′,θ)sinθ′dθ′∂Ixθ∂x=−κλ+σsλIxθ+σsλ2∫0πIxθ′Φλθ′θsinθ′dθ′ |
(5) |
where I is the intensity in direction θ at location x. The absorption and scattering coefficients are denoted by κλ and σsλ, respectively. The subscript λ stands for the wavelength of the incident laser. The scattering phase function Φλ(θ′, θ) represents the probability that radiation, which propagates from the incoming direction θ′, will be scattered into the direction θ. According to Mie theory, when the sphere particle is irradiated by a laser beam, the extinction, scattering, and absorption efficiency can be calculated by the following equation [41]:
|
Qext=2χ2Re[∑n=1∞(2n+1)(an+bn)]Qext=2χ2Re∑n=1∞2n+1an+bn |
(6) |
|
Qsca=2χ2∑n=1∞(2n+1)(|an|2+|bn|2)Qsca=2χ2∑n=1∞2n+1an2+bn2 |
(7) |
|
(8) |
where the extinction, scattering, and absorption efficiency are denoted by Qext, Qsca, and Qabs, respectively. Re is the real part of complex number, χ represents the scale parameter that can be expressed as πD/λ for spherical particles. D represents the diameter of particles. λ stands for the wavelength of incident laser. Mie scattering coefficients are denoted by an and bn, which can be written as follows [41]:
|
an=ψ′n(mχ)ψn(χ)−mψn(mχ)ψ′n(χ)ψ′n(mχ)ξn(χ)−mψn(mχ)ξ′n(χ)an=ψn′mχψnχ−mψnmχψn′χψn′mχξnχ−mψnmχξn′χ |
(9) |
|
bn=mψ′n(mχ)ψn(χ)−ψn(mχ)ψ′n(χ)mψ′n(mχ)ξn(χ)−ψn(mχ)ξ′n(χ)bn=mψn′mχψnχ−ψnmχψn′χmψn′mχξnχ−ψnmχξn′χ |
(10) |
where m (= n + ik) represents the optical constant of particle system. ξn(χ) = ψn + iχn, ψn(χ), and χn(χ) are Ricatti–Bessel functions. The details of Mie theory are well demonstrated in Ref. [41]. The absorption coefficient κλ and scattering coefficient σsλ of particle system can be calculated by the following equations:
|
κλ=∫0∞Cabs,λ(D)N(D)dDκλ=∫0∞Cabs,λDNDdD |
(11) |
|
σsλ=∫0∞Csca,λ(D)N(D)dDσsλ=∫0∞Csca,λDNDdD |
(12) |
where Cabs, λ and Csca, λ stand for absorption and scattering cross section, respectively. For sphere particles, absorption and scattering cross section can be calculated by using Mie theory [41]; N(D) represents the number density of particles with diameter D.
It is known that the PSDs can be described by certain distribution function, among which the commonly used parameter distribution function is as follows [42]:
|
fR−R(D)=σD×(DD)σ−1×exp[−(DD)σ]fR−RD=σD×DDσ−1×exp−DDσ |
(13) |
|
fS−N(D)=12π√σ×exp[−(D−D)22σ2]fS−ND=12πσ×exp−D−D22σ2 |
(14) |
|
fL−N(D)=12π√Dlnσ×exp[−(lnD−lnD)22(lnσ)2]fL−ND=12πDlnσ×exp−lnD−lnD22lnσ2 |
(15) |
where D¯¯¯D¯ represents the characteristic diameter parameter, and σ is the dispersion ratio. The volume frequency distribution is denoted by f(D).
Results and discussion
To demonstrate the validity of the proposed MSPSO-BP algorithm in the inverse radiative analysis, several different test cases are considered in this section. Firstly, the PSDs are retrieved using multi-wavelength method. Secondly, the optical constant is obtained for L–N distribution. The influence of measurement errors is considered. Finally, the PSDs and optical constant are retrieved simultaneously using a secondary optimization method. Each case was implemented using FORTRAN code, and the developed program was executed on an Intel Core i7 PC.
THE RETRIEVAL OF PSDS
Multi-wavelength method [43] has been applied to obtain the PSDs accurately and effectively. Therefore, to test the influence of wavelength number on the retrieval accuracy and efficiency, different amounts of wavelengths were used in the present work, which, to be specific, are one, two, and three wavelengths. For the single wavelength, the wavelength was set as 0.55 μm. For the two wavelengths, 0.55 and 0.63 μm were applied. And 0.55, 0.633, 0.649 μm were used for the three-wavelength method. The optical constants (n, k) corresponding to different wavelengths were set as (1.381, 0.00426), (1.377, 0.00162), and (1,376, 0.00504). The retrieval of the PSDs was solved by minimizing the objective function, which was defined as the sum of the square residuals between the predicted and the expected reflectance and transmittance. Therefore, the objective function can be expressed as follows:
|
Fobj=∑i=1n[(Rpre(i)−Rexp(i)Rexp(i))2+(Tpre(i)−Texp(i)Texp(i))2]Fobj=∑i=1nRprei−RexpiRexpi2+Tprei−TexpiTexpi2 |
(16) |
The number density of particle system was set to a fixed value N0 = 5.0 cm−3. The range of parameters in monomodal distributions are set as D¯¯¯∈(0.1,1.1)D¯∈0.11.1 and σ ∈ (1.4, 3.4). To obtain adequate training samples for training network, 1000 sets of data were generated randomly within the initial range. Nine hundred sets were randomly selected as training samples, and others were chosen as prediction samples to test the performance of trained network. Afterwards, the well-trained and well-tested network will be used to retrieve the PSD of the particle system. The control parameters of MSPSO-BP are shown in Table 1. Also, the relative errors were applied to demonstrate accuracy retrieval results, which can be expressed as follows:
|
ε=∣∣PRE−EXPEXP∣∣×100%ε=PRE−EXPEXP×100% |
(17) |
where PRE and EXP stand for the predicted and expected values of the results, respectively.
In the present work, three commonly used PSDs were applied, i.e. R–R, S–N, and L–N (see Equations 13–15). The training results, which are the mean relative errors of all the prediction samples, for different PSDs and different number of incident wavelengths are shown in Table 2.
It can be seen that MSPSO-BP has excellent prediction performance for the three kinds of PSD functions. Besides, with the number of wavelengths increasing, the error between predicted values and expected values of prediction samples gradually decreases, which means that the multiple-wavelength method can improve the prediction accuracy effectively. Also, it is noticeable that the calculation time is also increasing with the number of wavelengths used, which is caused by the increasing calculation times of the direct problem.
After the test of the network is finished, the network can be used to retrieve the PSD of the particle system. The results of the 100 sets of test data are illustrated in Figure 3. It demonstrates clearly that increasing the incident wavelength will improve the accuracy the prediction results. To show the performance of the MSPSO-BP algorithm, a case was presented to retrieve the PSD parameters of the L–N distribution. The true value of PSD parameters is set as σ = 2.0 and D¯¯¯=0.6D¯=0.6. Each case was run for 10 times. The inversion results are shown in Table 3. It is evident that the retrieval results are accurate when using two- and three-wavelength methods. However, the results of single wavelength are not acceptable. This can be contributed to the insufficient information supplied by the single wavelength method. Therefore, it is suggested that more than one incident wavelength should be used in the retrieval of PSDs.
The number of training sample is a very important parameter for the BP neural network. It should be dealt with very carefully. Therefore, different numbers of training samples are applied to the inversion of PSDs using three-wavelength method. The numbers of training samples are set as 600, 700, 800, and 900, respectively. The retrieval results are listed in Table 4.
It can be seen that by increasing the amount of the training samples, more accurate retrieval results can be obtained, which can be seen more clearly in Figure 4. Network training time is obviously related to the number of training data. Training the network with more data tends to spend more time. Therefore, to save computation time, less training samples should be applied. However, it will reduce the accuracy of the retrieval results at the same time (see Figure 4). Therefore, when using the neural network-based algorithm for an inverse problem, the number of the training samples should be decided very carefully
THE ESTIMATION OF OPTICAL CONSTANT
In this section, MSPSO-BP algorithm is applied to retrieve the optical constant, i.e. n and k. It is known that the optical constant is different for different wavelengths. Multi-wavelength method is no longer effective in this situation as more wavelengths mean more unknown optical constants. Therefore, other methods should be applied to obtain more information about the particle system. In the present work, different thicknesses of particle system were applied. To test the effectiveness of the multi-thickness method, the optical constant of a 1-D sample with thickness of 0.1, 0.2, and 0.3 m was retrieved.
First, the MSPSO-BP should be trained. The interval of the optical constant was set as n ∈ (1.3, 1.6) and k ∈ (0.001, 0.01). The PSD was set as L–N distribution and D¯¯¯=1.0D¯=1.0 and σ = 2.0. The incident wavelength was set as 0.55 μm, which means n = 1.381 and k = 0.00426. Other parameters were set the same as last section. One thousand sets of n and k were generated within the interval randomly, 900 of which were used as training, and others were used to test performance of the trained network. The test results are shown in Figure 5.
It is evident that the predicted values and expected values of n and k are in a good agreement, which indicates that the trained network has an excellent ability for retrieving optical constant. Therefore, the trained network will be used to conduct the inversion of optical constant. The optical constant that needs to be estimated was set as (n, k) = (1.4, 0.006). Considering the fact that the inevitable measurement errors will occur in practical measurement, random errors were added to the original reflectance and transmittance. The results are shown in Table 5. It can be seen that even when the measurement noise is 5%, the retrieval results of the optical constants are still acceptable.
THE RETRIEVAL OF OPTICAL CONSTANT AND PSDS
In this section, a secondary optimization method [44] was applied to obtain the optical constant and PSDs simultaneously by combining the multi-wavelength and multi-thickness methods. The diffuse transmittance, diffuse reflectance, and collimated transmittance were applied. The details of this method are described in Ref. [44] and will not be repeated here. The wavelengths of two incident beams were set as 0.55 and 0.65 μm. And its corresponding optical constants were 1.432–0.00753i and 1.426–0.00794i, respectively. Two samples, which were basically the same except for their thicknesses, are used. The original PSDs that need to be retrieved are shown in Table 6. Particle number density was set as 5.0 cm−3. The searching range of the optical constant and PSDs was set as n ∈ (1.3, 1.6) and k ∈ (0.001, 0.01), D¯¯¯∈(0.1,1.1)D¯∈0.1,1.1, and σ ∈ [1.4, 3.4]. Ten thousand sets of optical constant and PSDs were generated randomly within the parameter range. After network was trained, it was applied to obtain the optical constant and PSDs. Each case was run for 10 times. The retrieval results are shown in Table 7.
It can be seen that for the same distribution, the retrieval results of optical constant are much smaller than those of the PSDs. The optical constant can be retrieved with a tolerable error. However, the average retrieved results of the PSDs, i.e. D¯¯¯D¯ and σ, are not tolerable, and the maximum retrieval error can reach more than 20%. The reason the relative errors of the PSDs are so much larger than those of nand k is that the sensitivity coefficients of D¯¯¯D¯ and σ are much lower than those of n and k [44]. In other words, the radiative signals are insensitive to the PSDs, which mean even when the relative errors of the retrieved values of D¯¯¯D¯ and σ are very large, the deviation in the calculated reflectance and transmittance signals can be ignored.
It can be concluded that the influences of the PSDs on the diffused reflectance and transmittance signals are less significant than that of the optical constant. Even though the retrieved results of the PSDs were not accurate, those of the complex refractive index were relatively satisfactory. Hence, the accuracy of the retrieved PSDs needed to be improved which meant much more information, which was sensitive to the PSDs, should be included in the inverse procedure. It was found that the collimated transmittance was more sensitive to D¯¯¯D¯ and σ [44]. Therefore, the collimated transmittance signal was applied to improve the inversion accuracy of the PSDs. Based on the above analysis, a secondary optimization was performed in which D¯¯¯D¯ and σ were the unknown parameters that needed to be retrieved, and the optical constant retrieved by the first optimization was applied as the original values. The inverse procedure is shown in Figure 6. The results of the secondary optimization are shown in Table 8. It is evident that the estimated results of the secondary optimization were much more accurate than those of the first optimization, which proved that the applied inversion method is suitable for the simultaneous estimation of the complex refractive indexes and PSDs.
Conclusion
A multi-strategy PSO was applied to improve the performance of the BP multi-layer feed-forward ANN algorithm. The proposed MSPSO-based ANN method was introduced to solve the inverse radiation problems for the first time. Three commonly used PSD functions in a 1-D particle system were retrieved using the proposed algorithm. In addition, the optical constant was estimated, and the influence of measurement errors upon the precision of the estimated results was also investigated. Finally, the proposed algorithm was applied to the simultaneous estimation of the PSDs and optical constant using the multi-wavelength and multi-thickness method. All the retrieval results show that the proposed MSPSO-BP ANN can be applied to estimate the PSDs and optical constant in 1-D absorbing, scattering, and non-emitting system accurately even with noise data. Meanwhile, the secondary optimization can be recognized as an efficient way to obtain the PSDs and optical constant simultaneously by applying the multi-wavelength and multi-thickness method.