Merge branch 'feature/hyperelasticity' into 'develop'

Hyperelastic Lagrangian solver (Eucclhyd or Glace-type).

See merge request !163

src/language/modules/SchemeModule.cpp

@@ -32,6 +32,7 @@
 #include <scheme/FixedBoundaryConditionDescriptor.hpp>
 #include <scheme/FourierBoundaryConditionDescriptor.hpp>
 #include <scheme/FreeBoundaryConditionDescriptor.hpp>
+#include <scheme/HyperelasticSolver.hpp>
 #include <scheme/IBoundaryConditionDescriptor.hpp>
 #include <scheme/IDiscreteFunctionDescriptor.hpp>
 #include <scheme/NeumannBoundaryConditionDescriptor.hpp>
@@ -276,6 +277,18 @@ SchemeModule::SchemeModule()
 
                                           ));
 
+  this->_addBuiltinFunction("normalstress",
+                            std::function(
+
+                              [](std::shared_ptr<const IBoundaryDescriptor> boundary,
+                                 const FunctionSymbolId& normal_stress_id)
+                                -> std::shared_ptr<const IBoundaryConditionDescriptor> {
+                                return std::make_shared<DirichletBoundaryConditionDescriptor>("normal-stress", boundary,
+                                                                                              normal_stress_id);
+                              }
+
+                              ));
+
   this->_addBuiltinFunction("velocity", std::function(
 
                                           [](std::shared_ptr<const IBoundaryDescriptor> boundary,
@@ -400,6 +413,117 @@ SchemeModule::SchemeModule()
 
                               ));
 
+  this->_addBuiltinFunction("hyperelastic_eucclhyd_fluxes",
+                            std::function(
+
+                              [](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aL,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aT,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& u,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list)
+                                -> std::tuple<std::shared_ptr<const ItemValueVariant>,
+                                              std::shared_ptr<const SubItemValuePerItemVariant>> {
+                                return HyperelasticSolverHandler{getCommonMesh({rho, aL, aT, u, sigma})}
+                                  .solver()
+                                  .compute_fluxes(HyperelasticSolverHandler::SolverType::Eucclhyd, rho, aL, aT, u,
+                                                  sigma, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("hyperelastic_eucclhyd_solver",
+                            std::function(
+
+                              [](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& u,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& E,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& CG,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aL,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aT,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const IMesh>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return HyperelasticSolverHandler{getCommonMesh({rho, u, E, CG, aL, aT, sigma})}
+                                  .solver()
+                                  .apply(HyperelasticSolverHandler::SolverType::Eucclhyd, dt, rho, u, E, CG, aL, aT,
+                                         sigma, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("hyperelastic_glace_fluxes",
+                            std::function(
+
+                              [](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aL,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aT,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& u,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list)
+                                -> std::tuple<std::shared_ptr<const ItemValueVariant>,
+                                              std::shared_ptr<const SubItemValuePerItemVariant>> {
+                                return HyperelasticSolverHandler{getCommonMesh({rho, aL, aT, u, sigma})}
+                                  .solver()
+                                  .compute_fluxes(HyperelasticSolverHandler::SolverType::Glace,   //
+                                                  rho, aL, aT, u, sigma, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("hyperelastic_glace_solver",
+                            std::function(
+
+                              [](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& u,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& E,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& CG,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aL,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& aT,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const IMesh>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return HyperelasticSolverHandler{getCommonMesh({rho, u, E, CG, aL, aT, sigma})}
+                                  .solver()
+                                  .apply(HyperelasticSolverHandler::SolverType::Glace, dt, rho, u, E, CG, aL, aT, sigma,
+                                         bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("apply_hyperelastic_fluxes",
+                            std::function(
+
+                              [](const std::shared_ptr<const DiscreteFunctionVariant>& rho,      //
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& u,        //
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& E,        //
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& CG,       //
+                                 const std::shared_ptr<const ItemValueVariant>& ur,              //
+                                 const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr,   //
+                                 const double& dt) -> std::tuple<std::shared_ptr<const IMesh>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return HyperelasticSolverHandler{getCommonMesh({rho, u, E, CG})}   //
+                                  .solver()
+                                  .apply_fluxes(dt, rho, u, E, CG, ur, Fjr);
+                              }
+
+                              ));
+
   this->_addBuiltinFunction("lagrangian",
                             std::function(
 
@@ -458,6 +582,14 @@ SchemeModule::SchemeModule()
 
                               ));
 
+  this->_addBuiltinFunction("hyperelastic_dt", std::function(
+
+                                                 [](const std::shared_ptr<const DiscreteFunctionVariant>& c) -> double {
+                                                   return hyperelastic_dt(c);
+                                                 }
+
+                                                 ));
+
   MathFunctionRegisterForVh{*this};
 }
 

src/scheme/AcousticSolver.cpp

@@ -49,7 +49,7 @@ acoustic_dt(const std::shared_ptr<const DiscreteFunctionVariant>& c)
     return acoustic_dt(c->get<DiscreteFunctionP0<2, const double>>());
   }
   case 3: {
-    return acoustic_dt(c->get<DiscreteFunctionP0<2, const double>>());
+    return acoustic_dt(c->get<DiscreteFunctionP0<3, const double>>());
   }
   default: {
     throw UnexpectedError("invalid mesh dimension");

src/scheme/CMakeLists.txt

@@ -3,6 +3,7 @@
 add_library(
   PugsScheme
   AcousticSolver.cpp
+  HyperelasticSolver.cpp
   DiscreteFunctionIntegrator.cpp
   DiscreteFunctionInterpoler.cpp
   DiscreteFunctionUtils.cpp

src/scheme/HyperelasticSolver.hpp

+#ifndef HYPERELASTIC_SOLVER_HPP
+#define HYPERELASTIC_SOLVER_HPP
+
+#include <memory>
+#include <tuple>
+#include <vector>
+
+class IBoundaryConditionDescriptor;
+class IMesh;
+class ItemValueVariant;
+class SubItemValuePerItemVariant;
+class DiscreteFunctionVariant;
+
+double hyperelastic_dt(const std::shared_ptr<const DiscreteFunctionVariant>& c);
+
+class HyperelasticSolverHandler
+{
+ public:
+  enum class SolverType
+  {
+    Glace,
+    Eucclhyd
+  };
+
+ private:
+  struct IHyperelasticSolver
+  {
+    virtual std::tuple<const std::shared_ptr<const ItemValueVariant>,
+                       const std::shared_ptr<const SubItemValuePerItemVariant>>
+    compute_fluxes(
+      const SolverType& solver_type,
+      const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+      const std::shared_ptr<const DiscreteFunctionVariant>& aL,
+      const std::shared_ptr<const DiscreteFunctionVariant>& aT,
+      const std::shared_ptr<const DiscreteFunctionVariant>& u,
+      const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
+      const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) const = 0;
+
+    virtual std::tuple<std::shared_ptr<const IMesh>,
+                       std::shared_ptr<const DiscreteFunctionVariant>,
+                       std::shared_ptr<const DiscreteFunctionVariant>,
+                       std::shared_ptr<const DiscreteFunctionVariant>,
+                       std::shared_ptr<const DiscreteFunctionVariant>>
+    apply_fluxes(const double& dt,
+                 const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+                 const std::shared_ptr<const DiscreteFunctionVariant>& u,
+                 const std::shared_ptr<const DiscreteFunctionVariant>& E,
+                 const std::shared_ptr<const DiscreteFunctionVariant>& CG,
+                 const std::shared_ptr<const ItemValueVariant>& ur,
+                 const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr) const = 0;
+
+    virtual std::tuple<std::shared_ptr<const IMesh>,
+                       std::shared_ptr<const DiscreteFunctionVariant>,
+                       std::shared_ptr<const DiscreteFunctionVariant>,
+                       std::shared_ptr<const DiscreteFunctionVariant>,
+                       std::shared_ptr<const DiscreteFunctionVariant>>
+    apply(const SolverType& solver_type,
+          const double& dt,
+          const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+          const std::shared_ptr<const DiscreteFunctionVariant>& u,
+          const std::shared_ptr<const DiscreteFunctionVariant>& E,
+          const std::shared_ptr<const DiscreteFunctionVariant>& CG,
+          const std::shared_ptr<const DiscreteFunctionVariant>& aL,
+          const std::shared_ptr<const DiscreteFunctionVariant>& aT,
+          const std::shared_ptr<const DiscreteFunctionVariant>& p,
+          const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) const = 0;
+
+    IHyperelasticSolver()                                 = default;
+    IHyperelasticSolver(IHyperelasticSolver&&)            = default;
+    IHyperelasticSolver& operator=(IHyperelasticSolver&&) = default;
+
+    virtual ~IHyperelasticSolver() = default;
+  };
+
+  template <size_t Dimension>
+  class HyperelasticSolver;
+
+  std::unique_ptr<IHyperelasticSolver> m_hyperelastic_solver;
+
+ public:
+  const IHyperelasticSolver&
+  solver() const
+  {
+    return *m_hyperelastic_solver;
+  }
+
+  HyperelasticSolverHandler(const std::shared_ptr<const IMesh>& mesh);
+};
+
+#endif   // HYPERELASTIC_SOLVER_HPP