#723 FSI/Newton: start the computation of the 'dH' term of Bruno's Freefem program.

Sources/ModelInstances/FSI_EI/ImplicitStepFluidVariationalFormulation.cpp

@@ -202,6 +202,7 @@ namespace HappyHeart
                         GlobalVectorWithCoefficient vector(rhs, 1.);
                         GetInletPressureOperator().Assemble(std::make_tuple(std::ref(vector)), time);
                     }
+                    break;
                 }
             }
             

Sources/ModelInstances/FSI_EI/Newton/Model.hpp

@@ -157,7 +157,7 @@ namespace HappyHeart
             void ImplicitFluidStep(differential is_differential);
             
             //! Update fluid Dirichlet condition.
-            void UpdateFluidDirichletCondition();
+            void UpdateFluidDirichletCondition(differential is_differential);
             
             //! Update solid displacements before SoF.
             void UpdateSolidDisplacementBeforeSoF();
@@ -219,8 +219,13 @@ namespace HappyHeart
             
             //! Evolution of solid displacement in the latest SoF computation.
             GlobalVector& GetNonCstSolidDisplacementEvolutionInSoF() noexcept;
-
             
+            //! Evolution of solid displacement in the latest SoF computation.
+            const GlobalVector& GetSolidDisplacementEvolutionInSoFOnInterface() const noexcept;
+            
+            //! Evolution of solid displacement in the latest SoF computation.
+            GlobalVector& GetNonCstSolidDisplacementEvolutionInSoFOnInterface() noexcept;
+
             //! Solid displacement before SoF is applied (in Aitken loop). See Bruno's note for more details about Aitken.
             const GlobalVector& GetSolidDisplacementBeforeSoFOnInterface() const noexcept;
             
@@ -251,7 +256,8 @@ namespace HappyHeart
             GlobalVector::unique_ptr solid_displacement_after_SoF_ = nullptr;
             
             GlobalVector::unique_ptr solid_displacement_evolution_in_SoF_ = nullptr;
-
+            GlobalVector::unique_ptr solid_displacement_evolution_in_SoF_on_interface_ = nullptr;
+            
             GlobalVector::unique_ptr solid_velocity_on_interface_ = nullptr;
             GlobalVector::unique_ptr p2_solid_velocity_on_interface_ = nullptr;
             GlobalVector::unique_ptr solid_displacement_previous_time_iteration_ = nullptr;

Sources/ModelInstances/FSI_EI/Newton/Model.hxx

@@ -127,6 +127,8 @@ namespace HappyHeart
                 {
                     auto& solid_displacement_before_SoF_on_interface = GetNonCstSolidDisplacementBeforeSoFOnInterface();
                     HappyHeart::AllocateGlobalVector(god_of_dof, solid_displacement_before_SoF_on_interface);
+                    
+                    solid_displacement_evolution_in_SoF_on_interface_ = std::make_unique<GlobalVector>(solid_displacement_before_SoF_on_interface);
                 }
                 
                 solid_displacement_after_SoF_ =

Sources/ModelInstances/FSI_EI/Newton/ModelAccessors.hxx

@@ -123,6 +123,21 @@ namespace HappyHeart
         }
         
         
+        template<class SolidVariationalFormulationPolicyT>
+        GlobalVector& Model<SolidVariationalFormulationPolicyT>::GetNonCstSolidDisplacementEvolutionInSoFOnInterface() noexcept
+        {
+            return const_cast<GlobalVector&>(GetSolidDisplacementEvolutionInSoFOnInterface());
+        }
+        
+        
+        template<class SolidVariationalFormulationPolicyT>
+        const GlobalVector& Model<SolidVariationalFormulationPolicyT>::GetSolidDisplacementEvolutionInSoFOnInterface() const noexcept
+        {
+            assert(!(!solid_displacement_evolution_in_SoF_on_interface_));
+            return *solid_displacement_evolution_in_SoF_on_interface_;
+        }
+        
+        
         template<class SolidVariationalFormulationPolicyT>
         GlobalVector& Model<SolidVariationalFormulationPolicyT>::GetNonCstSolidDisplacementEvolutionInSoF() noexcept
         {

Sources/ModelInstances/FSI_EI/Newton/ModelForward.hxx

@@ -36,8 +36,7 @@ namespace HappyHeart
             
             UpdateSolidDisplacementBeforeSoF();
             
-            UpdateFluidDirichletCondition();
-            
+            UpdateFluidDirichletCondition(differential::no);
             ImplicitFluidStep(differential::no);
             
             ComputeSolidResidual();
@@ -46,7 +45,12 @@ namespace HappyHeart
             
             UpdateSolidDisplacementAfterSoF();
             
+            UpdateFluidDirichletCondition(differential::yes);
+            ImplicitFluidStep(differential::yes);
             
+//            const auto& fluid_varf = GetImplicitStepFluidVariationalFormulation();
+//            const auto& sol = fluid_varf.GetSystemSolution(fluid_varf.GetNumberingSubset());
+//            sol.View(mpi, __FILE__, __LINE__);
             
             
 //            do
@@ -116,13 +120,13 @@ namespace HappyHeart
         
         template<class SolidVariationalFormulationPolicyT>
         void Model<SolidVariationalFormulationPolicyT>
-        ::UpdateFluidDirichletCondition()
+        ::UpdateFluidDirichletCondition(differential is_differential)
         {
             // Update the value of the Dirichlet condition.
             // Compute first the solid velocity from the previous aitken iteration.
             const double inv_time_step = 1. / this->GetTimeManager().GetTimeStep();
             
-            const auto& solid_displacement_previous_time_iteration_on_interface = *solid_displacement_previous_time_iteration_on_interface_;
+            
             const auto& p1_2_p2_interpolation_matrix = solid_to_fluid_velocity_on_fsi_interpolator_->GetInterpolationMatrix();
             
             auto& solid_velocity_on_interface = *solid_velocity_on_interface_; // \todo #608 Proper accessor!
@@ -130,13 +134,28 @@ namespace HappyHeart
             auto& implicit_fluid_boundary_condition =
             Private::DirichletBoundaryConditionManager::GetInstance().GetNonCstDirichletBoundaryCondition(EnumUnderlyingType(BoundaryConditionIndex::fluid_first));
             
-            solid_velocity_on_interface.Copy(GetSolidDisplacementBeforeSoFOnInterface(),
-                                             __FILE__, __LINE__);
-            
-            Wrappers::Petsc::AXPY(-1.,
-                                  solid_displacement_previous_time_iteration_on_interface,
-                                  solid_velocity_on_interface,
-                                  __FILE__, __LINE__);
+            switch(is_differential)
+            {
+                case differential::no:
+                {
+                    const auto& solid_displacement_previous_time_iteration_on_interface = *solid_displacement_previous_time_iteration_on_interface_;
+                    
+                    solid_velocity_on_interface.Copy(GetSolidDisplacementBeforeSoFOnInterface(),
+                                                     __FILE__, __LINE__);
+            
+                    Wrappers::Petsc::AXPY(-1.,
+                                          solid_displacement_previous_time_iteration_on_interface,
+                                          solid_velocity_on_interface,
+                                          __FILE__, __LINE__);
+                    break;
+                }
+                case differential::yes:
+                {
+                    solid_velocity_on_interface.Copy(GetSolidDisplacementEvolutionInSoFOnInterface(),
+                                                     __FILE__, __LINE__);
+                    break;
+                }
+            }
             
             solid_velocity_on_interface.Scale(inv_time_step, __FILE__, __LINE__);
                         
@@ -263,7 +282,7 @@ namespace HappyHeart
             solid_displacement_after_SoF.Copy(solid_variational_formulation.GetSystemSolution(solid_variational_formulation.GetNumberingSubset()),
                                               __FILE__, __LINE__);
             
-            auto& solid_displacement_evolution_in_SoF = GetNonCstSolidDisplacementEvolutionInSoF();            
+            auto& solid_displacement_evolution_in_SoF = GetNonCstSolidDisplacementEvolutionInSoF();
             
             solid_displacement_evolution_in_SoF.Copy(solid_displacement_after_SoF, __FILE__, __LINE__);
             
@@ -272,6 +291,12 @@ namespace HappyHeart
                                   solid_displacement_evolution_in_SoF,
                                   __FILE__, __LINE__);
             
+            auto& solid_disp_on_interface = GetNonCstSolidDisplacementEvolutionInSoFOnInterface();
+            InterpolateSolidDisplacementToInterface(solid_displacement_evolution_in_SoF,
+                                                    solid_disp_on_interface,
+                                                    __FILE__, __LINE__);
+
+            
             
 //            const auto& mpi = parent::MpiHappyHeart();