#1022 Compiles but doesn't run well: I need first to deal with #1040.

Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStep/ImplicitStepFluid/NewtonFixedPoint/GlobalVariationalOperatorInstances/T21.hpp

@@ -91,6 +91,7 @@ namespace HappyHeart
                  */
                 explicit T21(const FEltSpace& felt_space,
                              const Unknown& fluid_mass,
+                             const Unknown& fluid_velocity,
                              const ScalarParameter<>& fluid_density,
                              const HyperelasticLawT& law,
                              const cauchy_green_tensor_type& cauchy_green_tensor,

Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStep/ImplicitStepFluid/NewtonFixedPoint/GlobalVariationalOperatorInstances/T21.hxx

@@ -36,12 +36,14 @@ namespace HappyHeart
             template<class HyperelasticLawT>
             T21<HyperelasticLawT>::T21(const FEltSpace& felt_space,
                                        const Unknown& fluid_mass,
+                                       const Unknown& fluid_velocity,
                                        const ScalarParameter<>& fluid_density,
                                        const HyperelasticLawT& law,
                                        const cauchy_green_tensor_type& cauchy_green_tensor,
                                        const QuadratureRulePerTopology* const quadrature_rule_per_topology)
             : parent(felt_space,
                      fluid_mass,
+                     fluid_velocity,
                      quadrature_rule_per_topology,
                      AllocateGradientFEltPhi::yes,
                      DoComputeProcessorWiseLocal2Global::yes, // \todo #820 See if no ok

Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStep/ImplicitStepFluid/NewtonFixedPoint/LocalVariationalOperatorInstances/T21.hxx

@@ -44,7 +44,7 @@ namespace HappyHeart
             fluid_density_(fluid_density),
             hyperelastic_law_(law)
             {
-                assert(unknown_list.size() == 1ul);
+                assert(unknown_list.size() == 2ul);
 
                 decltype(auto) elementary_data = parent::GetElementaryData();
 

Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStep/ImplicitStepFluid/NewtonFixedPoint/VariationalFormulationInit.hxx

@@ -171,7 +171,7 @@ namespace HappyHeart
                     decltype(auto) god_of_dof = parent::GetGodOfDof();
 
                     decltype(auto) unknown_manager = UnknownManager::GetInstance();
-                    decltype(auto) fluid_velocity =
+                    decltype(auto) fluid_velocity_unknown =
                         unknown_manager.GetUnknown(EnumUnderlyingType(UnknownIndex::fluid_velocity));
                     decltype(auto) fluid_mass_unknown =
                         unknown_manager.GetUnknown(EnumUnderlyingType(UnknownIndex::fluid_mass));
@@ -187,7 +187,7 @@ namespace HappyHeart
 
                     inlet_pressure_source_operator_ =
                         std::make_unique<inlet_pressure_operator_type>(inlet_border_felt_space,
-                                                                       fluid_velocity,
+                                                                       fluid_velocity_unknown,
                                                                        porosity,
                                                                        GetInletPressureParam());
 
@@ -205,28 +205,30 @@ namespace HappyHeart
                     decltype(auto) solid_displacement_data = this->GetSolidDisplacementData();
                     decltype(auto) new_fluid_pressure_data = this->GetNonCstNewFluidPressureData();
                     
-                    {
-                        decltype(auto) god_of_dof_manager = GodOfDofManager::GetInstance();
-                        
-                        decltype(auto) solid_god_of_dof =
+                    decltype(auto) god_of_dof_manager = GodOfDofManager::GetInstance();
+                    decltype(auto) solid_god_of_dof =
                         god_of_dof_manager.GetGodOfDof(EnumUnderlyingType(MeshIndex::solid));
+
+                    {
+                        decltype(auto) hyperelastic_law = GetHyperelasticLaw();
                         
                         decltype(auto) solid_felt_space =
-                        solid_god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::solid));
-                        
-                        decltype(auto) hyperelastic_law = GetHyperelasticLaw();
+                            solid_god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::solid));
+
                         
                         t21_operator_ =
-                        std::make_unique<t21_type>(solid_felt_space,
-                                                   fluid_mass_unknown,
-                                                   variable_holder.GetFluidDensity(),
-                                                   hyperelastic_law,
-                                                   this->GetCauchyGreenTensor());
+                            std::make_unique<t21_type>(solid_felt_space, // \todo #1040 Check adequate FEltSpace to use here.
+                                                       fluid_mass_unknown,
+                                                       fluid_velocity_unknown,
+                                                       variable_holder.GetFluidDensity(),
+                                                       hyperelastic_law,
+                                                       this->GetCauchyGreenTensor());
+
                     }
                     
                     darcy_operator_ =
                         std::make_unique<darcy_operator_type>(monolithic_felt_space,
-                                                              fluid_velocity,
+                                                              fluid_velocity_unknown,
                                                               variable_holder.GetFluidDensity(),
                                                               porosity,
                                                               fluid_velocity_data.GetNewAsParam(),
@@ -246,7 +248,7 @@ namespace HappyHeart
                     hybrid_vector_operator_ =
                         std::make_unique<GVO::HybridVector<HyperelasticLawT>>(monolithic_felt_space,
                                                                               fluid_mass_unknown,
-                                                                              fluid_velocity,
+                                                                              fluid_velocity_unknown,
                                                                               variable_holder.GetFluidDensity(),
                                                                               porosity,
                                                                               new_fluid_pressure_data.GetNonCstUpdatePressureParamOnFluidCurrent(),
@@ -255,7 +257,7 @@ namespace HappyHeart
 
                     mass_div_velocity_
                         = std::make_unique<scalar_div_vectorial_type>(fluid_felt_space,
-                                                                      fluid_velocity,
+                                                                      fluid_velocity_unknown,
                                                                       fluid_mass_unknown,
                                                                       1.,
                                                                       0.,
@@ -263,7 +265,7 @@ namespace HappyHeart
 
                     pressure_div_velocity_
                         = std::make_unique<scalar_div_vectorial_type>(fluid_felt_space,
-                                                                      fluid_velocity,
+                                                                      fluid_velocity_unknown,
                                                                       fluid_pressure,
                                                                       0.,
                                                                       -1.,
@@ -271,11 +273,9 @@ namespace HappyHeart
 
 
                     {
-                        decltype(auto) god_of_dof_manager = GodOfDofManager::GetInstance();
-
-                        decltype(auto) solid_god_of_dof =
-                        god_of_dof_manager.GetGodOfDof(EnumUnderlyingType(MeshIndex::solid));
+                        
 
+                      
                         decltype(auto) felt_space =
                         solid_god_of_dof.GetFEltSpace(EnumUnderlyingType(FEltSpaceIndex::solid));
 
@@ -314,7 +314,7 @@ namespace HappyHeart
                         t22_operator_ =
                             std::make_unique<GlobalVariationalOperatorNS::T22>(parent::GetTimeManager(),
                                                                                monolithic_felt_space,
-                                                                               fluid_velocity,
+                                                                               fluid_velocity_unknown,
                                                                                porosity,
                                                                                variable_holder.GetFluidDensity(),
                                                                                GetInternalFriction());
@@ -322,7 +322,7 @@ namespace HappyHeart
                         t22_operator_vel_felt_space_ =
                             std::make_unique<GlobalVariationalOperatorNS::T22>(parent::GetTimeManager(),
                                                                                fluid_felt_space,
-                                                                               fluid_velocity,
+                                                                               fluid_velocity_unknown,
                                                                                porosity,
                                                                                variable_holder.GetFluidDensity(),
                                                                                GetInternalFriction());

Sources/Operators/LocalVariationalOperator/CauchyAndInvariant/CauchyGreenTensor.cpp

+//! \file 
+//
+//
+//  CauchyGreenTensor.cpp
+//  HappyHeart
+//
+//  Created by Sebastien Gilles on 26/01/15.
+//  Copyright (c) 2015 Inria. All rights reserved.
+//
+
+#include "Utilities/Numeric/Numeric.hpp"
+
+#include "Operators/LocalVariationalOperator/CauchyAndInvariant/CauchyGreenTensor.hpp"
+#include "Operators/LocalVariationalOperator/Advanced/InformationsAtQuadraturePoint.hpp"
+#include "Operators/LocalVariationalOperator/Advanced/GradientDisplacementMatrix.hpp"
+
+
+namespace HappyHeart
+{
+    
+    
+    CauchyGreenTensor::CauchyGreenTensor(unsigned int mesh_dimension)
+    : mesh_dimension_(mesh_dimension)
+    {
+        switch(mesh_dimension)
+        {
+            case 2u:
+                value_list_.resize(3ul);
+                break;
+            case 3u:
+                value_list_.resize(6ul);
+                break;
+            default:
+                assert(false && "Case not foreseen!");
+                exit(EXIT_FAILURE);
+        }
+        
+        gradient_displacement_.Resize(static_cast<int>(mesh_dimension),
+                                      static_cast<int>(mesh_dimension));
+    }
+    
+    
+    
+    namespace // anonymous
+    {
+        
+        
+
+        /*!
+         * \brief Update for a tridimensional mesh.
+         *
+         */
+        void Update3D(const LocalMatrix& gradient_component_disp,
+                      std::vector<double>& value_list);
+        
+        /*!
+         * \brief Update for a bidimensional mesh.
+         */
+        void Update2D(const LocalMatrix& gradient_component_disp,
+                      std::vector<double>& value_list);
+        
+        
+    } // namespace anonymous
+
+               
+   
+    void CauchyGreenTensor::Update(const LocalMatrix& gradient_displacement)
+    {
+        const auto mesh_dimension = GetMeshDimension();
+        
+        auto& value_list = GetNonCstValueList();
+        
+        switch (mesh_dimension)
+        {
+            case 2u:
+                Update2D(gradient_displacement,
+                         value_list);
+                break;
+            case 3u:
+                Update3D(gradient_displacement,
+                         value_list);
+
+                break;
+        }
+    }
+    
+    
+    void CauchyGreenTensor::Update(const std::vector<double>& felt_displacement,
+                                   const Advanced::LocalVariationalOperatorNS::InformationsAtQuadraturePoint& infos_at_quad_pt)
+    {
+        auto& gradient_displacement = GetNonCstGradientDisplacement();
+        
+        //ASSUMES THE OPERATOR DOESN'T DEAL WITH TWO UNKNOWNS THAT DO NOT SHARE THE SAME NODES (which T21 totally does...)
+        assert(felt_displacement.size()
+               == static_cast<std::size_t>(infos_at_quad_pt.Nnode() * infos_at_quad_pt.GetMeshDimension()));
+        
+        Advanced::OperatorNS::ComputeGradientDisplacementMatrix(infos_at_quad_pt,
+                                                                felt_displacement,
+                                                                gradient_displacement);
+        
+        Update(gradient_displacement);
+    }
+
+      
+    namespace // anonymous
+    {
+        
+        
+        void Update3D(const LocalMatrix& gradient_component_disp,
+                      std::vector<double>& value_list)
+        {                               
+            using NumericNS::Square;
+            
+            // Component Cxx
+            value_list[0] = 1. + 2. * gradient_component_disp(0, 0)
+            + Square(gradient_component_disp(0, 0))
+            + Square(gradient_component_disp(1, 0))
+            + Square(gradient_component_disp(2, 0));
+            
+            // Component Cyy
+            value_list[1] = 1. + 2. * gradient_component_disp(1, 1)
+            + Square(gradient_component_disp(0, 1))
+            + Square(gradient_component_disp(1, 1))
+            + Square(gradient_component_disp(2, 1));
+            
+            // Component Czz
+            value_list[2] = 1. + 2. * gradient_component_disp(2, 2)
+            + Square(gradient_component_disp(0, 2))
+            + Square(gradient_component_disp(1, 2))
+            + Square(gradient_component_disp(2, 2));
+            
+            // Component Cxy
+            value_list[3] = gradient_component_disp(0, 1) + gradient_component_disp(1, 0)
+            + gradient_component_disp(0, 0) * gradient_component_disp(0, 1)
+            + gradient_component_disp(1, 0) * gradient_component_disp(1, 1)
+            + gradient_component_disp(2, 0) * gradient_component_disp(2, 1);
+            
+            // Component Cyz
+            value_list[4] = gradient_component_disp(1, 2) + gradient_component_disp(2, 1)
+            + gradient_component_disp(0, 2) * gradient_component_disp(0, 1)
+            + gradient_component_disp(1, 2) * gradient_component_disp(1, 1)
+            + gradient_component_disp(2, 2) * gradient_component_disp(2, 1);
+            
+            // Component Cxz
+            value_list[5] = gradient_component_disp(2, 0) + gradient_component_disp(0, 2)
+            + gradient_component_disp(0, 2) * gradient_component_disp(0, 0)
+            + gradient_component_disp(1, 2) * gradient_component_disp(1, 0)
+            + gradient_component_disp(2, 2) * gradient_component_disp(2, 0);
+            
+            
+            
+        }
+        
+        
+        void Update2D(const LocalMatrix& gradient_component_disp,
+                      std::vector<double>& value_list)
+        {
+            using NumericNS::Square;
+            
+            // Component Cxx
+            value_list[0] = 1. + 2. * gradient_component_disp(0, 0)
+            + Square(gradient_component_disp(0, 0))
+            + Square(gradient_component_disp(1, 0));
+            
+            // Component Cyy
+            value_list[1] = 1. + 2. * gradient_component_disp(1, 1)
+            + Square(gradient_component_disp(0, 1))
+            + Square(gradient_component_disp(1, 1));
+            
+            // Component Cxy
+            value_list[2] = gradient_component_disp(0, 1) + gradient_component_disp(1, 0)
+            + gradient_component_disp(0, 0) * gradient_component_disp(0, 1)
+            + gradient_component_disp(1, 0) * gradient_component_disp(1, 1);
+        }
+
+        
+    } // namespace anonymous
+
+    
+} // namespace HappyHeart