From 3f96e465cb85dc130118a48a3ce8f12d74e59745 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Fri, 14 Mar 2025 15:50:13 +0100
Subject: [PATCH] Create VariationalSchemeModule

---
 src/language/modules/CMakeLists.txt           |   1 +
 src/language/modules/SchemeModule.cpp         | 223 ----------------
 .../modules/VariationalSchemeModule.cpp       | 247 ++++++++++++++++++
 .../modules/VariationalSchemeModule.hpp       |  23 ++
 4 files changed, 271 insertions(+), 223 deletions(-)
 create mode 100644 src/language/modules/VariationalSchemeModule.cpp
 create mode 100644 src/language/modules/VariationalSchemeModule.hpp

diff --git a/src/language/modules/CMakeLists.txt b/src/language/modules/CMakeLists.txt
index 1ecbe1db5..0cc79a381 100644
--- a/src/language/modules/CMakeLists.txt
+++ b/src/language/modules/CMakeLists.txt
@@ -14,6 +14,7 @@ add_library(
   SchemeModule.cpp
   SocketModule.cpp
   UnaryOperatorRegisterForVh.cpp
+  VariationalSchemeModule.cpp
   WriterModule.cpp
 )
 
diff --git a/src/language/modules/SchemeModule.cpp b/src/language/modules/SchemeModule.cpp
index ba19d0314..5287c41d8 100644
--- a/src/language/modules/SchemeModule.cpp
+++ b/src/language/modules/SchemeModule.cpp
@@ -49,13 +49,8 @@
 #include <scheme/InflowListBoundaryConditionDescriptor.hpp>
 #include <scheme/NeumannBoundaryConditionDescriptor.hpp>
 #include <scheme/OutflowBoundaryConditionDescriptor.hpp>
-#include <scheme/P1P0AnalyticVariationalSolver.hpp>
-#include <scheme/P1P0VariationalSolver.hpp>
 #include <scheme/SymmetryBoundaryConditionDescriptor.hpp>
 #include <scheme/VariableBCDescriptor.hpp>
-#include <scheme/VariationalSolver.hpp>
-#include <scheme/VariationalSolverDevelopedReconstruction.hpp>
-#include <scheme/VariationalSolverO1.hpp>
 #include <scheme/WallBoundaryConditionDescriptor.hpp>
 #include <utils/Socket.hpp>
 
@@ -529,215 +524,6 @@ SchemeModule::SchemeModule()
 
                               ));
 
-  this->_addBuiltinFunction("variational_solver",
-                            std::function(
-
-                              [](const size_t& degree, 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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return VariationalSolverHandler{getCommonMesh({rho, u, E, c, a})}
-                                  .solver()
-                                  .apply(degree, dt, VariationalSolverHandler::VelocityBCTreatment::elimination, rho, u,
-                                         E, c, a, bc_descriptor_list, {});
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("variational_solver_with_source",
-                            std::function(
-
-                              [](const size_t& degree, 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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt, const FunctionSymbolId& function_id)
-                                -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                              std::shared_ptr<const DiscreteFunctionVariant>,
-                                              std::shared_ptr<const DiscreteFunctionVariant>,
-                                              std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return VariationalSolverHandler{getCommonMesh({rho, u, E, c, a})}
-                                  .solver()
-                                  .apply(degree, dt, VariationalSolverHandler::VelocityBCTreatment::elimination, rho, u,
-                                         E, c, a, bc_descriptor_list, function_id);
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("variational_solver_dr",
-                            std::function(
-
-                              [](const size_t& degree, 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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return VariationalSolverDevelopedReconstructionHandler{getCommonMesh({rho, u, E, c, a})}
-                                  .solver()
-                                  .apply(degree, dt,
-                                         VariationalSolverDevelopedReconstructionHandler::VelocityBCTreatment::
-                                           elimination,
-                                         rho, u, E, c, a, bc_descriptor_list, {});
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("variational_solver_dr_with_source",
-                            std::function(
-
-                              [](const size_t& degree, 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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt, const FunctionSymbolId& function_id)
-                                -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                              std::shared_ptr<const DiscreteFunctionVariant>,
-                                              std::shared_ptr<const DiscreteFunctionVariant>,
-                                              std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return VariationalSolverDevelopedReconstructionHandler{getCommonMesh({rho, u, E, c, a})}
-                                  .solver()
-                                  .apply(degree, dt,
-                                         VariationalSolverDevelopedReconstructionHandler::VelocityBCTreatment::
-                                           elimination,
-                                         rho, u, E, c, a, bc_descriptor_list, function_id);
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("variational_solver_o1",
-                            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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return VariationalSolverO1Handler{getCommonMesh({rho, u, E, c, a, p})}
-                                  .solver()
-                                  .apply(dt, VariationalSolverO1Handler::VelocityBCTreatment::elimination, rho, u, E, c,
-                                         a, p, bc_descriptor_list);
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("penalized_variational_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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return VariationalSolverO1Handler{getCommonMesh({rho, u, E, c, a, p})}
-                                  .solver()
-                                  .apply(dt, VariationalSolverO1Handler::VelocityBCTreatment::penalty, rho, u, E, c, a,
-                                         p, bc_descriptor_list);
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("p1_p0_analytic_variational_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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return P1P0AnalyticVariationalSolverHandler{getCommonMesh({rho, u, E, c, a, p})}
-                                  .solver()
-                                  .apply(dt, P1P0AnalyticVariationalSolverHandler::VelocityBCTreatment::elimination,
-                                         rho, u, E, c, a, p, bc_descriptor_list);
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("p1_p0_variational_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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return P1P0VariationalSolverHandler{getCommonMesh({rho, u, E, c, a, p})}
-                                  .solver()
-                                  .apply(dt, P1P0VariationalSolverHandler::VelocityBCTreatment::elimination, rho, u, E,
-                                         c, a, p, bc_descriptor_list);
-                              }
-
-                              ));
-
-  this->_addBuiltinFunction("p1_p0_analytic_penalized_variational_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>& c,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
-                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
-                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
-                                   bc_descriptor_list,
-                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
-                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
-                                return P1P0AnalyticVariationalSolverHandler{getCommonMesh({rho, u, E, c, a, p})}
-                                  .solver()
-                                  .apply(dt, P1P0AnalyticVariationalSolverHandler::VelocityBCTreatment::penalty, rho, u,
-                                         E, c, a, p, bc_descriptor_list);
-                              }
-
-                              ));
-
   this->_addBuiltinFunction("apply_acoustic_fluxes",
                             std::function(
 
@@ -885,15 +671,6 @@ SchemeModule::SchemeModule()
 
                                              ));
 
-  this->_addBuiltinFunction("variational_acoustic_dt",
-                            std::function(
-
-                              [](const std::shared_ptr<const DiscreteFunctionVariant>& c) -> double {
-                                return variational_acoustic_dt(c);
-                              }
-
-                              ));
-
   this->_addBuiltinFunction("cell_volume",
                             std::function(
 
diff --git a/src/language/modules/VariationalSchemeModule.cpp b/src/language/modules/VariationalSchemeModule.cpp
new file mode 100644
index 000000000..b6a8cad4f
--- /dev/null
+++ b/src/language/modules/VariationalSchemeModule.cpp
@@ -0,0 +1,247 @@
+#include <language/modules/VariationalSchemeModule.hpp>
+
+#include <language/modules/MeshModuleTypes.hpp>
+#include <language/modules/ModuleRepository.hpp>
+#include <language/modules/SchemeModuleTypes.hpp>
+#include <language/utils/BuiltinFunctionEmbedder.hpp>
+#include <language/utils/TypeDescriptor.hpp>
+#include <scheme/DiscreteFunctionUtils.hpp>
+#include <scheme/P1P0AnalyticVariationalSolver.hpp>
+#include <scheme/P1P0VariationalSolver.hpp>
+#include <scheme/VariationalSolver.hpp>
+#include <scheme/VariationalSolverDevelopedReconstruction.hpp>
+#include <scheme/VariationalSolverO1.hpp>
+#include <scheme/WallBoundaryConditionDescriptor.hpp>
+
+#include <memory>
+
+VariationalSchemeModule::VariationalSchemeModule()
+{
+  this->_addBuiltinFunction("variational_solver",
+                            std::function(
+
+                              [](const size_t& degree, 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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return VariationalSolverHandler{getCommonMesh({rho, u, E, c, a})}
+                                  .solver()
+                                  .apply(degree, dt, VariationalSolverHandler::VelocityBCTreatment::elimination, rho, u,
+                                         E, c, a, bc_descriptor_list, {});
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("variational_solver_with_source",
+                            std::function(
+
+                              [](const size_t& degree, 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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt, const FunctionSymbolId& function_id)
+                                -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                              std::shared_ptr<const DiscreteFunctionVariant>,
+                                              std::shared_ptr<const DiscreteFunctionVariant>,
+                                              std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return VariationalSolverHandler{getCommonMesh({rho, u, E, c, a})}
+                                  .solver()
+                                  .apply(degree, dt, VariationalSolverHandler::VelocityBCTreatment::elimination, rho, u,
+                                         E, c, a, bc_descriptor_list, function_id);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("variational_solver_dr",
+                            std::function(
+
+                              [](const size_t& degree, 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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return VariationalSolverDevelopedReconstructionHandler{getCommonMesh({rho, u, E, c, a})}
+                                  .solver()
+                                  .apply(degree, dt,
+                                         VariationalSolverDevelopedReconstructionHandler::VelocityBCTreatment::
+                                           elimination,
+                                         rho, u, E, c, a, bc_descriptor_list, {});
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("variational_solver_dr_with_source",
+                            std::function(
+
+                              [](const size_t& degree, 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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt, const FunctionSymbolId& function_id)
+                                -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                              std::shared_ptr<const DiscreteFunctionVariant>,
+                                              std::shared_ptr<const DiscreteFunctionVariant>,
+                                              std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return VariationalSolverDevelopedReconstructionHandler{getCommonMesh({rho, u, E, c, a})}
+                                  .solver()
+                                  .apply(degree, dt,
+                                         VariationalSolverDevelopedReconstructionHandler::VelocityBCTreatment::
+                                           elimination,
+                                         rho, u, E, c, a, bc_descriptor_list, function_id);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("variational_solver_o1",
+                            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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return VariationalSolverO1Handler{getCommonMesh({rho, u, E, c, a, p})}
+                                  .solver()
+                                  .apply(dt, VariationalSolverO1Handler::VelocityBCTreatment::elimination, rho, u, E, c,
+                                         a, p, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("penalized_variational_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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return VariationalSolverO1Handler{getCommonMesh({rho, u, E, c, a, p})}
+                                  .solver()
+                                  .apply(dt, VariationalSolverO1Handler::VelocityBCTreatment::penalty, rho, u, E, c, a,
+                                         p, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("p1_p0_analytic_variational_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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return P1P0AnalyticVariationalSolverHandler{getCommonMesh({rho, u, E, c, a, p})}
+                                  .solver()
+                                  .apply(dt, P1P0AnalyticVariationalSolverHandler::VelocityBCTreatment::elimination,
+                                         rho, u, E, c, a, p, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("p1_p0_variational_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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return P1P0VariationalSolverHandler{getCommonMesh({rho, u, E, c, a, p})}
+                                  .solver()
+                                  .apply(dt, P1P0VariationalSolverHandler::VelocityBCTreatment::elimination, rho, u, E,
+                                         c, a, p, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("p1_p0_analytic_penalized_variational_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>& c,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& a,
+                                 const std::shared_ptr<const DiscreteFunctionVariant>& p,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const double& dt) -> std::tuple<std::shared_ptr<const MeshVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>,
+                                                                 std::shared_ptr<const DiscreteFunctionVariant>> {
+                                return P1P0AnalyticVariationalSolverHandler{getCommonMesh({rho, u, E, c, a, p})}
+                                  .solver()
+                                  .apply(dt, P1P0AnalyticVariationalSolverHandler::VelocityBCTreatment::penalty, rho, u,
+                                         E, c, a, p, bc_descriptor_list);
+                              }
+
+                              ));
+
+  this->_addBuiltinFunction("variational_acoustic_dt",
+                            std::function(
+
+                              [](const std::shared_ptr<const DiscreteFunctionVariant>& c) -> double {
+                                return variational_acoustic_dt(c);
+                              }
+
+                              ));
+}
+
+void
+VariationalSchemeModule::registerOperators() const
+{}
+
+void
+VariationalSchemeModule::registerCheckpointResume() const
+{}
+
+ModuleRepository::Subscribe<VariationalSchemeModule> variational_scheme_module;
diff --git a/src/language/modules/VariationalSchemeModule.hpp b/src/language/modules/VariationalSchemeModule.hpp
new file mode 100644
index 000000000..7aa2ab54a
--- /dev/null
+++ b/src/language/modules/VariationalSchemeModule.hpp
@@ -0,0 +1,23 @@
+#ifndef VARIATIONAL_SCHEME_MODULE_HPP
+#define VARIATIONAL_SCHEME_MODULE_HPP
+
+#include <language/modules/BuiltinModule.hpp>
+
+class VariationalSchemeModule : public BuiltinModule
+{
+ public:
+  std::string_view
+  name() const final
+  {
+    return "variational_scheme";
+  }
+
+  void registerOperators() const final;
+  void registerCheckpointResume() const final;
+
+  VariationalSchemeModule();
+
+  ~VariationalSchemeModule() = default;
+};
+
+#endif   // VARIATIONAL_SCHEME_MODULE_HPP
-- 
GitLab