OGS 6.1.0-1721-g6382411ad
MatrixTranslator.h
Go to the documentation of this file.
1 
10 #pragma once
11 
12 #include <memory>
13 
14 #include "TimeDiscretization.h"
15 #include "Types.h"
16 
17 #include "MathLib/LinAlg/LinAlg.h"
18 
19 namespace NumLib
20 {
23 
30 template <ODESystemTag ODETag>
32 
38 template <>
40 {
41 public:
43  virtual void computeA(GlobalMatrix const& M, GlobalMatrix const& K,
44  GlobalMatrix& A) const = 0;
45 
47  virtual void computeRhs(const GlobalMatrix& M, const GlobalMatrix& K,
48  const GlobalVector& b, GlobalVector& rhs) const = 0;
49 
54  virtual void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
55  GlobalVector const& b,
56  GlobalVector const& x_new_timestep,
57  GlobalVector const& xdot,
58  GlobalVector& res) const = 0;
59 
61  virtual void computeJacobian(GlobalMatrix const& Jac_in,
62  GlobalMatrix& Jac_out) const = 0;
63 
71  virtual void pushMatrices(GlobalMatrix const& /*M*/,
72  GlobalMatrix const& /*K*/,
73  GlobalVector const& /*b*/)
74  {
75  }
76 
77  virtual ~MatrixTranslator() = default;
78 };
79 
85 template <ODESystemTag ODETag>
87 
96 template <>
98  : public MatrixTranslator<ODESystemTag::FirstOrderImplicitQuasilinear>
99 {
100 public:
106  : _time_disc(timeDisc)
107  {
108  }
109 
111  void computeA(GlobalMatrix const& M, GlobalMatrix const& K,
112  GlobalMatrix& A) const override;
113 
115  void computeRhs(const GlobalMatrix& M, const GlobalMatrix& /*K*/,
116  const GlobalVector& b, GlobalVector& rhs) const override;
117 
119  void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
120  GlobalVector const& b,
121  GlobalVector const& x_new_timestep,
122  GlobalVector const& xdot,
123  GlobalVector& res) const override;
124 
127  void computeJacobian(GlobalMatrix const& Jac_in,
128  GlobalMatrix& Jac_out) const override;
129 
130 private:
131  TimeDiscretization const&
133 
135  mutable std::size_t _tmp_id = 0u;
136 };
137 
142 template <ODESystemTag ODETag>
144 
153 template <>
155  : public MatrixTranslator<ODESystemTag::FirstOrderImplicitQuasilinear>
156 {
157 public:
163  : _fwd_euler(timeDisc)
164  {
165  }
166 
168  void computeA(GlobalMatrix const& M, GlobalMatrix const& /*K*/,
169  GlobalMatrix& A) const override;
170 
172  void computeRhs(const GlobalMatrix& M, const GlobalMatrix& K,
173  const GlobalVector& b, GlobalVector& rhs) const override;
174 
176  void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
177  GlobalVector const& b,
178  GlobalVector const& x_new_timestep,
179  GlobalVector const& xdot,
180  GlobalVector& res) const override;
181 
184  void computeJacobian(GlobalMatrix const& Jac_in,
185  GlobalMatrix& Jac_out) const override;
186 
187 private:
189 
191  mutable std::size_t _tmp_id = 0u;
192 };
193 
198 template <ODESystemTag ODETag>
200 
209 template <>
211  : public MatrixTranslator<ODESystemTag::FirstOrderImplicitQuasilinear>
212 {
213 public:
219  : _crank_nicolson(timeDisc),
220  _M_bar(NumLib::GlobalMatrixProvider::provider
221  .getMatrix()),
222  _b_bar(
223  NumLib::GlobalVectorProvider::provider.getVector())
224  {
225  }
226 
228  {
230  _M_bar);
232  _b_bar);
233  }
234 
237  void computeA(GlobalMatrix const& M, GlobalMatrix const& K,
238  GlobalMatrix& A) const override;
239 
242  void computeRhs(const GlobalMatrix& M, const GlobalMatrix& /*K*/,
243  const GlobalVector& b, GlobalVector& rhs) const override;
244 
247  void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
248  GlobalVector const& b,
249  GlobalVector const& x_new_timestep,
250  GlobalVector const& xdot,
251  GlobalVector& res) const override;
252 
259  void computeJacobian(GlobalMatrix const& Jac_in,
260  GlobalMatrix& Jac_out) const override;
261 
273  void pushMatrices(GlobalMatrix const& M, GlobalMatrix const& K,
274  GlobalVector const& b) override;
275 
276 private:
278 
279  GlobalMatrix&
281  GlobalVector& _b_bar;
283 
286  mutable std::size_t _tmp_id = 0u;
287 };
288 
291 template <ODESystemTag ODETag>
292 std::unique_ptr<MatrixTranslator<ODETag>> createMatrixTranslator(
293  TimeDiscretization const& timeDisc)
294 {
295  if (auto* fwd_euler = dynamic_cast<ForwardEuler const*>(&timeDisc))
296  {
297  return std::unique_ptr<MatrixTranslator<ODETag>>(
298  new MatrixTranslatorForwardEuler<ODETag>(*fwd_euler));
299  }
300  if (auto* crank = dynamic_cast<CrankNicolson const*>(&timeDisc))
301  {
302  return std::unique_ptr<MatrixTranslator<ODETag>>(
304  }
305 
306  return std::unique_ptr<MatrixTranslator<ODETag>>(
307  new MatrixTranslatorGeneral<ODETag>(timeDisc));
308 }
309 
311 }
Generalized Crank-Nicolson scheme.
TimeDiscretization const & _time_disc
the time discretization used.
virtual void releaseMatrix(GlobalMatrix const &A)=0
std::unique_ptr< MatrixTranslator< ODETag > > createMatrixTranslator(TimeDiscretization const &timeDisc)
virtual void pushMatrices(GlobalMatrix const &, GlobalMatrix const &, GlobalVector const &)
static NUMLIB_EXPORT MatrixProvider & provider
static NUMLIB_EXPORT VectorProvider & provider
ODESystemTag
Tag used to specify the type of ODE.
Definition: Types.h:25
Forward Euler scheme.
virtual void releaseVector(GlobalVector const &x)=0