Merge pull request #8 from aurora-multiphysics/EdwardPalmer99/impleme…

…nt-fespaces-and-gridfunctions

Implements MFEMFESpaces, MFEMVariable, MFEMFECollections

include/fespaces/MFEMFECollection.h

@@ -0,0 +1,36 @@
+#pragma once
+#include "MFEMGeneralUserObject.h"
+#include "mfem.hpp"
+#include <string>
+
+/**
+ * Constructs and stores an mfem::FiniteElementCollection object. Access using
+ * the getFEC() accessor.
+ */
+class MFEMFECollection : public MFEMGeneralUserObject
+{
+public:
+  static InputParameters validParams();
+
+  MFEMFECollection(const InputParameters & parameters);
+
+  ~MFEMFECollection() override = default;
+
+  /// Returns a shared pointer to the constructed fec.
+  inline std::shared_ptr<mfem::FiniteElementCollection> getFEC() const { return _fec; }
+
+protected:
+  const int _fec_order;
+  const std::string _fec_type;
+  const std::string _fec_name;
+
+private:
+  /// Constructs the fec name from the order and fespace type.
+  const std::string buildFECName();
+
+  /// Constructs the fec from the fec name.
+  const std::shared_ptr<mfem::FiniteElementCollection> buildFEC();
+
+  /// Stores the constructed fec.
+  const std::shared_ptr<mfem::FiniteElementCollection> _fec{nullptr};
+};

include/fespaces/MFEMFESpace.h

@@ -1,17 +1,38 @@
 #pragma once
 #include "MFEMGeneralUserObject.h"
+#include "MFEMFECollection.h"
 
+/**
+ * Constructs and stores an mfem::ParFiniteElementSpace object. Access using the
+ * getFESpace() accessor.
+ */
 class MFEMFESpace : public MFEMGeneralUserObject
 {
 public:
   static InputParameters validParams();
 
   MFEMFESpace(const InputParameters & parameters);
-  virtual ~MFEMFESpace();
 
-  static const std::string createFECName(const std::string & fespace_type, const int order);
-  const int order;
-  const int vdim;
-  const std::string fespace_type;
-  const std::string fec_name;
+  /// Returns a shared pointer to the constructed fespace.
+  inline std::shared_ptr<mfem::ParFiniteElementSpace> getFESpace() const { return _fespace; }
+
+  /// Returns a shared pointer to the constructed fec.
+  inline std::shared_ptr<mfem::FiniteElementCollection> getFEC() const { return _fec.getFEC(); }
+
+protected:
+  /// Vector dimension (number of unknowns per degree of freedom).
+  const int _vdim;
+
+  /// Type of ordering of the vector dofs when _vdim > 1.
+  const int _ordering;
+
+  /// Constructs and stores the fec.
+  const MFEMFECollection _fec;
+
+private:
+  /// Constructs the fespace.
+  const std::shared_ptr<mfem::ParFiniteElementSpace> buildFESpace();
+
+  /// Stores the constructed fespace.
+  const std::shared_ptr<mfem::ParFiniteElementSpace> _fespace{nullptr};
 };

include/problem/MFEMProblem.h

@@ -179,6 +179,6 @@ class MFEMProblem : public ExternalProblem
 
   std::shared_ptr<platypus::ProblemBuilder> mfem_problem_builder{nullptr};
 
-  std::unique_ptr<platypus::Problem> mfem_problem{nullptr};
+  std::shared_ptr<platypus::Problem> mfem_problem{nullptr};
   std::unique_ptr<platypus::Executioner> executioner{nullptr};
 };

include/problem_builders/problem_builder_base.h

@@ -52,13 +52,7 @@ class ProblemBuilder
   ProblemBuilder() = delete;
 
   // Virtual destructor required to prevent leaks.
-  virtual ~ProblemBuilder()
-  {
-    if (_problem != nullptr)
-    {
-      delete _problem;
-    }
-  }
+  virtual ~ProblemBuilder() = default;
 
   void SetMesh(std::shared_ptr<mfem::ParMesh> pmesh);
   void SetFESpaces(platypus::FESpaces & fespaces);
@@ -99,9 +93,15 @@ class ProblemBuilder
   /// operator has already been constructed earlier to avoid rebuilding it.
   void FinalizeProblem(bool build_operator = true);
 
+  /// Returns a shared pointer to the problem.
+  std::shared_ptr<platypus::Problem> ReturnProblem() { return _problem; }
+
 protected:
   /// Protected constructor. Derived classes must call this constructor.
-  ProblemBuilder(platypus::Problem * problem) : _problem{problem} {}
+  ProblemBuilder(platypus::Problem * problem)
+    : _problem(std::shared_ptr<platypus::Problem>(problem))
+  {
+  }
 
   /// Supported Jacobian solver types.
   enum class SolverType
@@ -148,25 +148,14 @@ class ProblemBuilder
       MFEM_ABORT("platypus::Problem instance is NULL.");
     }
 
-    return static_cast<TDerivedProblem *>(_problem);
-  }
-
-  /// Helper template for returning a unique pointer of the correct class.
-  /// Sets _problem = nullptr to avoid double-free.
-  template <class TDerivedProblem>
-  [[nodiscard]] std::unique_ptr<TDerivedProblem> ReturnProblem()
-  {
-    auto * problem = GetProblem<TDerivedProblem>();
-    _problem = nullptr; // Avoid double-free.
-
-    return std::unique_ptr<TDerivedProblem>(problem);
+    return static_cast<TDerivedProblem *>(_problem.get());
   }
 
   /// Coefficient used in some derived classes.
   mfem::ConstantCoefficient _one_coef{1.0};
 
 private:
-  platypus::Problem * _problem{nullptr};
+  std::shared_ptr<platypus::Problem> _problem{nullptr};
 };
 
 /// Interface for problem builders that are constructing problems with an equation system.

include/problem_builders/steady_state_equation_system_problem_builder.h

@@ -56,8 +56,6 @@ class SteadyStateEquationSystemProblemBuilder : public SteadyStateProblemBuilder
   /// equation system is initialized.
   void InitializeKernels() final;
 
-  auto ReturnProblem() { return ProblemBuilder::ReturnProblem<SteadyStateEquationSystemProblem>(); }
-
 protected:
   [[nodiscard]] platypus::SteadyStateEquationSystemProblem * GetProblem() const override
   {

include/problem_builders/steady_state_problem_builder.h

@@ -44,8 +44,6 @@ class SteadyStateProblemBuilder : public ProblemBuilder
 
   ~SteadyStateProblemBuilder() override = default;
 
-  auto ReturnProblem() { return ProblemBuilder::ReturnProblem<SteadyStateProblem>(); }
-
   void RegisterFESpaces() override {}
 
   void RegisterGridFunctions() override {}

include/problem_builders/time_domain_equation_system_problem_builder.h

@@ -56,8 +56,6 @@ class TimeDomainEquationSystemProblemBuilder : public TimeDomainProblemBuilder,
   /// NB: - note use of final. Ensure that the equation system is initialized.
   void InitializeKernels() final;
 
-  auto ReturnProblem() { return ProblemBuilder::ReturnProblem<TimeDomainEquationSystemProblem>(); }
-
 protected:
   [[nodiscard]] platypus::TimeDomainEquationSystemProblem * GetProblem() const override
   {

include/problem_builders/time_domain_problem_builder.h

@@ -43,8 +43,6 @@ class TimeDomainProblemBuilder : public ProblemBuilder
 
   ~TimeDomainProblemBuilder() override = default;
 
-  auto ReturnProblem() { return ProblemBuilder::ReturnProblem<TimeDomainProblem>(); }
-
   static std::vector<mfem::ParGridFunction *>
   RegisterTimeDerivatives(std::vector<std::string> gridfunction_names,
                           platypus::GridFunctions & gridfunctions);

include/variables/MFEMVariable.h

@@ -3,14 +3,26 @@
 #include "MFEMFESpace.h"
 #include "MFEMGeneralUserObject.h"
 
+/**
+ * Constructs and stores an mfem::ParGridFunction object.
+ */
 class MFEMVariable : public MFEMGeneralUserObject
 {
 public:
   static InputParameters validParams();
 
   MFEMVariable(const InputParameters & parameters);
-  virtual ~MFEMVariable();
 
-  const MFEMFESpace & fespace;
-  unsigned int components;
+  /// Returns a shared pointer to the constructed gridfunction.
+  inline std::shared_ptr<mfem::ParGridFunction> getGridFunction() const { return _gridfunction; }
+
+protected:
+  const MFEMFESpace & _fespace;
+
+private:
+  /// Constructs the gridfunction.
+  const std::shared_ptr<mfem::ParGridFunction> buildGridFunction();
+
+  /// Stores the constructed gridfunction.
+  const std::shared_ptr<mfem::ParGridFunction> _gridfunction{nullptr};
 };

src/base/PlatypusApp.C

@@ -47,6 +47,7 @@ associateSyntaxInner(Syntax & syntax, ActionFactory & /*action_factory*/)
 
   // add FESpaces
   registerMooseObjectTask("add_mfem_fespaces", MFEMFESpace, false);
+  appendMooseObjectTask("add_mfem_fespaces", MFEMFECollection);
   registerSyntaxTask("AddFESpaceAction", "FESpaces/*", "add_mfem_fespaces");
   addTaskDependency("add_variable", "add_mfem_fespaces");
   addTaskDependency("add_aux_variable", "add_mfem_fespaces");

src/fespaces/MFEMFECollection.C

@@ -0,0 +1,49 @@
+#include "MFEMFECollection.h"
+
+registerMooseObject("PlatypusApp", MFEMFECollection);
+
+InputParameters
+MFEMFECollection::validParams()
+{
+  InputParameters params = MFEMGeneralUserObject::validParams();
+  params.registerBase("MFEMFECollection");
+
+  MooseEnum fec_order(
+      "CONSTANT FIRST SECOND THIRD FOURTH FIFTH SIXTH SEVENTH EIGHTH NINTH TENTH ELEVENTH TWELFTH"
+      "THIRTEENTH FOURTEENTH FIFTEENTH SIXTEENTH SEVENTEENTH EIGHTTEENTH NINETEENTH TWENTIETH "
+      "TWENTYFIRST TWENTYSECOND TWENTYTHIRD TWENTYFOURTH TWENTYFIFTH TWENTYSIXTH TWENTYSEVENTH "
+      "TWENTYEIGHTH TWENTYNINTH THIRTIETH THIRTYFIRST THIRTYSECOND THIRTYTHIRD THIRTYFOURTH "
+      "THIRTYFIFTH THIRTYSIXTH THIRTYSEVENTH THIRTYEIGHTH THIRTYNINTH FORTIETH FORTYFIRST"
+      "FORTYSECOND FORTYTHIRD",
+      "FIRST",
+      true);
+  params.addParam<MooseEnum>("fec_order", fec_order, "Order of the FE shape function to use.");
+
+  MooseEnum fec_types("H1 ND RT L2", "H1", true);
+  params.addParam<MooseEnum>(
+      "fec_type", fec_types, "Specifies the family of FE shape functions (FE space).");
+
+  return params;
+}
+
+MFEMFECollection::MFEMFECollection(const InputParameters & parameters)
+  : MFEMGeneralUserObject(parameters),
+    _fec_order(parameters.get<MooseEnum>("fec_order")),
+    _fec_type(parameters.get<MooseEnum>("fec_type")),
+    _fec_name(buildFECName()),
+    _fec(buildFEC())
+{
+}
+
+const std::string
+MFEMFECollection::buildFECName()
+{
+  return _fec_type + "_3D_P" + std::to_string(_fec_order);
+}
+
+const std::shared_ptr<mfem::FiniteElementCollection>
+MFEMFECollection::buildFEC()
+{
+  auto * fec_ptr = mfem::FiniteElementCollection::New(_fec_name.c_str());
+  return std::shared_ptr<mfem::FiniteElementCollection>(fec_ptr);
+}

src/fespaces/MFEMFESpace.C

@@ -1,4 +1,5 @@
 #include "MFEMFESpace.h"
+#include "MFEMProblem.h"
 
 registerMooseObject("PlatypusApp", MFEMFESpace);
 
@@ -7,42 +8,31 @@ MFEMFESpace::validParams()
 {
   InputParameters params = MFEMGeneralUserObject::validParams();
   params.registerBase("MFEMFESpace");
-  MooseEnum fespace_types("H1 ND RT L2", "H1", true);
-  params.addParam<MooseEnum>(
-      "fespace_type",
-      fespace_types,
-      "Specifies the family of FE shape functions (FE space) to use for this variable.");
-  MooseEnum order(
-      "CONSTANT FIRST SECOND THIRD FOURTH FIFTH SIXTH SEVENTH EIGHTH NINTH TENTH ELEVENTH TWELFTH"
-      "THIRTEENTH FOURTEENTH FIFTEENTH SIXTEENTH SEVENTEENTH EIGHTTEENTH NINETEENTH TWENTIETH "
-      "TWENTYFIRST TWENTYSECOND TWENTYTHIRD TWENTYFOURTH TWENTYFIFTH TWENTYSIXTH TWENTYSEVENTH "
-      "TWENTYEIGHTH TWENTYNINTH THIRTIETH THIRTYFIRST THIRTYSECOND THIRTYTHIRD THIRTYFOURTH "
-      "THIRTYFIFTH THIRTYSIXTH THIRTYSEVENTH THIRTYEIGHTH THIRTYNINTH FORTIETH FORTYFIRST"
-      "FORTYSECOND FORTYTHIRD",
-      "FIRST",
-      true);
-  params.addParam<MooseEnum>("order",
-                             order,
-                             "Order of the FE shape function to use for this variable (additional"
-                             "orders not listed here are allowed,"
-                             "depending on the family.");
+
   params.addParam<int>("vdim", 1, "Dimension of degrees of freedom");
+
+  MooseEnum ordering("NODES VDIM", "NODES", false);
+  params.addParam<MooseEnum>("ordering", ordering, "Type of ordering of the vector dofs.");
+
+  params += MFEMFECollection::validParams();
+
   return params;
 }
 
 MFEMFESpace::MFEMFESpace(const InputParameters & parameters)
   : MFEMGeneralUserObject(parameters),
-    order(parameters.get<MooseEnum>("order")),
-    vdim(parameters.get<int>("vdim")),
-    fespace_type(parameters.get<MooseEnum>("fespace_type")),
-    fec_name(createFECName(fespace_type, order))
+    _vdim(parameters.get<int>("vdim")),
+    _ordering(parameters.get<MooseEnum>("ordering")),
+    _fec(parameters),
+    _fespace(buildFESpace())
 {
 }
 
-const std::string
-MFEMFESpace::createFECName(const std::string & fespace_type, const int order)
+const std::shared_ptr<mfem::ParFiniteElementSpace>
+MFEMFESpace::buildFESpace()
 {
-  return fespace_type + "_3D_P" + std::to_string(order);
-}
+  mfem::ParMesh & pmesh = getMFEMProblem().mesh().getMFEMParMesh();
+  mfem::FiniteElementCollection & fec = *_fec.getFEC();
 
-MFEMFESpace::~MFEMFESpace() {}
+  return std::make_shared<mfem::ParFiniteElementSpace>(&pmesh, &fec, _vdim, _ordering);
+}