5th World Congress on Industrial Process Tomography, Bergen, Norway
Intel® MKL Trust-Region Solvers in Mathematical Problems of Geophysics
Nikolay I. Gorbenko, Nikita A. Shustrov, Alexander V. Avdeev
Intel Corporation R&D Lab, Novosibirsk, Russia, e-mail:
ABSTRACT
Trust-Region (TR) algorithms are relatively new iterative algorithms for solving nonlinear optimization problems.High efficiency of TR methods was demonstrated in a number of recent papers and books. They have global convergence and local super convergence, which differ them from line search methods, used quite often for solving Inverse Problems. TR techniques is implemented in a number of well-known SW libraries such as IMSL, TAO, GALAHAD, LANCELOT, etc.
Let us explain the main difference of TR method from classical Newton one. Assume we have a current guess of the solution for the optimization problem. An approximate model can be constructed near the current point. Solution of the approximate model can be taken as the next iterate point. The classical line search algorithms also solve approximate models to obtain search directions. However, in TR algorithms, the approximation model is only “trusted” in the region near current iterate. This seems reasonable, because for general nonlinear functions local approximate models (such as linear and quadratic approximations) can only fit the original function locally. The region that the approximate model is trusted is called “trust region”. The trust region is adjusted from iteration to iteration, i.e. if computations indicate the approximation model fit the original problem quite well, the trust region can be enlarged. Otherwise when the approximation works not good enough the trust region will be reduced.
TR solvers developed and implemented by the authors are based on Intel MKL (Math Kernel Library). There are versions of TR for Intel IA32, EM64T, IA64 platforms (Windows and Linux operating systems). In addition to Fortran-interface, TR solvers include C-language interface for all functions and routines. TR solvers are implemented with OpenMP support and can be used in multiprocessing mode.
Keywords Nonlinear optimization, Trust-Region solvers, Intel MKL, Oil & Gas Problems
1 Trust-Region algorithm
Trust-Region algorithms allow to solve (Conn et al. 2000):
– nonlinear least squares problem (without and with bound constraints);
– system of nonlinear equations (without and with bound constraints);
– minimization of functional (without constraints).
Let us explain some details of TR approach through the example of a nonlinear least squares problem:
(1)
where is continuously differentiable function.
The step for TR method is the solution of the subproblem
, where <k (2)
that is the approximation of the goal function in a neighborhood of current point х, and - Jacobian of .
*Trademark Information Intel, Itanium, Dual-Core Xeon and Xeon are trademarks of Intel Corporation in the U.S. and other countries. Intel's trademarks may be used publicly with permission only from Intel. Fair use of Intel's trademarks in advertising and promotion of Intel products requires proper acknowledgement.
*Other names and brands may be claimed as the property of others.
TR pseudo code for solving a nonlinear least squares problem:
Step 1: Choose initial guess
Step 2: Solution of (2) gives value of dk
If then stop;
Computes
Step 3: Choose
Define
Step 4: k = k +1, go to Step 2.
The main problem of TR method is to solve the subproblem (2).
The following lemma is well known (Conn et al. 2000):
Lemma: A vector is the solution of the problem
where is a symmetric matrix, and , if and only if there exits such that and is semi-definite, and
The case where has zero eigenvalues is called “hard” case and can be written by , where is a vector in the null space of . Unless in the hard case is the unique solution of the following equation:
For solving this equation we use Newton’s method to calculate.
Consider the TR usage for solving the nonlinear least square problem with simple bound constraints.
Let . The first-order necessary conditions of vector-minimizer of (1) are stated as
(3)
where is the diagonal scaling matrix (Coleman, 1996):
(4)
with
Let , and we write for the strict nonempty interior of . To update the current iterate we need the solution of an elliptical trust-region subproblem and computation of a search step at each iteration.
Let and trust-region size be given. We consider the elliptical trust-region subproblem:
(5)
where is quadratic model for at , i.e.
(6)
We let and to the solution and Cauchy point of the trust-region subproblem (5), respectively. Namely, is the minimizer of along the scaling steepest descent direction subject to the satisfaction of the trust-region bound, i.e.
(7)
where
The function is quadratic in , and the unconstrained minimizer has the form . is given by
(8)
The search step used in TR algorithm is defined by a linear combination of two vectors and :
(9)
where is suitable scalar. The vector in (9) is equal to q, where is given component-wise by
(10)
At the k-th iteration we consider the second degree polynomial . On the basis of the current information, models the restriction of to the points for varying from 0 to 1. Therefore, we compute the value that minimizes this model and choose according fixed lower and upper bounds, . This means that we set
.
Now we consider the choose of the vector . If the point we simply let . Otherwise, is a simple scaling of obtained by a step-back along it. In other words, we let to be the stepsize along to the boundary:
The step is given by
(11)
for some fixed
The scalar in (9) is chosen in order to satisfy the following conditions
(12)
for the given constant and
(13)
To force conditions (12), (13), we compute the scalar in (9) according to the following strategy. If , we simply set , i.e. we choose as the candidate step.
Otherwise, we fix so that is the point which lies on the segment from to and satisfies
.
In fact, is given by (9) setting equal to the following scalar
(14)
where
A good agreement between the model function and objective function is ensured for the following standard condition
(15)
with the given constant . If this condition is not met, we reject and the trust-region size is adjusted be successive reduction so that satisfies the accuracy requirement (15).
Algorithm of k-th iteration:
Let be given
Compute the matrix by (3). Set , set .
Repeat
Set .
Find the solution for problem (5).
Compute the Cauchy step by (7).
Compute and by (9), (10), (11) respectively.
Find such that satisfies (12), (13).
Set .
2 Intel Trust-Region solvers performance
Intel MKL TR solvers are tested carefully on standard test sets for the unconstrained and constrained optimization from CUTEr collection (CUTEr is testing environment for optimization and linear algebra solvers. The package contains a collection of test problems, along with Fortran77, Fortran90/95 and Matlab tools).
TR solvers are implemented with OpenMP support and can be used in multiprocessing mode (fig.1):
Fig.1. Intel TR solvers scaling
TR solvers show excellent performance vs. VNI IMSL Fortran library v5.0 and TAO v1.8.1. On benchmark tests TR give up to 10 times speed-up vs. IMSL (Fig.2) and up to 5 times speed-up vs. TAO.
Fig.2. Intel TR solvers vs. IMSL
3 Intel Trust-Region solvers usage in mathematical problems of Geophysics
TR solvers were applied for solving the inverse problem of pore pressure evaluation (Madatov et al. 1995), namely – it was required to reconstruct some lythological values of strata (initial porosity, constants of compaction, “source” term and specific surface) using real porosity and pore pressure data.
The numerical solutions of this inverse problem (Lavrentiev, et al, 2003) were obtained with the required precision and rapid convergence (see Fig.3, 4).
Fig. 3. Inverse problem of pore pressure evaluation: Porosity curve reconstruction
Fig. 4. Inverse problem of pore pressure evaluation: Pore pressure curve reconstruction
4 REFERENCES
CONN A. R., GOULD N. I.M., TOINT P. L., (2000), Trust-Region Methods. SIAM Society for Industrial & Applied Mathematics, Englewood Cliffs, New Jersey.
LAVRENTIEV M.M., AVDEEV A.V., LAVRENTIEV Jr. M.M., PRIIMENKO V.I., (2003) Inverse Problems of Mathematical Physics. VSP Publ., The Netherlands, 272 p.
COLEMAN T.F., LI Y., (1996) An interior trust-region approach for minimization subject to bounds, SIAM J. Optim. 6 (3) 418–445.
CUTEr, testing environment for optimization and linear algebra solvers, http://hsl.rl.ac.uk/cuter-www/
MADATOV A.G., HELLE H.B., and SEREDA V.I. (1995) , A modelling and inversion technique applied to pore pressure evaluation, 57th EAEG Conference, Glasgow, May 29-June 2
INTEL MATH KERNEL LIBRARY 10.0, http://www.intel.com/cd/software/products/asmo-na/eng/perflib/mkl/index.htm
7