20 #ifndef GEOS_PHYSICSSOLVERS_PHYSICSSOLVERBASE_HPP_
21 #define GEOS_PHYSICSSOLVERS_PHYSICSSOLVERBASE_HPP_
23 #include "codingUtilities/traits.hpp"
32 #include "physicsSolvers/NonlinearSolverParameters.hpp"
34 #include "physicsSolvers/SolverStatistics.hpp"
41 class DomainPartition;
61 Group *
const parent );
122 real64 const eventProgress,
132 real64 const eventProgress,
403 bool const setSparsity =
true );
410 virtual std::unique_ptr< PreconditionerBase< LAInterface > >
470 integer const nonlinearIteration,
484 integer const nonlinearIteration,
554 real64 const scalingFactor );
592 real64 const scalingFactor,
848 template<
typename LAMBDA >
851 for(
auto const & target: m_meshTargets )
853 string const meshBodyName = target.first.first;
854 string const meshLevelName = target.first.second;
859 if( meshLevelPtr==
nullptr )
863 lambda( meshBodyName, *meshLevelPtr, regionNames );
874 template<
typename LAMBDA >
877 for(
auto const & target: m_meshTargets )
879 string const meshBodyName = target.first.first;
880 string const meshLevelName = target.first.second;
885 if( meshLevelPtr==
nullptr )
889 lambda( meshBodyName, *meshLevelPtr, regionNames );
907 virtual bool registerCallback(
void * func,
const std::type_info & funcType )
final override;
937 #if defined(GEOS_USE_PYGEOSX)
942 virtual PyTypeObject * getPythonType()
const override;
951 return m_meshTargets;
983 real64 const oldNewtonNorm,
993 template<
typename CONSTITUTIVE_BASE_TYPE >
1003 template<
typename CONSTITUTIVE_BASE_TYPE >
1014 template<
typename CONSTITUTIVE >
1034 template<
typename BASETYPE = constitutive::ConstitutiveBase,
typename LOOKUP_TYPE >
1038 return constitutiveModels.
getGroup< BASETYPE >( key );
1049 template<
typename BASETYPE = constitutive::ConstitutiveBase,
typename LOOKUP_TYPE >
1053 return constitutiveModels.
getGroup< BASETYPE >( key );
1062 template<
typename CONSTITUTIVE_TYPE >
1065 return getConstitutiveModel< CONSTITUTIVE_TYPE >( subRegion, getConstitutiveName< CONSTITUTIVE_TYPE >( subRegion ) );
1108 std::unique_ptr< PreconditionerBase< LAInterface > >
m_precond;
1136 std::map< std::string, std::chrono::system_clock::duration >
m_timers;
1161 bool solveNonlinearSystem(
real64 const & time_n,
1172 void logEndOfCycleInformation(
integer const cycleNumber,
1178 template<
typename CONSTITUTIVE_BASE_TYPE >
1186 GEOS_ERROR_IF( !validName.empty(),
"A valid constitutive model was already found." );
1187 validName = model.getName();
1193 template<
typename CONSTITUTIVE_BASE_TYPE >
1201 GEOS_ERROR_IF( !validName.empty(),
"A valid constitutive model was already found." );
1202 validName = model.getName();
1207 template<
typename CONSTITUTIVE >
1215 string & constitutiveName = subRegion.
getReference<
string >( wrapperName );
1216 constitutiveName = getConstitutiveName< CONSTITUTIVE >( subRegion );
1217 GEOS_ERROR_IF( constitutiveName.empty(), GEOS_FMT(
"{}: {} constitutive model not found on subregion {}",
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
#define GEOS_ERROR_IF(EXP, msg)
Conditionally raise a hard error and terminate the program.
Class containing convergence information given a time-step.
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Partition of the decomposed physical domain. It also manages the connexion information to its neighbo...
dataRepository::Group const & getConstitutiveModels() const
Get the group in which the constitutive models of this subregion are registered.
Class containing solver iterations data for a time-step.
The class is used to manage mesh body.
Group & getMeshLevels()
Get the meshLevels group.
Class facilitating the representation of a multi-level discretization of a MeshBody.
dataRepository::Group const & getConstitutiveModels() const
Get the group in which the constitutive models of this subregion are registered.
Base class for all physics solvers.
SolverStatistics const & getSolverStatistics() const
const accessor for the solver statistics.
void setSystemSetupTimestamp(Timestamp timestamp)
set the timestamp of the system setup
map< std::pair< string, string >, string_array > const & getMeshTargets() const
accessor for m_meshTargets
integer m_allowNonConvergedLinearSolverSolution
behavior in case of linear solver failure
virtual void registerDataOnMesh(Group &MeshBodies) override
Register wrappers that contain data on the mesh objects.
PhysicsSolverBase(PhysicsSolverBase const &)=delete
Deleted copy constructor.
static BASETYPE const & getConstitutiveModel(dataRepository::Group const &dataGroup, LOOKUP_TYPE const &key)
Get the Constitutive Model object.
virtual string getCatalogName() const =0
DofManager & getDofManager()
Getter for degree-of-freedom manager.
SolverStatistics m_solverStatistics
Solver statistics.
real64 m_cflFactor
Courant–Friedrichs–Lewy factor for the timestep.
string m_discretizationName
name of the FV discretization object in the data repository
IterationsStatistics & getIterationStats()
string getDiscretizationName() const
return the name of the discretization object
ParallelVector const & getSystemRhs() const
Getter for system rhs vector.
ParallelMatrix m_matrix
System matrix.
std::map< std::string, std::chrono::system_clock::duration > m_timers
Timers for the aggregate profiling of the solver.
LinearSolverParameters & getLinearSolverParameters()
accessor for the linear solver parameters.
static BASETYPE & getConstitutiveModel(dataRepository::Group &dataGroup, LOOKUP_TYPE const &key)
Get the Constitutive Model object.
void forDiscretizationOnMeshTargets(Group const &meshBodies, LAMBDA &&lambda) const
Loop over the target discretization on all mesh targets and apply callback.
virtual bool checkSequentialSolutionIncrements(DomainPartition &domain) const
Check if the solution increments are ok to use.
PhysicsSolverBase(string const &name, Group *const parent)
Constructor for PhysicsSolverBase.
ParallelMatrix & getSystemMatrix()
Getter for system matrix.
static string getConstitutiveName(ElementSubRegionBase const &subRegion)
Get the Constitutive Name object.
integer m_numTimestepsSinceLastDtCut
Number of cycles since last timestep cut.
virtual void setConstitutiveNamesCallSuper(ElementSubRegionBase &subRegion) const
This function sets constitutive name fields on an ElementSubRegionBase, and calls the base function i...
ParallelVector & getSystemRhs()
Getter for system rhs vector.
virtual void cleanup(real64 const time_n, integer const cycleNumber, integer const eventCounter, real64 const eventProgress, DomainPartition &domain) override
Called as the code exits the main run loop.
std::unique_ptr< PreconditionerBase< LAInterface > > m_precond
Custom preconditioner for the "native" iterative solver.
Timestamp getMeshModificationTimestamp(DomainPartition &domain) const
getter for the timestamp of the mesh modification on the mesh levels
ParallelVector m_scaling
Diagonal scaling vector D (Ahat = D * A * D, bhat = D * b, x = D * xhat)
virtual void postInputInitialization() override
void forDiscretizationOnMeshTargets(Group &meshBodies, LAMBDA &&lambda) const
Loop over the target discretization on all mesh targets and apply callback.
localIndex targetRegionIndex(string const ®ionName) const
Get position of a given region within solver's target region list.
ParallelMatrix const & getSystemMatrix() const
Getter for system rhs vector.
integer m_writeLinearSystem
flag for debug output of matrix, rhs, and solution
PhysicsSolverBase()=delete
Deleted constructor.
LinearSolverParameters const & getLinearSolverParameters() const
const accessor for the linear solver parameters.
LinearSolverResult m_linearSolverResult
Result of the last linear solver.
virtual ~PhysicsSolverBase() override
Destructor for PhysicsSolverBase.
CRSMatrix< real64, globalIndex > m_localMatrix
Local system matrix and rhs.
ParallelVector const & getSystemSolution() const
Getter for system solution vector.
virtual void synchronizeNonlinearSolverParameters()
synchronize the nonlinear solver parameters.
virtual bool registerCallback(void *func, const std::type_info &funcType) final override
function to set the value of m_assemblyCallback
SolverStatistics & getSolverStatistics()
accessor for the solver statistics.
static CONSTITUTIVE_TYPE & getConstitutiveModel(ElementSubRegionBase &subRegion)
Get the Constitutive Model object.
virtual bool execute(real64 const time_n, real64 const dt, integer const cycleNumber, integer const eventCounter, real64 const eventProgress, DomainPartition &domain) override
Main extension point of executable targets.
virtual Group * createChild(string const &childKey, string const &childName) override
creates a child group of of this PhysicsSolverBase instantiation
virtual void saveSequentialIterationState(DomainPartition &domain)
Save the state of the solver for sequential iteration.
PhysicsSolverBase(PhysicsSolverBase &&)=default
Move constructor for PhysicsSolverBase.
real64 getTimestepRequest()
getter for the next timestep size
void setConstitutiveName(ElementSubRegionBase &subRegion, string const &wrapperName, string const &constitutiveType) const
Register wrapper with given name and store constitutive model name on the subregion.
ConvergenceStatistics & getConvergenceStats()
PhysicsSolverBase & operator=(PhysicsSolverBase const &)=delete
Deleted copy assignment operator.
integer m_usePhysicsScaling
Flag to decide whether to apply physics-based scaling to the linear system.
static CatalogInterface::CatalogType & getCatalog()
Get the singleton catalog for PhysicsSolverBase.
PhysicsSolverBase & operator=(PhysicsSolverBase &&)=delete
Deleted move assignment operator.
NonlinearSolverParameters const & getNonlinearSolverParameters() const
const accessor for the nonlinear solver parameters.
DofManager m_dofManager
Data structure to handle degrees of freedom.
R1Tensor const gravityVector() const
return the value of the gravity vector specified in PhysicsSolverManager
integer m_writeStatisticsCSV
virtual void initialize_postMeshGeneration() override
Initialization tasks after mesh generation is completed.
real64 eisenstatWalker(real64 const newNewtonNorm, real64 const oldNewtonNorm, LinearSolverParameters::Krylov const &krylovParams)
Eisenstat-Walker adaptive tolerance.
ParallelVector m_solution
System solution vector.
ParallelVector m_rhs
System right-hand side vector.
std::function< void(CRSMatrix< real64, globalIndex >, array1d< real64 >) > m_assemblyCallback
Callback function for assembly step.
CRSMatrixView< real64 const, globalIndex const > getLocalMatrix() const
Getter for local matrix.
Timestamp getSystemSetupTimestamp() const
getter for the timestamp of the system setup
ParallelVector & getSystemSolution()
Getter for system solution vector.
NonlinearSolverParameters & getNonlinearSolverParameters()
accessor for the nonlinear solver parameters.
Timestamp m_systemSetupTimestamp
Timestamp of the last call to setup system.
LinearSolverParametersInput m_linearSolverParameters
Linear solver parameters.
CRSMatrix< real64, globalIndex > & getLocalMatrix()
Getter for local matrix.
std::unique_ptr< LinearSolverBase< LAInterface > > m_linearSolver
Custom linear solver for the "native" solver type.
NonlinearSolverParameters m_nonlinearSolverParameters
Nonlinear solver parameters.
real64 m_nextDt
timestep of the next cycle
void generateMeshTargetsFromTargetRegions(Group const &meshBodies)
Generate mesh targets from target regions.
DofManager const & getDofManager() const
Getter for degree-of-freedom manager.
This class records solver statistics for each time step.
IterationsStatistics m_iterationsStats
Contain iteration data given a time step.
ConvergenceStatistics m_convergenceStats
Contain convergence data given a time step.
This class provides the base class/interface for the catalog value objects.
std::unordered_map< std::string, std::unique_ptr< CatalogInterface< BASETYPE, ARGS... > > > CatalogType
This is the type that will be used for the catalog. The catalog is actually instantiated in the BASET...
Wrapper< TBASE > & registerWrapper(string const &name, wrapperMap::KeyIndex::index_type *const rkey=nullptr)
Create and register a Wrapper around a new object.
T * getGroupPointer(KEY const &key)
Return a pointer to a sub-group of the current Group.
DataContext const & getDataContext() const
void setRestartFlags(RestartFlags flags)
Set flags that control restart output of this group.
string const & getName() const
Get group name.
GEOS_DECLTYPE_AUTO_RETURN getReference(LOOKUP_TYPE const &lookup) const
Look up a wrapper and get reference to wrapped object.
T & getGroup(KEY const &key)
Return a reference to a sub-group of the current Group.
void forSubGroups(LAMBDA &&lambda)
Apply the given functor to subgroups that can be casted to one of specified types.
Group & setSizedFromParent(int val)
Set whether this wrapper is resized when its parent is resized.
Base template for ordered and unordered maps.
virtual void updateState(DomainPartition &domain)
Recompute all dependent quantities from primary variables (including constitutive models)
virtual void applySystemSolution(DofManager const &dofManager, arrayView1d< real64 const > const &localSolution, real64 const scalingFactor, real64 const dt, DomainPartition &domain)
Function to apply the solution vector to the state.
virtual real64 setNextDtBasedOnStateChange(real64 const ¤tDt, DomainPartition &domain)
function to set the next dt based on state change
virtual void resetStateToBeginningOfStep(DomainPartition &domain)
reset state of physics back to the beginning of the step.
virtual bool lineSearchWithParabolicInterpolation(real64 const &time_n, real64 const &dt, integer const cycleNumber, integer const newtonIter, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, ParallelVector &rhs, ParallelVector &solution, real64 const scaleFactor, real64 &lastResidual, real64 &residualNormT)
Function to perform line search using a parabolic interpolation to find the scaling factor.
virtual void assembleSystem(real64 const time, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs)
function to assemble the linear system matrix and rhs
virtual bool resetConfigurationToDefault(DomainPartition &domain) const
resets the configuration to the default value.
virtual real64 linearImplicitStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain)
Function for a linear implicit integration step.
virtual bool updateConfiguration(DomainPartition &domain)
updates the configuration (if needed) based on the state after a converged Newton loop.
virtual void applyBoundaryConditions(real64 const time, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs)
apply boundary condition to system
virtual real64 nonlinearImplicitStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain)
Function for a nonlinear implicit integration step.
virtual real64 scalingForSystemSolution(DomainPartition &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localSolution)
Function to determine if the solution vector should be scaled back in order to maintain a known const...
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain)
function to perform setup for implicit timestep
virtual std::unique_ptr< PreconditionerBase< LAInterface > > createPreconditioner(DomainPartition &domain) const
Create a preconditioner for this solver's linear system.
void debugOutputSolution(real64 const &time, integer const cycleNumber, integer const nonlinearIteration, ParallelVector const &solution) const
Output the linear system solution for debug purposes.
virtual bool lineSearch(real64 const &time_n, real64 const &dt, integer const cycleNumber, integer const newtonIter, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, ParallelVector &rhs, ParallelVector &solution, real64 const scaleFactor, real64 &lastResidual)
Function to perform line search.
virtual void setupSystem(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, ParallelVector &rhs, ParallelVector &solution, bool const setSparsity=true)
Set up the linear system (DOF indices and sparsity patterns)
virtual real64 setNextDtBasedOnIterNumber(real64 const ¤tDt)
function to set the next time step size based on convergence
virtual bool checkSystemSolution(DomainPartition &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localSolution, real64 const scalingFactor)
Function to check system solution for physical consistency and constraint violation.
virtual real64 calculateResidualNorm(real64 const &time, real64 const &dt, DomainPartition const &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localRhs)
calculate the norm of the global system residual
virtual real64 solverStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain)
entry function to perform a solver step
virtual void setupDofs(DomainPartition const &domain, DofManager &dofManager) const
Populate degree-of-freedom manager with fields relevant to this solver.
virtual void implicitStepComplete(real64 const &time, real64 const &dt, DomainPartition &domain)
perform cleanup for implicit timestep
virtual real64 setNextDt(real64 const ¤tTime, real64 const ¤tDt, DomainPartition &domain)
function to set the next time step size
virtual void updateAndWriteConvergenceStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, integer const iteration)
Update the convergence information and write then into a CSV file.
virtual void outputConfigurationStatistics(DomainPartition const &domain) const
void debugOutputSystem(real64 const &time, integer const cycleNumber, integer const nonlinearIteration, ParallelMatrix const &matrix, ParallelVector const &rhs) const
Output the assembled linear system for debug purposes.
virtual void solveLinearSystem(DofManager const &dofManager, ParallelMatrix &matrix, ParallelVector &rhs, ParallelVector &solution)
function to apply a linear system solver to the assembled system.
virtual real64 getTimestepRequest(real64 const) override
getter for the next timestep size
virtual void resetConfigurationToBeginningOfStep(DomainPartition &domain)
resets the configuration to the beginning of the time-step.
virtual real64 explicitStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain)
Entry function for an explicit time integration step.
@ NOPLOT
Do not ever write to plot file.
@ NO_WRITE
Do not write into restart.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
stdVector< string > string_array
A 1-dimensional array of geos::string types.
unsigned long long int Timestamp
Timestamp type (used to perform actions such a sparsity pattern computation after mesh modifications)
LvArray::CRSMatrix< T, COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > CRSMatrix
Alias for CRS Matrix class.
double real64
64-bit floating point type.
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
LvArray::CRSMatrixView< T, COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
LAInterface::ParallelMatrix ParallelMatrix
Alias for ParallelMatrix.
int integer
Signed integer type.
Array< T, 1 > array1d
Alias for 1D array.
LAInterface::ParallelVector ParallelVector
Alias for ParallelVector.
internal::StdVectorWrapper< T, Allocator, USE_STD_CONTAINER_BOUNDS_CHECKING > stdVector
static constexpr auto constitutiveModelsString()
Krylov-method parameters.
Set of parameters for a linear solver or preconditioner.
Results/stats of a linear solve.
static constexpr char const * baseDiscretizationString()
Structure to hold scoped key names.
static constexpr char const * nonlinearSolverParametersString()
static constexpr char const * solverStatisticsString()
static constexpr char const * linearSolverParametersString()
Structure to hold scoped key names.
static constexpr char const * cflFactorString()
static constexpr char const * discretizationString()
static constexpr char const * writeLinearSystemString()
static constexpr char const * targetRegionsString()
static constexpr char const * usePhysicsScalingString()
static constexpr char const * minDtIncreaseIntervalString()
static constexpr char const * initialDtString()
static constexpr char const * writeStatisticsCSVString()
static constexpr char const * numTimestepsSinceLastDtCutString()
static constexpr char const * allowNonConvergedLinearSolverSolutionString()