15 #include "../config/config.hpp"
24 #ifdef MFEM_USE_SUITESPARSE
89 bool final=
false)
const;
121 {
return dot_prod_type == 0 ? MPI_COMM_NULL : comm; }
149 const double damping=1.0);
162 const double damping=1.0);
184 const double damping;
192 const bool allow_updates;
213 int order,
double max_eig_estimate);
219 int order,
double max_eig_estimate);
231 int order, MPI_Comm comm = MPI_COMM_NULL,
232 int power_iterations = 10,
233 double power_tolerance = 1e-8);
239 int order, MPI_Comm comm = MPI_COMM_NULL,
240 int power_iterations = 10,
241 double power_tolerance = 1e-8);
245 int order,
int power_iterations = 10,
246 double power_tolerance = 1e-8);
252 int order,
int power_iterations = 10,
253 double power_tolerance = 1e-8);
271 double max_eig_estimate;
278 mutable Vector helperVector;
305 void SLI(
const Operator &A,
const Vector &
b, Vector &x,
306 int print_iter = 0,
int max_num_iter = 1000,
307 double RTOLERANCE = 1e-12,
double ATOLERANCE = 1e-24);
310 void SLI(
const Operator &A, Solver &B,
const Vector &
b, Vector &x,
311 int print_iter = 0,
int max_num_iter = 1000,
312 double RTOLERANCE = 1e-12,
double ATOLERANCE = 1e-24);
337 void CG(
const Operator &A,
const Vector &
b, Vector &x,
338 int print_iter = 0,
int max_num_iter = 1000,
339 double RTOLERANCE = 1e-12,
double ATOLERANCE = 1e-24);
342 void PCG(
const Operator &A, Solver &B,
const Vector &
b, Vector &x,
343 int print_iter = 0,
int max_num_iter = 1000,
344 double RTOLERANCE = 1e-12,
double ATOLERANCE = 1e-24);
385 int GMRES(
const Operator &A, Vector &x,
const Vector &
b, Solver &M,
386 int &max_iter,
int m,
double &tol,
double atol,
int printit);
389 void GMRES(
const Operator &A, Solver &B,
const Vector &
b, Vector &x,
390 int print_iter = 0,
int max_num_iter = 1000,
int m = 50,
391 double rtol = 1e-12,
double atol = 1e-24);
416 int BiCGSTAB(
const Operator &A, Vector &x,
const Vector &
b, Solver &M,
417 int &max_iter,
double &tol,
double atol,
int printit);
420 void BiCGSTAB(
const Operator &A, Solver &B,
const Vector &
b, Vector &x,
421 int print_iter = 0,
int max_num_iter = 1000,
422 double rtol = 1e-12,
double atol = 1e-24);
451 void MINRES(
const Operator &A,
const Vector &
b, Vector &x,
int print_it = 0,
452 int max_it = 1000,
double rtol = 1e-12,
double atol = 1e-24);
455 void MINRES(
const Operator &A, Solver &B,
const Vector &
b, Vector &x,
456 int print_it = 0,
int max_it = 1000,
457 double rtol = 1e-12,
double atol = 1e-24);
493 const double fnorm)
const;
500 const double fnorm)
const;
541 const double rtol0 = 0.5,
542 const double rtol_max = 0.9,
543 const double alpha = 0.5 * (1.0 + sqrt(5.0)),
544 const double gamma = 1.0);
568 { MFEM_WARNING(
"L-BFGS won't use the given preconditioner."); }
570 { MFEM_WARNING(
"L-BFGS won't use the given solver."); }
577 int aGMRES(
const Operator &A, Vector &x,
const Vector &
b,
578 const Operator &M,
int &max_iter,
579 int m_max,
int m_min,
int m_step,
double cf,
580 double &tol,
double &atol,
int printit);
616 { MFEM_ABORT(
"The objective gradient is not implemented."); }
654 { MFEM_ABORT(
"Not meaningful for this solver."); }
656 { MFEM_ABORT(
"Not meaningful for this solver."); }
832 bool final)
override;
836 #ifdef MFEM_USE_SUITESPARSE
851 mutable double Info[UMFPACK_INFO];
908 #endif // MFEM_USE_SUITESPARSE
917 std::unique_ptr<DenseMatrixInverse[]> block_solvers;
935 bool ownA,
bool ownS0,
bool ownS1)
936 :
Solver(A_->
NumRows()), A(A_, ownA), S0(S0_, ownS0), S1(S1_, ownS1) { }
954 void Mult(
const Vector &x,
Vector &y,
bool transpose)
const;
957 bool op_is_symmetric =
true,
bool own_aux_map =
false);
963 #endif // MFEM_USE_MPI
967 #endif // MFEM_SOLVERS
const Array< int > * ess_dofs_list
Not owned.
MINRESSolver(MPI_Comm comm_)
const Vector * GetInequalityVec_Hi() const
GMRESSolver(MPI_Comm comm_)
Conjugate gradient method.
OpType * As() const
Return the Operator pointer statically cast to a specified OpType. Similar to the method Get()...
HypreSmoother & GetSmoother()
virtual void MultTranspose(const Vector &b, Vector &x) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
double Info[UMFPACK_INFO]
void Mult(const Vector &b, Vector &x) const
Solve the system LUx = b, where L and U are the block ILU factors.
Chebyshev accelerated smoothing with given vector, no matrix necessary.
int GetNumIterations() const
MPI_Comm GetComm() const
Return the associated MPI communicator, or MPI_COMM_NULL if no communicator is set.
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
void SetHistorySize(int dim)
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
void SetSize(int s)
Resize the vector to size s.
void Monitor(int it, double norm, const Vector &r, const Vector &x, bool final=false) const
virtual void ProcessNewState(const Vector &x) const
This method can be overloaded in derived classes to perform computations that need knowledge of the n...
OperatorJacobiSmoother(const double damping=1.0)
Default constructor: the diagonal will be computed by subsequent calls to SetOperator() using the Ope...
Pointer to an Operator of a specified type.
~OperatorChebyshevSmoother()
class IterativeSolver * iter_solver
The last IterativeSolver to which this monitor was attached.
T * GetData()
Returns the data.
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
Direct sparse solver using KLU.
virtual void Mult(const Vector &x, Vector &y) const =0
Operator application: y=A(x).
virtual void SetSolver(Solver &solver)
Set the linear solver for inverting the Jacobian.
virtual void MultTranspose(const Vector &b, Vector &x) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
void SetOperator(const Operator &op)
UMFPackSolver(SparseMatrix &A, bool use_long_ints_=false)
Factorize the given SparseMatrix using the defaults. For larger matrices, if the solver fails...
const Operator * GetC() const
Reordering
The reordering method used by the BlockILU factorization.
virtual ~IterativeSolverMonitor()
ResidualBCMonitor(const Array< int > &ess_dofs_list_)
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
virtual double CalcObjective(const Vector &x) const =0
Objective F(x). In parallel, the result should be reduced over tasks.
SLISolver(MPI_Comm comm_)
void add(const Vector &v1, const Vector &v2, Vector &v)
virtual ~OptimizationSolver()
void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
double Norm(const Vector &x) const
double GetFinalNorm() const
int BiCGSTAB(const Operator &A, Vector &x, const Vector &b, Solver &M, int &max_iter, double &tol, double atol, int printit)
BiCGSTAB method. (tolerances are squared)
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
Direct sparse solver using UMFPACK.
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
void SetPrintLevel(int print_lvl)
const Vector * GetBoundsVec_Lo() const
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
double solve(double l, const Vector &xt, Vector &x, int &nclip) const
Solve QP at fixed lambda.
int GetNumConstraints() const
void MINRES(const Operator &A, const Vector &b, Vector &x, int print_it, int max_it, double rtol, double atol)
MINRES method without preconditioner. (tolerances are squared)
Jacobi smoothing for a given bilinear form (no matrix necessary).
~OperatorJacobiSmoother()
void SetMaxIter(int max_it)
OperatorChebyshevSmoother(const Operator &oper_, const Vector &d, const Array< int > &ess_tdof_list, int order, double max_eig_estimate)
void SetLinearConstraint(const Vector &w_, double a_)
void median(const Vector &lo, const Vector &hi)
v = median(v,lo,hi) entrywise. Implementation assumes lo <= hi.
DirectSubBlockSolver(const SparseMatrix &A, const SparseMatrix &block_dof)
void SetKDim(int dim)
Set the number of iteration to perform between restarts, default is 50.
Newton's method for solving F(x)=b for a given operator F.
void CG(const Operator &A, const Vector &b, Vector &x, int print_iter, int max_num_iter, double RTOLERANCE, double ATOLERANCE)
Conjugate gradient method. (tolerances are squared)
void print_iteration(int it, double r, double l) const
void Setup(const Vector &diag)
Parallel smoothers in hypre.
void PCG(const Operator &A, Solver &B, const Vector &b, Vector &x, int print_iter, int max_num_iter, double RTOLERANCE, double ATOLERANCE)
Preconditioned conjugate gradient method. (tolerances are squared)
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
double Dot(const Vector &x, const Vector &y) const
void SetIterativeSolver(const IterativeSolver &solver)
This method is invoked by ItertiveSolver::SetMonitor, informing the monitor which IterativeSolver is ...
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
const Vector * GetBoundsVec_Hi() const
BlockILU(int block_size_, Reordering reordering_=Reordering::MINIMUM_DISCARDED_FILL, int k_fill_=0)
virtual void Mult(const Vector &b, Vector &x) const
Solve the nonlinear system with right-hand side b.
void SetPrintLevel(int print_lvl)
Set the print level field in the Control data member.
void AdaptiveLinRtolPostSolve(const Vector &x, const Vector &b, const int it, const double fnorm) const
Method for the adaptive linear solver rtol invoked after the linear solve.
Stationary linear iteration: x <- x + B (b - A x)
void SetInequalityConstraint(const Vector &dl, const Vector &dh)
void SetEqualityConstraint(const Vector &c)
Abstract base class for an iterative solver monitor.
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
void SetOperator(const Operator &op_)
Set/update the solver for the given operator.
ProductSolver(Operator *A_, Solver *S0_, Solver *S1_, bool ownA, bool ownS0, bool ownS1)
void SetAbsTol(double atol)
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
void SetRelTol(double rtol)
double Control[UMFPACK_CONTROL]
Monitor that checks whether the residual is zero at a given set of dofs.
Abstract base class for iterative solver.
UMFPackSolver(bool use_long_ints_=false)
For larger matrices, if the solver fails, set the parameter use_long_ints_ = true.
virtual void Mult(const Vector &b, Vector &x) const
Solve the nonlinear system with right-hand side b.
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
int NumRows() const
Get the number of rows (size of output) of the Operator. Synonym with Height().
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
void SetAdaptiveLinRtol(const int type=2, const double rtol0=0.5, const double rtol_max=0.9, const double alpha=0.5 *(1.0+sqrt(5.0)), const double gamma=1.0)
Enable adaptive linear solver relative tolerance algorithm.
LBFGSSolver(MPI_Comm comm_)
virtual void SetSolver(Solver &solver)
Set the linear solver for inverting the Jacobian.
const Vector * GetEqualityVec() const
FGMRESSolver(MPI_Comm comm_)
virtual void SetOptimizationProblem(const OptimizationProblem &prob)
void SetMonitor(IterativeSolverMonitor &m)
Set the iterative solver monitor.
const Operator * GetD() const
Abstract solver for OptimizationProblems.
virtual double ComputeScalingFactor(const Vector &x, const Vector &b) const
This method can be overloaded in derived classes to implement line search algorithms.
OptimizationProblem(int insize, const Operator *C_, const Operator *D_)
In parallel, insize is the number of the local true dofs.
void SetSolutionBounds(const Vector &xl, const Vector &xh)
const Vector * GetInequalityVec_Lo() const
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
AuxSpaceSmoother(const HypreParMatrix &op, HypreParMatrix *aux_map, bool op_is_symmetric=true, bool own_aux_map=false)
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
virtual void Mult(const Vector &xt, Vector &x) const =0
Operator application: y=A(x).
BiCGSTABSolver(MPI_Comm comm_)
void SetOperator(const Operator &op)
Recompute the diagonal using the method AssembleDiagonal of the given new Operator, op.
SLBQPOptimizer(MPI_Comm comm_)
void AdaptiveLinRtolPreSolve(const Vector &x, const int it, const double fnorm) const
Method for the adaptive linear solver rtol invoked before the linear solve.
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
virtual void Mult(const Vector &xt, Vector &x) const
void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
virtual void CalcObjectiveGrad(const Vector &x, Vector &grad) const
The result grad is expected to enter with the correct size.
void MonitorResidual(int it, double norm, const Vector &r, bool final) override
Monitor the residual vector r.
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
virtual void MonitorSolution(int it, double norm, const Vector &x, bool final)
Monitor the solution vector x.
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
KLUSolver(SparseMatrix &A)
int aGMRES(const Operator &A, Vector &x, const Vector &b, const Operator &M, int &max_iter, int m_max, int m_min, int m_step, double cf, double &tol, double &atol, int printit)
OptimizationSolver(MPI_Comm comm_)
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
const OptimizationProblem * problem
Wrapper for hypre's ParCSR matrix class.
Block diagonal solver for A, each block is inverted by direct solver.
void SetBounds(const Vector &lo_, const Vector &hi_)
Rank 3 tensor (array of matrices)
virtual void SetOptimizationProblem(const OptimizationProblem &prob)
NewtonSolver(MPI_Comm comm_)
IterativeSolverMonitor * monitor
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
virtual void MonitorResidual(int it, double norm, const Vector &r, bool final)
Monitor the residual vector r.
int width
Dimension of the input / number of columns in the matrix.
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
const Operator * C
Not owned, some can remain unused (NULL).
void SLI(const Operator &A, const Vector &b, Vector &x, int print_iter, int max_num_iter, double RTOLERANCE, double ATOLERANCE)
Stationary linear iteration. (tolerances are squared)
int GMRES(const Operator &A, Vector &x, const Vector &b, Solver &M, int &max_iter, int m, double &tol, double atol, int printit)
GMRES method. (tolerances are squared)
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
virtual void SetOperator(const Operator &op)
Factorize the given Operator op which must be a SparseMatrix.