This paper considers an application of metaheuristic algorithms for solving the problem of convection velocity identification in the convection-diffusion transport model. Algorithms based on numerical minimization of the parameter identification criterion for a discrete linear stochastic model using the simulated annealing and a genetic algorithm are proposed. The log-likelihood function is used as the identification criterion. Numerical experiments were conducted to compare the computational properties of the proposed algorithms.
Convection-diffusion transport models are an indispensable tool for describing various natural and anthropogenic processes [1], [2]. These models contain parameters that must be specified to uniquely determine the solution of boundary value problems, but in practice, situations often arise where some of these parameters are unknown or given approximately and they need to be determined or refined. Such problems belong to the class of inverse problems for the models of matter transfer.
In the simplest one-dimensional case, the convection-diffusion transport model can be described by
equation (
When solving a wide range of problems in ecology, geophysics, seismology, and other areas, the
problem of determining (identifying) the convection velocity in the convection-diffusion transport
model often arises. Depending on the equation under consideration and the boundary conditions,
various methods can be used to solve this problem. In [3], [4], to find the parameters of equation (
In this paper, we propose the use of metaheuristic algorithms for numerical optimization to find the
optimal estimate of a parameter of a discrete linear stochastic model with respect to a given criterion
of identification quality. As in [4], this model with noisy measurements is constructed from the
original model (
The choice of metaheuristic algorithms to optimize a parameter identification criterion is due to the fact that deterministic numerical methods applicability is guaranteed under the conditions of convergence theorems [5]. The convergence of numerical methods is influenced by various factors, including a good choice of the initial approximation.
If the initial approximation is not chosen correctly, the exact algorithm for finding the estimates of the parameters may diverge, which means that it is impossible to solve the identification problem. Metaheuristic algorithms can be used both directly for solving the problem of minimizing the identification criterion, and for finding a good initial approximation for exact methods. A similar approach was applied by the authors in [6–9] and proved its efficiency.
We want to move from the model (
In computational practice, for the numerical solution of non-stationary problems for
convectiondiffusion equations, two- and three-layer finite-difference schemes are most often used. To obtain a
discrete linear stochastic model, consider in the plane a regular grid (
Let us denote
difference scheme for (
. Let us write (
The desired discrete linear stochastic model can be represented in the following form:
Note that the first equation in model (12) is deterministic, the initial value of the state vector is (
Suppose that the diffusion coefficient and the characteristics of the noise in the measurer are known, and the steps of the space-time grid and are given, then the unknown parameter of model (12) to be determined is the convection velocity , on which coefficients and of equation (10) and, consequently, matrix depend. Briefly, the model (12) can be written in the form: where
.
Let us consider the problem of parameter identification of model (13) from noisy measurements data for estimating an unknown parameter. The parameter identification problem consists in finding an unknown parameter from known input signals and the output observation data in accordance with the chosen identification criterion . In this case, the problem of estimating an unknown parameter requires solving the nonlinear programming problem with constraints (12) (13) (14) where (the domain of ).
In this paper, we use metaheuristic algorithms to solve problem (14). Metaheuristic is a high-level search strategy for finding solutions, applicable to a wide range of optimization tasks. Metaheuristics have the following properties: they are based on fairly simple ideas, for example, imitating biological or physical processes, they are problem-independent, practically all of them are nondeterministic.
Most metaheuristic optimization algorithms can be divided into two large groups according to the method of obtaining a solution: trajectory and population-based ones. In trajectory algorithms, the process of finding a solution can be considered as a movement between individual solutions of a problem while in population-based algorithms a group of solutions called population changes in the process of finding the solution.
One of the most popular trajectory algorithms used in solving global optimization problems is the simulated annealing (SA) method. A key feature of the method is the use of a control parameter– a temperature, which allows controlling the nondeterministic process of solution search. As a rule, the temperature decreases during the operation of the algorithm according to a certain law, starting with some initial value. At each iteration of the algorithm, the randomly generated new solution from the neighborhood of the current solution is taken with probability 1 if it improves it, and with probability less than 1, if worsens, and the probability of making the worst decision decreases with a decreasing temperature. The quality of decisions is estimated using a cost function (integer or real).
Simulated Annealing (SA)
The genetic algorithm (GA) is a popular version of the so-called evolutionary optimization algorithms based on the simulation of natural selection processes. In evolutionary algorithms, the quality of solutions is estimated using the fitness function, and the main idea of algorithms is that solutions with the best values of a given function “survive” in the course of evolution.
In GA, at each iteration of the evolutionary process, a new population is obtained from the current population using one or more genetic operators successively. The most common genetic operators are the crossover (recombination) which is used to generate the descendant solutions from the parent solutions and the mutation – an accidental change in the solution.
Genetic Algorithm (GA) 1: Population InitialPopulation() 2: for all Population do 3: EvaluateFitness( ) 4: end for 5: while not StopCondition() do 6: Parents SelectParents(Population) 7: Offspring Crossover(Parents) 8: Offspring Mutation(Offspring) 9: for all Offspring do 10: EvaluateFitness( ) 11: end for 12: Population 13: end while 14: Solution ChooseBestOf(Population) Output: Solution.
UpdatePopulation(Population
Offspring) The algorithms SA and GA are discussed in more detail, for example, in [10].
Since (13) is a discrete linear stochastic model with Gaussian noise, it is advisable to select the identification criterion (14) in the form of a negative logarithmic likelihood function as the cost/fitness function for implementing metaheuristic algorithms [11] where the residual vector and its covariance matrix are calculated from the known Kalman filter equations [12]:
A. Time update.
For the Kalman filter computes extrapolated estimates through the temporal update from to as with and the covariance matrices for a given value of the parameter
(15)
in (
B. Measurement update.
For the Kalman filter computes the so-called filtered (i.e., measurement updated) estimates . They are obtained through the measurement update using with noise covariance , as , , with filter gain and filtered estimates covariance matrices
At present, when solving practical problems with the use of a computer, it is preferable to apply square-root and UD-implementations that are numerically stable against machine round-off errors instead of the conventional form of the Kalman algorithm [12, Chapter 6].
The maximum likelihood method consists in optimizing criterion (15) with respect to the system parameter . It is often used in practice to solve parameter identification problems of discrete linear stochastic systems [12, 13].
Consider the following problem: where is the concentration of the pollutant in one-dimensional flow, (17) is the initial distribution of the pollutant and boundary conditions (18) correspond to the case of two absorbing walls.
The exact solution of (16)–(18) has the form:
Suppose that the value of parameter in equation (16) is known and it is required to determine the value of parameter , provided that noisy measurements from sensors located at nodes of some regular (16) (17) (18) (19) grid are available. Let a grid with 10 nodes ( ) be given on the axis, and the time step . , . The plot of the exact solution for this case is shown in figure 1.
To simulate noisy measurements we add random errors (white noise) to the values of the exact solution (19) at the grid nodes. The covariance matrix of the noise is , where is the identity matrix and is the known variance. An example of a noisy solution is shown in figure 2.
Figure 3 demonstrates the averaged plot of criterion (15), obtained from the results of 100 experiments. To minimize it, we have used functions simulannealbnd() and ga() from the Global Optimization Toolbox of the MATLAB system. Numerical experiments were conducted on a hardwaresoftware platform: Intel Core 2 Quad Q6600 @ 2.40 GHz, 4 Gb RAM, Microsoft Windows 10 Pro x64, MATLAB R2017a.
For each of the algorithms, the corresponding cost/fitness and output functions were written. Basic settings of the algorithms are given in table 1. As a stopping criterion for both algorithms, a time limit of 5 seconds was used. The search for solutions was carried out on the interval [0; 5], the number of steps in time (i.e., the number of time series) is . Table 2 shows the results of computational experiments for different values of the noise variance . For each value of , a series of 100 experiments was conducted and for each series, the following values were calculated: mean value of the identified parameter (Mean), mean absolute percentage error (MAPE) and root mean squared error (RMSE). The obtained results show that, with the selected settings, both algorithms allow identifying the convection velocity with an acceptable accuracy. At high noise values, the quality of parameter identification is approximately the same for both algorithms, with a low noise level, the simulated annealing method yields smaller MAPE and RMSE errors.
This paper demonstrates the practical applicability of metaheuristic algorithms for solving the problem of parameter identification in the model of convection-diffusion transport. The convection velocity was considered as an unknown parameter of the model. The problem solution was obtained with the use of the maximum likelihood method, in which the simulated annealing method and the genetic algorithm were used to numerically minimize the negative logarithmic likelihood function.
Numerical experiments were conducted in MATLAB. The results obtained make it possible to conclude that the application of metaheuristic algorithms is expedient since it allows to obtain acceptable estimates of the parameter for different levels of noisy measurements.
Further research will focus on the construction and software implementation of new exact and hybrid algorithms for parameter identification of convection-diffusion transport models and their application to real life problems. [10] Hromkovič J. 2002 Algorithmics for Hard Problems: Introduction to Combinatorial Optimization, Randomization, Approximation and Heuristics ed 2 (Springer-Verlag Berlin, Heidelberg) p 557 [11] Åström K J 1980 Maximum Likelihood and Prediction Error Methods Automatica 16 5 pp 551574 [12] Grewal M S and Andrews A P 2001 Kalman Filtering: Theory and Practice Using MATLAB ed 2 (New Jersey: Prentice Hall) p 401 [13] Gibbs V P 2011 Advanced Kalman Filtering, Least-squares and Modeling: Practical Handbook (Hoboken, New Jersey: John Wiley & Sons, Inc.) p 632