Merge pull request idaholab#29309 from grmnptr/linearfv-timeintegral-…

…29419

Add transient capabilities to linear systems.

framework/doc/content/source/linearfvkernels/LinearFVTimeDerivative.md

@@ -0,0 +1,47 @@
+# 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.
+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].
+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 contribution to the right hand side becomes:
+
+!equation
+\frac{c_C}{\Delta t}V_C,
+
+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
+\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

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 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,
+respectively. 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

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;
 };

framework/include/linearfvkernels/LinearFVTimeDerivative.h

@@ -0,0 +1,52 @@
+//* 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 to a linear system
+ * populated using the finite volume method.
+ */
+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;
+
+  virtual void setCurrentElemInfo(const ElemInfo * elem_info) override;
+
+protected:
+  /// The functor for the material property multipler
+  const Moose::Functor<Real> & _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 and older values of the material property multiplier.
+  std::vector<Real> _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<Moose::StateArg> _state_args;
+};

framework/include/timeintegrators/CrankNicolson.h

@@ -45,7 +45,7 @@ class CrankNicolson : public TimeIntegrator
 
   virtual Real duDotDuCoeff() const override;
 
-  NumericVector<Number> & _residual_old;
+  NumericVector<Number> * _residual_old;
 };
 
 template <typename T, typename T2>

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<Number> & _residual_old;
+  NumericVector<Number> * _residual_old;
 
   /// The older solution
   const NumericVector<Number> & _solution_older;

framework/include/timeintegrators/ExplicitSSPRungeKutta.h

@@ -57,9 +57,9 @@ class ExplicitSSPRungeKutta : public ExplicitTimeIntegrator
   std::vector<const NumericVector<Number> *> _solution_stage;
 
   /// Solution vector for intermediate stage
-  NumericVector<Number> & _solution_intermediate_stage;
+  NumericVector<Number> * _solution_intermediate_stage;
   /// Temporary solution vector
-  NumericVector<Number> & _tmp_solution;
+  NumericVector<Number> * _tmp_solution;
   /// Temporary mass-matrix/solution vector product
-  NumericVector<Number> & _tmp_mass_solution_product;
+  NumericVector<Number> * _tmp_mass_solution_product;
 };

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<Number> & _residual_old;
+  NumericVector<Number> * _residual_old;
 
   /// The older solution
   const NumericVector<Number> & _solution_older;

framework/include/timeintegrators/ExplicitTimeIntegrator.h

@@ -66,13 +66,13 @@ class ExplicitTimeIntegrator : public TimeIntegrator, public MeshChangedInterfac
   MooseEnum _solve_type;
 
   /// Residual used for the RHS
-  NumericVector<Real> & _explicit_residual;
+  NumericVector<Real> * _explicit_residual;
 
   /// Solution vector for the linear solve
-  NumericVector<Real> & _solution_update;
+  NumericVector<Real> * _solution_update;
 
   /// Diagonal of the lumped mass matrix (and its inversion)
-  NumericVector<Real> & _mass_matrix_diag;
+  NumericVector<Real> * _mass_matrix_diag;
 
   /// Vector of 1's to help with creating the lumped mass matrix
   NumericVector<Real> * _ones;

framework/include/timeintegrators/ImplicitEuler.h

@@ -30,6 +30,10 @@ class ImplicitEuler : public TimeIntegrator
   virtual void postResidual(NumericVector<Number> & residual) override;
   virtual bool overridesSolve() const override { return false; }
 
+  virtual Real timeDerivativeRHSContribution(const dof_id_type dof_id,
+                                             const std::vector<Real> & factors) const override;
+  virtual Real timeDerivativeMatrixContribution(const Real factor) const override;
+
 protected:
   /**
    * Helper function that actually does the math for computing the time derivative

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<Number> & _residual_stage1;
+  NumericVector<Number> * _residual_stage1;
 };
 
 template <typename T, typename T2>

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<Number> & _residual_stage1;
+  NumericVector<Number> * _residual_stage1;
 
   //! Buffer to store non-time residual from second stage solve
-  NumericVector<Number> & _residual_stage2;
+  NumericVector<Number> * _residual_stage2;
 
   // The parameter of the method, set at construction time and cannot be changed.
   const Real _alpha;

framework/include/timeintegrators/LinearTimeIntegratorInterface.h

@@ -0,0 +1,52 @@
+//* 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
+#include "MooseTypes.h"
+
+// Forward declarations
+class SystemBase;
+class LinearSystem;
+
+namespace libMesh
+{
+template <typename T>
+class NumericVector;
+class LinearImplicitSystem;
+} // namespace libMesh
+
+/**
+ * Interface class for routines and member variables for time integrators
+ * relying on linear system assembly method.
+ */
+class LinearTimeIntegratorInterface
+{
+public:
+  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
+  /// @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<Real> & 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 future.
+  virtual Real timeDerivativeMatrixContribution(const Real factor) const;
+
+protected:
+  /// Pointer to the linear system, can happen that we dont have any
+  LinearSystem * _linear_system;
+
+  /// Nonlinear implicit system, if applicable; otherwise, nullptr
+  libMesh::LinearImplicitSystem * _linear_implicit_system;
+};

framework/include/timeintegrators/NonlinearTimeIntegratorInterface.h

@@ -0,0 +1,90 @@
+//* 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
+#include "MooseTypes.h"
+
+// Libmesh includes
+#include "libmesh/enum_parallel_type.h"
+
+// Forward declarations
+class SystemBase;
+class FEProblemBase;
+class NonlinearSystemBase;
+
+namespace libMesh
+{
+template <typename T>
+class NumericVector;
+class NonlinearImplicitSystem;
+} // namespace libMesh
+
+/**
+ * Interface class for routines and member variables for time integrators
+ * relying on Newton's method.
+ */
+class NonlinearTimeIntegratorInterface
+{
+public:
+  NonlinearTimeIntegratorInterface(FEProblemBase & problem, SystemBase & system);
+
+  /**
+   * 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
+   * vectors in a way that makes sense for the particular
+   * TimeIntegration method.
+   */
+  virtual void postResidual(NumericVector<Number> & /*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. 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 The parallel distribution of the vetor
+  NumericVector<Number> *
+  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;
+
+  /// libMesh nonlinear implicit system, if applicable; otherwise, nullptr
+  libMesh::NonlinearImplicitSystem * _nonlinear_implicit_system;
+
+  /// residual vector for time contributions
+  NumericVector<Number> * _Re_time;
+
+  /// residual vector for non-time contributions
+  NumericVector<Number> * _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;
+};