https://mooseframework.inl.gov
EquationSystem.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #ifdef MOOSE_MFEM_ENABLED
11 
12 #pragma once
13 #include "libmesh/ignore_warnings.h"
14 #include "mfem/miniapps/common/pfem_extras.hpp"
15 #include "libmesh/restore_warnings.h"
16 #include "MFEMIntegratedBC.h"
17 #include "MFEMEssentialBC.h"
18 #include "MFEMContainers.h"
19 #include "MFEMKernel.h"
21 #include "ScaleIntegrator.h"
22 
23 namespace Moose::MFEM
24 {
25 
30 class EquationSystem : public mfem::Operator
31 {
32 
33 public:
36 
37  EquationSystem() = default;
38  ~EquationSystem() override;
39 
41  virtual void AddTestVariableNameIfMissing(const std::string & test_var_name);
43  virtual void AddTrialVariableNameIfMissing(const std::string & trial_var_name);
44 
46  virtual void AddKernel(std::shared_ptr<MFEMKernel> kernel);
47  virtual void AddIntegratedBC(std::shared_ptr<MFEMIntegratedBC> kernel);
48  virtual void AddEssentialBC(std::shared_ptr<MFEMEssentialBC> bc);
49  virtual void ApplyEssentialBCs();
50 
52  virtual void Init(Moose::MFEM::GridFunctions & gridfunctions,
53  const Moose::MFEM::FESpaces & fespaces,
54  mfem::AssemblyLevel assembly_level);
55  virtual void BuildLinearForms();
56  virtual void BuildBilinearForms();
57  virtual void BuildMixedBilinearForms();
58  virtual void BuildEquationSystem();
59 
61  virtual void FormLinearSystem(mfem::OperatorHandle & op,
62  mfem::BlockVector & trueX,
63  mfem::BlockVector & trueRHS);
64 
65  virtual void
66  FormSystem(mfem::OperatorHandle & op, mfem::BlockVector & trueX, mfem::BlockVector & trueRHS);
67  virtual void FormLegacySystem(mfem::OperatorHandle & op,
68  mfem::BlockVector & trueX,
69  mfem::BlockVector & trueRHS);
70 
72  virtual void BuildJacobian(mfem::BlockVector & trueX, mfem::BlockVector & trueRHS);
73 
75  void Mult(const mfem::Vector & u, mfem::Vector & residual) const override;
76 
78  mfem::Operator & GetGradient(const mfem::Vector & u) const override;
79 
81  virtual void RecoverFEMSolution(mfem::BlockVector & trueX,
82  Moose::MFEM::GridFunctions & gridfunctions);
83 
84  std::vector<mfem::Array<int>> _ess_tdof_lists;
85 
86  const std::vector<std::string> & TrialVarNames() const { return _trial_var_names; }
87  const std::vector<std::string> & TestVarNames() const { return _test_var_names; }
88 
89 private:
91  using mfem::Operator::RecoverFEMSolution;
92 
93 protected:
97  void DeleteAllBlocks();
98 
99  bool VectorContainsName(const std::vector<std::string> & the_vector,
100  const std::string & name) const;
101 
102  // Test variables are associated with LinearForms,
103  // whereas trial variables are associated with gridfunctions.
104 
107  std::vector<std::string> _trial_var_names;
111  std::vector<std::string> _test_var_names;
113  std::vector<mfem::ParFiniteElementSpace *> _test_pfespaces;
114 
115  // Components of weak form. // Named according to test variable
120  _mblfs; // named according to trial variable
121 
126  template <class FormType>
128  const std::string & trial_var_name,
129  const std::string & test_var_name,
130  std::shared_ptr<FormType> form,
132  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMKernel>>>> & kernels_map);
133 
135  const std::string & test_var_name,
136  std::shared_ptr<mfem::ParLinearForm> form,
138  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMKernel>>>> & kernels_map);
139 
140  template <class FormType>
142  const std::string & trial_var_name,
143  const std::string & test_var_name,
144  std::shared_ptr<FormType> form,
146  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMIntegratedBC>>>> &
147  integrated_bc_map);
148 
150  const std::string & test_var_name,
151  std::shared_ptr<mfem::ParLinearForm> form,
153  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMIntegratedBC>>>> &
154  integrated_bc_map);
155 
157  std::vector<std::unique_ptr<mfem::ParGridFunction>> _xs;
158  std::vector<std::unique_ptr<mfem::ParGridFunction>> _dxdts;
159 
160  mfem::Array2D<const mfem::HypreParMatrix *> _h_blocks;
161 
174 
175  mutable mfem::OperatorHandle _jacobian;
176 
177  mfem::AssemblyLevel _assembly_level;
178 };
179 
180 template <class FormType>
181 void
183  const std::string & trial_var_name,
184  const std::string & test_var_name,
185  std::shared_ptr<FormType> form,
187  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMKernel>>>> & kernels_map)
188 {
189  if (kernels_map.Has(test_var_name) && kernels_map.Get(test_var_name)->Has(trial_var_name))
190  {
191  auto kernels = kernels_map.GetRef(test_var_name).GetRef(trial_var_name);
192  for (auto & kernel : kernels)
193  {
194  mfem::BilinearFormIntegrator * integ = kernel->createBFIntegrator();
195  if (integ != nullptr)
196  {
197  kernel->isSubdomainRestricted()
198  ? form->AddDomainIntegrator(std::move(integ), kernel->getSubdomainMarkers())
199  : form->AddDomainIntegrator(std::move(integ));
200  }
201  }
202  }
203 }
204 
205 inline void
207  const std::string & test_var_name,
208  std::shared_ptr<mfem::ParLinearForm> form,
210  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMKernel>>>> & kernels_map)
211 {
212  if (kernels_map.Has(test_var_name) && kernels_map.Get(test_var_name)->Has(test_var_name))
213  {
214  auto kernels = kernels_map.GetRef(test_var_name).GetRef(test_var_name);
215  for (auto & kernel : kernels)
216  {
217  mfem::LinearFormIntegrator * integ = kernel->createLFIntegrator();
218  if (integ != nullptr)
219  {
220  kernel->isSubdomainRestricted()
221  ? form->AddDomainIntegrator(std::move(integ), kernel->getSubdomainMarkers())
222  : form->AddDomainIntegrator(std::move(integ));
223  }
224  }
225  }
226 }
227 
228 template <class FormType>
229 void
231  const std::string & trial_var_name,
232  const std::string & test_var_name,
233  std::shared_ptr<FormType> form,
235  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMIntegratedBC>>>> &
236  integrated_bc_map)
237 {
238  if (integrated_bc_map.Has(test_var_name) &&
239  integrated_bc_map.Get(test_var_name)->Has(trial_var_name))
240  {
241  auto bcs = integrated_bc_map.GetRef(test_var_name).GetRef(trial_var_name);
242  for (auto & bc : bcs)
243  {
244  mfem::BilinearFormIntegrator * integ = bc->createBFIntegrator();
245  if (integ != nullptr)
246  {
247  bc->isBoundaryRestricted()
248  ? form->AddBoundaryIntegrator(std::move(integ), bc->getBoundaryMarkers())
249  : form->AddBoundaryIntegrator(std::move(integ));
250  }
251  }
252  }
253 }
254 
255 inline void
257  const std::string & test_var_name,
258  std::shared_ptr<mfem::ParLinearForm> form,
260  Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMIntegratedBC>>>> &
261  integrated_bc_map)
262 {
263  if (integrated_bc_map.Has(test_var_name) &&
264  integrated_bc_map.Get(test_var_name)->Has(test_var_name))
265  {
266  auto bcs = integrated_bc_map.GetRef(test_var_name).GetRef(test_var_name);
267  for (auto & bc : bcs)
268  {
269  mfem::LinearFormIntegrator * integ = bc->createLFIntegrator();
270  if (integ != nullptr)
271  {
272  bc->isBoundaryRestricted()
273  ? form->AddBoundaryIntegrator(std::move(integ), bc->getBoundaryMarkers())
274  : form->AddBoundaryIntegrator(std::move(integ));
275  }
276  }
277  }
278 }
279 
284 {
285 public:
287 
288  void AddTrialVariableNameIfMissing(const std::string & trial_var_name) override;
289 
290  virtual void SetTimeStep(double dt);
291  virtual void UpdateEquationSystem();
292 
293  virtual void AddKernel(std::shared_ptr<MFEMKernel> kernel) override;
294  virtual void BuildBilinearForms() override;
295  virtual void FormLegacySystem(mfem::OperatorHandle & op,
296  mfem::BlockVector & truedXdt,
297  mfem::BlockVector & trueRHS) override;
298  virtual void FormSystem(mfem::OperatorHandle & op,
299  mfem::BlockVector & truedXdt,
300  mfem::BlockVector & trueRHS) override;
301 
302  const std::vector<std::string> & TrialVarTimeDerivativeNames() const;
303 
304 protected:
306  mfem::ConstantCoefficient _dt_coef;
307  std::vector<std::string> _trial_var_time_derivative_names;
308 
313 };
314 
315 inline const std::vector<std::string> &
317 {
319 }
320 
321 } // namespace Moose::MFEM
322 
323 #endif
void ApplyDomainBLFIntegrators(const std::string &trial_var_name, const std::string &test_var_name, std::shared_ptr< FormType > form, Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel >>>> &kernels_map)
Template method for applying BilinearFormIntegrators on domains from kernels to a BilinearForm...
mfem::ConstantCoefficient _dt_coef
Coefficient for timestep scaling.
virtual void AddKernel(std::shared_ptr< MFEMKernel > kernel) override
Add kernels.
virtual void AddTestVariableNameIfMissing(const std::string &test_var_name)
Add test variable to EquationSystem.
virtual void AddKernel(std::shared_ptr< MFEMKernel > kernel)
Add kernels.
virtual void RecoverFEMSolution(mfem::BlockVector &trueX, Moose::MFEM::GridFunctions &gridfunctions)
Update variable from solution vector after solve.
virtual void BuildBilinearForms()
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel > > > > _td_kernels_map
const std::vector< std::string > & TestVarNames() const
Moose::MFEM::NamedFieldsMap< mfem::ParNonlinearForm > _nlfs
virtual void FormSystem(mfem::OperatorHandle &op, mfem::BlockVector &truedXdt, mfem::BlockVector &trueRHS) override
const std::vector< std::string > & TrialVarNames() const
void ApplyBoundaryLFIntegrators(const std::string &test_var_name, std::shared_ptr< mfem::ParLinearForm > form, Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMIntegratedBC >>>> &integrated_bc_map)
bool VectorContainsName(const std::vector< std::string > &the_vector, const std::string &name) const
std::vector< std::unique_ptr< mfem::ParGridFunction > > _xs
Gridfunctions for setting Dirichlet BCs.
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< mfem::ParMixedBilinearForm > > _mblfs
Class to store weak form components (bilinear and linear forms, and optionally mixed and nonlinear fo...
std::vector< std::unique_ptr< mfem::ParGridFunction > > _dxdts
virtual void AddEssentialBC(std::shared_ptr< MFEMEssentialBC > bc)
Problem operator for time-dependent problems with an equation system.
std::vector< mfem::Array< int > > _ess_tdof_lists
mfem::AssemblyLevel _assembly_level
void AddTrialVariableNameIfMissing(const std::string &trial_var_name) override
Add trial variable to EquationSystem.
Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMEssentialBC > > > _essential_bc_map
Arrays to store essential BCs to act on each component of weak form.
std::vector< mfem::ParFiniteElementSpace * > _test_pfespaces
Pointers to finite element spaces associated with test variables.
const std::vector< std::string > & TrialVarTimeDerivativeNames() const
Class to store weak form components for time dependent PDEs.
void Mult(const mfem::Vector &u, mfem::Vector &residual) const override
Compute residual y = Mu.
std::vector< std::string > _trial_var_names
Names of all variables corresponding to gridfunctions.
virtual void FormLegacySystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
mfem::OperatorHandle _jacobian
virtual void AddTrialVariableNameIfMissing(const std::string &trial_var_name)
Add trial variable to EquationSystem.
std::vector< std::string > _test_var_names
Names of all test variables corresponding to linear forms in this equation system.
Steady-state problem operator with an equation system.
virtual void BuildEquationSystem()
virtual void Init(Moose::MFEM::GridFunctions &gridfunctions, const Moose::MFEM::FESpaces &fespaces, mfem::AssemblyLevel assembly_level)
Build forms.
virtual void FormLegacySystem(mfem::OperatorHandle &op, mfem::BlockVector &truedXdt, mfem::BlockVector &trueRHS) override
virtual void BuildMixedBilinearForms()
mfem::Array2D< const mfem::HypreParMatrix * > _h_blocks
virtual void FormLinearSystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
Form linear system, with essential boundary conditions accounted for.
std::vector< std::string > _trial_var_time_derivative_names
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMIntegratedBC > > > > _integrated_bc_map
Arrays to store integrated BCs to act on each component of weak form.
void DeleteAllBlocks()
Deletes the HypreParMatrix associated with any pointer stored in _h_blocks, and then proceeds to dele...
Moose::MFEM::NamedFieldsMap< mfem::ParBilinearForm > _td_blfs
Container to store contributions to weak form of the form (F du/dt, v)
Moose::MFEM::GridFunctions _trial_variables
Pointers to trial variables.
virtual void BuildJacobian(mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
Build linear system, with essential boundary conditions accounted for.
void ApplyDomainLFIntegrators(const std::string &test_var_name, std::shared_ptr< mfem::ParLinearForm > form, Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel >>>> &kernels_map)
virtual void ApplyEssentialBCs()
Moose::MFEM::NamedFieldsMap< mfem::ParLinearForm > _lfs
virtual void BuildBilinearForms() override
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel > > > > _kernels_map
Arrays to store kernels to act on each component of weak form.
mfem::Operator & GetGradient(const mfem::Vector &u) const override
Compute J = M + grad_H(u)
virtual void FormSystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
Moose::MFEM::NamedFieldsMap< mfem::ParBilinearForm > _blfs
void ApplyBoundaryBLFIntegrators(const std::string &trial_var_name, const std::string &test_var_name, std::shared_ptr< FormType > form, Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMIntegratedBC >>>> &integrated_bc_map)
virtual void AddIntegratedBC(std::shared_ptr< MFEMIntegratedBC > kernel)