Gathers Solver functions. More...
Files | |
file | OptimAux.h |
File that contains auxiliary functions for optimization. | |
file | OptimSA.h |
Function to solve an optimization problem using the Simulated Annealing method. | |
file | OptimTN.h |
Function to solve an optimization problem using the Truncated Newton method. | |
Classes | |
class | Reconstruction |
To perform various reconstruction operations. More... | |
class | Eigen |
Class to find the smallest eigenvalues and corresponding eigenvectors in a generalized eigenproblem using the Bathe subspace iteration method. More... | |
class | Iter< T_ > |
Class to drive an iterative process. More... | |
class | Prec< T_ > |
To set a preconditioner. More... | |
class | Precond< T_ > |
An abstract class from which derive all preconditioning classes. More... | |
Macros | |
#define | MAX_NB_EQUATIONS 5 |
Maximum number of equations. More... | |
#define | MAX_NB_INPUT_FIELDS 3 |
Maximum number of fields for an equation. More... | |
#define | MAX_NB_MESHES 10 |
Maximum number of meshes. More... | |
#define | TIME_LOOP(ts, t, ft, n) |
A macro to loop on time steps to integrate on time ts : Time step t : Initial time value updated at each time step ft : Final time value n : Time step index. | |
#define | TimeLoop |
A macro to loop on time steps to integrate on time. More... | |
#define | IterationLoop while (++theIteration<MaxNbIterations && Converged==false) |
A macro to loop on iterations for an iterative procedure. More... | |
Functions | |
ostream & | operator<< (ostream &s, const Muscl3DT &m) |
Output mesh data as calculated in class Muscl3DT. | |
template<class T_ , class M_ , class P_ > | |
int | BiCG (const M_ &A, const P_ &P, const Vect< T_ > &b, Vect< T_ > &x, int max_it, double &toler, int verbose) |
Biconjugate gradient solver function. More... | |
template<class T_ , class M_ , class P_ > | |
int | BiCGStab (const M_ &A, const P_ &P, const Vect< T_ > &b, Vect< T_ > &x, int max_it, double &toler, int verbose) |
Biconjugate gradient stabilized solver function. More... | |
template<class T_ , class M_ , class P_ > | |
int | CG (const M_ &A, const P_ &P, const Vect< T_ > &b, Vect< T_ > &x, int max_it, double &toler, int verbose) |
Conjugate gradient solver function. More... | |
template<class T_ , class M_ , class P_ > | |
int | CGS (const M_ &A, const P_ &P, const Vect< T_ > &b, Vect< T_ > &x, int max_it, real_t &toler, int verbose) |
Conjugate Gradient Squared solver function. More... | |
template<class T_ , class M_ , class P_ > | |
int | GMRes (const M_ &A, const P_ &P, const Vect< T_ > &b, Vect< T_ > &x, size_t m, int max_it, double &toler, int verbose) |
GMRes solver function. More... | |
template<class T_ , class M_ > | |
int | GS (const M_ &A, const Vect< T_ > &b, Vect< T_ > &x, real_t omega, int max_it, real_t &toler, int verbose) |
Gauss-Seidel solver function. More... | |
template<class T_ , class M_ > | |
int | Jacobi (const M_ &A, const Vect< T_ > &b, Vect< T_ > &x, real_t omega, int max_it, real_t toler, int verbose) |
Jacobi solver function. More... | |
void | BCAsConstraint (const Mesh &m, const Vect< real_t > &bc, Vect< real_t > &up, Vect< real_t > &low) |
To impose Dirichlet boundary conditions in an optimization problem. If such conditions are to present, this function has to be invoked by giving on input bc(i) as the value to impose for the i -th optimization variable. More... | |
template<class OPT_ > | |
int | OptimSA (OPT_ &theOpt, Vect< real_t > &x, real_t &rt, real_t &eps, int &ns, int &nt, int &neps, int &maxevl, Vect< real_t > &lb, Vect< real_t > &ub, Vect< real_t > &c, int &msg_lvl, int &seed1, int &seed2, real_t &t, Vect< real_t > &vm, Vect< real_t > &xopt, real_t &fopt, int &nacc, int &nfcnev, int &nobds) |
Simulated annealing optimization solver. More... | |
template<class OPT_ > | |
int | OptimTN (OPT_ &theOpt, Vect< real_t > &x, Vect< real_t > &low, Vect< real_t > &up, Vect< int > &pivot, int max_it, real_t toler, int msg_lvl) |
Truncated Newton optimization solver. More... | |
template<class M_ , class P1_ , class P2_ > | |
int | QMR (const M_ &A, const P1_ &P1, const P2_ &P2, const Vect< real_t > &b, Vect< real_t > &x, int max_it, real_t &toler, int verbose) |
QMR solver function. More... | |
template<class T_ , class M_ > | |
int | Richardson (const M_ &A, const Vect< T_ > &b, Vect< T_ > &x, real_t omega, int max_it, real_t toler, int verbose) |
Richardson solver function. More... | |
template<class T_ , class MAT_A_ , class MAT_D_ , class MAT_U_ > | |
void | Schur (MAT_A_ &A, MAT_U_ &U, MAT_U_ &L, MAT_D_ &D, Vect< T_ > &b, Vect< T_ > &c) |
Solve a linear system of equations with a 2x2-block matrix. More... | |
template<class T_ , class MAT_A_ , class MAT_D_ , class MAT_U_ > | |
void | Schur (Vect< MAT_A_ > &A, Vect< MAT_U_ > &U, Vect< MAT_U_ > &L, MAT_D_ &D, const Vect< Vect< T_ > > &b, const Vect< T_ > &c, Vect< Vect< T_ > > &x, Vect< T_ > &y) |
Solve a linear arrow block system by block factorization. More... | |
template<class T_ , class M_ > | |
int | SSOR (const M_ &A, const Vect< T_ > &b, Vect< T_ > &x, int max_it, double toler, int verbose) |
SSOR solver function. More... | |
Detailed Description
Gathers Solver functions.
Macro Definition Documentation
#define MAX_NB_EQUATIONS 5 |
Maximum number of equations.
Useful for coupled problems
#define MAX_NB_INPUT_FIELDS 3 |
Maximum number of fields for an equation.
Useful for coupled problems
#define MAX_NB_MESHES 10 |
Maximum number of meshes.
Useful for coupled problems
#define TimeLoop |
A macro to loop on time steps to integrate on time.
It uses the following global variables defined in OFELI: theStep, theTime, theTimeStep, theFinalTime
#define IterationLoop while (++theIteration<MaxNbIterations && Converged==false) |
A macro to loop on iterations for an iterative procedure.
It uses the following global variables defined in OFELI: theIteration, MaxNbIterations, Converged
Function Documentation
int BiCG | ( | const M_ & | A, |
const P_ & | P, | ||
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
int | max_it, | ||
double & | toler, | ||
int | verbose | ||
) |
Biconjugate gradient solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] P Preconditioner (Instance of abstract class P_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] max_it Maximum number of iterations. toler [in,out] Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. verbose [in] Information output parameter - 0: No output
- 1: Output iteration information,
- 2 and greater: Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations,
- Template Parameters
-
<T_> Data type (double, float, complex<double>, ...) <M_> Matrix storage class <P_> Preconditioner class
int BiCGStab | ( | const M_ & | A, |
const P_ & | P, | ||
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
int | max_it, | ||
double & | toler, | ||
int | verbose | ||
) |
Biconjugate gradient stabilized solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] P Preconditioner (Instance of abstract class P_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] max_it Maximum number of iterations. [in,out] toler Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. [in] verbose Information output parameter - 0: No output
- 1: Output iteration information,
- 2 and greater: Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations,
- Template Parameters
-
<T_> Data type (double, float, complex<double>, ...) <M_> Matrix storage class <P_> Preconditioner class
int CG | ( | const M_ & | A, |
const P_ & | P, | ||
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
int | max_it, | ||
double & | toler, | ||
int | verbose | ||
) |
Conjugate gradient solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] P Preconditioner (Instance of abstract class P_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] max_it Maximum number of iterations. [in] toler Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. [in] verbose Information output parameter - 0: No output
- 1: Output iteration information,
- 2 and greater: Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations,
- Template Parameters
-
<T_> Data type (double, float, complex<double>, ...) <M_> Matrix storage class <P_> Preconditioner class
int CGS | ( | const M_ & | A, |
const P_ & | P, | ||
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
int | max_it, | ||
real_t & | toler, | ||
int | verbose | ||
) |
Conjugate Gradient Squared solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] P Preconditioner (Instance of abstract class P_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] max_it Maximum number of iterations. [in,out] toler Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. [in] verbose Information output parameter -
0
: No output -
1
: Output iteration information, -
2
and greater: Output iteration information and solution at each iteration.
-
- Returns
- Number of performed iterations,
- Template Parameters
-
<T_> Data type (real_t, float, complex<real_t>, ...) <M_> Matrix storage class
<P_> Preconditioner class
int GMRes | ( | const M_ & | A, |
const P_ & | P, | ||
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
size_t | m, | ||
int | max_it, | ||
double & | toler, | ||
int | verbose | ||
) |
GMRes solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] P Preconditioner (Instance of abstract class P_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] m Number of subspaces to generate for iterations. [in] max_it Maximum number of iterations. [in,out] toler Tolerance for convergence (measured in relative weighted 2-Norm) on input. Actual convergence tolerance on output [in] verbose Information output parameter (0: No output, 1: Output iteration information, 2 and greater: Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations,
- Template Parameters
-
<T_> Data type (double, float, complex<double>, ...) <M_> Matrix storage class <P_> Preconditioner class
int GS | ( | const M_ & | A, |
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
real_t | omega, | ||
int | max_it, | ||
real_t & | toler, | ||
int | verbose | ||
) |
Gauss-Seidel solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] omega Relaxation parameter. [in] max_it Maximum number of iterations. [in] toler Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. [in] verbose Information output parameter - 0: No output
- 1: Output iteration information
- 2 and greater: Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations
- Template Parameters
-
<T_> Data type (real_t, float, complex<real_t>, ...) <M_> Matrix storage class
int Jacobi | ( | const M_ & | A, |
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
real_t | omega, | ||
int | max_it, | ||
real_t | toler, | ||
int | verbose | ||
) |
Jacobi solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] omega Relaxation parameter. [in] max_it Maximum number of iterations. [in,out] toler Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. [in] verbose Information output parameter (0: No output, 1: Output iteration information, 2 and greater: Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations,
- Template Parameters
-
<T_> Data type (real_t, float, complex<real_t>, ...) <M_> Matrix storage class
int void BCAsConstraint | ( | const Mesh & | m, |
const Vect< real_t > & | bc, | ||
Vect< real_t > & | up, | ||
Vect< real_t > & | low | ||
) |
To impose Dirichlet boundary conditions in an optimization problem. If such conditions are to present, this function has to be invoked by giving on input bc(i)
as the value to impose for the i
-th optimization variable.
int OptimSA | ( | OPT_ & | theOpt, |
Vect< real_t > & | x, | ||
real_t & | rt, | ||
real_t & | eps, | ||
int & | ns, | ||
int & | nt, | ||
int & | neps, | ||
int & | maxevl, | ||
Vect< real_t > & | lb, | ||
Vect< real_t > & | ub, | ||
Vect< real_t > & | c, | ||
int & | msg_lvl, | ||
int & | seed1, | ||
int & | seed2, | ||
real_t & | t, | ||
Vect< real_t > & | vm, | ||
Vect< real_t > & | xopt, | ||
real_t & | fopt, | ||
int & | nacc, | ||
int & | nfcnev, | ||
int & | nobds | ||
) |
Simulated annealing optimization solver.
Simulated annealing is a global optimization method that distinguishes between different local optima. Starting from an initial point, the algorithm takes a step and the function is evaluated. When minimizing a function, any downhill step is accepted and the process repeats from this new point. An uphill step may be accepted. Thus, it can escape from local optima. This uphill decision is made by the Metropolis criteria. As the optimization process proceeds, the length of the stx decline and the algorithm closes in on the global optimum. Since the algorithm makes very few assumptions regarding the function to be optimized, it is quite robust with respect to non-quadratic surfaces. The degree of robustness can be adjusted by the user. In fact, simulated annealing can be used as a local optimizer for difficult functions.
This implementation of simulated annealing was used in "Global Optimization of Statistical Functions with Simulated Annealing," Goffe, Ferrier and Rogers, Journal of Econometrics, vol. 60, no. 1/2, Jan./Feb. 1994, pp. 65-100. Briefly, we found it competitive, if not superior, to multiple restarts of conventional optimization routines for difficult optimization problems.
For more information on this routine, contact its author: Bill Goffe, bgoff.nosp@m.e@wh.nosp@m.ale.s.nosp@m.t.us.nosp@m.m.edu
Synopsis: This function implements the continuous simulated annealing global optimization algorithm described in Corana et al.'s article "Minimizing Multimodal Functions of Continuous Variables with the "Simulated Annealing" Algorithm" in the September 1987 (vol. 13, no. 3, pp. 262-280) issue of the ACM Transactions on Mathematical Software.
A very quick (perhaps too quick) overview of OptimSA:
OptimSA tries to find the global optimum of an N dimensional function. It moves both up and downhill and as the optimization process proceeds, it focuses on the most promising area.
To start, it randomly chooses a trial point within the step length vm of the user selected starting point. The function is evaluated at this trial point and its value is compared to its value at the initial point.
In a maximization problem, all uphill moves are accepted and the algorithm continues from that trial point. Downhill moves may be accepted; the decision is made by the Metropolis criteria. It uses T (temperature) and the size of the downhill move in a probabilistic manner. The smaller t and the size of the downhill move are, the more likely that move will be accepted. If the trial is accepted, the algorithm moves on from that point. If it is rejected, another point is chosen instead for a trial evaluation.
Each element of vm periodically adjusted so that half of all function evaluations in that direction are accepted.
A fall in t is imposed upon the system with the rt variable by t(i+1) = rt*t(i) where i is the i-th iteration. Thus, as t declines, downhill moves are less likely to be accepted and the percentage of rejections rise. Given the scheme for the selection for vm, vm falls. Thus, as t declines, vm falls and OptimSA focuses upon the most promising area for optimization.
The importance of the parameter t:
The parameter t is crucial in using OptimSA successfully. It influences vm, the step length over which the algorithm searches for optima. For a small initial t, the step length may be too small; thus not enough of the function might be evaluated to find the global optima. The user should carefully examine vm in the intermediate output (set msg_lvl = 1) to make sure that vm is appropriate. The relationship between the initial temperature and the resulting step length is function dependent.
To determine the starting temperature that is consistent with optimizing a function, it is worthwhile to run a trial run first. Set rt = 1.5 and t = 1.0. With rt > 1.0, the temperature increases and vm rises as well. Then select the T that produces a large enough vm.
For modifications to the algorithm and many details on its use, (particularly for econometric applications) see Goffe, Ferrier and Rogers, "Global Optimization of Statistical Functions with
Simulated Annealing," Journal of Econometrics, vol. 60, no. 1/2, Jan./Feb. 1994, pp. 65-100.
For more information, contact
Bill Goffe
Department of Economics and International Business
University of Southern Mississippi
Hattiesburg, MS 39506-5072
(601) 266-4484 (office)
(601) 266-4920 (fax)
bgoff.nosp@m.e@wh.nosp@m.ale.s.nosp@m.t.us.nosp@m.m.edu (Internet)
As far as possible, the parameters here have the same name as in the description of the algorithm on pp. 266-8 of Corana et al.
Note: The suggested values generally come from Corana et al. To drastically reduce runtime, see Goffe et al., pp. 90-1 for suggestions on choosing the appropriate rt
and nt
.
- Parameters
-
[in] theOpt Instance of class OPT_ that is implemented by the user and that provides the objective function. [in] x The starting values for the variables of the function to be optimized. [in] rt The temperature reduction factor. The value suggested by Corana et al. is .85. See Goffe et al. for more advice. [in] eps Error tolerance for termination. If the final function values from the last neps temperatures differ from the corresponding value at the current temperature by less than eps and the final function value at the current temperature differs from the current optimal function value by less than eps, execution terminates and the value 0 is returned. [in] ns Number of cycles. After ns*n function evaluations, each element of vm is adjusted so that approximately half of all function evaluations are accepted. The suggested value is 20. [in] nt Number of iterations before temperature reduction. After nt*ns*n function evaluations, temperature (t) is changed by the factor rt. Value suggested by Corana et al. is max(100,5*n). See Goffe et al. for further advice. [in] neps Number of final function values used to decide upon termination. See eps
. Suggested value is 4maxevl [in] The maximum number of function evaluations. If it is exceeded, the return code=1. [in] lb The lower bound for the allowable solution variables. [in] ub The upper bound for the allowable solution variables. If the algorithm chooses x(i) < lb(i) or x(i) > ub(i) i=1,...,n, a point is from inside is randomly selected. This focuses the algorithm on the region inside ub and lb. Unless the user wishes to concentrate the search to a particular region, ub and lb should be set to very large positive and negative values, respectively. Note that the starting vector x should be inside this region. Also note that lb and ub are fixed in position, while vm is centered on the last accepted trial set of variables that optimizes the function. c [in] Vector that controls the step length adjustment. The suggested value for all elements is 2. [in] msg_lvl controls printing inside OptimSA. - 0 - Nothing printed.
- 1 - Function value for the starting value and summary results before each temperature reduction. This includes the optimal function value found so far, the total number of moves (broken up into uphill, downhill, accepted and rejected), the number of out of bounds trials, the number of new optima found at this temperature, the current optimal x and the step length vm. Note that there are n*ns*nt function evalutations before each temperature reduction. Finally, notice is is also given upon achieveing the termination criteria.
- 2 - Each new step length (vm), the current optimal x (xopt) and the current trial x (x). This gives the user some idea about how far x strays from xopt as well as how vm is adapting to the function.
- 3 - Each function evaluation, its acceptance or rejection and new optima. For many problems, this option will likely require a small tree if hard copy is used. This option is best used to learn about the algorithm. A small value for MAXEVL is thus recommended when using msg_lvl=3. Suggested value: 1
[in] seed1 The first seed for the random number generator ranmar_. 0 <= seed1 <= 31328. (integer) [in] seed2 The second seed for the random number generator ranmar_. 0 <= seed <= 30081. Different values for seed1 and seed2 will lead to an entirely different sequence of trial points and decisions on downhill moves (when maximizing). See Goffe et al. on how this can be used to test the results of OptimSA. [in] t On input, the initial temperature. See Goffe et al. for advice. On output, the final temperature. [in] vm The step length vector. On input it should encompass the region of interest given the starting value x. For point x[i], the next trial point is selected is from x[i]-vm[i] to x[i]+vm[i]. Since vm is adjusted so that about half of all points are accepted, the input value is not very important (i.e. is the value is off, OptimSA adjusts vm to the correct value). [out] xopt The variables that optimize the function. [out] fopt The optimal value of the function. [out] nacc The number of accepted function evaluations. [out] nfcnev The total number of function evaluations. In a minor point, note that the first evaluation is not used in the core of the algorithm; it simply initializes the algorithm. [out] nobds The total number of trial function evaluations that would have been out of bounds of lb and ub. Note that a trial point is randomly selected between lb and ub.
- Returns
-
0 - Normal return; termination criteria achieved.
-
1 - Number of function evaluations (nfcnev) is greater than the maximum number (maxevl).
-
2 - The starting value (x) is not inside the bounds (lb and ub).
-
3 - The initial temperature is not positive.
99 - Should not be seen; only used internally.
-
0 - Normal return; termination criteria achieved.
- Template Parameters
-
<OPT_> Class that defined the objective function and its gradient
int OptimTN | ( | OPT_ & | theOpt, |
Vect< real_t > & | x, | ||
Vect< real_t > & | low, | ||
Vect< real_t > & | up, | ||
Vect< int > & | pivot, | ||
int | max_it = -1 , |
||
real_t | toler = 100*OFELI_EPSMCH , |
||
int | msg_lvl = 5 |
||
) |
Truncated Newton optimization solver.
Solves a bounds-constrained optimization problem using the Nash's Truncated Newton Algorithm (See paper by S.G. Nash, Newton-type Minimization via the Lanczos method, SIAM J. Numer. Anal. 21 (1984) 770-778). All vector variables are instances of class Vect<real_t>.
- Parameters
-
[in] theOpt Instance of class OPT_ that is implemented by the user and that provides the objective function. [in,out] x Vector that contains an initial guess of the solution and as output the final optimization variables if optimization has succeeded [in] low Vector of the same size as x
that contains for each variable the lower bound to impose. Note that Dirichlet boundary conditions are treated as equality conditions (i.e. lower and upper bounds) and that these ones can be imposed via an auxiliary optimization file (BCAsConstraint)[in] up Vector of the same size as x
that contains for each variable the upper bound to impose. Note that Dirichlet boundary conditions are treated as equality conditions (i.e. lower and upper bounds) and that these ones can be imposed via an auxiliary optimization function BCAsConstraint)[in] pivot Vector of the same size as x
that contains on return for each variable an integer value that says if the corresponding constraint was reached (different from0
) or not (= 0
). Note that Dirichlet boundary conditions are treated as equality constraints[in] max_it Maximum number of iterations for convergence [in] toler Tolerance for convergence (measured in relative weighted 2-Norm of projected gradient) [in] msg_lvl Output message level. Must be between 0
and10
- Returns
- Number of performed iterations
- Template Parameters
-
<OPT_> Class that provides the objective function. This class is defined by the user.
The OPT_ class:
This class is defined by the user. It must have the member function :
void Objective(Vect<real_t> &x, real_t &f, Vect<real_t> &g)
Here above, x
is the optimization variable vector, f
is the value of the objective to calculate for the given x
and g
is the gradient vector for x
.
The function BCAsConstraint:
This function is defined by the user:
void BCAsConstraint(const Mesh &m, const Vect<real_t> &bc, Vect<real_t> &up, Vect<real_t> &low)
This function imposes Dirichlet boundary conditions in an optimization problem as optimization constraints. If such conditions are to be present, this function has to be invoked by giving on input bc(i)
as the value to impose for the i
-th optimization variable.
int QMR | ( | const M_ & | A, |
const P1_ & | P1, | ||
const P2_ & | P2, | ||
const Vect< real_t > & | b, | ||
Vect< real_t > & | x, | ||
int | max_it, | ||
real_t & | toler, | ||
int | verbose | ||
) |
QMR solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] P1 First preconditioner for the Conjugate Gradient Squared. [in] P2 Second preconditioner for the Conjugate Gradient Squared. [in] A Problem matrix (Instance of abstract class M_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] max_it Maximum number of iterations. [in] toler Tolerance for convergence (measured in relative weighted 2-Norm). [in] verbose Information output parameter (0: No output, 1: Output iteration information, 2 and greater : Output iteration information and solution at each iteration.
- Returns
- Number of iterations:
-
0
, Convergence within max_iter iterations -
-1
, No convergence after max_iter iterations -
-2
, Breakdown in rho -
-3
, Breakdown in beta -
-4
, Breakdown in gamma -
-5
, Breakdown in delta -
-6
, Breakdown in ep -
-7
, Breakdown in xi
-
Template Arguments:
- M_ Matrix storage class
- P1_ Class of first preconditioner
- P2_ Class of second preconditioner
int Richardson | ( | const M_ & | A, |
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
real_t | omega, | ||
int | max_it, | ||
real_t | toler, | ||
int | verbose | ||
) |
Richardson solver function.
- Parameters
-
[in] A Problem matrix problem (Instance of abstract class M_). [in] b Right-hand side vector (class Vect) x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] omega Relaxation parameter. [in] max_it Maximum number of iterations. [in,out] toler Tolerance for convergence (measured in relative weighted 2-Norm) in input, effective discrepancy in output. [in] verbose Information output parameter ( 0
: No output,1
: Output iteration information,2
and greater: Output iteration information and solution at each iteration.
- Returns
- nb_it Number of performed iterations,
- Template Parameters
-
<T_> Data type (real_t, float, complex<real_t>, ...) <M_> Matrix storage class
void Schur | ( | MAT_A_ & | A, |
MAT_U_ & | U, | ||
MAT_U_ & | L, | ||
MAT_D_ & | D, | ||
Vect< T_ > & | b, | ||
Vect< T_ > & | c | ||
) |
Solve a linear system of equations with a 2x2-block matrix.
The linear system is of the form
| A U | |x| |b| | | | | = | | | L D | |y| |c|
- Parameters
-
[in] A Instance of a matrix class for the first diagonal block. The matrix must be invertible and factorizable (Do not use SpMatrix class) where A
,U
,L
,D
are instances of matrix classes,[in] U Instance of a matrix class for the upper triangle block. The matrix can be rectangular [in] L Instance of a matrix class for the lower triangle block. The matrix can be rectangular [in] D Instance of a matrix class for the second diagonal block. The matrix must be factorizable (Do not use SpMatrix class) [in,out] b Vector (Instance of class Vect) that contains the first block of right-hand side on input and the first block of the solution on output. b
must have the same size as the dimension ofA
.[in,out] c Vector (Instance of class Vect) that contains the second block of right-hand side on output and the first block of the solution on output. c
must have the same size as the dimension ofD
.
Template Arguments:
- Template Parameters
-
<T_> data type (real_t, float, ...) <MAT_A_> Matrix class for A
<MAT_U_> Matrix class for D
<MAT_D_> Matrix class for U
andL
void Schur | ( | Vect< MAT_A_ > & | A, |
Vect< MAT_U_ > & | U, | ||
Vect< MAT_U_ > & | L, | ||
MAT_D_ & | D, | ||
const Vect< Vect< T_ > > & | b, | ||
const Vect< T_ > & | c, | ||
Vect< Vect< T_ > > & | x, | ||
Vect< T_ > & | y | ||
) |
Solve a linear arrow block system by block factorization.
The linear system is of the form
| A[0] 0 ... 0 U[0] | | x[0] | | b[0] | | 0 A[1] ... 0 U[1] | | x[1] | | b[1] | | ... ... | | ... | | ... | | | | | = | | | | | | | | | 0 ... 0 A[nb-1] U[nb-1] | | x[nb-1] | | b[nb-1] | | L[0] L[1] ... L[nb-1] D | | y | | c |
-
Diagonal blocks
A[0],...,A[nb-1]
are instances of template class MAT_A_ -
Block
D
is a dense matrix (instance of template class MAT_D_) - Upper rectangular blocks are instances of template class MAT_U_
- Lower rectangular blocks are instances of template class MAT_U_
All vector blocks are instances of template class Vect<T_>
- Parameters
-
[in] A Vect instance whose coefficients are diagonal blocks of the matrix ( A[0], ...,A[nb-1]
). Each block is an instance of class MAT_A_[in] U Vect instance whose coefficients are upper last blocks of the matrix ( U[0], ..., U[nb-1]
). Each block is an instance of class MAT_U_[in] L Vect instance whose coefficients are lower last blocks of the matrix ( L[0], ...,L[nb-1]
). Each block is an instance of class MAT_U_[in] D Instance of class MAT_D_ (A matrix class) that contains last diagonal block of the matrix [in] b Vect instance whose entries are right-hand side blocks (as instances of class Vect) [in] c Vect instance that contains the last right-hand side block. [out] x Vect instance whose entries are solution blocks (as instances of class Vect) [out] y Vect instance that contains the last solution block.
- Template Parameters
-
<T_> Data type (real_t, float, complex<real_t>, ...) <MAT_A_> Matrix class for A
<MAT_D_> Matrix class for D
<MAT_U_> Matrix class for U
andL
int SSOR | ( | const M_ & | A, |
const Vect< T_ > & | b, | ||
Vect< T_ > & | x, | ||
int | max_it, | ||
double | toler, | ||
int | verbose | ||
) |
SSOR solver function.
- Parameters
-
[in] A Problem matrix (Instance of abstract class M_). [in] b Right-hand side vector (class Vect) [in,out] x Vect instance containing initial solution guess in input and solution of the linear system in output (If iterations have succeeded). [in] max_it Maximum number of iterations. [in] toler Tolerance for convergence (measured in relative weighted 2-Norm). [in] verbose Information output parameter (0: No output, 1: Output iteration information, 2 and greater : Output iteration information and solution at each iteration.
- Returns
- Number of performed iterations,
Template Arguments:
- T_ data type (double, float, ...)
- M_ Matrix storage class