OGS 6.2.1-76-gbb689931b
NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton > Class Template Referencefinal

Detailed Description

template<>
class NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >

Time discretized first order implicit quasi-linear ODE; to be solved using the Newton-Raphson method for resolving nonlinearities.

See also
ODESystemTag::FirstOrderImplicitQuasilinear

Definition at line 57 of file TimeDiscretizedODESystem.h.

#include <TimeDiscretizedODESystem.h>

Inheritance diagram for NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >:
Collaboration diagram for NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >:

Public Types

using ODE = ODESystem< ODETag, NonlinearSolverTag::Newton >
 The type of ODE. More...
 
using MatTrans = MatrixTranslator< ODETag >
 
using TimeDisc = TimeDiscretization
 A shortcut for a general time discretization scheme. More...
 

Public Member Functions

 TimeDiscretizedODESystem (const int process_id, ODE &ode, TimeDisc &time_discretization)
 
 ~TimeDiscretizedODESystem () override
 
void assemble (const GlobalVector &x_new_timestep) override
 
void getResidual (GlobalVector const &x_new_timestep, GlobalVector &res) const override
 
void getJacobian (GlobalMatrix &Jac) const override
 
void computeKnownSolutions (GlobalVector const &x) override
 
void applyKnownSolutions (GlobalVector &x) const override
 
void applyKnownSolutionsNewton (GlobalMatrix &Jac, GlobalVector &res, GlobalVector &minus_delta_x) const override
 
bool isLinear () const override
 
void preIteration (const unsigned iter, GlobalVector const &x) override
 
IterationResult postIteration (GlobalVector const &x) override
 
void pushMatrices () const override
 
TimeDiscgetTimeDiscretization () override
 Exposes the used time discretization scheme. More...
 
MathLib::MatrixSpecifications getMatrixSpecifications (const int process_id) const override
 
- Public Member Functions inherited from NumLib::InternalMatrixStorage
virtual ~InternalMatrixStorage ()=default
 

Static Public Attributes

static const ODESystemTag ODETag
 A tag indicating the type of ODE. More...
 

Private Types

using Index = MathLib::MatrixVectorTraits< GlobalMatrix >::Index
 

Private Attributes

ODE_ode
 ode the ODE being wrapped More...
 
TimeDisc_time_disc
 the time discretization to being used More...
 
std::unique_ptr< MatTrans_mat_trans
 the object used to compute the matrix/vector for the nonlinear solver More...
 
std::vector< NumLib::IndexValueVector< Index > > const * _known_solutions
 stores precomputed values for known solutions More...
 
GlobalMatrix * _Jac
 the Jacobian of the residual More...
 
GlobalMatrix * _M
 Matrix $ M $. More...
 
GlobalMatrix * _K
 Matrix $ K $. More...
 
GlobalVector * _b
 Matrix $ b $. More...
 
std::size_t _Jac_id = 0u
 ID of the _Jac matrix. More...
 
std::size_t _M_id = 0u
 ID of the _M matrix. More...
 
std::size_t _K_id = 0u
 ID of the _K matrix. More...
 
std::size_t _b_id = 0u
 ID of the _b vector. More...
 
std::size_t _xdot_id = 0u
 ID of the vector storing xdot in intermediate computations. More...
 

Member Typedef Documentation

◆ Index

◆ MatTrans

The auxiliary class that computes the matrix/vector used by the nonlinear solver.

Definition at line 70 of file TimeDiscretizedODESystem.h.

◆ ODE

◆ TimeDisc

Constructor & Destructor Documentation

◆ TimeDiscretizedODESystem()

Constructs a new instance.

Parameters
process_idID of the ODE to be solved.
odethe ODE to be wrapped.
time_discretizationthe time discretization to be used.

Definition at line 45 of file TimeDiscretizedODESystem.cpp.

References NumLib::FirstOrderImplicitQuasilinear, NumLib::MatrixProvider::getMatrix(), NumLib::VectorProvider::getVector(), NumLib::GlobalVectorProvider::provider, and NumLib::GlobalMatrixProvider::provider.

47  : _ode(ode),
48  _time_disc(time_discretization),
49  _mat_trans(createMatrixTranslator<ODETag>(time_discretization))
50 {
52  _ode.getMatrixSpecifications(process_id), _Jac_id);
54  _ode.getMatrixSpecifications(process_id), _M_id);
56  _ode.getMatrixSpecifications(process_id), _K_id);
58  _ode.getMatrixSpecifications(process_id), _b_id);
59 }
virtual GlobalVector & getVector()=0
Get an uninitialized vector.
std::unique_ptr< MatTrans > _mat_trans
the object used to compute the matrix/vector for the nonlinear solver
static NUMLIB_EXPORT MatrixProvider & provider
static NUMLIB_EXPORT VectorProvider & provider
virtual GlobalMatrix & getMatrix()=0
Get an uninitialized matrix.

◆ ~TimeDiscretizedODESystem()

Definition at line 63 of file TimeDiscretizedODESystem.cpp.

References NumLib::FirstOrderImplicitQuasilinear, NumLib::GlobalVectorProvider::provider, NumLib::GlobalMatrixProvider::provider, NumLib::MatrixProvider::releaseMatrix(), and NumLib::VectorProvider::releaseVector().

Member Function Documentation

◆ applyKnownSolutions()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::applyKnownSolutions ( GlobalVector &  x) const
override

Definition at line 134 of file TimeDiscretizedODESystem.cpp.

References detail::applyKnownSolutions(), NumLib::FirstOrderImplicitQuasilinear, and NumLib::Newton.

135 {
137 }
std::vector< NumLib::IndexValueVector< Index > > const * _known_solutions
stores precomputed values for known solutions
void applyKnownSolutions(std::vector< Solutions > const *const known_solutions, Vector &x)
Applies known solutions to the solution vector x.

◆ applyKnownSolutionsNewton()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::applyKnownSolutionsNewton ( GlobalMatrix &  Jac,
GlobalVector &  res,
GlobalVector &  minus_delta_x 
) const
override

Definition at line 141 of file TimeDiscretizedODESystem.cpp.

References MathLib::applyKnownSolution(), MathLib::LinAlg::copy(), NumLib::FirstOrderImplicitQuasilinear, and NumLib::Picard.

143 {
144  if (!_known_solutions)
145  {
146  return;
147  }
148 
150  std::vector<IndexType> ids;
151  for (auto const& bc : *_known_solutions)
152  {
153  std::copy(bc.ids.cbegin(), bc.ids.cend(), std::back_inserter(ids));
154  }
155 
156  // For the Newton method the values must be zero
157  std::vector<double> values(ids.size(), 0);
158  MathLib::applyKnownSolution(Jac, res, minus_delta_x, ids, values);
159 }
std::vector< NumLib::IndexValueVector< Index > > const * _known_solutions
stores precomputed values for known solutions
void applyKnownSolution(EigenMatrix &A_, EigenVector &b_, EigenVector &, const std::vector< EigenMatrix::IndexType > &vec_knownX_id, const std::vector< double > &vec_knownX_x, double)
Definition: EigenTools.cpp:19
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:36

◆ assemble()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::assemble ( const GlobalVector &  x_new_timestep)
override

Definition at line 73 of file TimeDiscretizedODESystem.cpp.

References MathLib::LinAlg::finalizeAssembly(), NumLib::FirstOrderImplicitQuasilinear, NumLib::VectorProvider::getVector(), NumLib::GlobalVectorProvider::provider, and NumLib::VectorProvider::releaseVector().

74 {
75  namespace LinAlg = MathLib::LinAlg;
76 
77  auto const t = _time_disc.getCurrentTime();
78  auto const& x_curr = _time_disc.getCurrentX(x_new_timestep);
79  auto const dxdot_dx = _time_disc.getNewXWeight();
80  auto const dx_dx = _time_disc.getDxDx();
81 
83  _time_disc.getXdot(x_new_timestep, xdot);
84 
85  _M->setZero();
86  _K->setZero();
87  _b->setZero();
88  _Jac->setZero();
89 
90  _ode.preAssemble(t, x_curr);
91  _ode.assembleWithJacobian(t, x_curr, xdot, dxdot_dx, dx_dx, *_M, *_K, *_b,
92  *_Jac);
93 
98 
100 }
virtual double getNewXWeight() const =0
Returns .
virtual double getCurrentTime() const =0
void finalizeAssembly(Matrix &)
virtual GlobalVector & getVector()=0
Get an uninitialized vector.
static NUMLIB_EXPORT VectorProvider & provider
virtual void releaseVector(GlobalVector const &x)=0
virtual double getDxDx() const
virtual GlobalVector const & getCurrentX(GlobalVector const &x_at_new_timestep) const
void getXdot(GlobalVector const &x_at_new_timestep, GlobalVector &xdot) const

◆ computeKnownSolutions()

◆ getJacobian()

Definition at line 120 of file TimeDiscretizedODESystem.cpp.

References NumLib::FirstOrderImplicitQuasilinear.

121 {
122  _mat_trans->computeJacobian(*_Jac, Jac);
123 }
std::unique_ptr< MatTrans > _mat_trans
the object used to compute the matrix/vector for the nonlinear solver

◆ getMatrixSpecifications()

Definition at line 120 of file TimeDiscretizedODESystem.h.

122  {
123  return _ode.getMatrixSpecifications(process_id);
124  }

◆ getResidual()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::getResidual ( GlobalVector const &  x_new_timestep,
GlobalVector &  res 
) const
override

Definition at line 104 of file TimeDiscretizedODESystem.cpp.

References NumLib::FirstOrderImplicitQuasilinear, NumLib::VectorProvider::getVector(), NumLib::GlobalVectorProvider::provider, and NumLib::VectorProvider::releaseVector().

106 {
107  // TODO Maybe the duplicate calculation of xdot here and in assembleJacobian
108  // can be optimuized. However, that would make the interface a bit more
109  // fragile.
111  _time_disc.getXdot(x_new_timestep, xdot);
112 
113  _mat_trans->computeResidual(*_M, *_K, *_b, x_new_timestep, xdot, res);
114 
116 }
virtual GlobalVector & getVector()=0
Get an uninitialized vector.
std::unique_ptr< MatTrans > _mat_trans
the object used to compute the matrix/vector for the nonlinear solver
static NUMLIB_EXPORT VectorProvider & provider
virtual void releaseVector(GlobalVector const &x)=0
void getXdot(GlobalVector const &x_at_new_timestep, GlobalVector &xdot) const

◆ getTimeDiscretization()

◆ isLinear()

◆ postIteration()

◆ preIteration()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::preIteration ( const unsigned  iter,
GlobalVector const &  x 
)
inlineoverride

Definition at line 104 of file TimeDiscretizedODESystem.h.

◆ pushMatrices()

Triggers a refresh of the internal matrix/vector storage.

Remarks
This method is needed in particular to fully implement the interaction of the CrankNicolson scheme with other classes.
Attention
This method must be called (if it is called) from within TimeDiscretization::pushState() after the internal state of the TimeDiscretization has been set to the new solution. Otherwise the pushMatrices() method of MatrixTranslator's will break!

Implements NumLib::InternalMatrixStorage.

Definition at line 114 of file TimeDiscretizedODESystem.h.

Member Data Documentation

◆ _b

◆ _b_id

ID of the _b vector.

Definition at line 145 of file TimeDiscretizedODESystem.h.

◆ _Jac

the Jacobian of the residual

Definition at line 137 of file TimeDiscretizedODESystem.h.

◆ _Jac_id

ID of the _Jac matrix.

Definition at line 142 of file TimeDiscretizedODESystem.h.

◆ _K

◆ _K_id

ID of the _K matrix.

Definition at line 144 of file TimeDiscretizedODESystem.h.

◆ _known_solutions

Initial value:
=
nullptr

stores precomputed values for known solutions

Definition at line 134 of file TimeDiscretizedODESystem.h.

◆ _M

◆ _M_id

ID of the _M matrix.

Definition at line 143 of file TimeDiscretizedODESystem.h.

◆ _mat_trans

the object used to compute the matrix/vector for the nonlinear solver

Definition at line 131 of file TimeDiscretizedODESystem.h.

◆ _ode

◆ _time_disc

the time discretization to being used

Definition at line 128 of file TimeDiscretizedODESystem.h.

◆ _xdot_id

ID of the vector storing xdot in intermediate computations.

Definition at line 148 of file TimeDiscretizedODESystem.h.

◆ ODETag


The documentation for this class was generated from the following files: