#1049 Fixes in several models.

Data/Lua/demo_input_ondomatic_2D.lua

@@ -351,6 +351,41 @@ Domain8 = {
     
 } -- Domain8
 
+
+Domain9 = {
+    
+    -- Index of the geometric mesh upon which the domain is defined (as defined in the present file). Might be
+    -- left empty if domain not limited to one mesh; at most one value is expected here.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_index = { 1 },
+    
+    -- List of dimensions encompassed by the domain. Might be left empty if no restriction at all upon
+    -- dimensions.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    -- Constraint: ops_in(v, {0, 1, 2, 3})
+    dimension_list = { },
+    
+    -- List of mesh labels encompassed by the domain. Might be left empty if no restriction at all upon mesh
+    -- labels. This parameter does not make sense if no mesh is defined for the domain.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_label_list = { },
+    
+    -- List of geometric element types considered in the domain. Might be left empty if no restriction upon
+    -- these. No constraint is applied at Ops level, as some geometric element types could be added after
+    -- generation of current input parameter file. Current list is below; if an incorrect value is put there it
+    -- will be detected a bit later when the domain object is built.
+    -- The known types when this file was generated are:
+    -- . Point1
+    -- . Segment2, Segment3
+    -- . Triangle3, Triangle6
+    -- . Quadrangle4, Quadrangle8, Quadrangle9
+    -- . Tetrahedron4, Tetrahedron10
+    -- . Hexahedron8, Hexahedron20, Hexahedron27.
+    -- Expected format: {"VALUE1", "VALUE2", ...}
+    geometric_element_type_list = { }
+    
+} -- Domain9
+
 EssentialBoundaryCondition1 = {
     
     -- Name of the boundary condition (must be unique).

Data/Lua/demo_input_ondomatic_3D.lua

@@ -351,6 +351,42 @@ Domain8 = {
     
 } -- Domain8
 
+
+Domain9 = {
+    
+    -- Index of the geometric mesh upon which the domain is defined (as defined in the present file). Might be
+    -- left empty if domain not limited to one mesh; at most one value is expected here.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_index = { 1 },
+    
+    -- List of dimensions encompassed by the domain. Might be left empty if no restriction at all upon
+    -- dimensions.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    -- Constraint: ops_in(v, {0, 1, 2, 3})
+    dimension_list = { },
+    
+    -- List of mesh labels encompassed by the domain. Might be left empty if no restriction at all upon mesh
+    -- labels. This parameter does not make sense if no mesh is defined for the domain.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_label_list = { },
+    
+    -- List of geometric element types considered in the domain. Might be left empty if no restriction upon
+    -- these. No constraint is applied at Ops level, as some geometric element types could be added after
+    -- generation of current input parameter file. Current list is below; if an incorrect value is put there it
+    -- will be detected a bit later when the domain object is built.
+    -- The known types when this file was generated are:
+    -- . Point1
+    -- . Segment2, Segment3
+    -- . Triangle3, Triangle6
+    -- . Quadrangle4, Quadrangle8, Quadrangle9
+    -- . Tetrahedron4, Tetrahedron10
+    -- . Hexahedron8, Hexahedron20, Hexahedron27.
+    -- Expected format: {"VALUE1", "VALUE2", ...}
+    geometric_element_type_list = { }
+    
+} -- Domain9
+
+
 EssentialBoundaryCondition1 = {
     
     -- Name of the boundary condition (must be unique).

Data/Lua/demo_input_surfacic_bidomain_CRN.lua

@@ -392,6 +392,42 @@ Domain8 = {
     
 } -- Domain8
 
+
+Domain9 = {
+    
+    -- Index of the geometric mesh upon which the domain is defined (as defined in the present file). Might be
+    -- left empty if domain not limited to one mesh; at most one value is expected here.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_index = { 1 },
+    
+    -- List of dimensions encompassed by the domain. Might be left empty if no restriction at all upon
+    -- dimensions.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    -- Constraint: ops_in(v, {0, 1, 2, 3})
+    dimension_list = { },
+    
+    -- List of mesh labels encompassed by the domain. Might be left empty if no restriction at all upon mesh
+    -- labels. This parameter does not make sense if no mesh is defined for the domain.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_label_list = { },
+    
+    -- List of geometric element types considered in the domain. Might be left empty if no restriction upon
+    -- these. No constraint is applied at Ops level, as some geometric element types could be added after
+    -- generation of current input parameter file. Current list is below; if an incorrect value is put there it
+    -- will be detected a bit later when the domain object is built.
+    -- The known types when this file was generated are:
+    -- . Point1
+    -- . Segment2, Segment3
+    -- . Triangle3, Triangle6
+    -- . Quadrangle4, Quadrangle8, Quadrangle9
+    -- . Tetrahedron4, Tetrahedron10
+    -- . Hexahedron8, Hexahedron20, Hexahedron27.
+    -- Expected format: {"VALUE1", "VALUE2", ...}
+    geometric_element_type_list = { }
+    
+} -- Domain9
+
+
 FiniteElementSpace1 = {
     
     -- Index of the god of dof into which the finite element space is defined.
@@ -709,7 +745,7 @@ Fiber_vector_1 = {
     
     -- Index of the domain upon which parameter is defined.
     -- Expected format: VALUE
-    mesh_index = 1,
+    domain_index = 9,
     
     -- Index of the finite element space upon which parameter is defined.
     -- Expected format: VALUE
@@ -730,7 +766,7 @@ Fiber_scalar_2 = {
 
 	-- Index of the domain upon which parameter is defined.
 	-- Expected format: VALUE
-	mesh_index = 1,
+	domain_index = 9,
 
 	-- Index of the finite element space upon which parameter is defined.
 	-- Expected format: VALUE

Data/Lua/demo_input_surfacic_bidomain_MS.lua

@@ -392,6 +392,41 @@ Domain8 = {
     
 } -- Domain8
 
+Domain9 = {
+    
+    -- Index of the geometric mesh upon which the domain is defined (as defined in the present file). Might be
+    -- left empty if domain not limited to one mesh; at most one value is expected here.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_index = { 1 },
+    
+    -- List of dimensions encompassed by the domain. Might be left empty if no restriction at all upon
+    -- dimensions.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    -- Constraint: ops_in(v, {0, 1, 2, 3})
+    dimension_list = { },
+    
+    -- List of mesh labels encompassed by the domain. Might be left empty if no restriction at all upon mesh
+    -- labels. This parameter does not make sense if no mesh is defined for the domain.
+    -- Expected format: {VALUE1, VALUE2, ...}
+    mesh_label_list = { },
+    
+    -- List of geometric element types considered in the domain. Might be left empty if no restriction upon
+    -- these. No constraint is applied at Ops level, as some geometric element types could be added after
+    -- generation of current input parameter file. Current list is below; if an incorrect value is put there it
+    -- will be detected a bit later when the domain object is built.
+    -- The known types when this file was generated are:
+    -- . Point1
+    -- . Segment2, Segment3
+    -- . Triangle3, Triangle6
+    -- . Quadrangle4, Quadrangle8, Quadrangle9
+    -- . Tetrahedron4, Tetrahedron10
+    -- . Hexahedron8, Hexahedron20, Hexahedron27.
+    -- Expected format: {"VALUE1", "VALUE2", ...}
+    geometric_element_type_list = { }
+    
+} -- Domain9
+
+
 FiniteElementSpace1 = {
     
     -- Index of the god of dof into which the finite element space is defined.
@@ -709,7 +744,7 @@ Fiber_vector_1 = {
     
     -- Index of the domain upon which parameter is defined.
     -- Expected format: VALUE
-    mesh_index = 1,
+    domain_index = 9,
     
     -- Index of the finite element space upon which parameter is defined.
     -- Expected format: VALUE
@@ -730,7 +765,7 @@ Fiber_scalar_2 = {
     
     -- Index of the domain upon which parameter is defined.
     -- Expected format: VALUE
-    mesh_index = 1,
+    domain_index = 9,
     
     -- Index of the finite element space upon which parameter is defined.
     -- Expected format: VALUE

Sources/ModelInstances/Ondomatic/InputParameterList.hpp

@@ -46,7 +46,8 @@ namespace HappyHeart
             face4 = 5,
             face5 = 6,
             face6 = 7,
-            face123 = 8
+            face123 = 8,
+            full_mesh = 9
         };
 
         
@@ -113,6 +114,7 @@ namespace HappyHeart
             InputParameter::Domain<EnumUnderlyingType(DomainIndex::face5)>,
             InputParameter::Domain<EnumUnderlyingType(DomainIndex::face6)>,
             InputParameter::Domain<EnumUnderlyingType(DomainIndex::face123)>,
+            InputParameter::Domain<EnumUnderlyingType(DomainIndex::full_mesh)>,
         
             InputParameter::DirichletBoundaryCondition<EnumUnderlyingType(BoundaryConditionIndex::face123)>,
         

Sources/ModelInstances/Ondomatic/VariationalFormulation.cpp

@@ -44,8 +44,6 @@ namespace HappyHeart
         {
             const auto& god_of_dof = GetGodOfDof();
             
-            const GeometricMeshRegion& mesh = god_of_dof.GetGeometricMeshRegion();
-            
             const FEltSpace& felt_space_highest_dimension = god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::highest_dimension));
             const FEltSpace& felt_space_surface_pressure_4 = god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::surface_pressure_4));
             const FEltSpace& felt_space_surface_pressure_5 = god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::surface_pressure_5));
@@ -66,13 +64,13 @@ namespace HappyHeart
             }
             
             namespace GVO = GlobalVariationalOperatorNS;
-            
+            decltype(auto) domain = DomainManager::GetInstance().GetDomain(EnumUnderlyingType(DomainIndex::full_mesh));
             {
                 using parameter_type = InputParameter::TransientSource<EnumUnderlyingType(ForceIndexList::surface_pressure_4)>;
                 
                 surface_pressure_4_parameter_ =
                     InitThreeDimensionalParameterFromInputData<parameter_type>("surface_pressure_4",
-                                                                  mesh,
+                                                                  domain,
                                                                   input_parameter_data);
                 
                 if (surface_pressure_4_parameter_ != nullptr)
@@ -88,7 +86,7 @@ namespace HappyHeart
                 
                 surface_pressure_5_parameter_ =
                 InitThreeDimensionalParameterFromInputData<parameter_type>("surface_pressure_5",
-                                                              mesh,
+                                                              domain,
                                                               input_parameter_data);
                 
                 if (surface_pressure_5_parameter_ != nullptr)
@@ -104,7 +102,7 @@ namespace HappyHeart
                 
                 surface_pressure_6_parameter_ =
                 InitThreeDimensionalParameterFromInputData<parameter_type>("surface_pressure_6",
-                                                              mesh,
+                                                              domain,
                                                               input_parameter_data);
                 
                 if (surface_pressure_6_parameter_ != nullptr)

Sources/ModelInstances/SurfacicBidomain/InputParameterList.hpp

@@ -52,7 +52,8 @@ namespace HappyHeart
             rai = 5,
             pm = 6,
             fo = 7,
-            regular_atria = 8
+            regular_atria = 8,
+            full_mesh = 9
         };
 
 
@@ -139,6 +140,7 @@ namespace HappyHeart
             InputParameter::Domain<EnumUnderlyingType(DomainIndex::pm)>,
             InputParameter::Domain<EnumUnderlyingType(DomainIndex::fo)>,
             InputParameter::Domain<EnumUnderlyingType(DomainIndex::regular_atria)>,
+            InputParameter::Domain<EnumUnderlyingType(DomainIndex::full_mesh)>,
             
             InputParameter::FEltSpace<EnumUnderlyingType(FEltSpaceIndex::highest_dimension)>,
             

Sources/ModelInstances/SurfacicBidomain/SurfacicBidomainVariationalFormulation.hxx

@@ -58,16 +58,15 @@ namespace HappyHeart
         template <Advanced::ReactionLawNS::ReactionLawName ReactionLawNameT>
         void SurfacicBidomainVariationalFormulation<ReactionLawNameT>::SupplInit(const InputParameterList& input_parameter_data)
         {
-            const auto& god_of_dof = this->GetGodOfDof();
-            const auto& mesh = god_of_dof.GetGeometricMeshRegion();
+            decltype(auto) domain = DomainManager::GetInstance().GetDomain(EnumUnderlyingType(DomainIndex::full_mesh));
             
             using Diffusion = InputParameter::Diffusion;
-            intracellular_trans_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::intracellular_trans_diffusion_tensor)>>("Intracellular Trans Diffusion tensor", mesh, input_parameter_data);
-            extracellular_trans_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::extracellular_trans_diffusion_tensor)>>("Extracellular Trans Diffusion tensor", mesh, input_parameter_data);
-            intracellular_fiber_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::intracellular_fiber_diffusion_tensor)>>("Intracellular Fiber Diffusion tensor", mesh, input_parameter_data);
-            extracellular_fiber_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::extracellular_fiber_diffusion_tensor)>>("Extracellular Fiber Diffusion tensor", mesh, input_parameter_data);
-            heterogeneous_conductivity_coefficient_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::heterogeneous_conductivity_coefficient)>>("Heterogeneous conductivity coefficient", mesh, input_parameter_data);
-            transcellular_density_ = InitScalarParameterFromInputData<Diffusion::Density>("Transcellular Density", mesh, input_parameter_data);
+            intracellular_trans_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::intracellular_trans_diffusion_tensor)>>("Intracellular Trans Diffusion tensor", domain, input_parameter_data);
+            extracellular_trans_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::extracellular_trans_diffusion_tensor)>>("Extracellular Trans Diffusion tensor", domain, input_parameter_data);
+            intracellular_fiber_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::intracellular_fiber_diffusion_tensor)>>("Intracellular Fiber Diffusion tensor", domain, input_parameter_data);
+            extracellular_fiber_diffusion_tensor_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::extracellular_fiber_diffusion_tensor)>>("Extracellular Fiber Diffusion tensor", domain, input_parameter_data);
+            heterogeneous_conductivity_coefficient_ = InitScalarParameterFromInputData<Diffusion::Tensor<EnumUnderlyingType(TensorIndex::heterogeneous_conductivity_coefficient)>>("Heterogeneous conductivity coefficient", domain, input_parameter_data);
+            transcellular_density_ = InitScalarParameterFromInputData<Diffusion::Density>("Transcellular Density", domain, input_parameter_data);
             
             if (!GetTranscellularDiffusionDensity().IsConstant())
                 throw Exception("Current Bidomain  model is restricted to a constant diffusion density.",
@@ -112,7 +111,6 @@ namespace HappyHeart
         void SurfacicBidomainVariationalFormulation<ReactionLawNameT>::DefineOperators(const InputParameterList& input_parameter_data)
         {
             const auto& god_of_dof = this->GetGodOfDof();
-            const auto& mesh = god_of_dof.GetGeometricMeshRegion();
             const auto& felt_space_highest_dimension = god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::highest_dimension));
             
             const auto& transcellular_potential = UnknownManager::GetInstance().GetUnknown(EnumUnderlyingType(UnknownIndex::transcellular_potential));
@@ -140,6 +138,10 @@ namespace HappyHeart
                                                                                   fibers,
                                                                                   angles);
             
+            decltype(auto) domain = DomainManager::GetInstance().GetDomain(EnumUnderlyingType(DomainIndex::full_mesh));
+            
+
+            
             capacity_operator_ = std::make_unique<GVO::Mass>(felt_space_highest_dimension,
                                                              transcellular_potential);
             
@@ -147,7 +149,7 @@ namespace HappyHeart
                                                                                      extracellular_potential);
             
             reaction_law_ = std::make_unique<reaction_law_type>(input_parameter_data,
-                                                                mesh,
+                                                                domain,
                                                                 this->GetTimeManager(),
                                                                 felt_space_highest_dimension.GetQuadratureRulePerTopology());
             
@@ -156,9 +158,10 @@ namespace HappyHeart
                                                                                                    GetNonCstReactionLaw());
             
             using parameter_type = InputParameter::TransientSource<EnumUnderlyingType(ForceIndexList::transcellular_current_applied)>;
+            
             transcellular_current_applied_parameter_ =
             InitThreeDimensionalParameterFromInputData<parameter_type>("Current Applied",
-                                                          mesh,
+                                                          domain,
                                                           input_parameter_data);
             
             if (transcellular_current_applied_parameter_ != nullptr)

Sources/OperatorInstances/VariationalOperator/NonlinearForm/Local/NonLinearSource/ReactionLaw/Instantiations/FitzHughNagumo.hpp

@@ -87,7 +87,7 @@ namespace HappyHeart
              */
             template <class InputParameterDataT>
             explicit ReactionLaw(const InputParameterDataT& input_parameter_data,
-                                 const GeometricMeshRegion& mesh,
+                                 const Domain& domain,
                                  const TimeManager& time_manager,
                                  const QuadratureRulePerTopology& default_quadrature_rule_set);
 

Sources/OperatorInstances/VariationalOperator/NonlinearForm/Local/NonLinearSource/ReactionLaw/Instantiations/FitzHughNagumo.hxx

@@ -26,7 +26,7 @@ namespace HappyHeart
 
         template <class InputParameterDataT>
         ReactionLaw<ReactionLawName::FitzHughNagumo>::ReactionLaw(const InputParameterDataT& input_parameter_data,
-                                                                  const GeometricMeshRegion& mesh,
+                                                                  const Domain& domain,
                                                                   const TimeManager& time_manager,
                                                                   const QuadratureRulePerTopology& default_quadrature_rule_set)
         : a_(Utilities::InputParameterListNS::Extract<input_param_fhn::ACoefficient>::Value(input_parameter_data)),
@@ -41,7 +41,7 @@ namespace HappyHeart
             const double initial_condition_gate = ipl::Extract<InitialConditionGate::Value>::Value(input_parameter_data);
 
             gate_ = std::make_unique<ParameterAtQuadraturePoint<ParameterNS::Type::scalar>>("Gate",
-                                                                                           mesh,
+                                                                                           domain,
                                                                                            default_quadrature_rule_set,
                                                                                            initial_condition_gate,
                                                                                            this->GetTimeManager());

Sources/OperatorInstances/VariationalOperator/NonlinearForm/Local/NonLinearSource/ReactionLaw/Instantiations/MitchellSchaeffer.hpp

@@ -85,7 +85,7 @@ namespace HappyHeart
              */
             template <class InputParameterDataT>
             explicit ReactionLaw(const InputParameterDataT& input_parameter_data,
-                                 const GeometricMeshRegion& mesh,
+                                 const Domain& domain,
                                  const TimeManager& time_manager,
                                  const QuadratureRulePerTopology& default_quadrature_rule_set);
 

Sources/OperatorInstances/VariationalOperator/NonlinearForm/Local/NonLinearSource/ReactionLaw/Instantiations/MitchellSchaeffer.hxx

@@ -29,7 +29,7 @@ namespace HappyHeart
 
         template <class InputParameterDataT>
         ReactionLaw<ReactionLawName::MitchellSchaeffer>::ReactionLaw(const InputParameterDataT& input_parameter_data,
-                                                                     const GeometricMeshRegion& mesh,
+                                                                     const Domain& domain,
                                                                      const TimeManager& time_manager,
                                                                      const QuadratureRulePerTopology& default_quadrature_rule_set)
         : tau_in_(Utilities::InputParameterListNS::Extract<input_param_ms::TauIn>::Value(input_parameter_data)),
@@ -47,13 +47,13 @@ namespace HappyHeart
             const double initial_condition_gate = ipl::Extract<InitialConditionGate::Value>::Value(input_parameter_data);
 
             gate_ = std::make_unique<ScalarParameterAtQuadPt>("Gate",
-                                                              mesh,
+                                                              domain,
                                                               default_quadrature_rule_set,
                                                               initial_condition_gate,
                                                               this->GetTimeManager());
 
             tau_close_ = InitScalarParameterFromInputData<input_param_ms::TauClose>("Tau Close",
-                                                                       mesh,
+                                                                       domain,
                                                                        input_parameter_data);
         }