GEOS
CoupledSolver.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 
21 #ifndef GEOS_PHYSICSSOLVERS_MULTIPHYSICS_COUPLEDSOLVER_HPP_
22 #define GEOS_PHYSICSSOLVERS_MULTIPHYSICS_COUPLEDSOLVER_HPP_
23 
26 
27 #include <tuple>
28 
29 namespace geos
30 {
31 
32 template< typename ... SOLVERS >
34 {
35 
36 public:
37 
43  CoupledSolver( const string & name,
44  Group * const parent )
45  : PhysicsSolverBase( name, parent )
46  {
47  forEachArgInTuple( m_solvers, [&]( auto solver, auto idx )
48  {
49  using SolverType = TYPEOFPTR( solver );
50  string const key = SolverType::coupledSolverAttributePrefix() + "SolverName";
51  registerWrapper( key, &m_names[idx()] ).
52  setRTTypeName( rtTypes::CustomTypes::groupNameRef ).
54  setDescription( "Name of the " + SolverType::coupledSolverAttributePrefix() + " solver used by the coupled solver" );
55  } );
56 
57  this->getWrapper< string >( PhysicsSolverBase::viewKeyStruct::discretizationString() ).
58  setInputFlag( dataRepository::InputFlags::FALSE );
59 
60  addLogLevel< logInfo::Coupling >();
61  }
62 
64  CoupledSolver( CoupledSolver const & ) = delete;
65 
67  CoupledSolver( CoupledSolver && ) = default;
68 
70  CoupledSolver & operator=( CoupledSolver const & ) = delete;
71 
74 
75 
79  void
81  {
82  forEachArgInTuple( m_solvers, [&]( auto & solver, auto idx )
83  {
84  using SolverPtr = TYPEOFREF( solver );
85  using SolverType = TYPEOFPTR( SolverPtr {} );
86  auto const & solverName = m_names[idx()];
87  auto const & solverType = LvArray::system::demangleType< SolverType >();
88  solver = this->getParent().template getGroupPointer< SolverType >( solverName );
89  GEOS_THROW_IF( solver == nullptr,
90  GEOS_FMT( "{}: Could not find solver '{}' of type {}",
92  solverName, solverType ),
93  InputError );
94  GEOS_LOG_LEVEL_RANK_0( logInfo::Coupling,
95  GEOS_FMT( "{}: found {} solver named {}",
96  getName(), solver->getCatalogName(), solverName ) );
97  } );
98  }
99 
100 
106  virtual void
108  DofManager & dofManager ) const
109  { GEOS_UNUSED_VAR( domain, dofManager ); }
110 
120  virtual void
122  real64 const dt,
123  DomainPartition const & domain,
124  DofManager const & dofManager,
125  CRSMatrixView< real64, globalIndex const > const & localMatrix,
126  arrayView1d< real64 > const & localRhs )
127  { GEOS_UNUSED_VAR( time_n, dt, domain, dofManager, localMatrix, localRhs ); }
128 
136  void
137  setupDofs( DomainPartition const & domain,
138  DofManager & dofManager ) const override
139  {
140  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
141  {
142  solver->setupDofs( domain, dofManager );
143  } );
144 
145  setupCoupling( domain, dofManager );
146  }
147 
148  virtual void
149  implicitStepSetup( real64 const & time_n,
150  real64 const & dt,
151  DomainPartition & domain ) override
152  {
153  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
154  {
155  solver->implicitStepSetup( time_n, dt, domain );
156  } );
157  }
158 
159  virtual void
160  implicitStepComplete( real64 const & time_n,
161  real64 const & dt,
162  DomainPartition & domain ) override
163  {
164  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
165  {
166  solver->implicitStepComplete( time_n, dt, domain );
167  } );
168  }
169 
170  // general version of assembleSystem function, keep in mind many solvers will override it
171  virtual void
172  assembleSystem( real64 const time_n,
173  real64 const dt,
174  DomainPartition & domain,
175  DofManager const & dofManager,
176  CRSMatrixView< real64, globalIndex const > const & localMatrix,
177  arrayView1d< real64 > const & localRhs ) override
178  {
180 
181  // 1. Assemble matrix blocks of each individual solver
182  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
183  {
184  solver->assembleSystem( time_n, dt, domain, dofManager, localMatrix, localRhs );
185  } );
186 
187  // 2. Assemble coupling blocks
188  assembleCouplingTerms( time_n, dt, domain, dofManager, localMatrix, localRhs );
189  }
190 
191  virtual void
192  applySystemSolution( DofManager const & dofManager,
193  arrayView1d< real64 const > const & localSolution,
194  real64 const scalingFactor,
195  real64 const dt,
196  DomainPartition & domain ) override
197  {
198  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
199  {
200  solver->applySystemSolution( dofManager, localSolution, scalingFactor, dt, domain );
201  } );
202  }
203 
204  virtual void
205  updateState( DomainPartition & domain ) override
206  {
207  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
208  {
209  solver->updateState( domain );
210  } );
211  }
212 
213  virtual void
215  {
216  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
217  {
218  solver->resetStateToBeginningOfStep( domain );
219  } );
220  }
221 
224  real64
225  solverStep( real64 const & time_n,
226  real64 const & dt,
227  int const cycleNumber,
228  DomainPartition & domain ) override final
229  {
231 
233  {
234  return fullyCoupledSolverStep( time_n, dt, cycleNumber, domain );
235  }
237  {
238  return sequentiallyCoupledSolverStep( time_n, dt, cycleNumber, domain );
239  }
240  else
241  {
242  GEOS_ERROR( getDataContext() << ": Invalid coupling type option." );
243  return 0;
244  }
245 
246  }
247 
248 
249  virtual real64
250  calculateResidualNorm( real64 const & time_n,
251  real64 const & dt,
252  DomainPartition const & domain,
253  DofManager const & dofManager,
254  arrayView1d< real64 const > const & localRhs ) override
255  {
256  real64 norm = 0.0;
257  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
258  {
259  real64 const singlePhysicsNorm = solver->calculateResidualNorm( time_n, dt, domain, dofManager, localRhs );
260  norm += singlePhysicsNorm * singlePhysicsNorm;
261  } );
262 
263  return sqrt( norm );
264  }
265 
266  virtual void
268  real64 const dt,
269  DomainPartition & domain,
270  DofManager const & dofManager,
271  CRSMatrixView< real64, globalIndex const > const & localMatrix,
272  arrayView1d< real64 > const & localRhs ) override
273  {
274  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
275  {
276  solver->applyBoundaryConditions( time_n, dt, domain, dofManager, localMatrix, localRhs );
277  } );
278  }
279 
280  virtual bool
282  DofManager const & dofManager,
283  arrayView1d< real64 const > const & localSolution,
284  real64 const scalingFactor ) override
285  {
286  bool validSolution = true;
287  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
288  {
289  bool const validSinglePhysicsSolution = solver->checkSystemSolution( domain, dofManager, localSolution, scalingFactor );
290  if( !validSinglePhysicsSolution )
291  {
292  GEOS_LOG_RANK_0( GEOS_FMT( " {}/{}: Solution check failed. Newton loop terminated.", getName(), solver->getName()) );
293  }
294  validSolution = validSolution && validSinglePhysicsSolution;
295  } );
296  return validSolution;
297  }
298 
299  virtual real64
301  DofManager const & dofManager,
302  arrayView1d< real64 const > const & localSolution ) override
303  {
304  real64 scalingFactor = 1e9;
305  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
306  {
307  real64 const singlePhysicsScalingFactor = solver->scalingForSystemSolution( domain, dofManager, localSolution );
308  scalingFactor = LvArray::math::min( scalingFactor, singlePhysicsScalingFactor );
309  } );
310  return scalingFactor;
311  }
312 
313  virtual real64
314  setNextDt( real64 const & currentTime,
315  real64 const & currentDt,
316  DomainPartition & domain ) override
317  {
318  real64 nextDt = PhysicsSolverBase::setNextDt( currentTime, currentDt, domain );
319  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
320  {
321  real64 const singlePhysicsNextDt =
322  solver->setNextDt( currentTime, currentDt, domain );
323  nextDt = LvArray::math::min( singlePhysicsNextDt, nextDt );
324  } );
325  return nextDt;
326  }
327 
328  virtual void cleanup( real64 const time_n,
329  integer const cycleNumber,
330  integer const eventCounter,
331  real64 const eventProgress,
332  DomainPartition & domain ) override
333  {
334  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
335  {
336  solver->cleanup( time_n, cycleNumber, eventCounter, eventProgress, domain );
337  } );
338  PhysicsSolverBase::cleanup( time_n, cycleNumber, eventCounter, eventProgress, domain );
339  }
340 
343  virtual bool checkSequentialSolutionIncrements( DomainPartition & domain ) const override
344  {
345  bool isConverged = true;
346  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
347  {
348  isConverged &= solver->checkSequentialSolutionIncrements( domain );
349  } );
350  return isConverged;
351  }
352 
353  virtual bool updateConfiguration( DomainPartition & domain ) override
354  {
355  bool result = true;
356  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
357  {
358  result &= solver->updateConfiguration( domain );
359  } );
360  return result;
361  }
362 
363  virtual void outputConfigurationStatistics( DomainPartition const & domain ) const override
364  {
365  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
366  {
367  solver->outputConfigurationStatistics( domain );
368  } );
369  }
370 
371  virtual void resetConfigurationToBeginningOfStep( DomainPartition & domain ) override
372  {
373  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
374  {
375  solver->resetConfigurationToBeginningOfStep( domain );
376  } );
377  }
378 
379  virtual bool resetConfigurationToDefault( DomainPartition & domain ) const override
380  {
381  bool result = true;
382  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
383  {
384  result &= solver->resetConfigurationToDefault( domain );
385  } );
386  return result;
387  }
388 
389 protected:
390 
400  virtual real64 fullyCoupledSolverStep( real64 const & time_n,
401  real64 const & dt,
402  int const cycleNumber,
403  DomainPartition & domain )
404  {
405  return PhysicsSolverBase::solverStep( time_n, dt, cycleNumber, domain );
406  }
407 
418  virtual real64 sequentiallyCoupledSolverStep( real64 const & time_n,
419  real64 const & dt,
420  int const cycleNumber,
421  DomainPartition & domain )
422  {
424 
425  // Only build the sparsity pattern if the mesh has changed
426  Timestamp const meshModificationTimestamp = getMeshModificationTimestamp( domain );
427  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
428  {
429  if( meshModificationTimestamp > solver->getSystemSetupTimestamp() )
430  {
431  solver->setupSystem( domain,
432  solver->getDofManager(),
433  solver->getLocalMatrix(),
434  solver->getSystemRhs(),
435  solver->getSystemSolution() );
436  solver->setSystemSetupTimestamp( meshModificationTimestamp );
437  }
438  } );
439 
440  implicitStepSetup( time_n, dt, domain );
441 
443  integer const maxNumberDtCuts = solverParams.m_maxTimeStepCuts;
444  real64 const dtCutFactor = solverParams.m_timeStepCutFactor;
445  integer & dtAttempt = solverParams.m_numTimeStepAttempts;
446 
447  bool isConverged = false;
448  // dt may be cut during the course of this step, so we are keeping a local
449  // value to track the achieved dt for this step.
450  real64 stepDt = dt;
451 
452  // outer loop attempts to apply full timestep, and managed the cutting of the timestep if
453  // required.
454  for( dtAttempt = 0; dtAttempt < maxNumberDtCuts; ++dtAttempt )
455  {
456  // TODO configuration loop
457 
458  // Reset the states of all solvers if any of them had to restart
459  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
460  {
461  solver->resetStateToBeginningOfStep( domain );
462  solver->getSolverStatistics().initializeTimeStepStatistics(); // initialize counters for subsolvers
463  } );
464  resetStateToBeginningOfStep( domain );
465 
466  integer & iter = solverParams.m_numNewtonIterations;
467 
469  for( iter = 0; iter < solverParams.m_maxIterNewton; iter++ )
470  {
471  // Increment the solver statistics for reporting purposes
472  // Pass a "0" as argument (0 linear iteration) to skip the output of linear iteration stats at the end
474 
475  startSequentialIteration( iter, domain );
476 
477  // Solve the subproblems nonlinearly
478  forEachArgInTuple( m_solvers, [&]( auto & solver, auto idx )
479  {
480  GEOS_LOG_LEVEL_RANK_0( logInfo::NonlinearSolver,
481  GEOS_FMT( " Iteration {:2}: {}", iter + 1, solver->getName() ) );
482  real64 solverDt = solver->nonlinearImplicitStep( time_n,
483  stepDt,
484  cycleNumber,
485  domain );
486 
487  // save fields (e.g. pressure and temperature) after inner solve
488  if( solver->getNonlinearSolverParameters().couplingType() == NonlinearSolverParameters::CouplingType::Sequential )
489  {
490  solver->saveSequentialIterationState( domain );
491  }
492 
493  mapSolutionBetweenSolvers( domain, idx() );
494 
495  if( solverDt < stepDt ) // subsolver had to cut the time step
496  {
497  iter = 0; // restart outer loop
498  stepDt = solverDt; // sync time step
500  }
501  } );
502 
503  // Check convergence of the outer loop
504  isConverged = checkSequentialConvergence( iter,
505  time_n,
506  stepDt,
507  domain );
508 
509  if( isConverged )
510  {
511  // we still want to count current iteration
512  ++iter;
513  // exit outer loop
514  break;
515  }
516  else
517  {
518  finishSequentialIteration( iter, domain );
519  }
520  }
521 
522  if( isConverged )
523  {
524  // Save time step statistics for the subsolvers
525  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
526  {
527  solver->getSolverStatistics().saveTimeStepStatistics();
528  } );
529  // get out of the time loop
530  break;
531  }
532  else
533  {
534  // cut timestep, go back to beginning of step and restart the Newton loop
535  stepDt *= dtCutFactor;
537  GEOS_LOG_LEVEL_RANK_0( logInfo::TimeStep, GEOS_FMT( "New dt = {}", stepDt ) );
538 
539  // notify the solver statistics counter that this is a time step cut
541  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
542  {
543  solver->getSolverStatistics().logTimeStepCut();
544  } );
545  }
546  }
547 
548  if( !isConverged )
549  {
550  GEOS_LOG_RANK_0( "Convergence not achieved." );
551 
553  {
554  GEOS_LOG_RANK_0( "The accepted solution may be inaccurate." );
555  }
556  else
557  {
558  GEOS_ERROR( "Nonconverged solutions not allowed. Terminating..." );
559  }
560  }
561 
562  implicitStepComplete( time_n, stepDt, domain );
563 
564  return stepDt;
565  }
566 
574  integer const solverType )
575  {
576  GEOS_UNUSED_VAR( domain, solverType );
577  }
578 
579  virtual bool checkSequentialConvergence( int const & iter,
580  real64 const & time_n,
581  real64 const & dt,
582  DomainPartition & domain )
583  {
585  bool isConverged = true;
586 
587  if( params.m_subcyclingOption == 0 )
588  {
589  GEOS_LOG_LEVEL_RANK_0( logInfo::Convergence, "***** Single Pass solver, no subcycling *****" );
590  }
591  else
592  {
593  GEOS_LOG_LEVEL_RANK_0( logInfo::Convergence, GEOS_FMT( " Iteration {:2}: outer-loop convergence check", iter + 1 ) );
594 
596  {
597  real64 residualNorm = 0;
598 
599  // loop over all the single-physics solvers
600  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
601  {
602 
603  solver->getLocalMatrix().toViewConstSizes().zero();
604  solver->getSystemRhs().zero();
605  arrayView1d< real64 > const localRhs = solver->getSystemRhs().open();
606 
607  // for each solver, we have to recompute the residual (and Jacobian, although not necessary)
608  solver->assembleSystem( time_n,
609  dt,
610  domain,
611  solver->getDofManager(),
612  solver->getLocalMatrix().toViewConstSizes(),
613  localRhs );
614  solver->applyBoundaryConditions( time_n,
615  dt,
616  domain,
617  solver->getDofManager(),
618  solver->getLocalMatrix().toViewConstSizes(),
619  localRhs );
620  solver->getSystemRhs().close();
621 
622  // once this is done, we recompute the single-physics residual
623  real64 const singlePhysicsNorm =
624  solver->calculateResidualNorm( time_n,
625  dt,
626  domain,
627  solver->getDofManager(),
628  solver->getSystemRhs().values() );
629  residualNorm += singlePhysicsNorm * singlePhysicsNorm;
630  } );
631 
632  // finally, we perform the convergence check on the multiphysics residual
633  residualNorm = sqrt( residualNorm );
634  GEOS_LOG_LEVEL_RANK_0( logInfo::Convergence,
635  GEOS_FMT( " ( R ) = ( {:4.2e} )", residualNorm ) );
636  isConverged = ( residualNorm < params.m_newtonTol );
637 
638  }
640  {
641  // TODO also make recursive?
642  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
643  {
644  NonlinearSolverParameters const & singlePhysicsParams = solver->getNonlinearSolverParameters();
645  if( singlePhysicsParams.m_numNewtonIterations > singlePhysicsParams.m_minIterNewton )
646  {
647  isConverged = false;
648  }
649  } );
650  }
652  {
653  isConverged = checkSequentialSolutionIncrements( domain );
654  }
655  else
656  {
657  GEOS_ERROR( getDataContext() << ": Invalid sequential convergence criterion." );
658  }
659 
660  if( isConverged )
661  {
662  GEOS_LOG_LEVEL_RANK_0( logInfo::Convergence,
663  GEOS_FMT( "***** The iterative coupling has converged in {} iteration(s) *****", iter + 1 ) );
664  }
665  }
666  return isConverged;
667  }
668 
669  virtual void
671  {
672  setSubSolvers();
673 
676  GEOS_THROW_IF( isSequential && usesLineSearch,
677  GEOS_FMT( "{}: line search is not supported by the coupled solver when {} is set to `{}`. Please set {} to `{}` to remove this error",
678  getNonlinearSolverParameters().getWrapperDataContext( NonlinearSolverParameters::viewKeysStruct::couplingTypeString() ),
679  NonlinearSolverParameters::viewKeysStruct::couplingTypeString(),
681  NonlinearSolverParameters::viewKeysStruct::lineSearchActionString(),
683  InputError );
684 
685  if( !isSequential )
686  {
688  }
689 
691  validateNonlinearAcceleration();
692  }
693 
694  virtual void validateNonlinearAcceleration()
695  {
696  GEOS_THROW ( GEOS_FMT( "{}: Nonlinear acceleration {} is not supported by {} solver '{}'",
697  getWrapperDataContext( NonlinearSolverParameters::viewKeysStruct::nonlinearAccelerationTypeString() ),
699  getCatalogName(), getName()),
700  InputError );
701  }
702 
703  virtual void
705  {
706  forEachArgInTuple( m_solvers, [&]( auto & solver, auto )
707  {
708  solver->getNonlinearSolverParameters() = getNonlinearSolverParameters();
709  } );
710  }
711 
712  virtual void startSequentialIteration( integer const & iter,
713  DomainPartition & domain )
714  {
715  GEOS_UNUSED_VAR( iter, domain );
716  }
717 
718  virtual void finishSequentialIteration( integer const & iter,
719  DomainPartition & domain )
720  {
721  GEOS_UNUSED_VAR( iter, domain );
722  }
723 
724 protected:
725 
727  std::tuple< SOLVERS *... > m_solvers;
728 
730  std::array< string, sizeof...( SOLVERS ) > m_names;
731 };
732 
733 } /* namespace geos */
734 
735 #endif /* GEOS_PHYSICSSOLVERS_MULTIPHYSICS_COUPLEDSOLVER_HPP_ */
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
Definition: GeosxMacros.hpp:84
#define GEOS_THROW(msg, TYPE)
Throw an exception.
Definition: Logger.hpp:164
#define GEOS_ERROR(msg)
Raise a hard error and terminate the program.
Definition: Logger.hpp:157
#define GEOS_LOG_RANK_0(msg)
Log a message on screen on rank 0.
Definition: Logger.hpp:101
#define GEOS_THROW_IF(EXP, msg, TYPE)
Conditionally throw an exception.
Definition: Logger.hpp:151
#define GEOS_MARK_FUNCTION
Mark function with both Caliper and NVTX if enabled.
virtual void assembleCouplingTerms(real64 const time_n, real64 const dt, DomainPartition const &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs)
Utility function to compute coupling terms.
CoupledSolver & operator=(CoupledSolver const &)=delete
deleted assignment operator
virtual real64 fullyCoupledSolverStep(real64 const &time_n, real64 const &dt, int const cycleNumber, DomainPartition &domain)
Fully coupled solution approach solution step.
virtual bool resetConfigurationToDefault(DomainPartition &domain) const override
resets the configuration to the default value.
CoupledSolver(const string &name, Group *const parent)
main constructor for CoupledSolver Objects
virtual void outputConfigurationStatistics(DomainPartition const &domain) const override
CoupledSolver(CoupledSolver &&)=default
default move constructor
virtual void resetConfigurationToBeginningOfStep(DomainPartition &domain) override
resets the configuration to the beginning of the time-step.
CoupledSolver(CoupledSolver const &)=delete
deleted copy constructor
virtual void synchronizeNonlinearSolverParameters() override
synchronize the nonlinear solver parameters.
virtual void setupCoupling(DomainPartition const &domain, DofManager &dofManager) const
Utility function to set the coupling between degrees of freedom.
virtual real64 sequentiallyCoupledSolverStep(real64 const &time_n, real64 const &dt, int const cycleNumber, DomainPartition &domain)
Sequentially coupled solver step. It solves a nonlinear system of equations using a sequential approa...
virtual bool checkSequentialSolutionIncrements(DomainPartition &domain) const override
Check if the solution increments are ok to use.
virtual void postInputInitialization() override
virtual bool updateConfiguration(DomainPartition &domain) override
updates the configuration (if needed) based on the state after a converged Newton loop.
std::array< string, sizeof...(SOLVERS) > m_names
Names of the single-physics solvers.
std::tuple< SOLVERS *... > m_solvers
Pointers of the single-physics solvers.
CoupledSolver & operator=(CoupledSolver &&)=delete
deleted move operator
void setSubSolvers()
Utility function to set the subsolvers pointers using the names provided by the user.
virtual void mapSolutionBetweenSolvers(DomainPartition &domain, integer const solverType)
Maps the solution obtained from one solver to the fields used by the other solver(s)
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...
SequentialConvergenceCriterion sequentialConvergenceCriterion() const
Getter for the sequential convergence criterion.
integer m_allowNonConverged
Flag to allow for a non-converged nonlinear solution and continue with the problem.
real64 m_newtonTol
The tolerance for the nonlinear convergence check.
NonlinearAccelerationType m_nonlinearAccelerationType
Type of nonlinear acceleration for sequential solver.
integer m_maxIterNewton
The maximum number of nonlinear iterations that are allowed.
real64 m_timeStepCutFactor
Factor by which the time step will be cut if a timestep cut is required.
integer m_numNewtonIterations
The number of nonlinear iterations that have been exectued.
integer m_numTimeStepAttempts
Number of times that the time-step had to be cut.
integer m_maxTimeStepCuts
Max number of time step cuts.
CouplingType couplingType() const
Getter for the coupling type.
LineSearchAction m_lineSearchAction
Flag to apply a line search.
@ ResidualNorm
convergence achieved when the residual drops below a given norm
@ NumberOfNonlinearIterations
convergence achieved when the subproblems convergence is achieved in less than minNewtonIteration
@ SolutionIncrements
convergence achieved when the solution increments are small enough
integer m_subcyclingOption
Flag to specify whether subcycling is allowed or not in sequential schemes.
Base class for all physics solvers.
virtual string getCatalogName() const =0
SolverStatistics m_solverStatistics
Solver statistics.
integer m_numTimestepsSinceLastDtCut
Number of cycles since last timestep cut.
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.
Timestamp getMeshModificationTimestamp(DomainPartition &domain) const
getter for the timestamp of the mesh modification on the mesh levels
NonlinearSolverParameters & getNonlinearSolverParameters()
accessor for the nonlinear solver parameters.
NonlinearSolverParameters m_nonlinearSolverParameters
Nonlinear solver parameters.
void logTimeStepCut()
Tell the solverStatistics that there is a time step cut.
void logNonlinearIteration(integer const numLinearIterations)
Tell the solverStatistics that we are doing a nonlinear iteration.
Wrapper< TBASE > & registerWrapper(string const &name, wrapperMap::KeyIndex::index_type *const rkey=nullptr)
Create and register a Wrapper around a new object.
DataContext const & getDataContext() const
Definition: Group.hpp:1345
string const & getName() const
Get group name.
Definition: Group.hpp:1331
Group & getParent()
Access the group's parent.
Definition: Group.hpp:1364
DataContext const & getWrapperDataContext(KEY key) const
Definition: Group.hpp:1356
#define GEOS_LOG_LEVEL_RANK_0(logInfoStruct, msg)
Output messages (only on rank 0) based on current Group's log level.
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override
function to perform setup for implicit timestep
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.
virtual real64 scalingForSystemSolution(DomainPartition &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localSolution) override
Function to determine if the solution vector should be scaled back in order to maintain a known const...
virtual void updateState(DomainPartition &domain) override
Recompute all dependent quantities from primary variables (including constitutive models)
virtual void implicitStepComplete(real64 const &time_n, real64 const &dt, DomainPartition &domain) override
perform cleanup for implicit timestep
virtual void applyBoundaryConditions(real64 const time_n, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) override
apply boundary condition to system
virtual real64 calculateResidualNorm(real64 const &time_n, real64 const &dt, DomainPartition const &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localRhs) override
calculate the norm of the global system residual
real64 solverStep(real64 const &time_n, real64 const &dt, int const cycleNumber, DomainPartition &domain) override final
virtual real64 setNextDt(real64 const &currentTime, real64 const &currentDt, DomainPartition &domain) override
function to set the next time step size
virtual real64 solverStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain)
entry function to perform a solver step
void setupDofs(DomainPartition const &domain, DofManager &dofManager) const override
Populate degree-of-freedom manager with fields relevant to this solver.
virtual void applySystemSolution(DofManager const &dofManager, arrayView1d< real64 const > const &localSolution, real64 const scalingFactor, real64 const dt, DomainPartition &domain) override
Function to apply the solution vector to the state.
virtual real64 setNextDt(real64 const &currentTime, real64 const &currentDt, DomainPartition &domain)
function to set the next time step size
virtual void assembleSystem(real64 const time_n, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) override
function to assemble the linear system matrix and rhs
virtual bool checkSystemSolution(DomainPartition &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localSolution, real64 const scalingFactor) override
Function to check system solution for physical consistency and constraint violation.
virtual void resetStateToBeginningOfStep(DomainPartition &domain) override
reset state of physics back to the beginning of the step.
@ FALSE
Not read from input.
@ REQUIRED
Required in input.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
Definition: DataTypes.hpp:179
unsigned long long int Timestamp
Timestamp type (used to perform actions such a sparsity pattern computation after mesh modifications)
Definition: DataTypes.hpp:126
std::string string
String type.
Definition: DataTypes.hpp:90
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
LvArray::CRSMatrixView< T, COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
Definition: DataTypes.hpp:309
int integer
Signed integer type.
Definition: DataTypes.hpp:81
Provides enum <-> string conversion facilities.
Exception class used to report errors in user input.
Definition: Logger.hpp:464
static constexpr char const * discretizationString()