OGS
NumLib::PETScNonlinearSolver Class Referencefinal

Detailed Description

Definition at line 24 of file PETScNonlinearSolver.h.

#include <PETScNonlinearSolver.h>

Inheritance diagram for NumLib::PETScNonlinearSolver:
[legend]
Collaboration diagram for NumLib::PETScNonlinearSolver:
[legend]

Public Types

using System = NonlinearSystem<NonlinearSolverTag::Newton>
 Type of the nonlinear equation system to be solved.
 

Public Member Functions

 PETScNonlinearSolver (GlobalLinearSolver &linear_solver, int maxiter, std::string prefix)
 
void setEquationSystem (System &eq, ConvergenceCriterion &conv_crit)
 Set the nonlinear equation system that will be solved.
 
void compensateNonEquilibriumInitialResiduum (bool const value)
 
void calculateNonEquilibriumInitialResiduum (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id) override
 
NonlinearSolverStatus solve (std::vector< GlobalVector * > &x, std::vector< GlobalVector * > const &x_prev, std::function< void(int, std::vector< GlobalVector * > const &)> const &postIterationCallback, int const process_id) override
 
- Public Member Functions inherited from NumLib::NonlinearSolverBase
virtual ~NonlinearSolverBase ()=default
 

Private Attributes

SNES _snes_solver
 
System_equation_system = nullptr
 
ConvergenceCriterion_convergence_criterion = nullptr
 
std::size_t _residual_id = 0u
 ID of the residual vector.
 
std::size_t _jacobian_id = 0u
 ID of the Jacobian matrix.
 
std::size_t _petsc_x_id = 0u
 
std::size_t _petsc_jacobian_id = 0u
 
std::size_t _petsc_residual_id = 0u
 
bool _compensate_non_equilibrium_initial_residuum = false
 

Member Typedef Documentation

◆ System

Type of the nonlinear equation system to be solved.

Definition at line 28 of file PETScNonlinearSolver.h.

Constructor & Destructor Documentation

◆ PETScNonlinearSolver()

NumLib::PETScNonlinearSolver::PETScNonlinearSolver ( GlobalLinearSolver & linear_solver,
int maxiter,
std::string prefix )

Constructs a new instance.

Parameters
linear_solverthe linear solver used by this nonlinear solver.
maxiterthe maximum number of iterations used to solve the equation.
prefixused to set the SNES options.

Definition at line 95 of file PETScNonlinearSolver.cpp.

98{
99 SNESCreate(PETSC_COMM_WORLD, &_snes_solver);
100 if (!prefix.empty())
101 {
102 prefix = prefix + "_";
103 SNESSetOptionsPrefix(_snes_solver, prefix.c_str());
104 }
105 // force SNESSolve() to take at least one iteration regardless of the
106 // initial residual norm
107 SNESSetForceIteration(_snes_solver, PETSC_TRUE);
108
109 // Set the maximum iterations.
110 PetscReal atol, rtol, stol;
111 PetscInt maxf;
112 SNESGetTolerances(_snes_solver, &atol, &rtol, &stol, nullptr, &maxf);
113 SNESSetTolerances(_snes_solver, atol, rtol, stol, maxiter, maxf);
114
115 SNESSetFromOptions(_snes_solver);
116#ifndef NDEBUG
117 PetscOptionsView(nullptr, PETSC_VIEWER_STDOUT_WORLD);
118#endif // NDEBUG
119}

References _snes_solver.

Member Function Documentation

◆ calculateNonEquilibriumInitialResiduum()

void NumLib::PETScNonlinearSolver::calculateNonEquilibriumInitialResiduum ( std::vector< GlobalVector * > const & x,
std::vector< GlobalVector * > const & x_prev,
int const process_id )
overridevirtual

Implements NumLib::NonlinearSolverBase.

Definition at line 134 of file PETScNonlinearSolver.cpp.

137{
139 {
140 return;
141 }
142
143 OGS_FATAL(
144 "Non-equilibrium initial residuum is not implemented for the "
145 "PETScNonlinearSolver.");
146}
#define OGS_FATAL(...)
Definition Error.h:26

References _compensate_non_equilibrium_initial_residuum, and OGS_FATAL.

◆ compensateNonEquilibriumInitialResiduum()

void NumLib::PETScNonlinearSolver::compensateNonEquilibriumInitialResiduum ( bool const value)

◆ setEquationSystem()

void NumLib::PETScNonlinearSolver::setEquationSystem ( System & eq,
ConvergenceCriterion & conv_crit )

Set the nonlinear equation system that will be solved.

Definition at line 121 of file PETScNonlinearSolver.cpp.

123{
124 _equation_system = &eq;
125 _convergence_criterion = &conv_crit;
126}
ConvergenceCriterion * _convergence_criterion

References _convergence_criterion, and _equation_system.

◆ solve()

NonlinearSolverStatus NumLib::PETScNonlinearSolver::solve ( std::vector< GlobalVector * > & x,
std::vector< GlobalVector * > const & x_prev,
std::function< void(int, std::vector< GlobalVector * > const &)> const & postIterationCallback,
int const process_id )
overridevirtual

Assemble and solve the equation system.

Parameters
xin: the initial guess, out: the solution.
x_prevprevious time step solution.
postIterationCallbackcalled after each iteration if set.
process_idusually used in staggered schemes.
Return values
trueif the equation system could be solved
falseotherwise

Implements NumLib::NonlinearSolverBase.

Definition at line 148 of file PETScNonlinearSolver.cpp.

155{
156 DBUG("PETScNonlinearSolver: solve()");
157 using TimeDiscretizedSystem =
158 TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
160
161 auto* system = static_cast<TimeDiscretizedSystem*>(_equation_system);
162
163 DBUG("PETScNonlinearSolver: create vectors");
164 // r and J on which the ogs assembly operates.
166 system->getMatrixSpecifications(process_id), _jacobian_id);
168 system->getMatrixSpecifications(process_id), _residual_id);
169
170 // Vectors and matrices used by SNES for solutions. These will be locked
171 // during the SNES' solve call.
173 system->getMatrixSpecifications(process_id), _petsc_jacobian_id);
174 BaseLib::RunTime timer_dirichlet;
175 double time_dirichlet = 0.0;
176
177 timer_dirichlet.start();
178 system->computeKnownSolutions(*x[process_id], process_id);
179 system->applyKnownSolutions(*x[process_id]);
180 time_dirichlet += timer_dirichlet.elapsed();
181 INFO("[time] Applying Dirichlet BCs took {} s.", time_dirichlet);
182
184 system->getMatrixSpecifications(process_id), _petsc_residual_id);
186 system->getMatrixSpecifications(process_id), _petsc_x_id);
187 MathLib::LinAlg::copy(*x[process_id], x_snes); // Initial guess.
188
189 PetscContext petsc_context{_equation_system, x, x_prev, &r, &J, process_id};
190
191 DBUG("PETScNonlinearSolver: set function");
192 SNESSetFunction(_snes_solver, r_snes.getRawVector(), updateResidual,
193 &petsc_context);
194
195 DBUG("PETScNonlinearSolver: set jacobian");
196 // The jacobian and the preconditioner matrices point to the same location.
197 SNESSetJacobian(_snes_solver, J_snes.getRawMatrix(), J_snes.getRawMatrix(),
198 updateJacobian, &petsc_context);
199
200 std::unique_ptr<GlobalVector> xl = nullptr;
201 std::unique_ptr<GlobalVector> xu = nullptr;
202
203 SNESType snes_type;
204 SNESGetType(_snes_solver, &snes_type);
205 if ((std::strcmp(snes_type, SNESVINEWTONRSLS) == 0) ||
206 (std::strcmp(snes_type, SNESVINEWTONSSLS) == 0))
207 {
208 // Set optional constraints via callback.
209 DBUG("PETScNonlinearSolver: set constraints");
211 system->getMatrixSpecifications(process_id));
213 system->getMatrixSpecifications(process_id));
214
215 system->updateConstraints(*xl, *xu, process_id);
218
219 SNESVISetVariableBounds(_snes_solver, xl->getRawVector(),
220 xu->getRawVector());
221 }
222
223 DBUG("PETScNonlinearSolver: call SNESSolve");
224 SNESSolve(_snes_solver, nullptr, x_snes.getRawVector());
225
226 SNESConvergedReason reason;
227 SNESGetConvergedReason(_snes_solver, &reason);
229
230 PetscInt iterations;
231 SNESGetIterationNumber(_snes_solver, &iterations);
232 INFO("PETScSNES used {} iterations.", iterations);
233
234 // Copy back the solution.
235 MathLib::LinAlg::copy(x_snes, *x[process_id]);
236
242
243 return {reason >= 0, iterations};
244}
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
Count the running time.
Definition RunTime.h:29
double elapsed() const
Get the elapsed time in seconds.
Definition RunTime.h:42
void start()
Start the timer.
Definition RunTime.h:32
virtual void releaseMatrix(GlobalMatrix const &A)=0
virtual GlobalMatrix & getMatrix(std::size_t &id)=0
Get an uninitialized matrix with the given id.
std::size_t _jacobian_id
ID of the Jacobian matrix.
std::size_t _residual_id
ID of the residual vector.
virtual GlobalVector & getVector(std::size_t &id)=0
Get an uninitialized vector with the given id.
virtual void releaseVector(GlobalVector const &x)=0
void copy(PETScVector const &x, PETScVector &y)
Definition LinAlg.cpp:37
void finalizeVectorAssembly(VEC_T &)
General function to finalize the vector assembly.
void printConvergenceReason(SNES const &snes_solver, SNESConvergedReason const &reason)
static NUMLIB_EXPORT MatrixProvider & provider
static NUMLIB_EXPORT VectorProvider & provider

References _equation_system, _jacobian_id, _petsc_jacobian_id, _petsc_residual_id, _petsc_x_id, _residual_id, _snes_solver, MathLib::LinAlg::copy(), DBUG(), BaseLib::RunTime::elapsed(), MathLib::finalizeVectorAssembly(), NumLib::FirstOrderImplicitQuasilinear, NumLib::MatrixProvider::getMatrix(), NumLib::VectorProvider::getVector(), INFO(), NumLib::Newton, NumLib::GlobalVectorProvider::provider, NumLib::GlobalMatrixProvider::provider, NumLib::MatrixProvider::releaseMatrix(), NumLib::VectorProvider::releaseVector(), and BaseLib::RunTime::start().

Member Data Documentation

◆ _compensate_non_equilibrium_initial_residuum

bool NumLib::PETScNonlinearSolver::_compensate_non_equilibrium_initial_residuum = false
private

◆ _convergence_criterion

ConvergenceCriterion* NumLib::PETScNonlinearSolver::_convergence_criterion = nullptr
private

Definition at line 62 of file PETScNonlinearSolver.h.

Referenced by setEquationSystem().

◆ _equation_system

System* NumLib::PETScNonlinearSolver::_equation_system = nullptr
private

Definition at line 61 of file PETScNonlinearSolver.h.

Referenced by setEquationSystem(), and solve().

◆ _jacobian_id

std::size_t NumLib::PETScNonlinearSolver::_jacobian_id = 0u
private

ID of the Jacobian matrix.

Definition at line 65 of file PETScNonlinearSolver.h.

Referenced by solve().

◆ _petsc_jacobian_id

std::size_t NumLib::PETScNonlinearSolver::_petsc_jacobian_id = 0u
private

Definition at line 68 of file PETScNonlinearSolver.h.

Referenced by solve().

◆ _petsc_residual_id

std::size_t NumLib::PETScNonlinearSolver::_petsc_residual_id = 0u
private

Definition at line 69 of file PETScNonlinearSolver.h.

Referenced by solve().

◆ _petsc_x_id

std::size_t NumLib::PETScNonlinearSolver::_petsc_x_id = 0u
private

Definition at line 67 of file PETScNonlinearSolver.h.

Referenced by solve().

◆ _residual_id

std::size_t NumLib::PETScNonlinearSolver::_residual_id = 0u
private

ID of the residual vector.

Definition at line 64 of file PETScNonlinearSolver.h.

Referenced by solve().

◆ _snes_solver

SNES NumLib::PETScNonlinearSolver::_snes_solver
private

Definition at line 59 of file PETScNonlinearSolver.h.

Referenced by PETScNonlinearSolver(), and solve().


The documentation for this class was generated from the following files: