OGS  v6.3.3
NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard > Class Referencefinal

Detailed Description

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

See also
ODESystemTag::FirstOrderImplicitQuasilinear

Definition at line 154 of file TimeDiscretizedODESystem.h.

#include <TimeDiscretizedODESystem.h>

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

Public Types

using ODE = ODESystem< ODETag, NonlinearSolverTag::Picard >
 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)
 Constructs a new instance. More...
 
 ~TimeDiscretizedODESystem () override
 
void assemble (std::vector< GlobalVector * > const &x_new_timestep, std::vector< GlobalVector * > const &x_prev, int const process_id) override
 
void getA (GlobalMatrix &A) const override
 
void getRhs (GlobalVector const &x_prev, GlobalVector &rhs) const override
 
void computeKnownSolutions (GlobalVector const &x, int const process_id) override
 
void applyKnownSolutions (GlobalVector &x) const override
 
void applyKnownSolutionsPicard (GlobalMatrix &A, GlobalVector &rhs, GlobalVector &x) const override
 
bool isLinear () const override
 
void preIteration (const unsigned iter, GlobalVector const &x) override
 
IterationResult postIteration (GlobalVector const &x) override
 
TimeDiscgetTimeDiscretization () override
 Exposes the used time discretization scheme. More...
 
MathLib::MatrixSpecifications getMatrixSpecifications (const int process_id) const override
 

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 * _M
 Matrix \( M \). More...
 
GlobalMatrix * _K
 Matrix \( K \). More...
 
GlobalVector * _b
 Matrix \( b \). 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...
 

Member Typedef Documentation

◆ Index

◆ MatTrans

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

Definition at line 167 of file TimeDiscretizedODESystem.h.

◆ ODE

◆ TimeDisc

Constructor & Destructor Documentation

◆ TimeDiscretizedODESystem()

Constructs a new instance.

Definition at line 187 of file TimeDiscretizedODESystem.cpp.

190  : _ode(ode),
191  _time_disc(time_discretization),
192  _mat_trans(createMatrixTranslator<ODETag>(time_discretization))
193 {
195  ode.getMatrixSpecifications(process_id), _M_id);
197  ode.getMatrixSpecifications(process_id), _K_id);
199  ode.getMatrixSpecifications(process_id), _b_id);
200 }
virtual GlobalMatrix & getMatrix()=0
Get an uninitialized matrix.
std::unique_ptr< MatTrans > _mat_trans
the object used to compute the matrix/vector for the nonlinear solver
virtual GlobalVector & getVector()=0
Get an uninitialized vector.
static NUMLIB_EXPORT MatrixProvider & provider
static NUMLIB_EXPORT VectorProvider & provider

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

◆ ~TimeDiscretizedODESystem()

Member Function Documentation

◆ applyKnownSolutions()

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

Definition at line 257 of file TimeDiscretizedODESystem.cpp.

258 {
260 }
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.

References detail::applyKnownSolutions().

◆ applyKnownSolutionsPicard()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >::applyKnownSolutionsPicard ( GlobalMatrix &  A,
GlobalVector &  rhs,
GlobalVector &  x 
) const
override

Definition at line 263 of file TimeDiscretizedODESystem.cpp.

267 {
268  if (!_known_solutions)
269  {
270  return;
271  }
272 
274  std::vector<IndexType> ids;
275  std::vector<double> values;
276  for (auto const& bc : *_known_solutions)
277  {
278  std::copy(bc.ids.cbegin(), bc.ids.cend(), std::back_inserter(ids));
279  std::copy(bc.values.cbegin(), bc.values.cend(),
280  std::back_inserter(values));
281  }
282  MathLib::applyKnownSolution(A, rhs, x, ids, values);
283 }
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:37
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

References MathLib::applyKnownSolution(), and MathLib::LinAlg::copy().

◆ assemble()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >::assemble ( std::vector< GlobalVector * > const &  x_new_timestep,
std::vector< GlobalVector * > const &  x_prev,
int const  process_id 
)
override

Definition at line 212 of file TimeDiscretizedODESystem.cpp.

216 {
217  namespace LinAlg = MathLib::LinAlg;
218 
219  auto const t = _time_disc.getCurrentTime();
220  auto const dt = _time_disc.getCurrentTimeIncrement();
221  auto const& x_curr = *x_new_timestep[process_id];
222  std::vector<GlobalVector*> xdot(x_new_timestep.size());
223  for (std::size_t i = 0; i < xdot.size(); i++)
224  {
226  _time_disc.getXdot(*x_new_timestep[i], *x_prev[i], *xdot[i]);
227  }
228 
229  _M->setZero();
230  _K->setZero();
231  _b->setZero();
232 
233  _ode.preAssemble(t, dt, x_curr);
234  _ode.assemble(t, dt, x_new_timestep, xdot, process_id, *_M, *_K, *_b);
235 
239 
240  for (auto& v : xdot)
241  {
243  }
244 }
void getXdot(GlobalVector const &x_at_new_timestep, GlobalVector const &x_old, GlobalVector &xdot) const
virtual double getCurrentTimeIncrement() const =0
virtual double getCurrentTime() const =0
void finalizeAssembly(Matrix &)

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

◆ computeKnownSolutions()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >::computeKnownSolutions ( GlobalVector const &  x,
int const  process_id 
)
override

Definition at line 248 of file TimeDiscretizedODESystem.cpp.

250 {
252  _ode.getKnownSolutions(_time_disc.getCurrentTime(), x, process_id);
253 }

◆ getA()

Definition at line 181 of file TimeDiscretizedODESystem.h.

182  {
183  _mat_trans->computeA(*_M, *_K, A);
184  }

◆ getMatrixSpecifications()

Definition at line 212 of file TimeDiscretizedODESystem.h.

214  {
215  return _ode.getMatrixSpecifications(process_id);
216  }

◆ getRhs()

void NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >::getRhs ( GlobalVector const &  x_prev,
GlobalVector &  rhs 
) const
inlineoverride

Definition at line 186 of file TimeDiscretizedODESystem.h.

187  {
188  _mat_trans->computeRhs(*_M, *_K, *_b, x_prev, rhs);
189  }

◆ getTimeDiscretization()

Exposes the used time discretization scheme.

Implements NumLib::TimeDiscretizedODESystemBase< NonlinearSolverTag::Picard >.

Definition at line 211 of file TimeDiscretizedODESystem.h.

211 { return _time_disc; }

◆ isLinear()

Definition at line 199 of file TimeDiscretizedODESystem.h.

199 { return _ode.isLinear(); }

◆ postIteration()

Definition at line 206 of file TimeDiscretizedODESystem.h.

207  {
208  return _ode.postIteration(x);
209  }

◆ preIteration()

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

Definition at line 201 of file TimeDiscretizedODESystem.h.

202  {
203  _ode.preIteration(iter, x);
204  }

Member Data Documentation

◆ _b

◆ _b_id

ID of the _b vector.

Definition at line 235 of file TimeDiscretizedODESystem.h.

◆ _K

◆ _K_id

ID of the _K matrix.

Definition at line 234 of file TimeDiscretizedODESystem.h.

◆ _known_solutions

Initial value:
=
nullptr

stores precomputed values for known solutions

Definition at line 226 of file TimeDiscretizedODESystem.h.

◆ _M

◆ _M_id

ID of the _M matrix.

Definition at line 233 of file TimeDiscretizedODESystem.h.

◆ _mat_trans

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

Definition at line 223 of file TimeDiscretizedODESystem.h.

◆ _ode

◆ _time_disc

the time discretization to being used

Definition at line 220 of file TimeDiscretizedODESystem.h.

◆ ODETag


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