GEOS
PhysicsSolverBase.hpp
Go to the documentation of this file.
1 /*
2  * ------------------------------------------------------------------------------------------------------------
3  * SPDX-License-Identifier: LGPL-2.1-only
4  *
5  * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
6  * Copyright (c) 2018-2024 TotalEnergies
7  * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
8  * Copyright (c) 2023-2024 Chevron
9  * Copyright (c) 2019- GEOS/GEOSX Contributors
10  * All rights reserved
11  *
12  * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
13  * ------------------------------------------------------------------------------------------------------------
14  */
15 
20 #ifndef GEOS_PHYSICSSOLVERS_PHYSICSSOLVERBASE_HPP_
21 #define GEOS_PHYSICSSOLVERS_PHYSICSSOLVERBASE_HPP_
22 
23 #include "codingUtilities/traits.hpp"
24 #include "common/DataTypes.hpp"
31 #include "mesh/MeshBody.hpp"
32 #include "physicsSolvers/NonlinearSolverParameters.hpp"
34 #include "physicsSolvers/SolverStatistics.hpp"
35 
36 #include <limits>
37 
38 namespace geos
39 {
40 
41 class DomainPartition;
42 
52 {
53 public:
54 
58  enum class StatsOutputType : integer
59  {
60  none, iteration, convergence, all
61  };
62 
68  explicit PhysicsSolverBase( string const & name,
69  Group * const parent );
70 
75 
79  virtual ~PhysicsSolverBase() override;
80 
84  PhysicsSolverBase() = delete;
85 
89  PhysicsSolverBase( PhysicsSolverBase const & ) = delete;
90 
95 
100 
104  virtual string getCatalogName() const = 0;
105 
106 
111  virtual void registerDataOnMesh( Group & MeshBodies ) override;
112 
116  virtual void initialize_postMeshGeneration() override;
117 
122  void generateMeshTargetsFromTargetRegions( Group const & meshBodies );
123 
127  virtual void cleanup( real64 const time_n,
128  integer const cycleNumber,
129  integer const eventCounter,
130  real64 const eventProgress,
131  DomainPartition & domain ) override;
132 
136  virtual bool execute( real64 const time_n,
137  real64 const dt,
138  integer const cycleNumber,
139  integer const eventCounter,
140  real64 const eventProgress,
141  DomainPartition & domain ) override;
142 
148 
153  ParallelMatrix const & getSystemMatrix() const { return m_matrix; }
154 
160 
165  ParallelVector const & getSystemRhs() const { return m_rhs; }
166 
172 
177  ParallelVector const & getSystemSolution() const { return m_solution; }
178 
184 
189  DofManager const & getDofManager() const { return m_dofManager; }
190 
196 
202 
221  virtual real64 solverStep( real64 const & time_n,
222  real64 const & dt,
223  integer const cycleNumber,
224  DomainPartition & domain );
225 
233  virtual real64 setNextDt( real64 const & currentTime,
234  real64 const & currentDt,
235  DomainPartition & domain );
236 
242  virtual real64 setNextDtBasedOnIterNumber( real64 const & currentDt );
243 
250  virtual real64 setNextDtBasedOnStateChange( real64 const & currentDt,
251  DomainPartition & domain );
252 
261  virtual real64 explicitStep( real64 const & time_n,
262  real64 const & dt,
263  integer const cycleNumber,
264  DomainPartition & domain );
265 
279  virtual real64 nonlinearImplicitStep( real64 const & time_n,
280  real64 const & dt,
281  integer const cycleNumber,
282  DomainPartition & domain );
283 
304  virtual bool
305  lineSearch( real64 const & time_n,
306  real64 const & dt,
307  integer const cycleNumber,
308  integer const newtonIter,
309  DomainPartition & domain,
310  DofManager const & dofManager,
311  CRSMatrixView< real64, globalIndex const > const & localMatrix,
312  ParallelVector & rhs,
313  ParallelVector & solution,
314  real64 const scaleFactor,
315  real64 & lastResidual );
316 
334  virtual bool
336  real64 const & dt,
337  integer const cycleNumber,
338  integer const newtonIter,
339  DomainPartition & domain,
340  DofManager const & dofManager,
341  CRSMatrixView< real64, globalIndex const > const & localMatrix,
342  ParallelVector & rhs,
343  ParallelVector & solution,
344  real64 const scaleFactor,
345  real64 & lastResidual,
346  real64 & residualNormT );
347 
362  virtual real64 linearImplicitStep( real64 const & time_n,
363  real64 const & dt,
364  integer const cycleNumber,
365  DomainPartition & domain );
366 
379  virtual void
380  implicitStepSetup( real64 const & time_n,
381  real64 const & dt,
382  DomainPartition & domain );
383 
389  virtual void
390  setupDofs( DomainPartition const & domain,
391  DofManager & dofManager ) const;
392 
405  virtual void
407  DofManager & dofManager,
408  CRSMatrix< real64, globalIndex > & localMatrix,
409  ParallelVector & rhs,
410  ParallelVector & solution,
411  bool const setSparsity = true );
412 
420  virtual void
422  DofManager & dofManager,
423  CRSMatrix< real64, globalIndex > & localMatrix,
424  SparsityPattern< globalIndex > & pattern );
425 
431  virtual std::unique_ptr< PreconditionerBase< LAInterface > >
433 
452  virtual void
453  assembleSystem( real64 const time,
454  real64 const dt,
455  DomainPartition & domain,
456  DofManager const & dofManager,
457  CRSMatrixView< real64, globalIndex const > const & localMatrix,
458  arrayView1d< real64 > const & localRhs );
459 
472  virtual void
474  real64 const dt,
475  DomainPartition & domain,
476  DofManager const & dofManager,
477  CRSMatrixView< real64, globalIndex const > const & localMatrix,
478  arrayView1d< real64 > const & localRhs );
479 
488  void
489  debugOutputSystem( real64 const & time,
490  integer const cycleNumber,
491  integer const nonlinearIteration,
492  ParallelMatrix const & matrix,
493  ParallelVector const & rhs ) const;
494 
502  void
504  integer const cycleNumber,
505  integer const nonlinearIteration,
506  ParallelVector const & solution ) const;
507 
515  virtual void
517  real64 const & dt,
518  integer const cycleNumber,
519  integer const iteration );
520 
533  virtual real64
535  real64 const & dt,
536  DomainPartition const & domain,
537  DofManager const & dofManager,
538  arrayView1d< real64 const > const & localRhs );
539 
553  virtual void
554  solveLinearSystem( DofManager const & dofManager,
555  ParallelMatrix & matrix,
556  ParallelVector & rhs,
557  ParallelVector & solution );
558 
571  virtual bool
573  DofManager const & dofManager,
574  arrayView1d< real64 const > const & localSolution,
575  real64 const scalingFactor );
576 
584  virtual real64
586  DofManager const & dofManager,
587  arrayView1d< real64 const > const & localSolution );
588 
610  virtual void
611  applySystemSolution( DofManager const & dofManager,
612  arrayView1d< real64 const > const & localSolution,
613  real64 const scalingFactor,
614  real64 const dt,
615  DomainPartition & domain );
616 
623  virtual bool updateConfiguration( DomainPartition & domain,
624  integer configurationLoopIter );
625 
630  virtual void outputConfigurationStatistics( DomainPartition const & domain ) const;
631 
637 
643  virtual bool resetConfigurationToDefault( DomainPartition & domain ) const;
644 
645 
650  virtual void updateState( DomainPartition & domain );
651 
663  virtual void
665 
679  virtual void
681  real64 const & dt,
682  DomainPartition & domain );
683 
684 
689  virtual real64 getTimestepRequest( real64 const ) override
690  {return m_nextDt;};
698  {return m_nextDt;};
699 
706  virtual Group * createChild( string const & childKey, string const & childName ) override;
707 
712 
718 
723  {
725  static constexpr char const * cflFactorString() { return "cflFactor"; }
726 
728  static constexpr char const * initialDtString() { return "initialDt"; }
729 
731  static constexpr char const * minDtIncreaseIntervalString() { return "minDtIncreaseInterval"; }
732 
734  static constexpr char const * discretizationString() { return "discretization"; }
735 
737  static constexpr char const * targetRegionsString() { return "targetRegions"; }
738 
740  static constexpr char const * writeLinearSystemString() { return "writeLinearSystem"; }
741 
743  static constexpr char const * usePhysicsScalingString() { return "usePhysicsScaling"; }
744 
746  static constexpr char const * allowNonConvergedLinearSolverSolutionString() { return "allowNonConvergedLinearSolverSolution"; }
747 
749  static constexpr char const * writeStatisticsCSVString() { return "writeStatistics"; }
750 
752  static constexpr char const * numTimestepsSinceLastDtCutString() { return "numTimestepsSinceLastDtCut"; }
753  };
754 
759  {
761  static constexpr char const * linearSolverParametersString() { return "LinearSolverParameters"; }
762 
764  static constexpr char const * nonlinearSolverParametersString() { return "NonlinearSolverParameters"; }
765 
767  static constexpr char const * solverStatisticsString() { return "SolverStatistics"; }
768  };
769 
775 
782 
788 
797  R1Tensor const gravityVector() const;
798 
804  virtual bool checkSequentialSolutionIncrements( DomainPartition & domain ) const;
805 
811 
817  {
818  return m_linearSolverParameters.get();
819  }
820 
826  {
827  return m_linearSolverParameters.get();
828  }
829 
835  {
837  }
838 
844  {
846  }
847 
851  virtual void
853  { /* empty here, overriden in CoupledSolver */ }
854 
860  localIndex targetRegionIndex( string const & regionName ) const;
861 
866  string_array const & getTargetRegionNames() const {return m_targetRegionNames;}
867 
868 
876  template< typename LAMBDA >
877  void forDiscretizationOnMeshTargets( Group const & meshBodies, LAMBDA && lambda ) const
878  {
879  for( auto const & target: m_meshTargets )
880  {
881  string const meshBodyName = target.first.first;
882  string const meshLevelName = target.first.second;
883  string_array const & regionNames = target.second;
884  MeshBody const & meshBody = meshBodies.getGroup< MeshBody >( meshBodyName );
885 
886  MeshLevel const * meshLevelPtr = meshBody.getMeshLevels().getGroupPointer< MeshLevel >( meshLevelName );
887  if( meshLevelPtr==nullptr )
888  {
890  }
891  lambda( meshBodyName, *meshLevelPtr, regionNames );
892  }
893  }
894 
902  template< typename LAMBDA >
903  void forDiscretizationOnMeshTargets( Group & meshBodies, LAMBDA && lambda ) const
904  {
905  for( auto const & target: m_meshTargets )
906  {
907  string const meshBodyName = target.first.first;
908  string const meshLevelName = target.first.second;
909  string_array const & regionNames = target.second;
910  MeshBody & meshBody = meshBodies.getGroup< MeshBody >( meshBodyName );
911 
912  MeshLevel * meshLevelPtr = meshBody.getMeshLevels().getGroupPointer< MeshLevel >( meshLevelName );
913  if( meshLevelPtr==nullptr )
914  {
916  }
917  lambda( meshBodyName, *meshLevelPtr, regionNames );
918  }
919  }
920 
926 
935  virtual bool registerCallback( void * func, const std::type_info & funcType ) final override;
936 
942  {
944  }
951  {
953  }
958  {
960  }
965  {
967  }
968 
974 
980 
981 #if defined(GEOS_USE_PYGEOSX)
986  virtual PyTypeObject * getPythonType() const override;
987 #endif
988 
994  {
995  return m_meshTargets;
996  }
997 
1002  bool detectOscillations() const;
1003 
1004 protected:
1005 
1006  virtual void postInputInitialization() override;
1007 
1033  real64 eisenstatWalker( real64 const newNewtonNorm,
1034  real64 const oldNewtonNorm,
1035  LinearSolverParameters::Krylov const & krylovParams );
1036 
1044  template< typename CONSTITUTIVE_BASE_TYPE >
1045  static string getConstitutiveName( ElementSubRegionBase const & subRegion );
1046 
1054  template< typename CONSTITUTIVE_BASE_TYPE >
1055  static string getConstitutiveName( ParticleSubRegionBase const & subRegion ); // particle overload
1056 
1065  template< typename CONSTITUTIVE >
1066  void setConstitutiveName( ElementSubRegionBase & subRegion, string const & wrapperName, string const & constitutiveType ) const;
1067 
1074  virtual void setConstitutiveNamesCallSuper( ElementSubRegionBase & subRegion ) const { GEOS_UNUSED_VAR( subRegion ); }
1075 
1076 
1085  template< typename BASETYPE = constitutive::ConstitutiveBase, typename LOOKUP_TYPE >
1086  static BASETYPE const & getConstitutiveModel( dataRepository::Group const & dataGroup, LOOKUP_TYPE const & key )
1087  {
1089  return constitutiveModels.getGroup< BASETYPE >( key );
1090  }
1091 
1100  template< typename BASETYPE = constitutive::ConstitutiveBase, typename LOOKUP_TYPE >
1101  static BASETYPE & getConstitutiveModel( dataRepository::Group & dataGroup, LOOKUP_TYPE const & key )
1102  {
1104  return constitutiveModels.getGroup< BASETYPE >( key );
1105  }
1106 
1113  template< typename CONSTITUTIVE_TYPE >
1114  static CONSTITUTIVE_TYPE & getConstitutiveModel( ElementSubRegionBase & subRegion )
1115  {
1116  return getConstitutiveModel< CONSTITUTIVE_TYPE >( subRegion, getConstitutiveName< CONSTITUTIVE_TYPE >( subRegion ) );
1117  }
1118 
1121 
1124 
1127 
1130 
1133 
1136 
1139 
1142 
1145 
1148 
1151 
1154 
1156  std::unique_ptr< LinearSolverBase< LAInterface > > m_linearSolver;
1157 
1159  std::unique_ptr< PreconditionerBase< LAInterface > > m_precond;
1160 
1163 
1167 
1170 
1173 
1176 
1179 
1182 
1185 
1187  std::map< std::string, std::chrono::system_clock::duration > m_timers;
1188 
1191 
1192 private:
1194  string_array m_targetRegionNames;
1195 
1198 
1205  virtual void setConstitutiveNames( ElementSubRegionBase & subRegion ) const { GEOS_UNUSED_VAR( subRegion ); }
1206 
1215  bool solveNonlinearSystem( real64 const & time_n,
1216  real64 const & dt,
1217  integer const cycleNumber,
1218  DomainPartition & domain );
1219 
1226  void logEndOfCycleInformation( integer const cycleNumber,
1227  integer const numOfSubSteps,
1228  stdVector< real64 > const & subStepDts ) const;
1229 };
1230 
1231 
1232 template< typename CONSTITUTIVE_BASE_TYPE >
1234 {
1235  string validName;
1236  dataRepository::Group const & constitutiveModels = subRegion.getConstitutiveModels();
1237 
1238  constitutiveModels.forSubGroups< CONSTITUTIVE_BASE_TYPE >( [&]( dataRepository::Group const & model )
1239  {
1240  GEOS_ERROR_IF( !validName.empty(), "A valid constitutive model was already found." );
1241  validName = model.getName();
1242  } );
1243 
1244  return validName;
1245 }
1246 
1247 template< typename CONSTITUTIVE_BASE_TYPE >
1248 string PhysicsSolverBase::getConstitutiveName( ParticleSubRegionBase const & subRegion ) // particle overload
1249 {
1250  string validName;
1251  dataRepository::Group const & constitutiveModels = subRegion.getConstitutiveModels();
1252 
1253  constitutiveModels.forSubGroups< CONSTITUTIVE_BASE_TYPE >( [&]( dataRepository::Group const & model )
1254  {
1255  GEOS_ERROR_IF( !validName.empty(), "A valid constitutive model was already found." );
1256  validName = model.getName();
1257  } );
1258  return validName;
1259 }
1260 
1261 template< typename CONSTITUTIVE >
1262 void PhysicsSolverBase::setConstitutiveName( ElementSubRegionBase & subRegion, string const & wrapperName, string const & constitutiveType ) const
1263 {
1264  subRegion.registerWrapper< string >( wrapperName ).
1265  setPlotLevel( dataRepository::PlotLevel::NOPLOT ).
1267  setSizedFromParent( 0 );
1268 
1269  string & constitutiveName = subRegion.getReference< string >( wrapperName );
1270  constitutiveName = getConstitutiveName< CONSTITUTIVE >( subRegion );
1271  GEOS_ERROR_IF( constitutiveName.empty(), GEOS_FMT( "{}: {} constitutive model not found on subregion {}",
1272  getDataContext(), constitutiveType, subRegion.getName() ) );
1273 }
1274 
1279  "none",
1280  "iteration",
1281  "convergence",
1282  "all" );
1283 
1284 } // namespace geos
1285 
1286 
1287 #endif /* GEOS_PHYSICSSOLVERS_PHYSICSSOLVERBASE_HPP_ */
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
Definition: GeosxMacros.hpp:84
#define GEOS_ERROR_IF(EXP, msg)
Conditionally raise a hard error and terminate the program.
Definition: Logger.hpp:142
Class containing convergence information given a time-step.
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Definition: DofManager.hpp:45
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.
Linear solver parameters with Group capabilities.
The class is used to manage mesh body.
Definition: MeshBody.hpp:36
Group & getMeshLevels()
Get the meshLevels group.
Definition: MeshBody.hpp:86
Class facilitating the representation of a multi-level discretization of a MeshBody.
Definition: MeshLevel.hpp:42
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
StatsOutputType m_writeStatisticsCSV
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.
bool detectOscillations() const
Detect oscillations in the solution.
string_array const & getTargetRegionNames() const
return the list of target regions
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...
StatsOutputType
Type of the stat output.
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.
ConvergenceStatistics const & getConvergenceStats() const
localIndex targetRegionIndex(string const &regionName) 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.
IterationsStatistics const & getIterationStats() const
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
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.
ArrayOfArrays< real64 > m_solutionHistory
History of the solution vector, used for oscillation detection.
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.
Definition: Group.hpp:299
DataContext const & getDataContext() const
Definition: Group.hpp:1345
void setRestartFlags(RestartFlags flags)
Set flags that control restart output of this group.
Definition: Group.hpp:1426
string const & getName() const
Get group name.
Definition: Group.hpp:1331
GEOS_DECLTYPE_AUTO_RETURN getReference(LOOKUP_TYPE const &lookup) const
Look up a wrapper and get reference to wrapped object.
Definition: Group.hpp:1275
T & getGroup(KEY const &key)
Return a reference to a sub-group of the current Group.
Definition: Group.hpp:318
void forSubGroups(LAMBDA &&lambda)
Apply the given functor to subgroups that can be casted to one of specified types.
Definition: Group.hpp:500
Group & setSizedFromParent(int val)
Set whether this wrapper is resized when its parent is resized.
Definition: Group.hpp:1411
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 &currentDt, DomainPartition &domain)
function to set the next dt based on state change
virtual bool updateConfiguration(DomainPartition &domain, integer configurationLoopIter)
updates the configuration (if needed) based on the state after a converged Newton loop.
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 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.
virtual void setSparsityPattern(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, SparsityPattern< globalIndex > &pattern)
Set the sparsity pattern of the linear system matrix.
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 &currentDt)
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 &currentTime, real64 const &currentDt, 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.
Definition: DataTypes.hpp:179
stdVector< string > string_array
A 1-dimensional array of geos::string types.
Definition: DataTypes.hpp:361
unsigned long long int Timestamp
Timestamp type (used to perform actions such a sparsity pattern computation after mesh modifications)
Definition: DataTypes.hpp:126
LvArray::CRSMatrix< T, COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > CRSMatrix
Alias for CRS Matrix class.
Definition: DataTypes.hpp:305
LvArray::SparsityPattern< COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > SparsityPattern
Alias for Sparsity pattern class.
Definition: DataTypes.hpp:297
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
Definition: DataTypes.hpp:84
LvArray::CRSMatrixView< T, COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
Definition: DataTypes.hpp:309
LAInterface::ParallelMatrix ParallelMatrix
Alias for ParallelMatrix.
int integer
Signed integer type.
Definition: DataTypes.hpp:81
Array< T, 1 > array1d
Alias for 1D array.
Definition: DataTypes.hpp:175
LvArray::ArrayOfArrays< T, INDEX_TYPE, LvArray::ChaiBuffer > ArrayOfArrays
Array of variable-sized arrays. See LvArray::ArrayOfArrays for details.
Definition: DataTypes.hpp:281
ENUM_STRINGS(LinearSolverParameters::SolverType, "direct", "cg", "gmres", "fgmres", "bicgstab", "richardson", "preconditioner")
Declare strings associated with enumeration values.
LAInterface::ParallelVector ParallelVector
Alias for ParallelVector.
internal::StdVectorWrapper< T, Allocator, USE_STD_CONTAINER_BOUNDS_CHECKING > stdVector
Set of parameters for a linear solver or preconditioner.
Results/stats of a linear solve.
static constexpr char const * baseDiscretizationString()
Definition: MeshBody.hpp:220
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()