10 #ifdef MOOSE_MFEM_ENABLED 13 #include "libmesh/ignore_warnings.h" 14 #include "mfem/miniapps/common/pfem_extras.hpp" 15 #include "libmesh/restore_warnings.h" 46 virtual void AddKernel(std::shared_ptr<MFEMKernel> kernel);
53 mfem::AssemblyLevel assembly_level);
67 mfem::BlockVector & trueX,
68 mfem::BlockVector & trueRHS);
71 FormSystem(mfem::OperatorHandle & op, mfem::BlockVector & trueX, mfem::BlockVector & trueRHS);
73 mfem::BlockVector & trueX,
74 mfem::BlockVector & trueRHS);
77 virtual void BuildJacobian(mfem::BlockVector & trueX, mfem::BlockVector & trueRHS);
80 void Mult(
const mfem::Vector & u, mfem::Vector & residual)
const override;
83 mfem::Operator &
GetGradient(
const mfem::Vector & u)
const override;
96 using mfem::Operator::RecoverFEMSolution;
108 const std::string & name)
const;
142 template <
class FormType>
144 const std::string & trial_var_name,
145 const std::string & test_var_name,
146 std::shared_ptr<FormType> form,
151 const std::string & test_var_name,
152 std::shared_ptr<mfem::ParLinearForm> form,
156 template <
class FormType>
158 const std::string & trial_var_name,
159 const std::string & test_var_name,
160 std::shared_ptr<FormType> form,
166 const std::string & test_var_name,
167 std::shared_ptr<mfem::ParLinearForm> form,
195 template <
class FormType>
198 const std::string & trial_var_name,
199 const std::string & test_var_name,
200 std::shared_ptr<FormType> form,
204 if (kernels_map.Has(test_var_name) && kernels_map.Get(test_var_name)->Has(trial_var_name))
206 auto kernels = kernels_map.GetRef(test_var_name).GetRef(trial_var_name);
207 for (
auto & kernel : kernels)
209 mfem::BilinearFormIntegrator * integ = kernel->createBFIntegrator();
210 if (integ !=
nullptr)
212 kernel->isSubdomainRestricted()
213 ? form->AddDomainIntegrator(std::move(integ), kernel->getSubdomainMarkers())
214 : form->AddDomainIntegrator(std::move(integ));
222 const std::string & test_var_name,
223 std::shared_ptr<mfem::ParLinearForm> form,
227 if (kernels_map.Has(test_var_name) && kernels_map.Get(test_var_name)->Has(test_var_name))
229 auto kernels = kernels_map.GetRef(test_var_name).GetRef(test_var_name);
230 for (
auto & kernel : kernels)
232 mfem::LinearFormIntegrator * integ = kernel->createLFIntegrator();
233 if (integ !=
nullptr)
235 kernel->isSubdomainRestricted()
236 ? form->AddDomainIntegrator(std::move(integ), kernel->getSubdomainMarkers())
237 : form->AddDomainIntegrator(std::move(integ));
243 template <
class FormType>
246 const std::string & trial_var_name,
247 const std::string & test_var_name,
248 std::shared_ptr<FormType> form,
253 if (integrated_bc_map.Has(test_var_name) &&
254 integrated_bc_map.Get(test_var_name)->Has(trial_var_name))
256 auto bcs = integrated_bc_map.GetRef(test_var_name).GetRef(trial_var_name);
257 for (
auto & bc : bcs)
259 mfem::BilinearFormIntegrator * integ = bc->createBFIntegrator();
260 if (integ !=
nullptr)
262 bc->isBoundaryRestricted()
263 ? form->AddBoundaryIntegrator(std::move(integ), bc->getBoundaryMarkers())
264 : form->AddBoundaryIntegrator(std::move(integ));
272 const std::string & test_var_name,
273 std::shared_ptr<mfem::ParLinearForm> form,
278 if (integrated_bc_map.Has(test_var_name) &&
279 integrated_bc_map.Get(test_var_name)->Has(test_var_name))
281 auto bcs = integrated_bc_map.GetRef(test_var_name).GetRef(test_var_name);
282 for (
auto & bc : bcs)
284 mfem::LinearFormIntegrator * integ = bc->createLFIntegrator();
285 if (integ !=
nullptr)
287 bc->isBoundaryRestricted()
288 ? form->AddBoundaryIntegrator(std::move(integ), bc->getBoundaryMarkers())
289 : form->AddBoundaryIntegrator(std::move(integ));
308 virtual void AddKernel(std::shared_ptr<MFEMKernel> kernel)
override;
311 mfem::BlockVector & truedXdt,
312 mfem::BlockVector & trueRHS)
override;
313 virtual void FormSystem(mfem::OperatorHandle & op,
314 mfem::BlockVector & truedXdt,
315 mfem::BlockVector & trueRHS)
override;
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 EliminateCoupledVariables()
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()
Build bilinear forms.
std::vector< mfem::ParFiniteElementSpace * > _coupled_pfespaces
Pointers to finite element spaces associated with coupled variables.
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::string > _eliminated_var_names
Names of all coupled variables without a corresponding test variable.
virtual void UpdateEquationSystem()
~EquationSystem() override
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...
TimeDependentEquationSystem()
Moose::MFEM::GridFunctions _eliminated_variables
Pointers to coupled variables not part of the reduced EquationSystem.
virtual void AddEssentialBC(std::shared_ptr< MFEMEssentialBC > bc)
Problem operator for time-dependent problems with an equation system.
void AddCoupledVariableNameIfMissing(const std::string &coupled_var_name) override
Add coupled variable to EquationSystem.
std::vector< mfem::Array< int > > _ess_tdof_lists
mfem::AssemblyLevel _assembly_level
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.
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
Subset of _coupled_var_names of all variables corresponding to gridfunctions with degrees of freedom ...
virtual void FormLegacySystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
virtual void SetTrialVariableNames() override
Set trial variable names from subset of coupled variables that have an associated test variable...
mfem::OperatorHandle _jacobian
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()
std::vector< std::unique_ptr< mfem::ParGridFunction > > _var_ess_constraints
Gridfunctions holding essential constraints from Dirichlet BCs.
virtual void Init(Moose::MFEM::GridFunctions &gridfunctions, const Moose::MFEM::FESpaces &fespaces, mfem::AssemblyLevel assembly_level)
Initialise.
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 SetTrialVariableNames()
Set trial variable names from subset of coupled variables that have an associated test variable...
virtual void AddCoupledVariableNameIfMissing(const std::string &coupled_var_name)
Add coupled variable to EquationSystem.
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 > _coupled_var_names
Names of all trial variables of kernels and boundary conditions added to this EquationSystem.
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)
virtual void BuildLinearForms()
Build linear forms and eliminate constrained DoFs.
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
Build bilinear forms.
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.
virtual void SetTimeStep(double dt)
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)