#820 Poromechanics / implicit: change the numbering subsets involved so that...

#820 Poromechanics / implicit: change the numbering subsets involved so that scalar div vectorial may correctly be assembled.

Data/Lua/demo_input_poromechanics.lua

@@ -593,7 +593,7 @@ FiniteElementSpace11 = {
     
     -- List of the numbering subset to use for each unknown;
     -- Expected format: {VALUE1, VALUE2, ...}
-    numbering_subset_list = { 10, 13 }
+    numbering_subset_list = { 15, 13 }
 }
 
 

Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStepFluidVariationalFormulation.cpp

@@ -46,7 +46,7 @@ namespace HappyHeart
                 decltype(auto) god_of_dof = parent::GetGodOfDof();
                 
                 decltype(auto) numbering_subset =
-                    god_of_dof.GetNumberingSubset(EnumUnderlyingType(NumberingSubsetIndex::fluid_velocity));
+                    god_of_dof.GetNumberingSubset(EnumUnderlyingType(NumberingSubsetIndex::fluid_mass_velocity_pressure));
             
                 auto& rhs = parent::GetNonCstSystemRhs(numbering_subset);
                 rhs.ZeroEntries(__FILE__, __LINE__);
@@ -62,8 +62,19 @@ namespace HappyHeart
 //                                         rhs,
 //                                         __FILE__, __LINE__);
                 
-                GlobalVectorWithCoefficient vector(rhs, -1.); // \todo #820 -1. empirical, to match Freefem's output. To check later!!!
-                GetInletPressureSourceOperator().Assemble(std::make_tuple(std::ref(vector)), time);
+                {
+                    GlobalVectorWithCoefficient vector(rhs, -1.); // \todo #820 -1. empirical, to match Freefem's output. To check later!!!
+                    GetInletPressureSourceOperator().Assemble(std::make_tuple(std::ref(vector)), time);
+                }
+                
+                auto& system_matrix = parent::GetNonCstSystemMatrix(numbering_subset, numbering_subset);
+                
+                {
+                    GlobalMatrixWithCoefficient system_matrix_with_coeff(system_matrix, 1.);
+                    
+                    GetScalarDivVectorial().Assemble(std::make_tuple(std::ref(system_matrix_with_coeff)));
+                    
+                }
 
                 
                 parent::ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset,
@@ -73,9 +84,6 @@ namespace HappyHeart
                                                     parent::GetOutputDirectory(numbering_subset) + "/implicit_rhs_time_" + std::to_string(time_manager.NtimeModified()) + ".hhdata",
                                                     __FILE__, __LINE__);
                 
-                
-                
-                
             
 //                parent::SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
 //                parent::WriteSolution(time_manager, numbering_subset);
@@ -88,7 +96,7 @@ namespace HappyHeart
                 decltype(auto) god_of_dof = parent::GetGodOfDof();
                 
                 decltype(auto) numbering_subset =
-                god_of_dof.GetNumberingSubset(EnumUnderlyingType(NumberingSubsetIndex::fluid_velocity));
+                god_of_dof.GetNumberingSubset(EnumUnderlyingType(NumberingSubsetIndex::fluid_mass_velocity_pressure));
 
                 
                 parent::AllocateSystemMatrix(numbering_subset, numbering_subset);