Detailed Description

Definition at line 24 of file PETScNonlinearSolver.h.

#include <PETScNonlinearSolver.h>

Inheritance diagram for NumLib::PETScNonlinearSolver:

Collaboration diagram for NumLib::PETScNonlinearSolver:

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

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. More...

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. More...

std::size_t	_jacobian_id = 0u
	ID of the Jacobian matrix. More...

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

using NumLib::PETScNonlinearSolver::System = NonlinearSystem<NonlinearSolverTag::Newton>

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_solver	the linear solver used by this nonlinear solver.
maxiter	the maximum number of iterations used to solve the equation.
prefix	used to set the SNES options.

Definition at line 80 of file PETScNonlinearSolver.cpp.

 {
     SNESCreate(PETSC_COMM_WORLD, &_snes_solver);
     if (!prefix.empty())
     {
         prefix = prefix + "_";
         SNESSetOptionsPrefix(_snes_solver, prefix.c_str());
     }
     // force SNESSolve() to take at least one iteration regardless of the
     // initial residual norm
     SNESSetForceIteration(_snes_solver, PETSC_TRUE);
  
     // Set the maximum iterations.
     PetscReal atol, rtol, stol;
     PetscInt maxf;
     SNESGetTolerances(_snes_solver, &atol, &rtol, &stol, nullptr, &maxf);
     SNESSetTolerances(_snes_solver, atol, rtol, stol, maxiter, maxf);
  
     SNESSetFromOptions(_snes_solver);
 #ifndef NDEBUG
     PetscOptionsView(nullptr, PETSC_VIEWER_STDOUT_WORLD);
 #endif  // NDEBUG
 }

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 119 of file PETScNonlinearSolver.cpp.

 {
     if (!_compensate_non_equilibrium_initial_residuum)
     {
         return;
     }
  
     OGS_FATAL(
         "Non-equilibrium initial residuum is not implemented for the "
         "PETScNonlinearSolver.");
 }

References _compensate_non_equilibrium_initial_residuum, and OGS_FATAL.

◆ compensateNonEquilibriumInitialResiduum()

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

Definition at line 113 of file PETScNonlinearSolver.cpp.

 {
     _compensate_non_equilibrium_initial_residuum = value;
 }

References _compensate_non_equilibrium_initial_residuum.

◆ setEquationSystem()

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

Set the nonlinear equation system that will be solved.

Definition at line 106 of file PETScNonlinearSolver.cpp.

 {
     _equation_system = &eq;
     _convergence_criterion = &conv_crit;
 }

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

x	in: the initial guess, out: the solution.
x_prev	previous time step solution.
postIterationCallback	called after each iteration if set.
process_id	usually used in staggered schemes.

Return values

true	if the equation system could be solved
false	otherwise

Implements NumLib::NonlinearSolverBase.

Definition at line 133 of file PETScNonlinearSolver.cpp.

 {
     DBUG("PETScNonlinearSolver: solve()");
     using TimeDiscretizedSystem =
         TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
                                  NonlinearSolverTag::Newton>;
  
     auto* system = static_cast<TimeDiscretizedSystem*>(_equation_system);
  
     DBUG("PETScNonlinearSolver: create vectors");
     // r and J on which the ogs assembly operates.
     auto& J = NumLib::GlobalMatrixProvider::provider.getMatrix(
         system->getMatrixSpecifications(process_id), _jacobian_id);
     auto& r = NumLib::GlobalVectorProvider::provider.getVector(
         system->getMatrixSpecifications(process_id), _residual_id);
  
     // Vectors and matrices used by SNES for solutions. These will be locked
     // during the SNES' solve call.
     auto& J_snes = NumLib::GlobalMatrixProvider::provider.getMatrix(
         system->getMatrixSpecifications(process_id), _petsc_jacobian_id);
     BaseLib::RunTime timer_dirichlet;
     double time_dirichlet = 0.0;
  
     timer_dirichlet.start();
     system->computeKnownSolutions(*x[process_id], process_id);
     system->applyKnownSolutions(*x[process_id]);
     time_dirichlet += timer_dirichlet.elapsed();
     INFO("[time] Applying Dirichlet BCs took {} s.", time_dirichlet);
  
     auto& r_snes = NumLib::GlobalVectorProvider::provider.getVector(
         system->getMatrixSpecifications(process_id), _petsc_residual_id);
     auto& x_snes = NumLib::GlobalVectorProvider::provider.getVector(
         system->getMatrixSpecifications(process_id), _petsc_x_id);
     MathLib::LinAlg::copy(*x[process_id], x_snes);  // Initial guess.
  
     PetscContext petsc_context{_equation_system, x, x_prev, &r, &J, process_id};
  
     DBUG("PETScNonlinearSolver: set function");
     SNESSetFunction(_snes_solver, r_snes.getRawVector(), updateResidual,
                     &petsc_context);
  
     DBUG("PETScNonlinearSolver: set jacobian");
     // The jacobian and the preconditioner matrices point to the same location.
     SNESSetJacobian(_snes_solver, J_snes.getRawMatrix(), J_snes.getRawMatrix(),
                     updateJacobian, &petsc_context);
  
     std::unique_ptr<GlobalVector> xl = nullptr;
     std::unique_ptr<GlobalVector> xu = nullptr;
  
     SNESType snes_type;
     SNESGetType(_snes_solver, &snes_type);
     if ((std::strcmp(snes_type, SNESVINEWTONRSLS) == 0) ||
         (std::strcmp(snes_type, SNESVINEWTONSSLS) == 0))
     {
         // Set optional constraints via callback.
         DBUG("PETScNonlinearSolver: set constraints");
         xl = MathLib::MatrixVectorTraits<GlobalVector>::newInstance(
             system->getMatrixSpecifications(process_id));
         xu = MathLib::MatrixVectorTraits<GlobalVector>::newInstance(
             system->getMatrixSpecifications(process_id));
  
         system->updateConstraints(*xl, *xu, process_id);
         MathLib::finalizeVectorAssembly(*xl);
         MathLib::finalizeVectorAssembly(*xu);
  
         SNESVISetVariableBounds(_snes_solver, xl->getRawVector(),
                                 xu->getRawVector());
     }
  
     DBUG("PETScNonlinearSolver: call SNESSolve");
     SNESSolve(_snes_solver, nullptr, x_snes.getRawVector());
  
     SNESConvergedReason reason = SNES_CONVERGED_ITERATING;
     SNESGetConvergedReason(_snes_solver, &reason);
     INFO("PETSsSNES convergence reason {}.", reason);
  
     PetscInt iterations;
     SNESGetIterationNumber(_snes_solver, &iterations);
     INFO("PETScSNES used {} iterations.", iterations);
  
     // Copy back the solution.
     MathLib::LinAlg::copy(x_snes, *x[process_id]);
  
     NumLib::GlobalVectorProvider::provider.releaseVector(x_snes);
     NumLib::GlobalVectorProvider::provider.releaseVector(r_snes);
     NumLib::GlobalMatrixProvider::provider.releaseMatrix(J_snes);
     NumLib::GlobalVectorProvider::provider.releaseVector(r);
     NumLib::GlobalMatrixProvider::provider.releaseMatrix(J);
  
     return {reason >= 0, iterations};
 }

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, MathLib::r, NumLib::MatrixProvider::releaseMatrix(), NumLib::VectorProvider::releaseVector(), BaseLib::RunTime::start(), anonymous_namespace{PETScNonlinearSolver.cpp}::updateJacobian(), and anonymous_namespace{PETScNonlinearSolver.cpp}::updateResidual().

Member Data Documentation

◆ _compensate_non_equilibrium_initial_residuum

bool NumLib::PETScNonlinearSolver::_compensate_non_equilibrium_initial_residuum = false

private

Enables computation of the non-equilibrium initial residuum \( r_{\rm neq} \) before the first time step. The forces are zero if the external forces are in equilibrium with the initial state/initial conditions. During the simulation the new residuum reads \( \tilde r = r - r_{\rm neq} \).

Definition at line 73 of file PETScNonlinearSolver.h.

Referenced by calculateNonEquilibriumInitialResiduum(), and compensateNonEquilibriumInitialResiduum().

◆ _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:

NumLib/ODESolver/PETScNonlinearSolver.h
NumLib/ODESolver/PETScNonlinearSolver.cpp

Detailed Description

Public Types

Public Member Functions

Private Attributes

Member Typedef Documentation

◆ System

Constructor & Destructor Documentation

◆ PETScNonlinearSolver()

Member Function Documentation

◆ calculateNonEquilibriumInitialResiduum()

◆ compensateNonEquilibriumInitialResiduum()

◆ setEquationSystem()

◆ solve()

Member Data Documentation

◆ _compensate_non_equilibrium_initial_residuum

◆ _convergence_criterion

◆ _equation_system

◆ _jacobian_id

◆ _petsc_jacobian_id

◆ _petsc_residual_id

◆ _petsc_x_id

◆ _residual_id

◆ _snes_solver