From 803a69f38c0319d80fedb407d9c9d6b7e4220a46 Mon Sep 17 00:00:00 2001 From: Peter German Date: Mon, 2 Dec 2024 17:02:52 -0700 Subject: [PATCH 01/18] Add inlude for linearFVtimederivative. (#29149) --- .../linearfvkernels/LinearFVElementalKernel.h | 8 ++- .../linearfvkernels/LinearFVTimeDerivative.h | 57 +++++++++++++++++++ .../include/timeintegrators/TimeIntegrator.h | 6 ++ 3 files changed, 69 insertions(+), 2 deletions(-) create mode 100644 framework/include/linearfvkernels/LinearFVTimeDerivative.h diff --git a/framework/include/linearfvkernels/LinearFVElementalKernel.h b/framework/include/linearfvkernels/LinearFVElementalKernel.h index b9af62b54dd8..0c1e08f490df 100644 --- a/framework/include/linearfvkernels/LinearFVElementalKernel.h +++ b/framework/include/linearfvkernels/LinearFVElementalKernel.h @@ -35,10 +35,11 @@ class LinearFVElementalKernel : public LinearFVKernel * Set the current ElemInfo object. * @param elem_info The new ElemInfo which will be used as the current */ - void setCurrentElemInfo(const ElemInfo * elem_info) { _current_elem_info = elem_info; } + virtual void setCurrentElemInfo(const ElemInfo * elem_info); /** - * Set the coordinate system specific volume + * Set the coordinate system specific volume, the multiplication with + * the transformation factor is done outside of the kernel. * @param volume the new coordinate specific volume */ void setCurrentElemVolume(const Real volume) { _current_elem_volume = volume; } @@ -55,4 +56,7 @@ class LinearFVElementalKernel : public LinearFVKernel /// The coordinate-specific element volume Real _current_elem_volume; + + /// The dof index for the current variable associated with the element + dof_id_type _dof_id; }; diff --git a/framework/include/linearfvkernels/LinearFVTimeDerivative.h b/framework/include/linearfvkernels/LinearFVTimeDerivative.h new file mode 100644 index 000000000000..54ef24cb96e9 --- /dev/null +++ b/framework/include/linearfvkernels/LinearFVTimeDerivative.h @@ -0,0 +1,57 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "LinearFVElementalKernel.h" +#include "TimeIntegrator.h" + +/** + * Kernel that adds contributions from a time derivative term discretized using the + * finite volume method to a linear system. + */ +class LinearFVTimeDerivative : public LinearFVElementalKernel +{ +public: + static InputParameters validParams(); + + /** + * Class constructor. + * @param params The InputParameters for the kernel. + */ + LinearFVTimeDerivative(const InputParameters & params); + + virtual Real computeMatrixContribution() override; + + virtual Real computeRightHandSideContribution() override; + + /** + * Set the current ElemInfo object. + * @param elem_info The new ElemInfo which will be used as the current + */ + virtual void setCurrentElemInfo(const ElemInfo * elem_info) override; + +protected: + /// The functor for the material property multipler + const Moose::Functor & _factor; + + /// The time integrator to use in this kernel, will provide information + /// on how many states are required in the history. + const TimeIntegrator & _time_integrator; + +private: + /// Current values of the material property multiplier. Used + /// for caching. + std::vector _factor_history; + + /// State args, the args which will help us fetch the different states of + /// the material property multiplier. 0th is the current, 1st is old + /// 2nd is the older. Might not need all, depends on the time integrator. + std::vector _state_args; +}; diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index 3a4c4662005d..f6ad3a1e2561 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -133,6 +133,12 @@ class TimeIntegrator : public MooseObject, public Restartable */ virtual bool isExplicit() const { return false; } + /** + * Return the number of states this requires in a linear + * system setting + */ + virtual unsigned int numStatesRequired() const { return 1; } + /** * Returns whether mass matrix is lumped */ From 200608c9023bb116914efeb43e0e18d812fbfb4c Mon Sep 17 00:00:00 2001 From: Peter German Date: Mon, 2 Dec 2024 17:03:26 -0700 Subject: [PATCH 02/18] Add source for LinearFVTimeDerivative. (#29149) --- .../linearfvkernels/LinearFVElementalKernel.C | 13 ++-- .../linearfvkernels/LinearFVTimeDerivative.C | 59 +++++++++++++++++++ 2 files changed, 68 insertions(+), 4 deletions(-) create mode 100644 framework/src/linearfvkernels/LinearFVTimeDerivative.C diff --git a/framework/src/linearfvkernels/LinearFVElementalKernel.C b/framework/src/linearfvkernels/LinearFVElementalKernel.C index 8965066af9e3..993ed4a7502c 100644 --- a/framework/src/linearfvkernels/LinearFVElementalKernel.C +++ b/framework/src/linearfvkernels/LinearFVElementalKernel.C @@ -29,8 +29,7 @@ LinearFVElementalKernel::addMatrixContribution() { // These only contribute to the diagonal of the matrix, so we just get // the contribution and insert it immediately. - const auto dof_id = _current_elem_info->dofIndices()[_sys_num][_var_num]; - (*_linear_system.matrix).add(dof_id, dof_id, computeMatrixContribution()); + (*_linear_system.matrix).add(_dof_id, _dof_id, computeMatrixContribution()); } void @@ -38,6 +37,12 @@ LinearFVElementalKernel::addRightHandSideContribution() { // These only contribute to one entry of the right hand side, so we just get // the contribution and insert it immediately. - const auto dof_id = _current_elem_info->dofIndices()[_sys_num][_var_num]; - (*_linear_system.rhs).add(dof_id, computeRightHandSideContribution()); + (*_linear_system.rhs).add(_dof_id, computeRightHandSideContribution()); +} + +void +LinearFVElementalKernel::setCurrentElemInfo(const ElemInfo * elem_info) +{ + _current_elem_info = elem_info; + _dof_id = _current_elem_info->dofIndices()[_sys_num][_var_num]; } diff --git a/framework/src/linearfvkernels/LinearFVTimeDerivative.C b/framework/src/linearfvkernels/LinearFVTimeDerivative.C new file mode 100644 index 000000000000..52617ac2d5d2 --- /dev/null +++ b/framework/src/linearfvkernels/LinearFVTimeDerivative.C @@ -0,0 +1,59 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "LinearFVTimeDerivative.h" + +registerMooseObject("MooseApp", LinearFVTimeDerivative); + +InputParameters +LinearFVTimeDerivative::validParams() +{ + InputParameters params = LinearFVElementalKernel::validParams(); + params.addClassDescription("Represents the matrix and right hand side contributions of a " + "time derivative term in a partial differential equation."); + params.addParam( + "factor", 1.0, "A multiplier on the variable within the time derivative."); + return params; +} + +LinearFVTimeDerivative::LinearFVTimeDerivative(const InputParameters & params) + : LinearFVElementalKernel(params), + _factor(getFunctor("factor")), + _time_integrator(_sys.getTimeIntegrator(_var_num)), + _factor_history(_time_integrator.numStatesRequired(), 0.0), + _state_args(_time_integrator.numStatesRequired(), determineState()) +{ + // In case we need older states + for (const auto i : index_range(_state_args)) + _state_args[i] = Moose::StateArg(i, Moose::SolutionIterationType::Time); +} + +Real +LinearFVTimeDerivative::computeMatrixContribution() +{ + const auto elem_arg = makeElemArg(_current_elem_info->elem()); + return _factor(elem_arg, _state_args[0]) * /*_time_integrator.matrixContribution() */ + _current_elem_volume; +} + +Real +LinearFVTimeDerivative::computeRightHandSideContribution() +{ + return /*_time_integrator.rhsContribution(_dof_id, _factor_history) */ _current_elem_volume; +} + +void +LinearFVTimeDerivative::setCurrentElemInfo(const ElemInfo * elem_info) +{ + LinearFVElementalKernel::setCurrentElemInfo(elem_info); + + const auto elem_arg = makeElemArg(_current_elem_info->elem()); + for (const auto i : index_range(_factor_history)) + _factor_history[i] = _factor(elem_arg, _state_args[i]); +} From 23c1457865e9a8672de4d5fa403ffa2c8141cde1 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 12:34:03 -0700 Subject: [PATCH 03/18] Add nonlinear TI interface. (#29149) --- .../timeintegrators/NonlinearTimeIntegrator.h | 91 +++++++++++++++++++ .../include/timeintegrators/TimeIntegrator.h | 58 +++--------- .../timeintegrators/NonlinearTimeIntegrator.C | 37 ++++++++ .../src/timeintegrators/TimeIntegrator.C | 9 +- 4 files changed, 145 insertions(+), 50 deletions(-) create mode 100644 framework/include/timeintegrators/NonlinearTimeIntegrator.h create mode 100644 framework/src/timeintegrators/NonlinearTimeIntegrator.C diff --git a/framework/include/timeintegrators/NonlinearTimeIntegrator.h b/framework/include/timeintegrators/NonlinearTimeIntegrator.h new file mode 100644 index 000000000000..76be3f72da15 --- /dev/null +++ b/framework/include/timeintegrators/NonlinearTimeIntegrator.h @@ -0,0 +1,91 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + + +// MOOSE includes +class SystemBase; +class NonlinearSystemBase; + +namespace libMesh +{ +template +class NumericVector; +class NonlinearImplicitSystem; +} // namespace libMesh + +/** + * Interface class for routines and member variables for time integrators + * relying on Newton's method. + */ +class NonlinearTimeIntegrator +{ +public: + NonlinearTimeIntegrator(FEProblemBase & problem, SystemBase & system); + + /** + * Callback to the NonlinearTimeIntegrator called immediately after the + * residuals are computed in NonlinearSystem::computeResidual(). + * The residual vector which is passed in to this function should + * be filled in by the user with the _Re_time and _Re_non_time + * vectors in a way that makes sense for the particular + * TimeIntegration method. + */ + virtual void postResidual(NumericVector & /*residual*/) {} + + /** + * Returns the tag for the nodal multiplication factor for the residual calculation of the udot + * term. + * + * By default, this tag will be associated with udot. + */ + TagID uDotFactorTag() const { return _u_dot_factor_tag; } + /** + * Returns the tag for the nodal multiplication factor for the residual calculation of the udotdot + * term. + * + * By default, this tag will be associated with udotdot. + */ + TagID uDotDotFactorTag() const { return _u_dotdot_factor_tag; } + +protected: + + /// Wrapper around vector addition for nonlinear time integrators. + /// @param name The name of the vector + /// @param project If the vector should be projected + /// @param type PThe parallel distribution of the vetor + NumericVector * addVectorForNonlinearTI(const std::string & name, + const bool project, + const ParallelType type); + + /// Pointer to the nonlinear system + NonlinearSystemBase * _nl; + + /// Boolean to check if this integrator belongs to a nonlinear system + const bool _integrates_nl; + + /// Nonlinear implicit system, if applicable; otherwise, nullptr + /// NonlinearImplicitSystem * _nonlinear_implicit_system; + NonlinearImplicitSystem * _nonlinear_implicit_system; + + /// residual vector for time contributions + NumericVector * _Re_time; + + /// residual vector for non-time contributions + NumericVector * _Re_non_time; + + /// The vector tag for the nodal multiplication factor for the residual calculation of the udot term + const TagID _u_dot_factor_tag; + + /// The vector tag for the nodal multiplication factor for the residual calculation of the udotdot term + const TagID _u_dotdot_factor_tag; + + +}; diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index f6ad3a1e2561..899ff40f850b 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -12,6 +12,7 @@ // MOOSE includes #include "MooseObject.h" #include "Restartable.h" +#include "NonlinearTimeIntegrator.h" class FEProblemBase; class SystemBase; @@ -35,7 +36,7 @@ class NonlinearImplicitSystem; * used * only by NonlinearSystem (AuxiliarySystem does not produce residual). */ -class TimeIntegrator : public MooseObject, public Restartable +class TimeIntegrator : public MooseObject, public Restartable, public NonlinearTimeIntegrator { public: static InputParameters validParams(); @@ -69,16 +70,6 @@ class TimeIntegrator : public MooseObject, public Restartable */ virtual void solve(); - /** - * Callback to the TimeIntegrator called immediately after the - * residuals are computed in NonlinearSystem::computeResidual(). - * The residual vector which is passed in to this function should - * be filled in by the user with the _Re_time and _Re_non_time - * vectors in a way that makes sense for the particular - * TimeIntegration method. - */ - virtual void postResidual(NumericVector & /*residual*/) {} - /** * Callback to the TimeIntegrator called immediately after * TimeIntegrator::solve() (so the name does make sense!). See @@ -144,21 +135,6 @@ class TimeIntegrator : public MooseObject, public Restartable */ virtual const bool & isLumped() const { return _is_lumped; } - /** - * Returns the tag for the nodal multiplication factor for the residual calculation of the udot - * term. - * - * By default, this tag will be associated with udot. - */ - TagID uDotFactorTag() const { return _u_dot_factor_tag; } - /** - * Returns the tag for the nodal multiplication factor for the residual calculation of the udotdot - * term. - * - * By default, this tag will be associated with udotdot. - */ - TagID uDotDotFactorTag() const { return _u_dotdot_factor_tag; } - /** * @returns whether this integrator integrates the given variable */ @@ -202,31 +178,32 @@ class TimeIntegrator : public MooseObject, public Restartable */ void computeDuDotDu(); + /// Reference to the problem FEProblemBase & _fe_problem; - SystemBase & _sys; - NonlinearSystemBase & _nl; - - /// Nonlinear implicit system, if applicable; otherwise, nullptr - libMesh::NonlinearImplicitSystem * _nonlinear_implicit_system; - /// residual vector for time contributions - NumericVector & _Re_time; - /// residual vector for non-time contributions - NumericVector & _Re_non_time; + /// Reference to the system this time integrator operates on + SystemBase & _sys; /// Derivative of time derivative with respect to current solution: \f$ {du^dot}\over{du} \f$ for /// the different variables. We will only modify the elements in this vector corresponding to the /// variables that we integrate std::vector & _du_dot_du; - /// solution vectors + + /// @{ + /// Solution vectors for different states and variable restrictions const NumericVector * const & _solution; const NumericVector & _solution_old; std::unique_ptr> & _solution_sub; std::unique_ptr> & _solution_old_sub; - // + ///@} + + /// The current time step number int & _t_step; - // + + /// The current time step size Real & _dt; + + /// The previous time step size Real & _dt_old; /// Total number of nonlinear iterations over all stages of the time step @@ -237,11 +214,6 @@ class TimeIntegrator : public MooseObject, public Restartable /// Boolean flag that is set to true if lumped mass matrix is used bool _is_lumped; - /// The vector tag for the nodal multiplication factor for the residual calculation of the udot term - const TagID _u_dot_factor_tag; - /// The vector tag for the nodal multiplication factor for the residual calculation of the udotdot term - const TagID _u_dotdot_factor_tag; - /// Whether the user has requested that the time integrator be applied to a subset of variables bool & _var_restriction; diff --git a/framework/src/timeintegrators/NonlinearTimeIntegrator.C b/framework/src/timeintegrators/NonlinearTimeIntegrator.C new file mode 100644 index 000000000000..1124da334aba --- /dev/null +++ b/framework/src/timeintegrators/NonlinearTimeIntegrator.C @@ -0,0 +1,37 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "NonlinearTimeIntegrator.h" +#include "FEProblem.h" +#include "NonlinearSystemBase.h" + +#include "libmesh/nonlinear_implicit_system.h" +#include "libmesh/nonlinear_solver.h" +#include "libmesh/dof_map.h" + +NonlinearTimeIntegrator::NonlinearTimeIntegrator(FEProblemBase & problem, SystemBase & system) + : _nl(dynamic_cast(&system)), + _integrates_nl(_nl), + _nonlinear_implicit_system(_integrates_nl ? dynamic_cast(&_nl->system()) : nullptr), + _Re_time(_integrates_nl ? &_nl->getResidualTimeVector() : nullptr), + _Re_non_time(_integrates_nl ? &_nl->getResidualNonTimeVector() : nullptr), + _u_dot_factor_tag(problem.addVectorTag("u_dot_factor", Moose::VECTOR_TAG_SOLUTION)), + _u_dotdot_factor_tag(problem.addVectorTag("u_dotdot_factor", Moose::VECTOR_TAG_SOLUTION)) +{ +} + +NumericVector * NonlinearTimeIntegrator::addVectorForNonlinearTI(const std::string & name, + const bool project, + const ParallelType type) +{ + if (_integrates_nl) + return &_nl->addVector(name, project, type); + else + return nullptr; +} diff --git a/framework/src/timeintegrators/TimeIntegrator.C b/framework/src/timeintegrators/TimeIntegrator.C index 75ace47ad6ac..1645f5226230 100644 --- a/framework/src/timeintegrators/TimeIntegrator.C +++ b/framework/src/timeintegrators/TimeIntegrator.C @@ -30,13 +30,10 @@ TimeIntegrator::validParams() TimeIntegrator::TimeIntegrator(const InputParameters & parameters) : MooseObject(parameters), Restartable(this, "TimeIntegrators"), + NonlinearTimeIntegrator(*getCheckedPointerParam("_fe_problem_base"), + *getCheckedPointerParam("_sys")), _fe_problem(*getCheckedPointerParam("_fe_problem_base")), _sys(*getCheckedPointerParam("_sys")), - _nl(_fe_problem.getNonlinearSystemBase( - dynamic_cast(&_sys) ? _sys.number() : 0)), - _nonlinear_implicit_system(dynamic_cast(&_sys.system())), - _Re_time(_nl.getResidualTimeVector()), - _Re_non_time(_nl.getResidualNonTimeVector()), _du_dot_du(_sys.duDotDus()), _solution(_sys.currentSolution()), _solution_old(_sys.solutionState(1)), @@ -50,8 +47,6 @@ TimeIntegrator::TimeIntegrator(const InputParameters & parameters) _n_nonlinear_iterations(0), _n_linear_iterations(0), _is_lumped(false), - _u_dot_factor_tag(_fe_problem.addVectorTag("u_dot_factor", Moose::VECTOR_TAG_SOLUTION)), - _u_dotdot_factor_tag(_fe_problem.addVectorTag("u_dotdot_factor", Moose::VECTOR_TAG_SOLUTION)), _var_restriction(declareRestartableData( "var_restriction", !getParam>("variables").empty())), _local_indices(declareRestartableData>("local_indices")), From 164d8c069f064f30a3fc63105f777d3760fad444 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 13:08:06 -0700 Subject: [PATCH 04/18] Adapt to the new notation in other time integrators. (#29149) --- .../include/timeintegrators/CrankNicolson.h | 2 +- .../include/timeintegrators/ExplicitRK2.h | 2 +- .../timeintegrators/ExplicitSSPRungeKutta.h | 6 +- .../include/timeintegrators/ExplicitTVDRK2.h | 2 +- .../timeintegrators/ExplicitTimeIntegrator.h | 6 +- .../timeintegrators/ImplicitMidpoint.h | 2 +- .../include/timeintegrators/LStableDirk2.h | 4 +- .../timeintegrators/NonlinearTimeIntegrator.h | 5 +- framework/src/timeintegrators/AStableDirk4.C | 12 ++-- .../timeintegrators/ActuallyExplicitEuler.C | 18 +++--- framework/src/timeintegrators/BDF2.C | 4 +- .../src/timeintegrators/CentralDifference.C | 8 +-- framework/src/timeintegrators/CrankNicolson.C | 56 ++++++++++--------- framework/src/timeintegrators/ExplicitEuler.C | 4 +- framework/src/timeintegrators/ExplicitRK2.C | 22 ++++---- .../timeintegrators/ExplicitSSPRungeKutta.C | 38 ++++++------- .../src/timeintegrators/ExplicitTVDRK2.C | 20 +++---- .../timeintegrators/ExplicitTimeIntegrator.C | 26 ++++----- framework/src/timeintegrators/ImplicitEuler.C | 12 ++-- .../src/timeintegrators/ImplicitMidpoint.C | 20 +++---- framework/src/timeintegrators/LStableDirk2.C | 32 +++++------ framework/src/timeintegrators/LStableDirk3.C | 14 ++--- framework/src/timeintegrators/LStableDirk4.C | 14 ++--- framework/src/timeintegrators/NewmarkBeta.C | 4 +- .../timeintegrators/DirectCentralDifference.C | 30 +++++----- 25 files changed, 183 insertions(+), 180 deletions(-) diff --git a/framework/include/timeintegrators/CrankNicolson.h b/framework/include/timeintegrators/CrankNicolson.h index 90e403707191..1a4eda505670 100644 --- a/framework/include/timeintegrators/CrankNicolson.h +++ b/framework/include/timeintegrators/CrankNicolson.h @@ -45,7 +45,7 @@ class CrankNicolson : public TimeIntegrator virtual Real duDotDuCoeff() const override; - NumericVector & _residual_old; + NumericVector * _residual_old; }; template diff --git a/framework/include/timeintegrators/ExplicitRK2.h b/framework/include/timeintegrators/ExplicitRK2.h index d4ed2aaed04c..73cdc3eaf5aa 100644 --- a/framework/include/timeintegrators/ExplicitRK2.h +++ b/framework/include/timeintegrators/ExplicitRK2.h @@ -79,7 +79,7 @@ class ExplicitRK2 : public TimeIntegrator unsigned int _stage; /// Buffer to store non-time residual from the first stage. - NumericVector & _residual_old; + NumericVector * _residual_old; /// The older solution const NumericVector & _solution_older; diff --git a/framework/include/timeintegrators/ExplicitSSPRungeKutta.h b/framework/include/timeintegrators/ExplicitSSPRungeKutta.h index 5176f9130d90..312c5b292f3c 100644 --- a/framework/include/timeintegrators/ExplicitSSPRungeKutta.h +++ b/framework/include/timeintegrators/ExplicitSSPRungeKutta.h @@ -57,9 +57,9 @@ class ExplicitSSPRungeKutta : public ExplicitTimeIntegrator std::vector *> _solution_stage; /// Solution vector for intermediate stage - NumericVector & _solution_intermediate_stage; + NumericVector * _solution_intermediate_stage; /// Temporary solution vector - NumericVector & _tmp_solution; + NumericVector * _tmp_solution; /// Temporary mass-matrix/solution vector product - NumericVector & _tmp_mass_solution_product; + NumericVector * _tmp_mass_solution_product; }; diff --git a/framework/include/timeintegrators/ExplicitTVDRK2.h b/framework/include/timeintegrators/ExplicitTVDRK2.h index b8cd92be8e2a..1fda1fed88d2 100644 --- a/framework/include/timeintegrators/ExplicitTVDRK2.h +++ b/framework/include/timeintegrators/ExplicitTVDRK2.h @@ -70,7 +70,7 @@ class ExplicitTVDRK2 : public TimeIntegrator unsigned int _stage; /// Buffer to store non-time residual from the first stage. - NumericVector & _residual_old; + NumericVector * _residual_old; /// The older solution const NumericVector & _solution_older; diff --git a/framework/include/timeintegrators/ExplicitTimeIntegrator.h b/framework/include/timeintegrators/ExplicitTimeIntegrator.h index fd3775bb2cc3..c1eb653c62c9 100644 --- a/framework/include/timeintegrators/ExplicitTimeIntegrator.h +++ b/framework/include/timeintegrators/ExplicitTimeIntegrator.h @@ -66,13 +66,13 @@ class ExplicitTimeIntegrator : public TimeIntegrator, public MeshChangedInterfac MooseEnum _solve_type; /// Residual used for the RHS - NumericVector & _explicit_residual; + NumericVector * _explicit_residual; /// Solution vector for the linear solve - NumericVector & _solution_update; + NumericVector * _solution_update; /// Diagonal of the lumped mass matrix (and its inversion) - NumericVector & _mass_matrix_diag; + NumericVector * _mass_matrix_diag; /// Vector of 1's to help with creating the lumped mass matrix NumericVector * _ones; diff --git a/framework/include/timeintegrators/ImplicitMidpoint.h b/framework/include/timeintegrators/ImplicitMidpoint.h index bdba146909c3..6869d6afebe4 100644 --- a/framework/include/timeintegrators/ImplicitMidpoint.h +++ b/framework/include/timeintegrators/ImplicitMidpoint.h @@ -65,7 +65,7 @@ class ImplicitMidpoint : public TimeIntegrator unsigned int _stage; /// Buffer to store non-time residual from the first stage. - NumericVector & _residual_stage1; + NumericVector * _residual_stage1; }; template diff --git a/framework/include/timeintegrators/LStableDirk2.h b/framework/include/timeintegrators/LStableDirk2.h index cc268a50f0bc..e705dbbeb3a9 100644 --- a/framework/include/timeintegrators/LStableDirk2.h +++ b/framework/include/timeintegrators/LStableDirk2.h @@ -62,10 +62,10 @@ class LStableDirk2 : public TimeIntegrator unsigned int _stage; //! Buffer to store non-time residual from first stage solve. - NumericVector & _residual_stage1; + NumericVector * _residual_stage1; //! Buffer to store non-time residual from second stage solve - NumericVector & _residual_stage2; + NumericVector * _residual_stage2; // The parameter of the method, set at construction time and cannot be changed. const Real _alpha; diff --git a/framework/include/timeintegrators/NonlinearTimeIntegrator.h b/framework/include/timeintegrators/NonlinearTimeIntegrator.h index 76be3f72da15..6d0b52c5837c 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegrator.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegrator.h @@ -57,7 +57,8 @@ class NonlinearTimeIntegrator protected: - /// Wrapper around vector addition for nonlinear time integrators. + /// Wrapper around vector addition for nonlinear time integrators. If we don't + /// operate on a nonlinear system we don't need to add the vector. /// @param name The name of the vector /// @param project If the vector should be projected /// @param type PThe parallel distribution of the vetor @@ -65,7 +66,7 @@ class NonlinearTimeIntegrator const bool project, const ParallelType type); - /// Pointer to the nonlinear system + /// Pointer to the nonlinear system, can happen that we dont have any NonlinearSystemBase * _nl; /// Boolean to check if this integrator belongs to a nonlinear system diff --git a/framework/src/timeintegrators/AStableDirk4.C b/framework/src/timeintegrators/AStableDirk4.C index e233631eef00..de24fb33e263 100644 --- a/framework/src/timeintegrators/AStableDirk4.C +++ b/framework/src/timeintegrators/AStableDirk4.C @@ -42,7 +42,7 @@ AStableDirk4::AStableDirk4(const InputParameters & parameters) { std::ostringstream oss; oss << "residual_stage" << stage + 1; - _stage_residuals[stage] = &(_nl.addVector(oss.str(), false, GHOSTED)); + _stage_residuals[stage] = addVectorForNonlinearTI(oss.str(), false, GHOSTED); } // Initialize parameters @@ -148,14 +148,14 @@ AStableDirk4::solve() } // Do the solve - _nl.system().solve(); + _nl->system().solve(); // Update the iteration counts _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; } } @@ -188,10 +188,10 @@ AStableDirk4::postResidual(NumericVector & residual) // Store this stage's non-time residual. We are calling operator= // here, and that calls close(). - *_stage_residuals[_stage - 1] = _Re_non_time; + *_stage_residuals[_stage - 1] = *_Re_non_time; // Build up the residual for this stage. - residual.add(1., _Re_time); + residual.add(1., *_Re_time); for (unsigned int j = 0; j < _stage; ++j) residual.add(_a[_stage - 1][j], *_stage_residuals[j]); residual.close(); @@ -206,7 +206,7 @@ AStableDirk4::postResidual(NumericVector & residual) // just do one more stage, but I think handling it separately like // this is easier to understand and doesn't create too much code // repitition. - residual.add(1., _Re_time); + residual.add(1., *_Re_time); for (unsigned int j = 0; j < 3; ++j) residual.add(_b[j], *_stage_residuals[j]); residual.close(); diff --git a/framework/src/timeintegrators/ActuallyExplicitEuler.C b/framework/src/timeintegrators/ActuallyExplicitEuler.C index c6488e3abb61..861b90235ea2 100644 --- a/framework/src/timeintegrators/ActuallyExplicitEuler.C +++ b/framework/src/timeintegrators/ActuallyExplicitEuler.C @@ -76,12 +76,12 @@ ActuallyExplicitEuler::solve() _nonlinear_implicit_system->update(); // Compute the residual - _explicit_residual.zero(); + _explicit_residual->zero(); _fe_problem.computeResidual( - *_nonlinear_implicit_system->current_local_solution, _explicit_residual, _nl.number()); + *_nonlinear_implicit_system->current_local_solution, *_explicit_residual, _nl->number()); // Move the residual to the RHS - _explicit_residual *= -1.0; + *_explicit_residual *= -1.0; // Compute the mass matrix auto & mass_matrix = _nonlinear_implicit_system->get_system_matrix(); @@ -93,14 +93,14 @@ ActuallyExplicitEuler::solve() bool converged = performExplicitSolve(mass_matrix); // Update the solution - *_nonlinear_implicit_system->solution = _nl.solutionOld(); - *_nonlinear_implicit_system->solution += _solution_update; + *_nonlinear_implicit_system->solution = _nl->solutionOld(); + *_nonlinear_implicit_system->solution += *_solution_update; // Constraints may be solved in an uncoupled way. For example, momentum-balance equations may be // solved node-wise and then the solution (e.g. velocities or positions)can be applied to those // nodes without solving for such constraints on a system level. This strategy is being used for // node-face contact in explicit dynamics. - _nl.overwriteNodeFace(*_nonlinear_implicit_system->solution); + _nl->overwriteNodeFace(*_nonlinear_implicit_system->solution); // Enforce contraints on the solution DofMap & dof_map = _nonlinear_implicit_system->get_dof_map(); @@ -108,7 +108,7 @@ ActuallyExplicitEuler::solve() _nonlinear_implicit_system->solution.get()); _nonlinear_implicit_system->update(); - _nl.setSolution(*_nonlinear_implicit_system->current_local_solution); + _nl->setSolution(*_nonlinear_implicit_system->current_local_solution); _nonlinear_implicit_system->nonlinear_solver->converged = converged; } @@ -116,8 +116,8 @@ ActuallyExplicitEuler::solve() void ActuallyExplicitEuler::postResidual(NumericVector & residual) { - residual += _Re_time; - residual += _Re_non_time; + residual += *_Re_time; + residual += *_Re_non_time; residual.close(); // Reset time to the time at which to evaluate nodal BCs, which comes next diff --git a/framework/src/timeintegrators/BDF2.C b/framework/src/timeintegrators/BDF2.C index 87f5123df0a1..29487f3558f7 100644 --- a/framework/src/timeintegrators/BDF2.C +++ b/framework/src/timeintegrators/BDF2.C @@ -72,8 +72,8 @@ BDF2::computeADTimeDerivatives(ADReal & ad_u_dot, void BDF2::postResidual(NumericVector & residual) { - residual += _Re_time; - residual += _Re_non_time; + residual += *_Re_time; + residual += *_Re_non_time; residual.close(); } diff --git a/framework/src/timeintegrators/CentralDifference.C b/framework/src/timeintegrators/CentralDifference.C index b6b07393357b..d8d971bf8b05 100644 --- a/framework/src/timeintegrators/CentralDifference.C +++ b/framework/src/timeintegrators/CentralDifference.C @@ -61,10 +61,10 @@ CentralDifference::initialSetup() ActuallyExplicitEuler::initialSetup(); // _nl here so that we don't create this vector in the aux system time integrator - _nl.disassociateVectorFromTag(*_nl.solutionUDot(), _u_dot_factor_tag); - _nl.addVector(_u_dot_factor_tag, true, GHOSTED); - _nl.disassociateVectorFromTag(*_nl.solutionUDotDot(), _u_dotdot_factor_tag); - _nl.addVector(_u_dotdot_factor_tag, true, GHOSTED); + _nl->disassociateVectorFromTag(*_nl->solutionUDot(), _u_dot_factor_tag); + _nl->addVector(_u_dot_factor_tag, true, GHOSTED); + _nl->disassociateVectorFromTag(*_nl->solutionUDotDot(), _u_dotdot_factor_tag); + _nl->addVector(_u_dotdot_factor_tag, true, GHOSTED); } void diff --git a/framework/src/timeintegrators/CrankNicolson.C b/framework/src/timeintegrators/CrankNicolson.C index 103b4e97e239..89ba0e2464cb 100644 --- a/framework/src/timeintegrators/CrankNicolson.C +++ b/framework/src/timeintegrators/CrankNicolson.C @@ -23,7 +23,7 @@ CrankNicolson::validParams() CrankNicolson::CrankNicolson(const InputParameters & parameters) : TimeIntegrator(parameters), - _residual_old(_nl.addVector("residual_old", false, libMesh::GHOSTED)) + _residual_old(addVectorForNonlinearTI("residual_old", false, libMesh::GHOSTED)) { } @@ -74,22 +74,26 @@ CrankNicolson::init() mooseError("CrankNicolson: Time derivative of solution (`u_dot`) is not stored. Please set " "uDotRequested() to true in FEProblemBase befor requesting `u_dot`."); + // time derivative is assumed to be zero on initial NumericVector & u_dot = *_sys.solutionUDot(); u_dot.zero(); computeDuDotDu(); - // compute residual for the initial time step - // Note: we can not directly pass _residual_old in computeResidualTag because - // the function will call postResidual, which will cause _residual_old - // to be added on top of itself prohibited by PETSc. - // Objects executed on initial have been executed by FEProblem, - // so we can and should directly call NonlinearSystem residual evaluation. - _fe_problem.setCurrentResidualVectorTags({_nl.nonTimeVectorTag()}); - _nl.computeResidualTag(_nl.RHS(), _nl.nonTimeVectorTag()); - _fe_problem.clearCurrentResidualVectorTags(); - - copyVector(_nl.RHS(), _residual_old); + if (_integrates_nl) + { + // compute residual for the initial time step + // Note: we can not directly pass _residual_old in computeResidualTag because + // the function will call postResidual, which will cause _residual_old + // to be added on top of itself prohibited by PETSc. + // Objects executed on initial have been executed by FEProblem, + // so we can and should directly call NonlinearSystem residual evaluation. + _fe_problem.setCurrentResidualVectorTags({_nl->nonTimeVectorTag()}); + _nl->computeResidualTag(_nl->RHS(), _nl->nonTimeVectorTag()); + _fe_problem.clearCurrentResidualVectorTags(); + + copyVector(_nl->RHS(), *_residual_old); + } } void @@ -99,7 +103,7 @@ CrankNicolson::postResidual(NumericVector & residual) // and that's probably a good idea with earlier versions too, but // we don't always get here with _Re_time closed. std::vector inputs_closed = { - _Re_time.closed(), _Re_non_time.closed(), _residual_old.closed()}; + _Re_time->closed(), _Re_non_time->closed(), _residual_old->closed()}; // We might have done work on one processor but not all processors, // so we have to sync our closed() checks. Congrats to the BISON @@ -107,38 +111,38 @@ CrankNicolson::postResidual(NumericVector & residual) comm().min(inputs_closed); if (!inputs_closed[0]) - _Re_time.close(); + _Re_time->close(); if (!inputs_closed[1]) - _Re_non_time.close(); + _Re_non_time->close(); if (!inputs_closed[2]) - _residual_old.close(); + _residual_old->close(); if (!_var_restriction) { - residual += _Re_time; - residual += _Re_non_time; - residual += _residual_old; + residual += *_Re_time; + residual += *_Re_non_time; + residual += *_residual_old; } else { auto residual_sub = residual.get_subvector(_local_indices); - auto re_time_sub = _Re_time.get_subvector(_local_indices); - auto re_non_time_sub = _Re_non_time.get_subvector(_local_indices); - auto residual_old_sub = _residual_old.get_subvector(_local_indices); + auto re_time_sub = _Re_time->get_subvector(_local_indices); + auto re_non_time_sub = _Re_non_time->get_subvector(_local_indices); + auto residual_old_sub = _residual_old->get_subvector(_local_indices); *residual_sub += *re_time_sub; *residual_sub += *re_non_time_sub; *residual_sub += *residual_old_sub; residual.restore_subvector(std::move(residual_sub), _local_indices); - _Re_time.restore_subvector(std::move(re_time_sub), _local_indices); - _Re_non_time.restore_subvector(std::move(re_non_time_sub), _local_indices); - _residual_old.restore_subvector(std::move(residual_old_sub), _local_indices); + _Re_time->restore_subvector(std::move(re_time_sub), _local_indices); + _Re_non_time->restore_subvector(std::move(re_non_time_sub), _local_indices); + _residual_old->restore_subvector(std::move(residual_old_sub), _local_indices); } } void CrankNicolson::postStep() { - copyVector(_Re_non_time, _residual_old); + copyVector(*_Re_non_time, *_residual_old); } Real diff --git a/framework/src/timeintegrators/ExplicitEuler.C b/framework/src/timeintegrators/ExplicitEuler.C index 26f0cd76d98e..76ac3d149d59 100644 --- a/framework/src/timeintegrators/ExplicitEuler.C +++ b/framework/src/timeintegrators/ExplicitEuler.C @@ -57,7 +57,7 @@ ExplicitEuler::computeADTimeDerivatives(ADReal & ad_u_dot, void ExplicitEuler::postResidual(NumericVector & residual) { - residual += _Re_time; - residual += _Re_non_time; + residual += *_Re_time; + residual += *_Re_non_time; residual.close(); } diff --git a/framework/src/timeintegrators/ExplicitRK2.C b/framework/src/timeintegrators/ExplicitRK2.C index bb5ed9e1920e..7bd463d8deb8 100644 --- a/framework/src/timeintegrators/ExplicitRK2.C +++ b/framework/src/timeintegrators/ExplicitRK2.C @@ -25,7 +25,7 @@ ExplicitRK2::validParams() ExplicitRK2::ExplicitRK2(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_old(_nl.addVector("residual_old", false, GHOSTED)), + _residual_old(addVectorForNonlinearTI("residual_old", false, GHOSTED)), _solution_older(_sys.solutionState(2)) { mooseInfo("ExplicitRK2-derived TimeIntegrators (ExplicitMidpoint, Heun, Ralston) and other " @@ -86,12 +86,12 @@ ExplicitRK2::solve() _stage = 2; _fe_problem.timeOld() = time_old; _fe_problem.time() = time_stage2; - _nl.system().solve(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; // Advance solutions old->older, current->old. Also moves Material @@ -105,7 +105,7 @@ ExplicitRK2::solve() _stage = 3; _fe_problem.timeOld() = time_stage2; _fe_problem.time() = time_new; - _nl.system().solve(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); @@ -136,11 +136,11 @@ ExplicitRK2::postResidual(NumericVector & residual) // .) The minus signs are "baked in" to the non-time residuals, so // they do not appear here. // .) The current non-time residual is saved for the next stage. - _residual_old = _Re_non_time; - _residual_old.close(); + *_residual_old = *_Re_non_time; + _residual_old->close(); - residual.add(1., _Re_time); - residual.add(a(), _residual_old); + residual.add(1., *_Re_time); + residual.add(a(), *_residual_old); residual.close(); } else if (_stage == 3) @@ -156,9 +156,9 @@ ExplicitRK2::postResidual(NumericVector & residual) // residuals, so it does not appear here. // .) Although this is an update step, we have to do a "solve" // using the mass matrix. - residual.add(1., _Re_time); - residual.add(b1(), _residual_old); - residual.add(b2(), _Re_non_time); + residual.add(1., *_Re_time); + residual.add(b1(), *_residual_old); + residual.add(b2(), *_Re_non_time); residual.close(); } else diff --git a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C index c1da27e2c998..02abad54e86d 100644 --- a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C +++ b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C @@ -32,13 +32,11 @@ ExplicitSSPRungeKutta::validParams() ExplicitSSPRungeKutta::ExplicitSSPRungeKutta(const InputParameters & parameters) : ExplicitTimeIntegrator(parameters), - _order(getParam("order")), - _stage(0), - _solution_intermediate_stage(_nl.addVector("solution_intermediate_stage", false, GHOSTED)), - _tmp_solution(_nl.addVector("tmp_solution", false, GHOSTED)), - _tmp_mass_solution_product(_nl.addVector("tmp_mass_solution_product", false, GHOSTED)) + _solution_intermediate_stage(addVectorForNonlinearTI("solution_intermediate_stage", false, GHOSTED)), + _tmp_solution(addVectorForNonlinearTI("tmp_solution", false, GHOSTED)), + _tmp_mass_solution_product(addVectorForNonlinearTI("tmp_mass_solution_product", false, GHOSTED)) { // For SSPRK methods up to order 3, the number of stages equals the order _n_stages = _order; @@ -125,9 +123,9 @@ ExplicitSSPRungeKutta::solve() else { // Else must be the intermediate stage of the 3-stage method - _solution_intermediate_stage = *_solution; - _solution_intermediate_stage.close(); - _solution_stage[_stage] = &_solution_intermediate_stage; + *_solution_intermediate_stage = *_solution; + _solution_intermediate_stage->close(); + _solution_stage[_stage] = _solution_intermediate_stage; } // Set stage time for residual evaluation @@ -156,12 +154,12 @@ ExplicitSSPRungeKutta::solveStage() *_nonlinear_implicit_system->current_local_solution, mass_matrix, _Ke_time_tag); // Compute RHS vector using previous stage solution in steady-state residual - _explicit_residual.zero(); + _explicit_residual->zero(); _fe_problem.computeResidual( - *_nonlinear_implicit_system->current_local_solution, _explicit_residual, _nl.number()); + *_nonlinear_implicit_system->current_local_solution, *_explicit_residual, _nl->number()); // Move the residual to the RHS - _explicit_residual *= -1.0; + *_explicit_residual *= -1.0; // Perform the linear solve bool converged = performExplicitSolve(mass_matrix); @@ -169,7 +167,7 @@ ExplicitSSPRungeKutta::solveStage() // Update the solution: u^(s) = u^(s-1) + du^(s) (*_nonlinear_implicit_system->solution).zero(); *_nonlinear_implicit_system->solution = *(_solution_stage[_stage]); - *_nonlinear_implicit_system->solution += _solution_update; + *_nonlinear_implicit_system->solution += *_solution_update; // Enforce contraints on the solution DofMap & dof_map = _nonlinear_implicit_system->get_dof_map(); @@ -177,7 +175,7 @@ ExplicitSSPRungeKutta::solveStage() _nonlinear_implicit_system->solution.get()); _nonlinear_implicit_system->update(); - _nl.setSolution(*_nonlinear_implicit_system->current_local_solution); + _nl->setSolution(*_nonlinear_implicit_system->current_local_solution); _nonlinear_implicit_system->nonlinear_solver->converged = converged; @@ -188,21 +186,21 @@ void ExplicitSSPRungeKutta::postResidual(NumericVector & residual) { // The time residual is not included in the steady-state residual - residual += _Re_non_time; + residual += *_Re_non_time; // Compute \sum_{k=0}^{s-1} a_{s,k} u^(k) - u^(s-1) - _tmp_solution.zero(); + _tmp_solution->zero(); for (unsigned int k = 0; k <= _stage; k++) - _tmp_solution.add(_a[_stage][k], *(_solution_stage[k])); - _tmp_solution.add(-1.0, *(_solution_stage[_stage])); - _tmp_solution.close(); + _tmp_solution->add(_a[_stage][k], *(_solution_stage[k])); + _tmp_solution->add(-1.0, *(_solution_stage[_stage])); + _tmp_solution->close(); // Perform mass matrix product with the above vector auto & mass_matrix = _nonlinear_implicit_system->get_system_matrix(); - mass_matrix.vector_mult(_tmp_mass_solution_product, _tmp_solution); + mass_matrix.vector_mult(*_tmp_mass_solution_product, *_tmp_solution); // Finish computing residual vector (before modification by nodal BCs) - residual -= _tmp_mass_solution_product; + residual -= *_tmp_mass_solution_product; residual.close(); diff --git a/framework/src/timeintegrators/ExplicitTVDRK2.C b/framework/src/timeintegrators/ExplicitTVDRK2.C index 61d4da4a66bb..3fed86ba9141 100644 --- a/framework/src/timeintegrators/ExplicitTVDRK2.C +++ b/framework/src/timeintegrators/ExplicitTVDRK2.C @@ -26,7 +26,7 @@ ExplicitTVDRK2::validParams() ExplicitTVDRK2::ExplicitTVDRK2(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_old(_nl.addVector("residual_old", false, libMesh::GHOSTED)), + _residual_old(addVectorForNonlinearTI("residual_old", false, libMesh::GHOSTED)), _solution_older(_sys.solutionState(2)) { mooseInfo("ExplicitTVDRK2 and other multistage TimeIntegrators are known not to work with " @@ -88,12 +88,12 @@ ExplicitTVDRK2::solve() _stage = 2; _fe_problem.timeOld() = time_old; _fe_problem.time() = time_stage2; - _nl.system().solve(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; // Advance solutions old->older, current->old. Also moves Material @@ -107,7 +107,7 @@ ExplicitTVDRK2::solve() _stage = 3; _fe_problem.timeOld() = time_stage2; _fe_problem.time() = time_new; - _nl.system().solve(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); @@ -138,11 +138,11 @@ ExplicitTVDRK2::postResidual(NumericVector & residual) // .) The minus signs are "baked in" to the non-time residuals, so // they do not appear here. // .) The current non-time residual is saved for the next stage. - _residual_old = _Re_non_time; - _residual_old.close(); + *_residual_old = *_Re_non_time; + _residual_old->close(); - residual.add(1.0, _Re_time); - residual.add(1.0, _residual_old); + residual.add(1.0, *_Re_time); + residual.add(1.0, *_residual_old); residual.close(); } else if (_stage == 3) @@ -158,8 +158,8 @@ ExplicitTVDRK2::postResidual(NumericVector & residual) // residuals, so it does not appear here. // .) Although this is an update step, we have to do a "solve" // using the mass matrix. - residual.add(1.0, _Re_time); - residual.add(0.5, _Re_non_time); + residual.add(1.0, *_Re_time); + residual.add(0.5, *_Re_non_time); residual.close(); } else diff --git a/framework/src/timeintegrators/ExplicitTimeIntegrator.C b/framework/src/timeintegrators/ExplicitTimeIntegrator.C index 53a2a9844fa7..0af2a1e2674f 100644 --- a/framework/src/timeintegrators/ExplicitTimeIntegrator.C +++ b/framework/src/timeintegrators/ExplicitTimeIntegrator.C @@ -41,9 +41,9 @@ ExplicitTimeIntegrator::ExplicitTimeIntegrator(const InputParameters & parameter MeshChangedInterface(parameters), _solve_type(getParam("solve_type")), - _explicit_residual(_nl.addVector("explicit_residual", false, PARALLEL)), - _solution_update(_nl.addVector("solution_update", true, PARALLEL)), - _mass_matrix_diag(_nl.addVector("mass_matrix_diag", false, PARALLEL)) + _explicit_residual(addVectorForNonlinearTI("explicit_residual", false, PARALLEL)), + _solution_update(addVectorForNonlinearTI("solution_update", true, PARALLEL)), + _mass_matrix_diag(addVectorForNonlinearTI("mass_matrix_diag", false, PARALLEL)) { _Ke_time_tag = _fe_problem.getMatrixTagID("TIME"); @@ -52,7 +52,7 @@ ExplicitTimeIntegrator::ExplicitTimeIntegrator(const InputParameters & parameter _fe_problem.solverParams()._type = Moose::ST_LINEAR; if (_solve_type == LUMPED || _solve_type == LUMP_PRECONDITIONED) - _ones = &_nl.addVector("ones", false, PARALLEL); + _ones = addVectorForNonlinearTI("ones", false, PARALLEL); } void @@ -86,7 +86,7 @@ ExplicitTimeIntegrator::meshChanged() if (_solve_type == LUMP_PRECONDITIONED) { - _preconditioner = std::make_unique(_mass_matrix_diag); + _preconditioner = std::make_unique(*_mass_matrix_diag); _linear_solver->attach_preconditioner(_preconditioner.get()); _linear_solver->init(); } @@ -113,16 +113,16 @@ ExplicitTimeIntegrator::performExplicitSolve(SparseMatrix & mass_matrix) // Computes the sum of each row (lumping) // Note: This is actually how PETSc does it // It's not "perfectly optimal" - but it will be fast (and universal) - mass_matrix.vector_mult(_mass_matrix_diag, *_ones); + mass_matrix.vector_mult(*_mass_matrix_diag, *_ones); // "Invert" the diagonal mass matrix - _mass_matrix_diag.reciprocal(); + _mass_matrix_diag->reciprocal(); // Multiply the inversion by the RHS - _solution_update.pointwise_mult(_mass_matrix_diag, _explicit_residual); + _solution_update->pointwise_mult(*_mass_matrix_diag, *_explicit_residual); // Check for convergence by seeing if there is a nan or inf - auto sum = _solution_update.sum(); + auto sum = _solution_update->sum(); converged = std::isfinite(sum); // The linear iteration count remains zero @@ -132,8 +132,8 @@ ExplicitTimeIntegrator::performExplicitSolve(SparseMatrix & mass_matrix) } case LUMP_PRECONDITIONED: { - mass_matrix.vector_mult(_mass_matrix_diag, *_ones); - _mass_matrix_diag.reciprocal(); + mass_matrix.vector_mult(*_mass_matrix_diag, *_ones); + _mass_matrix_diag->reciprocal(); converged = solveLinearSystem(mass_matrix); @@ -153,8 +153,8 @@ ExplicitTimeIntegrator::solveLinearSystem(SparseMatrix & mass_matrix) const auto num_its_and_final_tol = _linear_solver->solve(mass_matrix, - _solution_update, - _explicit_residual, + *_solution_update, + *_explicit_residual, es.parameters.get("linear solver tolerance"), es.parameters.get("linear solver maximum iterations")); diff --git a/framework/src/timeintegrators/ImplicitEuler.C b/framework/src/timeintegrators/ImplicitEuler.C index 5a2a6a447781..0ca502e6efbd 100644 --- a/framework/src/timeintegrators/ImplicitEuler.C +++ b/framework/src/timeintegrators/ImplicitEuler.C @@ -65,19 +65,19 @@ ImplicitEuler::postResidual(NumericVector & residual) { if (!_var_restriction) { - residual += _Re_time; - residual += _Re_non_time; + residual += *_Re_time; + residual += *_Re_non_time; residual.close(); } else { auto residual_sub = residual.get_subvector(_local_indices); - auto re_time_sub = _Re_time.get_subvector(_local_indices); - auto re_non_time_sub = _Re_non_time.get_subvector(_local_indices); + auto re_time_sub = _Re_time->get_subvector(_local_indices); + auto re_non_time_sub = _Re_non_time->get_subvector(_local_indices); *residual_sub += *re_time_sub; *residual_sub += *re_non_time_sub; residual.restore_subvector(std::move(residual_sub), _local_indices); - _Re_time.restore_subvector(std::move(re_time_sub), _local_indices); - _Re_non_time.restore_subvector(std::move(re_non_time_sub), _local_indices); + _Re_time->restore_subvector(std::move(re_time_sub), _local_indices); + _Re_non_time->restore_subvector(std::move(re_non_time_sub), _local_indices); } } diff --git a/framework/src/timeintegrators/ImplicitMidpoint.C b/framework/src/timeintegrators/ImplicitMidpoint.C index 5bde2899fdec..01bebab105db 100644 --- a/framework/src/timeintegrators/ImplicitMidpoint.C +++ b/framework/src/timeintegrators/ImplicitMidpoint.C @@ -25,7 +25,7 @@ ImplicitMidpoint::validParams() ImplicitMidpoint::ImplicitMidpoint(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_stage1(_nl.addVector("residual_stage1", false, libMesh::GHOSTED)) + _residual_stage1(addVectorForNonlinearTI("residual_stage1", false, libMesh::GHOSTED)) { mooseInfo("ImplicitMidpoint and other multistage TimeIntegrators are known not to work with " "Materials/AuxKernels that accumulate 'state' and should be used with caution."); @@ -72,12 +72,12 @@ ImplicitMidpoint::solve() _console << "1st stage" << std::endl; _stage = 1; _fe_problem.time() = time_half; - _nl.system().solve(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; // Compute second stage @@ -85,7 +85,7 @@ ImplicitMidpoint::solve() _console << "2nd stage" << std::endl; _stage = 2; _fe_problem.time() = time_new; - _nl.system().solve(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); } @@ -103,11 +103,11 @@ ImplicitMidpoint::postResidual(NumericVector & residual) // .) M is the mass matrix // .) f(t_n + dt/2, Y_1) is saved in _residual_stage1 // .) The minus sign is baked in to the non-time residuals, so it does not appear here. - _residual_stage1 = _Re_non_time; - _residual_stage1.close(); + *_residual_stage1 = *_Re_non_time; + _residual_stage1->close(); - residual.add(1., _Re_time); - residual.add(0.5, _Re_non_time); + residual.add(1., *_Re_time); + residual.add(0.5, *_Re_non_time); residual.close(); } else if (_stage == 2) @@ -123,8 +123,8 @@ ImplicitMidpoint::postResidual(NumericVector & residual) // been saved as _residual_stage1. // .) The minus signs are "baked in" to the non-time residuals, so // they do not appear here. - residual.add(1., _Re_time); - residual.add(1., _residual_stage1); + residual.add(1., *_Re_time); + residual.add(1., *_residual_stage1); residual.close(); } else diff --git a/framework/src/timeintegrators/LStableDirk2.C b/framework/src/timeintegrators/LStableDirk2.C index f84f1e4aae89..37c39bc5a09d 100644 --- a/framework/src/timeintegrators/LStableDirk2.C +++ b/framework/src/timeintegrators/LStableDirk2.C @@ -28,8 +28,8 @@ LStableDirk2::validParams() LStableDirk2::LStableDirk2(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_stage1(_nl.addVector("residual_stage1", false, GHOSTED)), - _residual_stage2(_nl.addVector("residual_stage2", false, GHOSTED)), + _residual_stage1(addVectorForNonlinearTI("residual_stage1", false, GHOSTED)), + _residual_stage2(addVectorForNonlinearTI("residual_stage2", false, GHOSTED)), _alpha(1. - 0.5 * std::sqrt(2)) { mooseInfo("LStableDirk2 and other multistage TimeIntegrators are known not to work with " @@ -82,13 +82,13 @@ LStableDirk2::solve() _console << "1st stage" << std::endl; _stage = 1; _fe_problem.time() = time_stage1; - _nl.system().solve(); - _nl.destroyColoring(); + _nl->system().solve(); + _nl->destroyColoring(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; // Compute second stage @@ -97,8 +97,8 @@ LStableDirk2::solve() _stage = 2; _fe_problem.timeOld() = time_stage1; _fe_problem.time() = time_new; - _nl.potentiallySetupFiniteDifferencing(); - _nl.system().solve(); + _nl->potentiallySetupFiniteDifferencing(); + _nl->system().solve(); _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); @@ -122,11 +122,11 @@ LStableDirk2::postResidual(NumericVector & residual) // .) (Y_1 - y_n)/dt corresponds to the residual of the time kernels. // .) The minus sign in front of alpha is already "baked in" to // the non-time residuals, so it does not appear here. - _residual_stage1 = _Re_non_time; - _residual_stage1.close(); + *_residual_stage1 = *_Re_non_time; + _residual_stage1->close(); - residual.add(1., _Re_time); - residual.add(_alpha, _residual_stage1); + residual.add(1., *_Re_time); + residual.add(_alpha, *_residual_stage1); residual.close(); } else if (_stage == 2) @@ -143,12 +143,12 @@ LStableDirk2::postResidual(NumericVector & residual) // residuals, so they do not appear here. // // The solution at the end of stage 2, i.e. Y_2, is also the final solution. - _residual_stage2 = _Re_non_time; - _residual_stage2.close(); + *_residual_stage2 = *_Re_non_time; + _residual_stage2->close(); - residual.add(1., _Re_time); - residual.add(1. - _alpha, _residual_stage1); - residual.add(_alpha, _residual_stage2); + residual.add(1., *_Re_time); + residual.add(1. - _alpha, *_residual_stage1); + residual.add(_alpha, *_residual_stage2); residual.close(); } else diff --git a/framework/src/timeintegrators/LStableDirk3.C b/framework/src/timeintegrators/LStableDirk3.C index 4b6771778c6a..80b1769fabf4 100644 --- a/framework/src/timeintegrators/LStableDirk3.C +++ b/framework/src/timeintegrators/LStableDirk3.C @@ -38,7 +38,7 @@ LStableDirk3::LStableDirk3(const InputParameters & parameters) { std::ostringstream oss; oss << "residual_stage" << stage + 1; - _stage_residuals[stage] = &(_nl.addVector(oss.str(), false, GHOSTED)); + _stage_residuals[stage] = addVectorForNonlinearTI(oss.str(), false, GHOSTED); } // Initialize parameters @@ -108,20 +108,20 @@ LStableDirk3::solve() // If we previously used coloring, destroy the old object so it doesn't leak when we allocate a // new object in the following lines - _nl.destroyColoring(); + _nl->destroyColoring(); // Potentially setup finite differencing contexts for the solve - _nl.potentiallySetupFiniteDifferencing(); + _nl->potentiallySetupFiniteDifferencing(); // Do the solve - _nl.system().solve(); + _nl->system().solve(); // Update the iteration counts _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; } } @@ -147,10 +147,10 @@ LStableDirk3::postResidual(NumericVector & residual) // Store this stage's non-time residual. We are calling operator= // here, and that calls close(). - *_stage_residuals[_stage - 1] = _Re_non_time; + *_stage_residuals[_stage - 1] = *_Re_non_time; // Build up the residual for this stage. - residual.add(1., _Re_time); + residual.add(1., *_Re_time); for (unsigned int j = 0; j < _stage; ++j) residual.add(_a[_stage - 1][j], *_stage_residuals[j]); residual.close(); diff --git a/framework/src/timeintegrators/LStableDirk4.C b/framework/src/timeintegrators/LStableDirk4.C index 202a993f6ae3..dcfb23829797 100644 --- a/framework/src/timeintegrators/LStableDirk4.C +++ b/framework/src/timeintegrators/LStableDirk4.C @@ -46,7 +46,7 @@ LStableDirk4::LStableDirk4(const InputParameters & parameters) { std::ostringstream oss; oss << "residual_stage" << stage + 1; - _stage_residuals[stage] = &(_nl.addVector(oss.str(), false, GHOSTED)); + _stage_residuals[stage] = addVectorForNonlinearTI(oss.str(), false, GHOSTED); } } @@ -104,20 +104,20 @@ LStableDirk4::solve() // If we previously used coloring, destroy the old object so it doesn't leak when we allocate a // new object in the following lines - _nl.destroyColoring(); + _nl->destroyColoring(); // Potentially setup finite differencing contexts for the solve - _nl.potentiallySetupFiniteDifferencing(); + _nl->potentiallySetupFiniteDifferencing(); // Do the solve - _nl.system().solve(); + _nl->system().solve(); // Update the iteration counts _n_nonlinear_iterations += getNumNonlinearIterationsLastSolve(); _n_linear_iterations += getNumLinearIterationsLastSolve(); // Abort time step immediately on stage failure - see TimeIntegrator doc page - if (!_fe_problem.converged(_nl.number())) + if (!_fe_problem.converged(_nl->number())) return; } } @@ -146,10 +146,10 @@ LStableDirk4::postResidual(NumericVector & residual) // Store this stage's non-time residual. We are calling operator= // here, and that calls close(). - *_stage_residuals[_stage - 1] = _Re_non_time; + *_stage_residuals[_stage - 1] = *_Re_non_time; // Build up the residual for this stage. - residual.add(1., _Re_time); + residual.add(1., *_Re_time); for (unsigned int j = 0; j < _stage; ++j) residual.add(_a[_stage - 1][j], *_stage_residuals[j]); residual.close(); diff --git a/framework/src/timeintegrators/NewmarkBeta.C b/framework/src/timeintegrators/NewmarkBeta.C index 4ff7475fe254..d02b8302619f 100644 --- a/framework/src/timeintegrators/NewmarkBeta.C +++ b/framework/src/timeintegrators/NewmarkBeta.C @@ -114,8 +114,8 @@ NewmarkBeta::computeADTimeDerivatives(ADReal & ad_u_dot, void NewmarkBeta::postResidual(NumericVector & residual) { - residual += _Re_time; - residual += _Re_non_time; + residual += *_Re_time; + residual += *_Re_non_time; residual.close(); } diff --git a/modules/solid_mechanics/src/timeintegrators/DirectCentralDifference.C b/modules/solid_mechanics/src/timeintegrators/DirectCentralDifference.C index 228eed5398ed..40cbff69613f 100644 --- a/modules/solid_mechanics/src/timeintegrators/DirectCentralDifference.C +++ b/modules/solid_mechanics/src/timeintegrators/DirectCentralDifference.C @@ -90,35 +90,35 @@ DirectCentralDifference::solve() _nonlinear_implicit_system->update(); // Calculating the lumped mass matrix for use in residual calculation - mass_matrix.vector_mult(_mass_matrix_diag, *_ones); + mass_matrix.vector_mult(*_mass_matrix_diag, *_ones); // Compute the residual - _explicit_residual.zero(); + _explicit_residual->zero(); _fe_problem.computeResidual( - *_nonlinear_implicit_system->current_local_solution, _explicit_residual, _nl.number()); + *_nonlinear_implicit_system->current_local_solution, *_explicit_residual, _nl->number()); // Move the residual to the RHS - _explicit_residual *= -1.0; + *_explicit_residual *= -1.0; // Perform the linear solve bool converged = performExplicitSolve(mass_matrix); - _nl.overwriteNodeFace(*_nonlinear_implicit_system->solution); + _nl->overwriteNodeFace(*_nonlinear_implicit_system->solution); // Update the solution - *_nonlinear_implicit_system->solution = _nl.solutionOld(); - *_nonlinear_implicit_system->solution += _solution_update; + *_nonlinear_implicit_system->solution = _nl->solutionOld(); + *_nonlinear_implicit_system->solution += *_solution_update; _nonlinear_implicit_system->update(); - _nl.setSolution(*_nonlinear_implicit_system->current_local_solution); + _nl->setSolution(*_nonlinear_implicit_system->current_local_solution); _nonlinear_implicit_system->nonlinear_solver->converged = converged; } void DirectCentralDifference::postResidual(NumericVector & residual) { - residual += _Re_time; - residual += _Re_non_time; + residual += *_Re_time; + residual += *_Re_non_time; residual.close(); // Reset time to the time at which to evaluate nodal BCs, which comes next @@ -131,10 +131,10 @@ DirectCentralDifference::performExplicitSolve(SparseMatrix &) bool converged = false; // "Invert" the diagonal mass matrix - _mass_matrix_diag.reciprocal(); + _mass_matrix_diag->reciprocal(); // Calculate acceleration auto & accel = *_sys.solutionUDotDot(); - accel.pointwise_mult(_mass_matrix_diag, _explicit_residual); + accel.pointwise_mult(*_mass_matrix_diag, *_explicit_residual); // Scaling the acceleration auto accel_scaled = accel.clone(); @@ -147,11 +147,11 @@ DirectCentralDifference::performExplicitSolve(SparseMatrix &) vel = *old_vel->clone(); vel += *accel_scaled; - _solution_update = vel; - _solution_update.scale(_dt); + *_solution_update = vel; + _solution_update->scale(_dt); // Check for convergence by seeing if there is a nan or inf - auto sum = _solution_update.sum(); + auto sum = _solution_update->sum(); converged = std::isfinite(sum); // The linear iteration count remains zero From 0861fed8ea4c5c5145c49aeee50b3f210320606d Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 15:02:01 -0700 Subject: [PATCH 05/18] Add linear time integrator interface. (#29149) --- .../timeintegrators/LinearTimeIntegrator.h | 53 +++++++++++++++++++ .../timeintegrators/LinearTimeIntegrator.C | 32 +++++++++++ 2 files changed, 85 insertions(+) create mode 100644 framework/include/timeintegrators/LinearTimeIntegrator.h create mode 100644 framework/src/timeintegrators/LinearTimeIntegrator.C diff --git a/framework/include/timeintegrators/LinearTimeIntegrator.h b/framework/include/timeintegrators/LinearTimeIntegrator.h new file mode 100644 index 000000000000..900ce0955436 --- /dev/null +++ b/framework/include/timeintegrators/LinearTimeIntegrator.h @@ -0,0 +1,53 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +// MOOSE includes +class SystemBase; +class LinearSystem; + +namespace libMesh +{ +template +class NumericVector; +class LinearImplicitSystem; +} // namespace libMesh + +/** + * Interface class for routines and member variables for time integrators + * relying on Newton's method. + */ +class LinearTimeIntegrator +{ +public: + LinearTimeIntegrator(SystemBase & system); + + /// The time derivative's contribution to the right hand side of a linear system + /// @param dof_id The dof index at which this contribution should be fetched at + /// @param factors Multiplicative factor (e.g. a material property) at multiple + /// states (old, older, etc) + virtual Real timeDerivativeRHSContribution(dof_id_type dof_id, + const std::vector & factors = {}) const; + + /// The time derivative's contribution to the right hand side of a linear system. + /// For now, this does not depend of the DoF index, might change in the furutre. + virtual Real timeDerivativeMatrixContribution(const Real factor) const; + +protected: + + /// Pointer to the nonlinear system, can happen that we dont have any + LinearSystem * _linear_system; + + /// Boolean to check if it integrates a linear system + const bool _integrates_linear_system; + + /// Nonlinear implicit system, if applicable; otherwise, nullptr + LinearImplicitSystem * _linear_implicit_system; +}; diff --git a/framework/src/timeintegrators/LinearTimeIntegrator.C b/framework/src/timeintegrators/LinearTimeIntegrator.C new file mode 100644 index 000000000000..73b466178ca6 --- /dev/null +++ b/framework/src/timeintegrators/LinearTimeIntegrator.C @@ -0,0 +1,32 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "LinearTimeIntegrator.h" +#include "FEProblem.h" +#include "LinearSystem.h" + +#include "libmesh/linear_implicit_system.h" + +LinearTimeIntegrator::LinearTimeIntegrator(SystemBase & system) + : _linear_system(dynamic_cast(&system)), + _integrates_linear_system(_linear_system), + _linear_implicit_system(_linear_system ? dynamic_cast(&_linear_system->system()) : nullptr) +{ +} + +Real LinearTimeIntegrator::timeDerivativeRHSContribution(const dof_id_type /*dof_id*/, + const std::vector & /*factors*/) const +{ + mooseError("The time derivative right hand side contribution has not been implemented yet!"); +} + +Real LinearTimeIntegrator::timeDerivativeMatrixContribution(const Real /*factor*/) const +{ + mooseError("The time derivative matrix contribution has not been implemented yet!"); +} From 9b64c3dfc269e7e99671224292bafd47071a1462 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 15:02:58 -0700 Subject: [PATCH 06/18] Inherit from the linear interface and adopt the changes in ImplicitEuler. (#29149) --- framework/include/timeintegrators/ImplicitEuler.h | 4 ++++ .../timeintegrators/NonlinearTimeIntegrator.h | 1 - framework/include/timeintegrators/TimeIntegrator.h | 3 ++- .../src/linearfvkernels/LinearFVTimeDerivative.C | 4 ++-- framework/src/timeintegrators/ImplicitEuler.C | 14 ++++++++++++++ framework/src/timeintegrators/TimeIntegrator.C | 1 + 6 files changed, 23 insertions(+), 4 deletions(-) diff --git a/framework/include/timeintegrators/ImplicitEuler.h b/framework/include/timeintegrators/ImplicitEuler.h index b71317e5db1b..1b5d246ceea8 100644 --- a/framework/include/timeintegrators/ImplicitEuler.h +++ b/framework/include/timeintegrators/ImplicitEuler.h @@ -30,6 +30,10 @@ class ImplicitEuler : public TimeIntegrator virtual void postResidual(NumericVector & residual) override; virtual bool overridesSolve() const override { return false; } + virtual Real timeDerivativeRHSContribution(const dof_id_type dof_id, + const std::vector & factors) const override; + virtual Real timeDerivativeMatrixContribution(const Real factor) const override; + protected: /** * Helper function that actually does the math for computing the time derivative diff --git a/framework/include/timeintegrators/NonlinearTimeIntegrator.h b/framework/include/timeintegrators/NonlinearTimeIntegrator.h index 6d0b52c5837c..dfeab044b5c8 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegrator.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegrator.h @@ -73,7 +73,6 @@ class NonlinearTimeIntegrator const bool _integrates_nl; /// Nonlinear implicit system, if applicable; otherwise, nullptr - /// NonlinearImplicitSystem * _nonlinear_implicit_system; NonlinearImplicitSystem * _nonlinear_implicit_system; /// residual vector for time contributions diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index 899ff40f850b..0cf3f616adeb 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -13,6 +13,7 @@ #include "MooseObject.h" #include "Restartable.h" #include "NonlinearTimeIntegrator.h" +#include "LinearTimeIntegrator.h" class FEProblemBase; class SystemBase; @@ -36,7 +37,7 @@ class NonlinearImplicitSystem; * used * only by NonlinearSystem (AuxiliarySystem does not produce residual). */ -class TimeIntegrator : public MooseObject, public Restartable, public NonlinearTimeIntegrator +class TimeIntegrator : public MooseObject, public Restartable, public NonlinearTimeIntegrator, public LinearTimeIntegrator { public: static InputParameters validParams(); diff --git a/framework/src/linearfvkernels/LinearFVTimeDerivative.C b/framework/src/linearfvkernels/LinearFVTimeDerivative.C index 52617ac2d5d2..765889f00bd2 100644 --- a/framework/src/linearfvkernels/LinearFVTimeDerivative.C +++ b/framework/src/linearfvkernels/LinearFVTimeDerivative.C @@ -38,14 +38,14 @@ Real LinearFVTimeDerivative::computeMatrixContribution() { const auto elem_arg = makeElemArg(_current_elem_info->elem()); - return _factor(elem_arg, _state_args[0]) * /*_time_integrator.matrixContribution() */ + return _time_integrator.timeDerivativeMatrixContribution(_factor(elem_arg, _state_args[0])) * _current_elem_volume; } Real LinearFVTimeDerivative::computeRightHandSideContribution() { - return /*_time_integrator.rhsContribution(_dof_id, _factor_history) */ _current_elem_volume; + return _time_integrator.timeDerivativeRHSContribution(_dof_id, _factor_history) * _current_elem_volume; } void diff --git a/framework/src/timeintegrators/ImplicitEuler.C b/framework/src/timeintegrators/ImplicitEuler.C index 0ca502e6efbd..2d27949b3768 100644 --- a/framework/src/timeintegrators/ImplicitEuler.C +++ b/framework/src/timeintegrators/ImplicitEuler.C @@ -81,3 +81,17 @@ ImplicitEuler::postResidual(NumericVector & residual) _Re_non_time->restore_subvector(std::move(re_non_time_sub), _local_indices); } } + +Real +ImplicitEuler::timeDerivativeRHSContribution(dof_id_type dof_id, + const std::vector & factors) const +{ + mooseAssert(factors.size() == numStatesRequired(), "Either too many or too few states are given!"); + return factors[0] * _solution_old(dof_id) / _dt; +} + +Real +ImplicitEuler::timeDerivativeMatrixContribution(const Real factor) const +{ + return factor / _dt; +} diff --git a/framework/src/timeintegrators/TimeIntegrator.C b/framework/src/timeintegrators/TimeIntegrator.C index 1645f5226230..7010e4df386f 100644 --- a/framework/src/timeintegrators/TimeIntegrator.C +++ b/framework/src/timeintegrators/TimeIntegrator.C @@ -32,6 +32,7 @@ TimeIntegrator::TimeIntegrator(const InputParameters & parameters) Restartable(this, "TimeIntegrators"), NonlinearTimeIntegrator(*getCheckedPointerParam("_fe_problem_base"), *getCheckedPointerParam("_sys")), + LinearTimeIntegrator(*getCheckedPointerParam("_sys")), _fe_problem(*getCheckedPointerParam("_fe_problem_base")), _sys(*getCheckedPointerParam("_sys")), _du_dot_du(_sys.duDotDus()), From 273e7a0ba336b36ee1bde8438532d9fbab526d4e Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 15:21:16 -0700 Subject: [PATCH 07/18] Add test file for implicit euler. (#29149) --- .../implicit-euler/ie-linearfv.i | 91 +++++++++++++++++++ .../time_integrators/implicit-euler/tests | 6 +- 2 files changed, 94 insertions(+), 3 deletions(-) create mode 100644 test/tests/time_integrators/implicit-euler/ie-linearfv.i diff --git a/test/tests/time_integrators/implicit-euler/ie-linearfv.i b/test/tests/time_integrators/implicit-euler/ie-linearfv.i new file mode 100644 index 000000000000..5fc548763e6c --- /dev/null +++ b/test/tests/time_integrators/implicit-euler/ie-linearfv.i @@ -0,0 +1,91 @@ +########################################################### +# This is a simple test with a time-dependent problem +# demonstrating the use of the TimeIntegrator system. +# +# Testing a solution that is second order in space +# and first order in time +# +# @Requirement F1.30 +########################################################### + +[Mesh] + type = GeneratedMesh + dim = 2 + xmin = -1 + xmax = 1 + ymin = -1 + ymax = 1 + nx = 10 + ny = 10 +[] + +[Problem] + linear_sys_names = 'u_sys' +[] + +[Variables] + [u] + type = MooseLinearVariableFVReal + solver_sys = 'u_sys' + initial_condition = 0.0 + [] +[] + +[Functions] + [forcing_fn] + type = ParsedFunction + expression = ((x*x)+(y*y))-(4*t) + [] + [exact_fn] + type = ParsedFunction + expression = t*((x*x)+(y*y)) + [] +[] + +[LinearFVKernels] + [ie] + type = LinearFVTimeDerivative + variable = u + [] + [diff] + type = LinearFVDiffusion + variable = u + [] + [source] + type = LinearFVSource + variable = u + source_density = forcing_fn + [] +[] + +[LinearFVBCs] + [all] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + variable = u + boundary = '0 1 2 3' + functor = exact_fn + [] +[] + +[Postprocessors] + [l2_err] + type = ElementL2Error + variable = u + function = exact_fn + [] +[] + +[Executioner] + type = Transient + + # Test of the TimeIntegrator System + scheme = 'implicit-euler' + + start_time = 0.0 + num_steps = 5 + dt = 0.25 +[] + +[Outputs] + exodus = true +[] diff --git a/test/tests/time_integrators/implicit-euler/tests b/test/tests/time_integrators/implicit-euler/tests index 2ead2c8c068c..fb3cda87c8f1 100644 --- a/test/tests/time_integrators/implicit-euler/tests +++ b/test/tests/time_integrators/implicit-euler/tests @@ -3,7 +3,7 @@ issues = '#1953' [group] - requirement = "The system shall support the use of an implicit Euler solver with" + requirement = "The system shall support the use of an implicit Euler solver " [test] type = 'Exodiff' @@ -11,7 +11,7 @@ exodiff = 'ie_out.e' use_old_floor = true abs_zero = 1e-9 - detail = "with and" + detail = "without and" [] [adapt] @@ -19,7 +19,7 @@ input = 'ie_adapt.i' exodiff = 'ie_adapt_out.e-s005' group = 'adaptive' - detail = "without mesh adaptivity." + detail = "with mesh adaptivity." [] [] From de0309d2f616d3a277af53b405045537fedc8fb5 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 16:35:45 -0700 Subject: [PATCH 08/18] Add mods to enable advanceState in transients. (#29149) --- framework/src/problems/FEProblemBase.C | 6 +++--- test/tests/time_integrators/implicit-euler/ie-linearfv.i | 6 ++++++ 2 files changed, 9 insertions(+), 3 deletions(-) diff --git a/framework/src/problems/FEProblemBase.C b/framework/src/problems/FEProblemBase.C index 6bb47530c5b5..8e6d293c08db 100644 --- a/framework/src/problems/FEProblemBase.C +++ b/framework/src/problems/FEProblemBase.C @@ -6457,13 +6457,13 @@ FEProblemBase::advanceState() { TIME_SECTION("advanceState", 5, "Advancing State"); - for (auto & nl : _nl) - nl->copyOldSolutions(); + for (auto & sys : _solver_systems) + sys->copyOldSolutions(); _aux->copyOldSolutions(); if (_displaced_problem) { - for (const auto i : index_range(_nl)) + for (const auto i : index_range(_solver_systems)) _displaced_problem->solverSys(i).copyOldSolutions(); _displaced_problem->auxSys().copyOldSolutions(); } diff --git a/test/tests/time_integrators/implicit-euler/ie-linearfv.i b/test/tests/time_integrators/implicit-euler/ie-linearfv.i index 5fc548763e6c..9e2d36b86ab0 100644 --- a/test/tests/time_integrators/implicit-euler/ie-linearfv.i +++ b/test/tests/time_integrators/implicit-euler/ie-linearfv.i @@ -78,6 +78,12 @@ [Executioner] type = Transient + system_names = u_sys + l_tol = 1e-10 + + petsc_options_iname = '-pc_type -pc_hypre_type' + petsc_options_value = 'hypre boomeramg' + # Test of the TimeIntegrator System scheme = 'implicit-euler' From 1474984f043624a7a8588873e69c53ccbb844152 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 16:50:03 -0700 Subject: [PATCH 09/18] Add test and gold files. (#29149) --- .../implicit-euler/gold/ie-linearfv_out.e | Bin 0 -> 66400 bytes .../tests/time_integrators/implicit-euler/tests | 10 +++++++++- 2 files changed, 9 insertions(+), 1 deletion(-) create mode 100644 test/tests/time_integrators/implicit-euler/gold/ie-linearfv_out.e diff --git a/test/tests/time_integrators/implicit-euler/gold/ie-linearfv_out.e b/test/tests/time_integrators/implicit-euler/gold/ie-linearfv_out.e new file mode 100644 index 0000000000000000000000000000000000000000..d1b96f4eac95ea7567c06e865fdb3aa113079c23 GIT binary patch literal 66400 zcmeHw37i|nmG^KN+Z-V#vB3$XIc*LbAL9!c+GCGx48E|94KW}zQqPP)Bgs0(GY-Lo zE0E<5B*Y|~Sr&3FNeBTJFbNRM2LZzRgOCFV5Uzz_65xX@$NvBIs#`6orBU@LU+~)f zn>STms`~eHq)ykq(yyicX$PELOBgEN?6PSRndIB2NcN7|F8sz?Oi$xRo@IeDFSTQ@wH zh~v2>;X7-LM2axu6%MgKrEww|CyQsQABeLf!v!bDK<3dOHX$HvVnA544%@+fPD;Tw z?SOD!7~#sp2v;6PI2ERNp^PFu(JaC!6NT5AOzR*Orev8&b)+Hgm58ItD49s1jGB%C z?FiZ)v_mwLi1awwTq2#qZPRw3ty;JG<;8HfuM)mnM=~AFCsLh}OgfRuM^>e?amZ{n zu_;6HG#E6bBauXT$Pd2X`1`Ki(cVNibpU!3fdX~ryAbT9x8UCq;FFu7NYQ+OMeSn4 z+q7Likx%l0xTsIV63!TAFPiP>iKGkpvHCa4dvDNCrIV|=(mD}&Ah{RynLOj#I353n zIPJv;bxjqw)@lrJ*Ig z3Z63AM)~9)c;9Bdc!##ZmAGsXM1ED`?#$;??i4TegL>eWk;~iA$m~Qixh$G2h(a(9 zpU(?^*Zb)WxTj40Sg!o*#G=S}4O3WP!i+}&~;12e+ z-L>N0!nhijjWFY~t#B_t?r#z3;`S4F6w|x7Y{-dAIk}e)cNF<{ar=wAh;cQ)6f1Fe z0=btTm$l5r?Jw@Nj4S%!$#-WL6fSWOmitlb;kt2MIebcfQ~&*O-}c3I^R`QH+|iaV zuFLPItHhm$@SRw()Ih=Q@kB3&=`^lh)>-#QQ(m9$i6(O%=6I)_%8rY=Oj!@+jlS4k zS>oS#BwW#F^v98Wv>o?M!*K55-QKX?EqWxD5F};54RD{Fma9 zZx@reY|nH){?-LQ`qz$hB~Mt;0Q6(p_o0uHZJ? zuB-X5U03sAyZq>IkMDB6^z&Q}DxO>ZoGrL4z_nZCOEGN1a&Yz7b=`X6y5-^j?0dJ~ zHrjn!%5hzzUEQv@?OdZ=x1DOV>$Q8D?}purT#M~@qwd{yxzVoHf3H7j)V;1hZoAni zx6yjAso}M!^ZB$gO>JtNH-_JAqO&XC7^qDWC7;f0l6*Ic>7>)KITwG^mioV?{x_;F zZHhRXd;Q-O%WMnJ4V3_GRIcrcb~@VcXfLDPo;I@~AljX2Bc^?hwpZFoY5&{>M0+Le zeYCf^HZ;zSjsWcj+8s0!^f}PyL3@Dq1bqQyf%XFJ4cZ5^FK9o|{vaBZM}rOk9SAxI zbTDWPXe{Uu&^Qpym9)uD1T}*ufhL2dfTn_`fer;7209$n0y+X@gQkO8K{G%zL2aN} zpxK}~pd&$ZK}Ug(2F(M_2Q2_C1T6v`16mAP0$K_>7PJhs9CRG$c+eL?D?lfJP6V9< zihxcAMM3SL7$^>MKpmh?5N(wS&?%r(K}irt>nTtg#KGH2P!>cxCdV=bP!DJoi1t>_ zCT#`U8ng|FWlX|cVf2#^>KDW4Cm+-`hS5(xsCNvbpL|dU8Ad<( zus{r>pS0F>EY|e7Niod>o;}=qDex4x^uZOb0QHe)7=@ny&riZ3c*8(?H~HCa6vO$=fU)Mn8F* zt;6UiZ*z1Q{p9UP9Y#NSo2$dPYZ;N#p{p4+l4x^vEE!AQ4laFI{82#j9nGU0$d@R>t^pm&abQpOe zZ^whasQu(^g$|>iyq%!K=qGO{>M;7r+etc%e)1O4Vf2%?lR*rlpS(puCu=`>YX>om ze)1LrF^oKsw>Zere)86#!{{e(ojQzu^46up=qGOp9Y#O-I7Nrie=_J)5X0yvA4w3y z=qDfDAcoOT-clfjktgz&24%FLysgw>^pm%&4x^vE<#ZVR}dG1=<_54`^S|exUt9qd=oU2Y?O) z9RxZUGzK&lbO>l1Xgp{FXdaVM()J^KBtB0<>xq21FZ|xxJ4E1Cy z(AJ>-%7gLhA(kD>g)*l+DJ#l>c_jbjjnCyX>!Et{LH_F@=7G9Jy`fI9yjhkkCzb)* z6Fzg>I?0tE^_V(L{Utv2mO4v)rLOXMuCDF{KXr>bGXumjcH0=YuA~u`0j&gOK{-$! zQ~M?QETkj zL*%bJoz6Ltypzk1NoVrbn1tn9Lpe49cx#y~47KKBp-8po&CX&q4witPs$ylku`)EB zVmA`@l5)_GKGZ5+or$~^gDl`o9M?Rp-I*9~k6Y=KHDTO@<}u?Zj2XvXt@4_eXrJ%o zy4>6pAtyCYojAF4NzL(z<0o`Xm^yXhr1to<=JDgEOrF5DuH;%Y zYw40$C7fkV96MpcSa!39*OGa&7Bv)bW+9t}l2}Wz61$jRh;lh}z0HYPRh_!)5cEe;@IE?3uC_bM#U^2Z5 z+aB7}1@tNbm$?k1NkB$YPA3$j$BC?5<#0;bEASS}tWkv86mpo7$DWOJ>{Jh?$f`Mq z8W2kz(We!GT%WlNxfmBsY(yctFhv(N<=CYh!Pa?MRuv!@(i~7Ng|3cj4#*|mE;|8k zi#4b9*dn#KtM~QoT-G@4iapoviAP(e>4hS`L^= z#G=^56$(pbL+fEQ-V;s5f*J|)29}D=w~1IJ7mH$VmtU(POmo0kmlHcR5-s%li=yZ^ z5TDL6UPdrH#1=Q<59M8wb#p!BUx!tZ6{VewcDKi)5sbof$UALQ z{a!NqDhEvF)2=dxkfHR;(e_-V2;-OGcxTL{y4A4B>*zBw*eW7vUuX`q4#u2bMvJ&8 z6lcshnmg;eEGwQs3TGT3c`ZDd6#@tLblJdxvzR7ow)I_Rwn7=@Mc2e;OUnI(r^u2kOVfw>1K=diWc{k_K--`rPXeC3Vx?BR7>`Z49 z`L6C_{#z_j(0UceviyJ@>1+(U7EwLC;k$X&^8>OYO%&fz$mhUh2#stkU_cwD3dN}T zN_0YI4tL5WsNuN?1c^>CpLTmyvnJNEl92#q@t3k`7m&n0)QX}g0UU>`Uae^&he0jB zhp1#Y-UzVRAA2Z1Pm87P4R!Je4vT=1Of+^XY$!fF=wSw@7{uy?V*^VG_IBV{fI=n? zgSH=x{uZlw!nEdTlc!9WR#&?Xz>W^IUsYfDTC7!>y4gcS0Ky)Vh~)ZUSfQE48^CH- z&FgYCzz9ZBIJuzH;b^RW<+iMp(!*3b#nQoU&J4z!nY1_VgxVNrA+wC(e?FSc7cy|2 zjMr@G9D}#|b&#?Kpwzv1q9;(psn`b@wb&Yv>;)l)=CE4*Q&rED!)n>3995Q$>N2vd zx)Q$aR1^sTLlg>YUq>Mo3&10XhRqi^tkod4SOx-g)rMV3LrXcJZ*uio6^hvyyIAS` zmDdMoY9;140WXcK@fwv*7L+YJp6KY1;dOf-#!%#vwW_WQRS(A{bY>IrpduDMtP9AE zB!Gda>-2PWl6W|a#;Dk*J6!Q7*=~a-3_x# zcV{q=3k^UqtW+B@nYV`;Fkm#Sa&F#yX~1eNj!7+Hw!JJ@YampXjYhCQiuNiK7Yu2! zeA+jc!aSREr~qM7@y5hTblyMn$pyS|J9)M#!IycVs8{$%PcCvw)H8o!pa_jqPY4M} z`vMjc5H#JMh^ElxcRKklCmW)10EZmg$b?yYG*?pbSSY5`a9O8-k(du=A9*SjGIp*J35QGs#3Okssp}lGP_! zi>*!W)I=uY^wLz!*;_Fq6*0|am`gdi?hc$0LQw{la{XPXG!9}9ulimQLWYn~o`;bY zv+4Jl>md$q5<_wD7}mV`0$?p;Rfs}S6|bsoDGZ$kOsC=_`QRY!Dr87+ERtO8$E-I~>}Ygx=64--tJ;}Mwt%7-`RitPtf zA}oa*5MU{kWIG(7`%Y0-AcRzfVq}BObCxaQR={FuXrH;-C_oIAy47r1E{_xJVoolX z)++$a8wWpNMF$Y3*AsJ%m~NtZWiVDEZFf3p-!{)O8){`4iDJT{azqjiki+IB^B1BY zkXZ^}(U=3mHq(xkhdzuFuPNXru!>+dOatnjry*+y`fZHNBUsjov4=RloGI)`HG)Z& zXw9W$Y9W~`msTYm_b%sg;>!pDIXhZ)DVcfm)z4Bg^Tq*ywH79rzo(#AscL|QJX)+p z(QM9%yMptHE6__B)|}=F!=clQ#v#Y!V+mh&wOIKBN0073u=jw2M)w`qcR+vZ@t8c< z*WHfwg26p=m{Qh&YVx+C;|l}GqxQ-&TD%PjSUIokLWTnjTH~e>%>h?n!FVHDg+Q%r zK%i>AXc56{KtmcZN-BKwRpi9W8v(a8Uz9@e+7z0vlE!Py=;i}!E-7Z-wlrVG5~k0* zwXrOWtQ*uMQ^|M&r z=5%zF1Da05ZJxcJT0!iNA23hKEzPds%GB)466owo=df8L+84B?Fifxa$?ryU4mOZi zyTZtj^A*^?kV!^kn54nJ+b$n^SU@Z%!~UzTv|0;=8iUR7GNl=&P}C{{^`Wa2ceP?? zm+$2?42PT9P^4z-8|FvVFDINv^~%zXHl|-r@T{6G{j#qLd){KdycARAZI~(VQJcl9 zX#f*`^84k@@#PWIgw@_1Iq6@mKm)^B%3HcoN}+g#=~Mn8-Kc!taZ9fksJZ{Y168k= zOBb*Ut}#6y6|OP8o^aJfz%6g1ltS?;_IgWH(eb~O{SY;C#8=y z2wHR>s-tiC%cC?rSh9-4Ir8N*0xed+M0}0<$Evrb8>JMA*A_c=4TX`PR`6`W#`MSn z*iWy_DRns`sLeH-EXt;_v%6bN#bl?H>3(`;0-A2ta4j}+5@k_ew9mZ7ZuqHgHES}n zZ2F#(NX3!`d^Z`EU2LPl_Bz!gcpDGQ2ko-Gt_-FqN^?MUSXv}e%|2&dC<0i%tIlIF z1RzhG!`y8ZCe8KGGv)>}X%2|-Yc8KHU?HqJJuR7lfma<2sJ@ZRTc}R`JB0PiEi9Xw z(}d$Q9zv4l94fUd0faJk$cc;!&~#U3%|)wWYD!L>MON)#JR|wBA;X#~ZflC;z&*9A zDeGV(TwWb@7KwJ?^eZe&@4&i*x`8#IhvL;?SF=q@SaZl-NNO*fi&5g5mvYU8zbuqF z4$h?G(Z!YReQE2j?l7`RCAbm7l~oCQxG~Q%8(5g_aTsslWeymP$5o}{R?2{V*5<8^ zA247!aY(>IdmP7-@u-*bO+Et!tl%*OpgM;WNBqV29u4p8i&dU*V z666=3D|NQ{8-LXI1i8=NW+B^DyO)cZBK;hA5PplOJJ5 z4xLe)bl=FJQPq?78*chL?T z*pZq8@*sz#6YZ$rP54ADsekwi<-tKE*rT3Hlnz^cc{Z0}s(Ncp+z=dAZ>>WDR&T8} zpjWNx@!(7h9Nmyi@<__S4P7l(+pO8G$Ih$wkZ}blL*SA_63Cpo%ugf>)@K(Fh0<)O zHUg;kx8j`&W7SSyIkiBquoJW%@u@hLm0gXM9BgBwy+PX8V)6Sy^@^dCJr^c}GZ~&F9CWBB8%Fcdmk*E!A!4INSp#eO z@)%1et1bsp92&K+O=`9RdxlEJ)qWPuyMUUD5-3B{=?@K9-XFVMnione*W)G*a?Jtb z>5lGbl{3}x1am-HE$!g^=^zFHIBM>Kxl8A^&a;-zosZAu%vuzgwV7 zY7!)erfQ=;-BWo+mCLUAP^7cF!G+S9DrA!B()#AI7phRY67j0lL;;#AcEhLPPRUhK z7|FM}9%3Yt=*C%`>Ug9?3g&-oEG!Dq9I)Z9w|ZVL(V-4Z<7vcsn2SA4H`ha)D%Ano zMKt4I8b6RM7Xk65iOMf&nOD6sKrG|$1~Na?r!D{*o}63`s7kwtXnr_q2F7stt_BWv zMvlE}TCXO?6_55A%9yx17xtTL3Fh*7b62 zmTA}=RF;c?Ww|zkQuyzjLyi^@n;Co4V zU_VD+g(>A=OR1N4EP}#P4ZJ3lvr^?AuND0KlUma5f%rY4Y$ol4rI{XZC3i}6qR%HEbEs?w6Q#c)Gxk}fDNBQ;>$7{n+o0}z+8s?e#gi3 zec3Hq{Q*f7wp&E12{%_87cikdR$Gzr0RmRO@r96~+G@+M_r%MyOLoI^n860z;-Q`% zt6|HuZ`rtGbsPzc37eh-&cJm`+2q?wP4ej?D6|Jy z);^jy4&@w{ofmBm$Qkf#dR6d%SgB3DXTB?y#HWj}+A9aEgDkSFFg=)eqbX+OJP@y^j`o<>D`B02-EZ0HD91%?l*}&|lEzfH1Qp6IfNG4m-v= z#IlKS)Z}qvr&<%LS|H1$iSMUqK)0zb8!(x7QrZy2k+iC7W&h08MguB4h&t+0pDV0^ zL&*?dcu}7l^?#wVI5Zhs|Y2*A=|92RtX8Sl0Q*r!AA14yoM;`hIo4tbr4J5`m|To2!lDu*%R2TtbEzebwJ;m+cz$m#_HE zWhn0_rJY`+q*XM+#Oa5P{pNyY6>VREeE?)=)O`60yfI#NT7pb5H)eGQi?XKC)ZAp+ z&uS6%iG+n}g^DZ`k?d5J4>{Be6Dm#}Rnfd2)bNSLLT8B`To*$)&`O)IQ8qGv5if@jU44A@ji2cD-#YC!*u3Y002C!%>MKh?(h!ljr zLIKwL-n)4>$^#W^7@orvivNAv&8g>x0@UWTqAY!DGa3WvzB|MA5=w@e!$3>+fB{QQ zy3!t+86=A)NwH-6LFX7){KZCopC*>FVg2a;v)9nW7d0|5-Q zwGTvwfwuMmaaepGDB9Ww!eO8#+qCXcH{`$^-867-D~yv7F*t=&>}qB&rHEl<=*>ur zDwqBb16XtQocY3uO4GRC1j{lpw82&amUIbWvRu>F9tN=TE1^`&vExU+e`jp_|zeWtj zA>JoeWfS?5*NGyB=6YDu*4|=GEPI=)eM?3lJhkG+WvavtqlcxXQEIaZRrO@EUe$A! z%w?ztNdd1>$b6Vm&N(X!j{j@Cej$dDVa-LCVOkq8-Urs?aAO;F=s>l1+)y&CH5-7s zY-CKM{A(KTLlsJ`k4}dItYy3pB}316?^i(W52ZF9V8CoZm@F&9bG4PQVOs2(mTUkT z+OhkwWQQpf?k`5y$>Zs;?xbvh6Gny&HH}j_YFe_xR6W(>MY`f0lLGb7w8h^1J`A*| zWHz`l=KI=10mEO_n#K~qG!yOf+|D%{ormV_>qv>IaM+vH7gtc9DSUs0n(4t}yLk z%>~q9`V`mPSN0w%{I-H#7?$!!CRhPuk74M5fGEhG4v{EB{*AN zb`1a!o53pa>Mq>?4zo@oR>mG~LvdKT(yJ8o6Y;@5c>%Ao<}!?8^MyOo>E!YL9my3s zIuib(4mF(nM|`@d^sIJYxr8a@nrGIA1yu9%nKag}Vn-pWjQ2Jk2GE_~@ny8=HUEm@ zW;nKyv2bIC4XUHBIc{KUBZoyV@}Ci{Z`>x!@Aflp2=l(&@J;0nKe4uAJzr6&>BFqq z{Dljb%(ArH`dd?sHA{iVe0iuX*6dj|%D`V9%CGhXWO;CM6;93DLu@oqfH;pK*r&jKE!Zh25&F0t^ z3Fvwj4uV7X1aBW8rmmgPfa>v)LOO!6S4sCIVl1b80$v?|iomH$_{rS8ZdTe=AU2~Zl z0)YJ%RAZCH7eV_qoLbE{=&!XLBj9WdJ- z&1TiMI2?l##hIaUtctk|Gg+r6fyJ9>TD&$J-O0)VPW!|cM$7>*oQvY{)=u@^0v>P> zb299fa-1YpJy^~J700*AAq48H3+ki_B|7h4y}g}V_k8EUHtK*`HSbv71ajcXb1NuA z)0(1=Z}qKDc157~Z3@0CE?SR|XXx^DN`w;8YH_J{%qp)IykZ zNFz(|A`Hc?Ifr_0g||DbtPe2f&;`VyEJZ_UQIyzd5KZRlJ@=$(Q#fsPSrlC!<_L2Q zvx~U3xU;-6yN6g6!$W(di|VZWZpPnWm>k%MqY9WL(tsCOX@iqLQ0%z zxPZs)*Lh!ZnfH4A7^_;W`DwJi^VBgCt_IFtwiqv#l}d;Za^n_?`f!$7w#$VIdM7E> zW+@H`4NM&n7!24`T;&%g%NQg@tLz#L)9dwPxe3EzEZv<+r|`yC)}2rsAw~{(r%zz$ zrBVm{3sX~SnlZxEl#W`*t%d|LZ3_MaOnrr8UMQR0d`duj zu_b7W1pXIEzZ?!kv^!N+6r{1$4f8x(4zVxHG{m6Oa`{rMMyZzAn z3+|rx(HL=`{HcWhe5}N~^D0UAtt$lQyTX2^4w z|4eYsy3f9`FzU_+mXEiuJoC&O|9aCh<^JW-SrUF*MB<&Xlcc+PqTt-OlgMS=PgH(x z`e=mWIES;f0{$Sx|~&C6dDoGX4I?zvMX z{I7Ea@8MTey8l|=qHw0(HB-5*Bkoe+FI{q$#CvF-@bUZgg7f)r3!k^`CwzY}Pvo=u zT$Qf}KU}NQzuEh|ayPDdNy0DxP~siFT+(%ZU2raH7dd|SDUs`0yNjH!|DD8t?PSHr z4G-=l?jDy)_}i~byw~24bSJD9oVznp&gZp>zP!6se!pOjJon`jRerudNV`KWm+)KO zl6a>tlyp}dA~w{`K8UZ@JBW?(uuR*?a5R_8-pr;Wa-z zJZ?Yz_8t39zDc{6|4PCa4U>3xT`%d@_t^^P)(_?Ww(IR@N4(jx=gDWFf!(55idb6avX|3Qa{fo%woMn=qA>UAZ zEPQ>PO8@9_!|g9Wyz9tO%m1L6gDT5!hiEPO6JMfm>X`;w0h_o{pyI(A!? z{)z3c5%=tORrs`5S7K;-!2^Q1i9ds4~y)&0&;@oyV>yyD}S#q)%p z_Yc(JTT8r|QXXA@+E#Gx?vip|w@AwU%I5_4wAPZOLU{mEW)`ysdd zedz1%F#_dMaSG4`ud*8KPM=F#?_PyYD#KmOU#_6vnOZn$X4 zbK)MghlHgTy=HbxHT=O@edWMVg;W)UP^zV>gO%DrmIJQe=&_IVZWq@|Bby4hb7oMjt@{~Jz_=N|s9;C%gl``+&y zbwS@tlk6YAa>ozOI(L3&7ItzSJ}@K&slcxPTMe7y9y;EZ@rm)CB>_tmRKJ`vz5Nw^1A%q_L@BRs@FsgU;9A2S4;Sp50`lMcO~674i=o*uZbQ$ z@@0Lm`E5H-_`mNKk;};$?fy{0?|)19X}d+@A}1$?Txy=ynkQk z2g6>!+J5_AFTQZkp`G@|#|FLpPULO-^);`q9Xa<}`_+xVIe+yT-x9ayNeN$ar^Jg- zlyt+_*a~N(E%!S`>^FxEo&WMTcd=jYxb+LqfAuH!^K;jJ@6{8}v!5+|+;!Qu;-u=REl5YLgg0r$+_@8pKJa^GzjdPFvt6{T#bHT97?4RHIeEat2 z|3Bppx@2bwpRiuyO85xX}Y7{LddrynAP>bla`}oWglG_jl#~{i3+wMQ)RL zcU&gvu3RrTZQs}Bb%yZ$m1UBTv)4$zp1WVAAN1#Ql)L>t(H87#oUTd2|8ZkQE*CD; z?)egad_u>YuKC$oe(SLn9}dd=%-%TZ&sVR1bfLHpbxU~RK8bh4U`coOo`N&yk8;2I1N(2MoH%BO zd-ky3XnAYgzy0_L`=zHR-yb>Wc>9GLAHMFlD~}Pk=O76`|89x*Ut=X5#@~vcbMkKev_no9rKk6*MeUNz>J zv$*$sQ^MOZo@V}bUaZm$e`tZidF_;!l>3)+X9!+&n8e%u`;zXk8wBT;<+{9X6uxgB zFY>XbsC?y?+^*7(czF-y?zZnAB>a#EB;L@yCEc3q1SkEN$Z_gLBG)(8ikv?;T;jL= zR`|I76>$siN%%~RPf%`a_Wqot`_lxiKS?R)HQ$hOf8<^H{jHch_omU3pNsYp_rfS!#dd Mjx}%XbLV&d2jee8PXGV_ literal 0 HcmV?d00001 diff --git a/test/tests/time_integrators/implicit-euler/tests b/test/tests/time_integrators/implicit-euler/tests index fb3cda87c8f1..234525b6669e 100644 --- a/test/tests/time_integrators/implicit-euler/tests +++ b/test/tests/time_integrators/implicit-euler/tests @@ -1,6 +1,6 @@ [Tests] design = 'ImplicitEuler.md' - issues = '#1953' + issues = '#1953 #29149' [group] requirement = "The system shall support the use of an implicit Euler solver " @@ -23,6 +23,14 @@ [] [] + [linearfv] + type = 'Exodiff' + input = 'ie-linearfv.i' + exodiff = 'ie-linearfv_out.e' + abs_zero = 1e-9 + requirement = "The system shall be able to use implicit euler time integration with linear FV spatial discretization." + [] + [monomials] requirement = "The system shall support the use of an implicit Euler solver with discontinuous " "(first-order Monomial) shape functions." From 9b1eacb97da75b750145349df35717acd6860807 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 16:51:15 -0700 Subject: [PATCH 10/18] Polish source and include files. (#29149) --- framework/include/timeintegrators/ImplicitEuler.h | 2 +- .../include/timeintegrators/LinearTimeIntegrator.h | 3 +-- .../include/timeintegrators/NonlinearTimeIntegrator.h | 9 ++------- framework/include/timeintegrators/TimeIntegrator.h | 5 ++++- .../src/linearfvkernels/LinearFVTimeDerivative.C | 3 ++- framework/src/timeintegrators/CrankNicolson.C | 1 - framework/src/timeintegrators/ExplicitSSPRungeKutta.C | 3 ++- framework/src/timeintegrators/ImplicitEuler.C | 5 +++-- framework/src/timeintegrators/LinearTimeIntegrator.C | 11 +++++++---- .../src/timeintegrators/NonlinearTimeIntegrator.C | 10 ++++++---- 10 files changed, 28 insertions(+), 24 deletions(-) diff --git a/framework/include/timeintegrators/ImplicitEuler.h b/framework/include/timeintegrators/ImplicitEuler.h index 1b5d246ceea8..2e19f8c38e7c 100644 --- a/framework/include/timeintegrators/ImplicitEuler.h +++ b/framework/include/timeintegrators/ImplicitEuler.h @@ -31,7 +31,7 @@ class ImplicitEuler : public TimeIntegrator virtual bool overridesSolve() const override { return false; } virtual Real timeDerivativeRHSContribution(const dof_id_type dof_id, - const std::vector & factors) const override; + const std::vector & factors) const override; virtual Real timeDerivativeMatrixContribution(const Real factor) const override; protected: diff --git a/framework/include/timeintegrators/LinearTimeIntegrator.h b/framework/include/timeintegrators/LinearTimeIntegrator.h index 900ce0955436..d1ffc3f651b6 100644 --- a/framework/include/timeintegrators/LinearTimeIntegrator.h +++ b/framework/include/timeintegrators/LinearTimeIntegrator.h @@ -34,14 +34,13 @@ class LinearTimeIntegrator /// @param factors Multiplicative factor (e.g. a material property) at multiple /// states (old, older, etc) virtual Real timeDerivativeRHSContribution(dof_id_type dof_id, - const std::vector & factors = {}) const; + const std::vector & factors = {}) const; /// The time derivative's contribution to the right hand side of a linear system. /// For now, this does not depend of the DoF index, might change in the furutre. virtual Real timeDerivativeMatrixContribution(const Real factor) const; protected: - /// Pointer to the nonlinear system, can happen that we dont have any LinearSystem * _linear_system; diff --git a/framework/include/timeintegrators/NonlinearTimeIntegrator.h b/framework/include/timeintegrators/NonlinearTimeIntegrator.h index dfeab044b5c8..c02c8f41f998 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegrator.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegrator.h @@ -9,7 +9,6 @@ #pragma once - // MOOSE includes class SystemBase; class NonlinearSystemBase; @@ -56,15 +55,13 @@ class NonlinearTimeIntegrator TagID uDotDotFactorTag() const { return _u_dotdot_factor_tag; } protected: - /// Wrapper around vector addition for nonlinear time integrators. If we don't /// operate on a nonlinear system we don't need to add the vector. /// @param name The name of the vector /// @param project If the vector should be projected /// @param type PThe parallel distribution of the vetor - NumericVector * addVectorForNonlinearTI(const std::string & name, - const bool project, - const ParallelType type); + NumericVector * + addVectorForNonlinearTI(const std::string & name, const bool project, const ParallelType type); /// Pointer to the nonlinear system, can happen that we dont have any NonlinearSystemBase * _nl; @@ -86,6 +83,4 @@ class NonlinearTimeIntegrator /// The vector tag for the nodal multiplication factor for the residual calculation of the udotdot term const TagID _u_dotdot_factor_tag; - - }; diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index 0cf3f616adeb..39880567e8ca 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -37,7 +37,10 @@ class NonlinearImplicitSystem; * used * only by NonlinearSystem (AuxiliarySystem does not produce residual). */ -class TimeIntegrator : public MooseObject, public Restartable, public NonlinearTimeIntegrator, public LinearTimeIntegrator +class TimeIntegrator : public MooseObject, + public Restartable, + public NonlinearTimeIntegrator, + public LinearTimeIntegrator { public: static InputParameters validParams(); diff --git a/framework/src/linearfvkernels/LinearFVTimeDerivative.C b/framework/src/linearfvkernels/LinearFVTimeDerivative.C index 765889f00bd2..9a95afa2590c 100644 --- a/framework/src/linearfvkernels/LinearFVTimeDerivative.C +++ b/framework/src/linearfvkernels/LinearFVTimeDerivative.C @@ -45,7 +45,8 @@ LinearFVTimeDerivative::computeMatrixContribution() Real LinearFVTimeDerivative::computeRightHandSideContribution() { - return _time_integrator.timeDerivativeRHSContribution(_dof_id, _factor_history) * _current_elem_volume; + return _time_integrator.timeDerivativeRHSContribution(_dof_id, _factor_history) * + _current_elem_volume; } void diff --git a/framework/src/timeintegrators/CrankNicolson.C b/framework/src/timeintegrators/CrankNicolson.C index 89ba0e2464cb..671027862e09 100644 --- a/framework/src/timeintegrators/CrankNicolson.C +++ b/framework/src/timeintegrators/CrankNicolson.C @@ -74,7 +74,6 @@ CrankNicolson::init() mooseError("CrankNicolson: Time derivative of solution (`u_dot`) is not stored. Please set " "uDotRequested() to true in FEProblemBase befor requesting `u_dot`."); - // time derivative is assumed to be zero on initial NumericVector & u_dot = *_sys.solutionUDot(); u_dot.zero(); diff --git a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C index 02abad54e86d..113a998a9108 100644 --- a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C +++ b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C @@ -34,7 +34,8 @@ ExplicitSSPRungeKutta::ExplicitSSPRungeKutta(const InputParameters & parameters) : ExplicitTimeIntegrator(parameters), _order(getParam("order")), _stage(0), - _solution_intermediate_stage(addVectorForNonlinearTI("solution_intermediate_stage", false, GHOSTED)), + _solution_intermediate_stage( + addVectorForNonlinearTI("solution_intermediate_stage", false, GHOSTED)), _tmp_solution(addVectorForNonlinearTI("tmp_solution", false, GHOSTED)), _tmp_mass_solution_product(addVectorForNonlinearTI("tmp_mass_solution_product", false, GHOSTED)) { diff --git a/framework/src/timeintegrators/ImplicitEuler.C b/framework/src/timeintegrators/ImplicitEuler.C index 2d27949b3768..a9f22b9e3d95 100644 --- a/framework/src/timeintegrators/ImplicitEuler.C +++ b/framework/src/timeintegrators/ImplicitEuler.C @@ -84,9 +84,10 @@ ImplicitEuler::postResidual(NumericVector & residual) Real ImplicitEuler::timeDerivativeRHSContribution(dof_id_type dof_id, - const std::vector & factors) const + const std::vector & factors) const { - mooseAssert(factors.size() == numStatesRequired(), "Either too many or too few states are given!"); + mooseAssert(factors.size() == numStatesRequired(), + "Either too many or too few states are given!"); return factors[0] * _solution_old(dof_id) / _dt; } diff --git a/framework/src/timeintegrators/LinearTimeIntegrator.C b/framework/src/timeintegrators/LinearTimeIntegrator.C index 73b466178ca6..85de3f7ed853 100644 --- a/framework/src/timeintegrators/LinearTimeIntegrator.C +++ b/framework/src/timeintegrators/LinearTimeIntegrator.C @@ -16,17 +16,20 @@ LinearTimeIntegrator::LinearTimeIntegrator(SystemBase & system) : _linear_system(dynamic_cast(&system)), _integrates_linear_system(_linear_system), - _linear_implicit_system(_linear_system ? dynamic_cast(&_linear_system->system()) : nullptr) + _linear_implicit_system( + _linear_system ? dynamic_cast(&_linear_system->system()) : nullptr) { } -Real LinearTimeIntegrator::timeDerivativeRHSContribution(const dof_id_type /*dof_id*/, - const std::vector & /*factors*/) const +Real +LinearTimeIntegrator::timeDerivativeRHSContribution(const dof_id_type /*dof_id*/, + const std::vector & /*factors*/) const { mooseError("The time derivative right hand side contribution has not been implemented yet!"); } -Real LinearTimeIntegrator::timeDerivativeMatrixContribution(const Real /*factor*/) const +Real +LinearTimeIntegrator::timeDerivativeMatrixContribution(const Real /*factor*/) const { mooseError("The time derivative matrix contribution has not been implemented yet!"); } diff --git a/framework/src/timeintegrators/NonlinearTimeIntegrator.C b/framework/src/timeintegrators/NonlinearTimeIntegrator.C index 1124da334aba..a231b2bb836e 100644 --- a/framework/src/timeintegrators/NonlinearTimeIntegrator.C +++ b/framework/src/timeintegrators/NonlinearTimeIntegrator.C @@ -18,7 +18,8 @@ NonlinearTimeIntegrator::NonlinearTimeIntegrator(FEProblemBase & problem, SystemBase & system) : _nl(dynamic_cast(&system)), _integrates_nl(_nl), - _nonlinear_implicit_system(_integrates_nl ? dynamic_cast(&_nl->system()) : nullptr), + _nonlinear_implicit_system( + _integrates_nl ? dynamic_cast(&_nl->system()) : nullptr), _Re_time(_integrates_nl ? &_nl->getResidualTimeVector() : nullptr), _Re_non_time(_integrates_nl ? &_nl->getResidualNonTimeVector() : nullptr), _u_dot_factor_tag(problem.addVectorTag("u_dot_factor", Moose::VECTOR_TAG_SOLUTION)), @@ -26,9 +27,10 @@ NonlinearTimeIntegrator::NonlinearTimeIntegrator(FEProblemBase & problem, System { } -NumericVector * NonlinearTimeIntegrator::addVectorForNonlinearTI(const std::string & name, - const bool project, - const ParallelType type) +NumericVector * +NonlinearTimeIntegrator::addVectorForNonlinearTI(const std::string & name, + const bool project, + const ParallelType type) { if (_integrates_nl) return &_nl->addVector(name, project, type); From f292071b28ca9a9348c80197a090b8bd97e0e761 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 17:11:44 -0700 Subject: [PATCH 11/18] Add documentation for LinearFVTimeDerivative. (#29149) --- .../linearfvkernels/LinearFVTimeDerivative.md | 45 +++++++++++++++++++ .../linearfvkernels/LinearFVTimeDerivative.h | 4 +- 2 files changed, 47 insertions(+), 2 deletions(-) create mode 100644 framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md diff --git a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md new file mode 100644 index 000000000000..b1cfd5a9c00f --- /dev/null +++ b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md @@ -0,0 +1,45 @@ +# LinearFVTimeDerivative + +## Description + +This kernel represents a time derivative term in a partial differential equation +discretized using the finite volume method: + +!equation +\int\limits_{V_C} \frac{\partial cu}{\partial t} dV \approx \left(\frac{\partial cu}{\partial t}\right)_C V_C~, + +where $\left(\frac{\partial cu}{\partial t}\right)_C$ and $V_C$ are the time derivative of the +field value at the cell center and the cell volume, respectively. +Note that we added a multiplier, $c$ which often represents a material property. +This can be defined through parameter [!param](/LinearFVKernels/LinearFVTimeDerivative/factor) +that accepts anything that supports functor-based evaluations. For more information on functors in +MOOSE, see [Functors/index.md]. +This kernel adds to the matrix diagonal and right hand side of a +linear system and the contributions depend on the +method chosen for time integration. For more information on available methods, see +the [TimeIntegrators](Executioner/TimeIntegrators/index.md) page. +For example, with an implicit euler scheme the controbution to the right hand side becomes: + +!equation +\frac{c_C}{\delta t}V_C, + +where $\delta t$ and $\frac{c_C}$ are the time step size and multiplier at the cell center, + repectively. With these, the contribution to the right hand side becomes: + +!equation +\frac{c_C u_{old,C}}{\delta t}V_C, + +where $u_{old,C}$ represents the solution at the previous time step. + +## Example Syntax + +The case below demonstrates the use of `LinearFVTimeDerivative` used in a simple +linear time-dependent diffusion problem: + +!listing test/tests/time_integrators/implicit-euler/ie-linearfv.i block=LinearFVKernels + +!syntax parameters /LinearFVKernels/LinearFVTimeDerivative + +!syntax inputs /LinearFVKernels/LinearFVTimeDerivative + +!syntax children /LinearFVKernels/LinearFVTimeDerivative diff --git a/framework/include/linearfvkernels/LinearFVTimeDerivative.h b/framework/include/linearfvkernels/LinearFVTimeDerivative.h index 54ef24cb96e9..f9e93712822b 100644 --- a/framework/include/linearfvkernels/LinearFVTimeDerivative.h +++ b/framework/include/linearfvkernels/LinearFVTimeDerivative.h @@ -13,8 +13,8 @@ #include "TimeIntegrator.h" /** - * Kernel that adds contributions from a time derivative term discretized using the - * finite volume method to a linear system. + * Kernel that adds contributions from a time derivative term to a linear system + * populated using the finite volume method. */ class LinearFVTimeDerivative : public LinearFVElementalKernel { From 7820818feddf0c8c7549168c8ede9fbbd1c2c9d9 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 17:23:21 -0700 Subject: [PATCH 12/18] Extend documentation for implicit euler. (#29149) --- .../linearfvkernels/LinearFVTimeDerivative.md | 6 +++--- .../source/timeintegrators/ImplicitEuler.md | 16 ++++++++++++++++ 2 files changed, 19 insertions(+), 3 deletions(-) diff --git a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md index b1cfd5a9c00f..8ecede86bea5 100644 --- a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md +++ b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md @@ -21,13 +21,13 @@ the [TimeIntegrators](Executioner/TimeIntegrators/index.md) page. For example, with an implicit euler scheme the controbution to the right hand side becomes: !equation -\frac{c_C}{\delta t}V_C, +\frac{c_C}{\Delta t}V_C, -where $\delta t$ and $\frac{c_C}$ are the time step size and multiplier at the cell center, +where $\Delta t$ and $\frac{c_C}$ are the time step size and multiplier at the cell center, repectively. With these, the contribution to the right hand side becomes: !equation -\frac{c_C u_{old,C}}{\delta t}V_C, +\frac{c_C u_{old,C}}{\Delta t}V_C, where $u_{old,C}$ represents the solution at the previous time step. diff --git a/framework/doc/content/source/timeintegrators/ImplicitEuler.md b/framework/doc/content/source/timeintegrators/ImplicitEuler.md index 6648a2cfeb90..e4922e5324ed 100644 --- a/framework/doc/content/source/timeintegrators/ImplicitEuler.md +++ b/framework/doc/content/source/timeintegrators/ImplicitEuler.md @@ -21,6 +21,22 @@ This is an implicit system with $U(t+\Delta t)$, the variable to solve for, appe equation. We solve this system iteratively, usually with a Newton or Newton-Krylov method as described in the non linear system solve [documentation](systems/NonlinearSystem.md). +## Contributions to linear systems + +For [linear systems](systems/LinearSystem.md), on top of creating the +time derivatives of the degrees of freedom, this provides contributions +to the matrix and the right hand side. Taking a finite volume system for example, +the contributions to the matrix will be: + +!equation +\frac{1}{\Delta t}V_C, + +where $\Delta t$ and $\frac{V_C}$ are the time step size and cell volume, + repectively. The contribution to the right hand side is: + +!equation +\frac{u_{old,C}}{\Delta t}V_C, + !syntax parameters /Executioner/TimeIntegrator/ImplicitEuler !syntax inputs /Executioner/TimeIntegrator/ImplicitEuler From 9913a6d39a8ced5e7ef10bad20244fb92326c6b8 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 17:28:27 -0700 Subject: [PATCH 13/18] Rename utility classes to interface. (#29149) --- ...meIntegrator.h => LinearTimeIntegratorInterface.h} | 4 ++-- ...ntegrator.h => NonlinearTimeIntegratorInterface.h} | 6 +++--- framework/include/timeintegrators/TimeIntegrator.h | 8 ++++---- ...meIntegrator.C => LinearTimeIntegratorInterface.C} | 10 +++++----- ...ntegrator.C => NonlinearTimeIntegratorInterface.C} | 11 ++++++----- framework/src/timeintegrators/TimeIntegrator.C | 7 ++++--- 6 files changed, 24 insertions(+), 22 deletions(-) rename framework/include/timeintegrators/{LinearTimeIntegrator.h => LinearTimeIntegratorInterface.h} (94%) rename framework/include/timeintegrators/{NonlinearTimeIntegrator.h => NonlinearTimeIntegratorInterface.h} (92%) rename framework/src/timeintegrators/{LinearTimeIntegrator.C => LinearTimeIntegratorInterface.C} (70%) rename framework/src/timeintegrators/{NonlinearTimeIntegrator.C => NonlinearTimeIntegratorInterface.C} (71%) diff --git a/framework/include/timeintegrators/LinearTimeIntegrator.h b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h similarity index 94% rename from framework/include/timeintegrators/LinearTimeIntegrator.h rename to framework/include/timeintegrators/LinearTimeIntegratorInterface.h index d1ffc3f651b6..30fc5ce05053 100644 --- a/framework/include/timeintegrators/LinearTimeIntegrator.h +++ b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h @@ -24,10 +24,10 @@ class LinearImplicitSystem; * Interface class for routines and member variables for time integrators * relying on Newton's method. */ -class LinearTimeIntegrator +class LinearTimeIntegratorInterface { public: - LinearTimeIntegrator(SystemBase & system); + LinearTimeIntegratorInterface(SystemBase & system); /// The time derivative's contribution to the right hand side of a linear system /// @param dof_id The dof index at which this contribution should be fetched at diff --git a/framework/include/timeintegrators/NonlinearTimeIntegrator.h b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h similarity index 92% rename from framework/include/timeintegrators/NonlinearTimeIntegrator.h rename to framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h index c02c8f41f998..fffdcaeef3a9 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegrator.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h @@ -24,13 +24,13 @@ class NonlinearImplicitSystem; * Interface class for routines and member variables for time integrators * relying on Newton's method. */ -class NonlinearTimeIntegrator +class NonLinearTimeIntegratorInterfaceInterface { public: - NonlinearTimeIntegrator(FEProblemBase & problem, SystemBase & system); + NonLinearTimeIntegratorInterfaceInterface(FEProblemBase & problem, SystemBase & system); /** - * Callback to the NonlinearTimeIntegrator called immediately after the + * Callback to the NonLinearTimeIntegratorInterfaceInterface called immediately after the * residuals are computed in NonlinearSystem::computeResidual(). * The residual vector which is passed in to this function should * be filled in by the user with the _Re_time and _Re_non_time diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index 39880567e8ca..30657ec45a43 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -12,8 +12,8 @@ // MOOSE includes #include "MooseObject.h" #include "Restartable.h" -#include "NonlinearTimeIntegrator.h" -#include "LinearTimeIntegrator.h" +#include "NonLinearTimeIntegratorInterfaceInterface.h" +#include "LinearTimeIntegratorInterface.h" class FEProblemBase; class SystemBase; @@ -39,8 +39,8 @@ class NonlinearImplicitSystem; */ class TimeIntegrator : public MooseObject, public Restartable, - public NonlinearTimeIntegrator, - public LinearTimeIntegrator + public NonLinearTimeIntegratorInterfaceInterface, + public LinearTimeIntegratorInterface { public: static InputParameters validParams(); diff --git a/framework/src/timeintegrators/LinearTimeIntegrator.C b/framework/src/timeintegrators/LinearTimeIntegratorInterface.C similarity index 70% rename from framework/src/timeintegrators/LinearTimeIntegrator.C rename to framework/src/timeintegrators/LinearTimeIntegratorInterface.C index 85de3f7ed853..5f8a73ef0277 100644 --- a/framework/src/timeintegrators/LinearTimeIntegrator.C +++ b/framework/src/timeintegrators/LinearTimeIntegratorInterface.C @@ -7,13 +7,13 @@ //* Licensed under LGPL 2.1, please see LICENSE for details //* https://www.gnu.org/licenses/lgpl-2.1.html -#include "LinearTimeIntegrator.h" +#include "LinearTimeIntegratorInterface.h" #include "FEProblem.h" #include "LinearSystem.h" #include "libmesh/linear_implicit_system.h" -LinearTimeIntegrator::LinearTimeIntegrator(SystemBase & system) +LinearTimeIntegratorInterface::LinearTimeIntegratorInterface(SystemBase & system) : _linear_system(dynamic_cast(&system)), _integrates_linear_system(_linear_system), _linear_implicit_system( @@ -22,14 +22,14 @@ LinearTimeIntegrator::LinearTimeIntegrator(SystemBase & system) } Real -LinearTimeIntegrator::timeDerivativeRHSContribution(const dof_id_type /*dof_id*/, - const std::vector & /*factors*/) const +LinearTimeIntegratorInterface::timeDerivativeRHSContribution( + const dof_id_type /*dof_id*/, const std::vector & /*factors*/) const { mooseError("The time derivative right hand side contribution has not been implemented yet!"); } Real -LinearTimeIntegrator::timeDerivativeMatrixContribution(const Real /*factor*/) const +LinearTimeIntegratorInterface::timeDerivativeMatrixContribution(const Real /*factor*/) const { mooseError("The time derivative matrix contribution has not been implemented yet!"); } diff --git a/framework/src/timeintegrators/NonlinearTimeIntegrator.C b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C similarity index 71% rename from framework/src/timeintegrators/NonlinearTimeIntegrator.C rename to framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C index a231b2bb836e..7327f6f38a8c 100644 --- a/framework/src/timeintegrators/NonlinearTimeIntegrator.C +++ b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C @@ -7,7 +7,7 @@ //* Licensed under LGPL 2.1, please see LICENSE for details //* https://www.gnu.org/licenses/lgpl-2.1.html -#include "NonlinearTimeIntegrator.h" +#include "NonLinearTimeIntegratorInterfaceInterface.h" #include "FEProblem.h" #include "NonlinearSystemBase.h" @@ -15,7 +15,8 @@ #include "libmesh/nonlinear_solver.h" #include "libmesh/dof_map.h" -NonlinearTimeIntegrator::NonlinearTimeIntegrator(FEProblemBase & problem, SystemBase & system) +NonLinearTimeIntegratorInterfaceInterface::NonLinearTimeIntegratorInterfaceInterface( + FEProblemBase & problem, SystemBase & system) : _nl(dynamic_cast(&system)), _integrates_nl(_nl), _nonlinear_implicit_system( @@ -28,9 +29,9 @@ NonlinearTimeIntegrator::NonlinearTimeIntegrator(FEProblemBase & problem, System } NumericVector * -NonlinearTimeIntegrator::addVectorForNonlinearTI(const std::string & name, - const bool project, - const ParallelType type) +NonLinearTimeIntegratorInterfaceInterface::addVectorForNonlinearTI(const std::string & name, + const bool project, + const ParallelType type) { if (_integrates_nl) return &_nl->addVector(name, project, type); diff --git a/framework/src/timeintegrators/TimeIntegrator.C b/framework/src/timeintegrators/TimeIntegrator.C index 7010e4df386f..250c050c8ae6 100644 --- a/framework/src/timeintegrators/TimeIntegrator.C +++ b/framework/src/timeintegrators/TimeIntegrator.C @@ -30,9 +30,10 @@ TimeIntegrator::validParams() TimeIntegrator::TimeIntegrator(const InputParameters & parameters) : MooseObject(parameters), Restartable(this, "TimeIntegrators"), - NonlinearTimeIntegrator(*getCheckedPointerParam("_fe_problem_base"), - *getCheckedPointerParam("_sys")), - LinearTimeIntegrator(*getCheckedPointerParam("_sys")), + NonLinearTimeIntegratorInterfaceInterface( + *getCheckedPointerParam("_fe_problem_base"), + *getCheckedPointerParam("_sys")), + LinearTimeIntegratorInterface(*getCheckedPointerParam("_sys")), _fe_problem(*getCheckedPointerParam("_fe_problem_base")), _sys(*getCheckedPointerParam("_sys")), _du_dot_du(_sys.duDotDus()), From c09ff7a0e6124f6796b3a84894783a0403252103 Mon Sep 17 00:00:00 2001 From: Peter German Date: Wed, 4 Dec 2024 17:49:27 -0700 Subject: [PATCH 14/18] Make sure it compiles in nonunity. (#29149) --- .../LinearTimeIntegratorInterface.h | 7 +++++-- .../NonlinearTimeIntegratorInterface.h | 15 +++++++++++---- .../include/timeintegrators/TimeIntegrator.h | 4 ++-- .../NonlinearTimeIntegratorInterface.C | 12 ++++++------ framework/src/timeintegrators/TimeIntegrator.C | 5 ++--- 5 files changed, 26 insertions(+), 17 deletions(-) diff --git a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h index 30fc5ce05053..23aa76697c34 100644 --- a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h @@ -9,7 +9,10 @@ #pragma once -// MOOSE includes +// Moose includes +#include "MooseTypes.h" + +// Forward declarations class SystemBase; class LinearSystem; @@ -22,7 +25,7 @@ class LinearImplicitSystem; /** * Interface class for routines and member variables for time integrators - * relying on Newton's method. + * relying on linear system assembly method. */ class LinearTimeIntegratorInterface { diff --git a/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h index fffdcaeef3a9..eae0e9f6b8fc 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h @@ -9,8 +9,15 @@ #pragma once -// MOOSE includes +// Moose includes +#include "MooseTypes.h" + +// Libmesh includes +#include "libmesh/enum_parallel_type.h" + +// Forward declarations class SystemBase; +class FEProblemBase; class NonlinearSystemBase; namespace libMesh @@ -24,13 +31,13 @@ class NonlinearImplicitSystem; * Interface class for routines and member variables for time integrators * relying on Newton's method. */ -class NonLinearTimeIntegratorInterfaceInterface +class NonlinearTimeIntegratorInterface { public: - NonLinearTimeIntegratorInterfaceInterface(FEProblemBase & problem, SystemBase & system); + NonlinearTimeIntegratorInterface(FEProblemBase & problem, SystemBase & system); /** - * Callback to the NonLinearTimeIntegratorInterfaceInterface called immediately after the + * Callback to the NonLinearTimeIntegratorInterface called immediately after the * residuals are computed in NonlinearSystem::computeResidual(). * The residual vector which is passed in to this function should * be filled in by the user with the _Re_time and _Re_non_time diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index 30657ec45a43..f696d46d3fb8 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -12,7 +12,7 @@ // MOOSE includes #include "MooseObject.h" #include "Restartable.h" -#include "NonLinearTimeIntegratorInterfaceInterface.h" +#include "NonlinearTimeIntegratorInterface.h" #include "LinearTimeIntegratorInterface.h" class FEProblemBase; @@ -39,7 +39,7 @@ class NonlinearImplicitSystem; */ class TimeIntegrator : public MooseObject, public Restartable, - public NonLinearTimeIntegratorInterfaceInterface, + public NonlinearTimeIntegratorInterface, public LinearTimeIntegratorInterface { public: diff --git a/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C index 7327f6f38a8c..c47998295475 100644 --- a/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C +++ b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C @@ -7,7 +7,7 @@ //* Licensed under LGPL 2.1, please see LICENSE for details //* https://www.gnu.org/licenses/lgpl-2.1.html -#include "NonLinearTimeIntegratorInterfaceInterface.h" +#include "NonlinearTimeIntegratorInterface.h" #include "FEProblem.h" #include "NonlinearSystemBase.h" @@ -15,8 +15,8 @@ #include "libmesh/nonlinear_solver.h" #include "libmesh/dof_map.h" -NonLinearTimeIntegratorInterfaceInterface::NonLinearTimeIntegratorInterfaceInterface( - FEProblemBase & problem, SystemBase & system) +NonlinearTimeIntegratorInterface::NonlinearTimeIntegratorInterface(FEProblemBase & problem, + SystemBase & system) : _nl(dynamic_cast(&system)), _integrates_nl(_nl), _nonlinear_implicit_system( @@ -29,9 +29,9 @@ NonLinearTimeIntegratorInterfaceInterface::NonLinearTimeIntegratorInterfaceInter } NumericVector * -NonLinearTimeIntegratorInterfaceInterface::addVectorForNonlinearTI(const std::string & name, - const bool project, - const ParallelType type) +NonlinearTimeIntegratorInterface::addVectorForNonlinearTI(const std::string & name, + const bool project, + const ParallelType type) { if (_integrates_nl) return &_nl->addVector(name, project, type); diff --git a/framework/src/timeintegrators/TimeIntegrator.C b/framework/src/timeintegrators/TimeIntegrator.C index 250c050c8ae6..64d6999e6edf 100644 --- a/framework/src/timeintegrators/TimeIntegrator.C +++ b/framework/src/timeintegrators/TimeIntegrator.C @@ -30,9 +30,8 @@ TimeIntegrator::validParams() TimeIntegrator::TimeIntegrator(const InputParameters & parameters) : MooseObject(parameters), Restartable(this, "TimeIntegrators"), - NonLinearTimeIntegratorInterfaceInterface( - *getCheckedPointerParam("_fe_problem_base"), - *getCheckedPointerParam("_sys")), + NonlinearTimeIntegratorInterface(*getCheckedPointerParam("_fe_problem_base"), + *getCheckedPointerParam("_sys")), LinearTimeIntegratorInterface(*getCheckedPointerParam("_sys")), _fe_problem(*getCheckedPointerParam("_fe_problem_base")), _sys(*getCheckedPointerParam("_sys")), From a90d7694567d13a6145db00e1761f13b36813b35 Mon Sep 17 00:00:00 2001 From: Peter German <31662443+grmnptr@users.noreply.github.com> Date: Thu, 5 Dec 2024 14:30:06 -0700 Subject: [PATCH 15/18] Address Guillaume's review comments Co-authored-by: Guillaume Giudicelli --- .../source/linearfvkernels/LinearFVTimeDerivative.md | 4 ++-- .../doc/content/source/timeintegrators/ImplicitEuler.md | 6 +++--- .../timeintegrators/LinearTimeIntegratorInterface.h | 2 +- .../src/timeintegrators/LinearTimeIntegratorInterface.C | 8 ++++++-- test/tests/time_integrators/implicit-euler/tests | 1 + 5 files changed, 13 insertions(+), 8 deletions(-) diff --git a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md index 8ecede86bea5..fe26bd6d0e6e 100644 --- a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md +++ b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md @@ -18,13 +18,13 @@ This kernel adds to the matrix diagonal and right hand side of a linear system and the contributions depend on the method chosen for time integration. For more information on available methods, see the [TimeIntegrators](Executioner/TimeIntegrators/index.md) page. -For example, with an implicit euler scheme the controbution to the right hand side becomes: +For example, with an implicit Euler scheme the contribution to the right hand side becomes: !equation \frac{c_C}{\Delta t}V_C, where $\Delta t$ and $\frac{c_C}$ are the time step size and multiplier at the cell center, - repectively. With these, the contribution to the right hand side becomes: +respectively. With these, the contribution to the right hand side becomes: !equation \frac{c_C u_{old,C}}{\Delta t}V_C, diff --git a/framework/doc/content/source/timeintegrators/ImplicitEuler.md b/framework/doc/content/source/timeintegrators/ImplicitEuler.md index e4922e5324ed..83a4378ede34 100644 --- a/framework/doc/content/source/timeintegrators/ImplicitEuler.md +++ b/framework/doc/content/source/timeintegrators/ImplicitEuler.md @@ -25,14 +25,14 @@ in the non linear system solve [documentation](systems/NonlinearSystem.md). For [linear systems](systems/LinearSystem.md), on top of creating the time derivatives of the degrees of freedom, this provides contributions -to the matrix and the right hand side. Taking a finite volume system for example, -the contributions to the matrix will be: +to the matrix diagonal and the right hand side. Taking a finite volume system for example, +the contributions to the matrix diagonal will be: !equation \frac{1}{\Delta t}V_C, where $\Delta t$ and $\frac{V_C}$ are the time step size and cell volume, - repectively. The contribution to the right hand side is: +respectively. The contribution to the right hand side is: !equation \frac{u_{old,C}}{\Delta t}V_C, diff --git a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h index 23aa76697c34..0f877bf7b8f0 100644 --- a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h @@ -40,7 +40,7 @@ class LinearTimeIntegratorInterface const std::vector & factors = {}) const; /// The time derivative's contribution to the right hand side of a linear system. - /// For now, this does not depend of the DoF index, might change in the furutre. + /// For now, this does not depend of the DoF index, might change in the future. virtual Real timeDerivativeMatrixContribution(const Real factor) const; protected: diff --git a/framework/src/timeintegrators/LinearTimeIntegratorInterface.C b/framework/src/timeintegrators/LinearTimeIntegratorInterface.C index 5f8a73ef0277..c9344dd2986b 100644 --- a/framework/src/timeintegrators/LinearTimeIntegratorInterface.C +++ b/framework/src/timeintegrators/LinearTimeIntegratorInterface.C @@ -25,11 +25,15 @@ Real LinearTimeIntegratorInterface::timeDerivativeRHSContribution( const dof_id_type /*dof_id*/, const std::vector & /*factors*/) const { - mooseError("The time derivative right hand side contribution has not been implemented yet!"); + mooseError("The time derivative right hand side contribution has not been implemented yet", + _linear_system ? " for time integrator of system " + _linear_system->name() : "", + "!"); } Real LinearTimeIntegratorInterface::timeDerivativeMatrixContribution(const Real /*factor*/) const { - mooseError("The time derivative matrix contribution has not been implemented yet!"); + mooseError("The time derivative matrix contribution has not been implemented yet", + _linear_system ? " for time integrator of system " + _linear_system->name() : "", + "!"); } diff --git a/test/tests/time_integrators/implicit-euler/tests b/test/tests/time_integrators/implicit-euler/tests index 234525b6669e..d897a7af61a9 100644 --- a/test/tests/time_integrators/implicit-euler/tests +++ b/test/tests/time_integrators/implicit-euler/tests @@ -29,6 +29,7 @@ exodiff = 'ie-linearfv_out.e' abs_zero = 1e-9 requirement = "The system shall be able to use implicit euler time integration with linear FV spatial discretization." + max_threads = 1 [] [monomials] From 6d33dbf1610662cac1d70bf8e8239f6276c1d29b Mon Sep 17 00:00:00 2001 From: Peter German <31662443+grmnptr@users.noreply.github.com> Date: Fri, 6 Dec 2024 16:39:29 -0700 Subject: [PATCH 16/18] Address Alex's code review. (#29149) Co-authored-by: Alex Lindsay --- .../source/linearfvkernels/LinearFVTimeDerivative.md | 2 +- .../include/linearfvkernels/LinearFVTimeDerivative.h | 7 +------ .../timeintegrators/LinearTimeIntegratorInterface.h | 4 ++-- .../timeintegrators/NonlinearTimeIntegratorInterface.h | 8 ++++---- framework/include/timeintegrators/TimeIntegrator.h | 1 - framework/src/timeintegrators/AStableDirk4.C | 2 +- framework/src/timeintegrators/CrankNicolson.C | 3 +-- framework/src/timeintegrators/ExplicitRK2.C | 2 +- framework/src/timeintegrators/ExplicitSSPRungeKutta.C | 7 +++---- framework/src/timeintegrators/ExplicitTVDRK2.C | 2 +- framework/src/timeintegrators/ExplicitTimeIntegrator.C | 8 ++++---- framework/src/timeintegrators/ImplicitMidpoint.C | 2 +- framework/src/timeintegrators/LStableDirk2.C | 4 ++-- framework/src/timeintegrators/LStableDirk3.C | 2 +- framework/src/timeintegrators/LStableDirk4.C | 2 +- .../timeintegrators/NonlinearTimeIntegratorInterface.C | 6 +++--- 16 files changed, 27 insertions(+), 35 deletions(-) diff --git a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md index fe26bd6d0e6e..a841dcf28154 100644 --- a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md +++ b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md @@ -23,7 +23,7 @@ For example, with an implicit Euler scheme the contribution to the right hand si !equation \frac{c_C}{\Delta t}V_C, -where $\Delta t$ and $\frac{c_C}$ are the time step size and multiplier at the cell center, +where $\Delta t$ and $c_C$ are the time step size and multiplier at the cell center, respectively. With these, the contribution to the right hand side becomes: !equation diff --git a/framework/include/linearfvkernels/LinearFVTimeDerivative.h b/framework/include/linearfvkernels/LinearFVTimeDerivative.h index f9e93712822b..95bc437d9963 100644 --- a/framework/include/linearfvkernels/LinearFVTimeDerivative.h +++ b/framework/include/linearfvkernels/LinearFVTimeDerivative.h @@ -31,10 +31,6 @@ class LinearFVTimeDerivative : public LinearFVElementalKernel virtual Real computeRightHandSideContribution() override; - /** - * Set the current ElemInfo object. - * @param elem_info The new ElemInfo which will be used as the current - */ virtual void setCurrentElemInfo(const ElemInfo * elem_info) override; protected: @@ -46,8 +42,7 @@ class LinearFVTimeDerivative : public LinearFVElementalKernel const TimeIntegrator & _time_integrator; private: - /// Current values of the material property multiplier. Used - /// for caching. + /// Current and older values of the material property multiplier. std::vector _factor_history; /// State args, the args which will help us fetch the different states of diff --git a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h index 0f877bf7b8f0..0f85e742ec44 100644 --- a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h @@ -44,12 +44,12 @@ class LinearTimeIntegratorInterface virtual Real timeDerivativeMatrixContribution(const Real factor) const; protected: - /// Pointer to the nonlinear system, can happen that we dont have any + /// Pointer to the linear system, can happen that we dont have any LinearSystem * _linear_system; /// Boolean to check if it integrates a linear system const bool _integrates_linear_system; /// Nonlinear implicit system, if applicable; otherwise, nullptr - LinearImplicitSystem * _linear_implicit_system; + libMesh::LinearImplicitSystem * _linear_implicit_system; }; diff --git a/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h index eae0e9f6b8fc..1430586ef497 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h @@ -66,9 +66,9 @@ class NonlinearTimeIntegratorInterface /// operate on a nonlinear system we don't need to add the vector. /// @param name The name of the vector /// @param project If the vector should be projected - /// @param type PThe parallel distribution of the vetor + /// @param type The parallel distribution of the vetor NumericVector * - addVectorForNonlinearTI(const std::string & name, const bool project, const ParallelType type); + addVector(const std::string & name, const bool project, const libMesh::ParallelType type); /// Pointer to the nonlinear system, can happen that we dont have any NonlinearSystemBase * _nl; @@ -76,8 +76,8 @@ class NonlinearTimeIntegratorInterface /// Boolean to check if this integrator belongs to a nonlinear system const bool _integrates_nl; - /// Nonlinear implicit system, if applicable; otherwise, nullptr - NonlinearImplicitSystem * _nonlinear_implicit_system; + /// libMesh nonlinear implicit system, if applicable; otherwise, nullptr + libMesh::NonlinearImplicitSystem * _nonlinear_implicit_system; /// residual vector for time contributions NumericVector * _Re_time; diff --git a/framework/include/timeintegrators/TimeIntegrator.h b/framework/include/timeintegrators/TimeIntegrator.h index f696d46d3fb8..e3410a83348a 100644 --- a/framework/include/timeintegrators/TimeIntegrator.h +++ b/framework/include/timeintegrators/TimeIntegrator.h @@ -23,7 +23,6 @@ namespace libMesh { template class NumericVector; -class NonlinearImplicitSystem; } // namespace libMesh /** diff --git a/framework/src/timeintegrators/AStableDirk4.C b/framework/src/timeintegrators/AStableDirk4.C index de24fb33e263..151f150cbc99 100644 --- a/framework/src/timeintegrators/AStableDirk4.C +++ b/framework/src/timeintegrators/AStableDirk4.C @@ -42,7 +42,7 @@ AStableDirk4::AStableDirk4(const InputParameters & parameters) { std::ostringstream oss; oss << "residual_stage" << stage + 1; - _stage_residuals[stage] = addVectorForNonlinearTI(oss.str(), false, GHOSTED); + _stage_residuals[stage] = addVector(oss.str(), false, GHOSTED); } // Initialize parameters diff --git a/framework/src/timeintegrators/CrankNicolson.C b/framework/src/timeintegrators/CrankNicolson.C index 671027862e09..c2dd27896aab 100644 --- a/framework/src/timeintegrators/CrankNicolson.C +++ b/framework/src/timeintegrators/CrankNicolson.C @@ -22,8 +22,7 @@ CrankNicolson::validParams() } CrankNicolson::CrankNicolson(const InputParameters & parameters) - : TimeIntegrator(parameters), - _residual_old(addVectorForNonlinearTI("residual_old", false, libMesh::GHOSTED)) + : TimeIntegrator(parameters), _residual_old(addVector("residual_old", false, libMesh::GHOSTED)) { } diff --git a/framework/src/timeintegrators/ExplicitRK2.C b/framework/src/timeintegrators/ExplicitRK2.C index 7bd463d8deb8..76107e1f4618 100644 --- a/framework/src/timeintegrators/ExplicitRK2.C +++ b/framework/src/timeintegrators/ExplicitRK2.C @@ -25,7 +25,7 @@ ExplicitRK2::validParams() ExplicitRK2::ExplicitRK2(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_old(addVectorForNonlinearTI("residual_old", false, GHOSTED)), + _residual_old(addVector("residual_old", false, GHOSTED)), _solution_older(_sys.solutionState(2)) { mooseInfo("ExplicitRK2-derived TimeIntegrators (ExplicitMidpoint, Heun, Ralston) and other " diff --git a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C index 113a998a9108..3c152f5058ed 100644 --- a/framework/src/timeintegrators/ExplicitSSPRungeKutta.C +++ b/framework/src/timeintegrators/ExplicitSSPRungeKutta.C @@ -34,10 +34,9 @@ ExplicitSSPRungeKutta::ExplicitSSPRungeKutta(const InputParameters & parameters) : ExplicitTimeIntegrator(parameters), _order(getParam("order")), _stage(0), - _solution_intermediate_stage( - addVectorForNonlinearTI("solution_intermediate_stage", false, GHOSTED)), - _tmp_solution(addVectorForNonlinearTI("tmp_solution", false, GHOSTED)), - _tmp_mass_solution_product(addVectorForNonlinearTI("tmp_mass_solution_product", false, GHOSTED)) + _solution_intermediate_stage(addVector("solution_intermediate_stage", false, GHOSTED)), + _tmp_solution(addVector("tmp_solution", false, GHOSTED)), + _tmp_mass_solution_product(addVector("tmp_mass_solution_product", false, GHOSTED)) { // For SSPRK methods up to order 3, the number of stages equals the order _n_stages = _order; diff --git a/framework/src/timeintegrators/ExplicitTVDRK2.C b/framework/src/timeintegrators/ExplicitTVDRK2.C index 3fed86ba9141..336486955685 100644 --- a/framework/src/timeintegrators/ExplicitTVDRK2.C +++ b/framework/src/timeintegrators/ExplicitTVDRK2.C @@ -26,7 +26,7 @@ ExplicitTVDRK2::validParams() ExplicitTVDRK2::ExplicitTVDRK2(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_old(addVectorForNonlinearTI("residual_old", false, libMesh::GHOSTED)), + _residual_old(addVector("residual_old", false, libMesh::GHOSTED)), _solution_older(_sys.solutionState(2)) { mooseInfo("ExplicitTVDRK2 and other multistage TimeIntegrators are known not to work with " diff --git a/framework/src/timeintegrators/ExplicitTimeIntegrator.C b/framework/src/timeintegrators/ExplicitTimeIntegrator.C index 0af2a1e2674f..616ddffbfcc2 100644 --- a/framework/src/timeintegrators/ExplicitTimeIntegrator.C +++ b/framework/src/timeintegrators/ExplicitTimeIntegrator.C @@ -41,9 +41,9 @@ ExplicitTimeIntegrator::ExplicitTimeIntegrator(const InputParameters & parameter MeshChangedInterface(parameters), _solve_type(getParam("solve_type")), - _explicit_residual(addVectorForNonlinearTI("explicit_residual", false, PARALLEL)), - _solution_update(addVectorForNonlinearTI("solution_update", true, PARALLEL)), - _mass_matrix_diag(addVectorForNonlinearTI("mass_matrix_diag", false, PARALLEL)) + _explicit_residual(addVector("explicit_residual", false, PARALLEL)), + _solution_update(addVector("solution_update", true, PARALLEL)), + _mass_matrix_diag(addVector("mass_matrix_diag", false, PARALLEL)) { _Ke_time_tag = _fe_problem.getMatrixTagID("TIME"); @@ -52,7 +52,7 @@ ExplicitTimeIntegrator::ExplicitTimeIntegrator(const InputParameters & parameter _fe_problem.solverParams()._type = Moose::ST_LINEAR; if (_solve_type == LUMPED || _solve_type == LUMP_PRECONDITIONED) - _ones = addVectorForNonlinearTI("ones", false, PARALLEL); + _ones = addVector("ones", false, PARALLEL); } void diff --git a/framework/src/timeintegrators/ImplicitMidpoint.C b/framework/src/timeintegrators/ImplicitMidpoint.C index 01bebab105db..e745b5ce0e03 100644 --- a/framework/src/timeintegrators/ImplicitMidpoint.C +++ b/framework/src/timeintegrators/ImplicitMidpoint.C @@ -25,7 +25,7 @@ ImplicitMidpoint::validParams() ImplicitMidpoint::ImplicitMidpoint(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_stage1(addVectorForNonlinearTI("residual_stage1", false, libMesh::GHOSTED)) + _residual_stage1(addVector("residual_stage1", false, libMesh::GHOSTED)) { mooseInfo("ImplicitMidpoint and other multistage TimeIntegrators are known not to work with " "Materials/AuxKernels that accumulate 'state' and should be used with caution."); diff --git a/framework/src/timeintegrators/LStableDirk2.C b/framework/src/timeintegrators/LStableDirk2.C index 37c39bc5a09d..536f083c6fff 100644 --- a/framework/src/timeintegrators/LStableDirk2.C +++ b/framework/src/timeintegrators/LStableDirk2.C @@ -28,8 +28,8 @@ LStableDirk2::validParams() LStableDirk2::LStableDirk2(const InputParameters & parameters) : TimeIntegrator(parameters), _stage(1), - _residual_stage1(addVectorForNonlinearTI("residual_stage1", false, GHOSTED)), - _residual_stage2(addVectorForNonlinearTI("residual_stage2", false, GHOSTED)), + _residual_stage1(addVector("residual_stage1", false, GHOSTED)), + _residual_stage2(addVector("residual_stage2", false, GHOSTED)), _alpha(1. - 0.5 * std::sqrt(2)) { mooseInfo("LStableDirk2 and other multistage TimeIntegrators are known not to work with " diff --git a/framework/src/timeintegrators/LStableDirk3.C b/framework/src/timeintegrators/LStableDirk3.C index 80b1769fabf4..4449399e7bdd 100644 --- a/framework/src/timeintegrators/LStableDirk3.C +++ b/framework/src/timeintegrators/LStableDirk3.C @@ -38,7 +38,7 @@ LStableDirk3::LStableDirk3(const InputParameters & parameters) { std::ostringstream oss; oss << "residual_stage" << stage + 1; - _stage_residuals[stage] = addVectorForNonlinearTI(oss.str(), false, GHOSTED); + _stage_residuals[stage] = addVector(oss.str(), false, GHOSTED); } // Initialize parameters diff --git a/framework/src/timeintegrators/LStableDirk4.C b/framework/src/timeintegrators/LStableDirk4.C index dcfb23829797..159653076d74 100644 --- a/framework/src/timeintegrators/LStableDirk4.C +++ b/framework/src/timeintegrators/LStableDirk4.C @@ -46,7 +46,7 @@ LStableDirk4::LStableDirk4(const InputParameters & parameters) { std::ostringstream oss; oss << "residual_stage" << stage + 1; - _stage_residuals[stage] = addVectorForNonlinearTI(oss.str(), false, GHOSTED); + _stage_residuals[stage] = addVector(oss.str(), false, GHOSTED); } } diff --git a/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C index c47998295475..ad74656ecd61 100644 --- a/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C +++ b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C @@ -29,9 +29,9 @@ NonlinearTimeIntegratorInterface::NonlinearTimeIntegratorInterface(FEProblemBase } NumericVector * -NonlinearTimeIntegratorInterface::addVectorForNonlinearTI(const std::string & name, - const bool project, - const ParallelType type) +NonlinearTimeIntegratorInterface::addVector(const std::string & name, + const bool project, + const libMesh::ParallelType type) { if (_integrates_nl) return &_nl->addVector(name, project, type); From 834dc42a546882486044f5bc3978a351f735947d Mon Sep 17 00:00:00 2001 From: Peter German Date: Tue, 10 Dec 2024 20:22:49 -0700 Subject: [PATCH 17/18] Address Mauricio's review. (#29149) --- .../content/source/linearfvkernels/LinearFVTimeDerivative.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md index a841dcf28154..62482e14afb0 100644 --- a/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md +++ b/framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md @@ -11,6 +11,8 @@ discretized using the finite volume method: where $\left(\frac{\partial cu}{\partial t}\right)_C$ and $V_C$ are the time derivative of the field value at the cell center and the cell volume, respectively. Note that we added a multiplier, $c$ which often represents a material property. +A good example for the multiplier can be the density in the momentum equation +in the Navier Stokes equation. This can be defined through parameter [!param](/LinearFVKernels/LinearFVTimeDerivative/factor) that accepts anything that supports functor-based evaluations. For more information on functors in MOOSE, see [Functors/index.md]. From fc547cd8de175e49e0f7d12781f1c81f43469e00 Mon Sep 17 00:00:00 2001 From: Peter German Date: Thu, 12 Dec 2024 12:39:48 -0700 Subject: [PATCH 18/18] Remove integrates_nl/linear and just use pointers. (#29149) --- .../timeintegrators/LinearTimeIntegratorInterface.h | 3 --- .../NonlinearTimeIntegratorInterface.h | 3 --- framework/src/timeintegrators/CrankNicolson.C | 2 +- .../timeintegrators/LinearTimeIntegratorInterface.C | 1 - .../NonlinearTimeIntegratorInterface.C | 11 +++++------ 5 files changed, 6 insertions(+), 14 deletions(-) diff --git a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h index 0f85e742ec44..376c65d95b0d 100644 --- a/framework/include/timeintegrators/LinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/LinearTimeIntegratorInterface.h @@ -47,9 +47,6 @@ class LinearTimeIntegratorInterface /// Pointer to the linear system, can happen that we dont have any LinearSystem * _linear_system; - /// Boolean to check if it integrates a linear system - const bool _integrates_linear_system; - /// Nonlinear implicit system, if applicable; otherwise, nullptr libMesh::LinearImplicitSystem * _linear_implicit_system; }; diff --git a/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h index 1430586ef497..9bc56bb7fd1d 100644 --- a/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h +++ b/framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h @@ -73,9 +73,6 @@ class NonlinearTimeIntegratorInterface /// Pointer to the nonlinear system, can happen that we dont have any NonlinearSystemBase * _nl; - /// Boolean to check if this integrator belongs to a nonlinear system - const bool _integrates_nl; - /// libMesh nonlinear implicit system, if applicable; otherwise, nullptr libMesh::NonlinearImplicitSystem * _nonlinear_implicit_system; diff --git a/framework/src/timeintegrators/CrankNicolson.C b/framework/src/timeintegrators/CrankNicolson.C index c2dd27896aab..64ee40182eb7 100644 --- a/framework/src/timeintegrators/CrankNicolson.C +++ b/framework/src/timeintegrators/CrankNicolson.C @@ -78,7 +78,7 @@ CrankNicolson::init() u_dot.zero(); computeDuDotDu(); - if (_integrates_nl) + if (_nl) { // compute residual for the initial time step // Note: we can not directly pass _residual_old in computeResidualTag because diff --git a/framework/src/timeintegrators/LinearTimeIntegratorInterface.C b/framework/src/timeintegrators/LinearTimeIntegratorInterface.C index c9344dd2986b..da94a504837e 100644 --- a/framework/src/timeintegrators/LinearTimeIntegratorInterface.C +++ b/framework/src/timeintegrators/LinearTimeIntegratorInterface.C @@ -15,7 +15,6 @@ LinearTimeIntegratorInterface::LinearTimeIntegratorInterface(SystemBase & system) : _linear_system(dynamic_cast(&system)), - _integrates_linear_system(_linear_system), _linear_implicit_system( _linear_system ? dynamic_cast(&_linear_system->system()) : nullptr) { diff --git a/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C index ad74656ecd61..7784e6c98df2 100644 --- a/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C +++ b/framework/src/timeintegrators/NonlinearTimeIntegratorInterface.C @@ -18,11 +18,10 @@ NonlinearTimeIntegratorInterface::NonlinearTimeIntegratorInterface(FEProblemBase & problem, SystemBase & system) : _nl(dynamic_cast(&system)), - _integrates_nl(_nl), - _nonlinear_implicit_system( - _integrates_nl ? dynamic_cast(&_nl->system()) : nullptr), - _Re_time(_integrates_nl ? &_nl->getResidualTimeVector() : nullptr), - _Re_non_time(_integrates_nl ? &_nl->getResidualNonTimeVector() : nullptr), + _nonlinear_implicit_system(_nl ? dynamic_cast(&_nl->system()) + : nullptr), + _Re_time(_nl ? &_nl->getResidualTimeVector() : nullptr), + _Re_non_time(_nl ? &_nl->getResidualNonTimeVector() : nullptr), _u_dot_factor_tag(problem.addVectorTag("u_dot_factor", Moose::VECTOR_TAG_SOLUTION)), _u_dotdot_factor_tag(problem.addVectorTag("u_dotdot_factor", Moose::VECTOR_TAG_SOLUTION)) { @@ -33,7 +32,7 @@ NonlinearTimeIntegratorInterface::addVector(const std::string & name, const bool project, const libMesh::ParallelType type) { - if (_integrates_nl) + if (_nl) return &_nl->addVector(name, project, type); else return nullptr;