Commit 1cc41e96 authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#718 FSI: introduce hyperelasticity back with the same format as for the...

#718 FSI: introduce hyperelasticity back with the same format as for the brand-new elasticity. Content (and noticeably update) have not been reviewed and the results aren't correct at the moment.
parent b7cad03e
......@@ -4467,6 +4467,10 @@
BEFC2B2D1BD7D6770067F358 /* VariationalFormulation.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = VariationalFormulation.cpp; sourceTree = "<group>"; };
BEFC2B2E1BD7D6770067F358 /* VariationalFormulation.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = VariationalFormulation.hpp; sourceTree = "<group>"; };
BEFC2B2F1BD7D6770067F358 /* VariationalFormulation.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = VariationalFormulation.hxx; sourceTree = "<group>"; };
BEFC2B341BD7EBFF0067F358 /* Policy.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = Policy.hpp; sourceTree = "<group>"; };
BEFC2B351BD7EBFF0067F358 /* Policy.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = Policy.hxx; sourceTree = "<group>"; };
BEFC2B381BD7EDAB0067F358 /* VariationalFormulation.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = VariationalFormulation.hpp; sourceTree = "<group>"; };
BEFC2B391BD7EDAB0067F358 /* VariationalFormulation.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = VariationalFormulation.hxx; sourceTree = "<group>"; };
BEFF48FC1859C519005BB48A /* TestInputOutputMesh */ = {isa = PBXFileReference; explicitFileType = "compiled.mach-o.executable"; includeInIndex = 0; path = TestInputOutputMesh; sourceTree = BUILT_PRODUCTS_DIR; };
/* End PBXFileReference section */
 
......@@ -6157,6 +6161,7 @@
isa = PBXGroup;
children = (
BEFC2B291BD7D6770067F358 /* Elasticity */,
BEFC2B321BD7E54F0067F358 /* Hyperelasticity */,
BE65A1DD1BBD2C8F00E99115 /* Hyperelasticity.hpp */,
BE65A1DE1BBD2C8F00E99115 /* Hyperelasticity.hxx */,
);
......@@ -7926,6 +7931,18 @@
path = SolidVariationalFormulationPolicy/Elasticity;
sourceTree = "<group>";
};
BEFC2B321BD7E54F0067F358 /* Hyperelasticity */ = {
isa = PBXGroup;
children = (
BEFC2B341BD7EBFF0067F358 /* Policy.hpp */,
BEFC2B351BD7EBFF0067F358 /* Policy.hxx */,
BEFC2B381BD7EDAB0067F358 /* VariationalFormulation.hpp */,
BEFC2B391BD7EDAB0067F358 /* VariationalFormulation.hxx */,
);
name = Hyperelasticity;
path = SolidVariationalFormulationPolicy/Hyperelasticity;
sourceTree = "<group>";
};
/* End PBXGroup section */
 
/* Begin PBXHeadersBuildPhase section */
......
......@@ -21,8 +21,6 @@ namespace HappyHeart
namespace Policy
{
inline const typename Elasticity::variational_formulation_type&
Elasticity::GetVariationalFormulation() const noexcept
......
//
// Policy.hpp
// HappyHeart
//
// Created by Sebastien Gilles on 21/10/15.
// Copyright © 2015 Inria. All rights reserved.
//
#ifndef HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_POLICY_HPP_
# define HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_POLICY_HPP_
# include <memory>
# include <vector>
# include "ModelInstances/Hyperelasticity/TimeSchemes/TimeScheme.hpp"
# include "ModelInstances/FSI_EI/SolidVariationalFormulationPolicy/Hyperelasticity/VariationalFormulation.hpp"
# include "ModelInstances/FSI_EI/InputParameterList.hpp"
namespace HappyHeart
{
namespace FSI_EINS
{
namespace Policy
{
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
class Hyperelasticity
{
public:
//! Alias to self.
using self = Hyperelasticity<LawPolicyT, TimeSchemeT>;
//! Alias to unique pointer.
using unique_ptr = std::unique_ptr<self>;
//! Alias to hyperelastic variational formulation.
using variational_formulation_type = VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>;
public:
/// \name Special members.
///@{
//! Constructor.
explicit Hyperelasticity() = default;
protected:
//! Destructor.
~Hyperelasticity() = default;
//! Copy constructor.
Hyperelasticity(const Hyperelasticity&) = delete;
//! Move constructor.
Hyperelasticity(Hyperelasticity&&) = delete;
//! Copy affectation.
Hyperelasticity& operator=(const Hyperelasticity&) = delete;
//! Move affectation.
Hyperelasticity& operator=(Hyperelasticity&&) = delete;
///@}
protected:
/*!
* \brief Constructor.
*
* \param[in] mpi Object that manages mpi informations.
* \param[in] time_manager Object in charge of time management.
* \param[in] god_of_dof God of dof into which the formulation works.
* \param[in] main_felt_space Finite element space upon which the variational formulation apply.
* \param[in] solid_displacement Unknown considered; it should be a solid displacement.
* \param[in] dof_source Contribution to the RHS that stems from the solution of the fluid variational
* formulations.
*
*/
void CreateSolidFormulation(const Wrappers::Mpi& mpi,
const InputParameterList& input_parameter_data,
const FEltSpace& main_felt_space,
const Unknown& solid_displacement,
const NumberingSubset& numbering_subset,
const TimeManager& time_manager,
const GodOfDof& god_of_dof,
DirichletBoundaryCondition::vector_shared_ptr&& boundary_condition_list,
const GlobalVector& dof_source);
//! Call the appropriate solver to deal with the equation.
void Solve();
/*!
* \brief Action to do in the variational formulation when Model::FinalizeStep() is called
*
* \param[in] newest_displacement Displacement that has just been computed at the end of the implicit
* loop.
* \param[in] former_displacement Actually the \a newest_displacement obtained one time iteration ago.
*/
void FinalizeStepForFormulation(const GlobalVector& newest_displacement,
const GlobalVector& former_displacement,
const std::string& output_dir);
protected:
//! Underlying variational formulation.
const variational_formulation_type& GetVariationalFormulation() const noexcept;
//! Underlying variational formulation.
variational_formulation_type& GetNonCstVariationalFormulation() noexcept;
private:
//! Underlying variational formulation for solid implicit step.
typename variational_formulation_type::unique_ptr variational_formulation_ = nullptr;
};
} //namespace Policy
} //namespace FSI_EINS
} // namespace HappyHeart
# include "ModelInstances/FSI_EI/SolidVariationalFormulationPolicy/Hyperelasticity/Policy.hxx"
#endif // HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_POLICY_HPP_
//
// Policy.hxx
// HappyHeart
//
// Created by Sebastien Gilles on 21/10/15.
// Copyright © 2015 Inria. All rights reserved.
//
#ifndef HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_POLICY_HXX_
# define HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_POLICY_HXX_
namespace HappyHeart
{
namespace FSI_EINS
{
namespace Policy
{
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
void Hyperelasticity<LawPolicyT, TimeSchemeT>
::CreateSolidFormulation(const Wrappers::Mpi& mpi,
const InputParameterList& input_parameter_data,
const FEltSpace& main_felt_space,
const Unknown& solid_displacement,
const NumberingSubset& numbering_subset,
const TimeManager& time_manager,
const GodOfDof& god_of_dof,
DirichletBoundaryCondition::vector_shared_ptr&& boundary_condition_list,
const GlobalVector& dof_source)
{
variational_formulation_ =
std::make_unique<variational_formulation_type>(mpi,
main_felt_space,
solid_displacement,
numbering_subset,
time_manager,
god_of_dof,
std::move(boundary_condition_list),
1., // contribution of DofSource must be removed from rhs in elastic case.
dof_source);
variational_formulation_->Init(input_parameter_data);
}
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
inline const typename Hyperelasticity<LawPolicyT, TimeSchemeT>::variational_formulation_type&
Hyperelasticity<LawPolicyT, TimeSchemeT>
::GetVariationalFormulation() const noexcept
{
assert(!(!variational_formulation_));
return *variational_formulation_;
}
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
inline typename Hyperelasticity<LawPolicyT, TimeSchemeT>::variational_formulation_type&
Hyperelasticity<LawPolicyT, TimeSchemeT>
::GetNonCstVariationalFormulation() noexcept
{
return const_cast<variational_formulation_type&>(GetVariationalFormulation());
}
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
void Hyperelasticity<LawPolicyT, TimeSchemeT>::Solve()
{
auto& variational_formulation = GetNonCstVariationalFormulation();
const auto& numbering_subset = variational_formulation.GetNumberingSubset();
static bool first = true;
if (first)
{
std::cout << "IMPLICIT SOLID MATRIX = " << std::endl;
variational_formulation.GetSystemMatrix(numbering_subset, numbering_subset).View(variational_formulation.MpiHappyHeart(), "/Users/sebastien/Desktop/ImplicitSolid.m",
__FILE__, __LINE__,
PETSC_VIEWER_ASCII_MATLAB
);
first = false;
}
variational_formulation.SolveNonLinear(numbering_subset, numbering_subset);
}
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
void Hyperelasticity<LawPolicyT, TimeSchemeT>::FinalizeStepForFormulation(const GlobalVector& newest_displacement,
const GlobalVector& former_displacement,
const std::string& output_dir)
{
// auto& variational_formulation = GetNonCstVariationalFormulation();
// variational_formulation.UpdatePreviousTimeIterationData(newest_displacement,
// former_displacement,
// output_dir);
}
} //namespace Policy
} //namespace FSI_EINS
} // namespace HappyHeart
#endif // HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_POLICY_HXX_
//
// VariationalFormulation.hpp
// HappyHeart
//
// Created by Sebastien Gilles on 21/10/15.
// Copyright © 2015 Inria. All rights reserved.
//
#ifndef HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_VARIATIONAL_FORMULATION_HPP_
# define HAPPY_HEART_x_MODEL_INSTANCES_x_F_S_I_x_E_I_x_SOLID_VARIATIONAL_FORMULATION_POLICY_x_HYPERELASTICITY_x_VARIATIONAL_FORMULATION_HPP_
# include "Core/TimeManager/TimeManagerInstance.hpp"
# include "Core/TimeManager/Policy/ConstantTimeStep.hpp"
# include "Core/InputParameter/Parameter/Solid/Solid.hpp"
# include "Operators/LocalVariationalOperatorInstances/NonlinearForm/Hyperelasticity/InvariantManager.hpp"
# include "Operators/LocalVariationalOperatorInstances/NonlinearForm/Hyperelasticity/Private/Helper.hpp"
# include "Operators/GlobalVariationalOperatorInstances/BilinearForm/Mass.hpp"
# include "Operators/GlobalVariationalOperatorInstances/LinearForm/TransientSource.hpp"
# include "Operators/GlobalVariationalOperatorInstances/NonlinearForm/GradOnGradientBasedHyperelasticityTensor.hpp"
# include "FormulationSolver/VariationalFormulation.hpp"
# include "FormulationSolver/DofSourcePolicy/DofSource.hpp"
# include "ModelInstances/Hyperelasticity/Impl/SnesFunctions.hpp"
# include "ModelInstances/Hyperelasticity/TimeSchemes/TimeScheme.hpp"
# include "ModelInstances/Hyperelasticity/Impl/SnesFunctions.hpp"
# include "ModelInstances/Hyperelasticity/TimeSchemes/HalfSum/HalfSum.hpp"
# include "ModelInstances/Hyperelasticity/TimeSchemes/Midpoint/Midpoint.hpp"
# include "ModelInstances/FSI_EI/InputParameterList.hpp"
namespace HappyHeart
{
namespace FSI_EINS
{
namespace Policy
{
/*!
* \brief Class in charge of hyperelasticity problem.
*
* This class is policy-based: the law actually used is specified by a template parameter
*
* For instance, to use CiarletGeymonat laws you have to call VariationalFormulationHyperElasticity<CiarletGeymonat> (and you must
* of course include the file in which CiarletGeymonat is defined).
*
* \tparam TimeSchemeT Approximation used to handle non-linear term.
*/
template
<
class LawPolicyT,
HyperelasticityNS::TimeScheme TimeSchemeT
>
class VariationalFormulationHyperElasticity final
: public VariationalFormulation<VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>>,
private ::HappyHeart::Crtp::Mutex<VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>>,
public VariationalFormulationNS::DofSourcePolicyNS::DofSource
{
private:
//! Convenient alias.
using self = VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>;
//! Alias to the parent class.
using parent = VariationalFormulation<self>;
//! Friendship to parent class, so this one can access private methods defined below through CRTP.
friend parent;
//! Friendship to the class which implements the prototyped functions required by Petsc Snes algorithm.
friend struct HyperelasticityNS::Impl::Snes<self>;
using dof_source_parent = VariationalFormulationNS::DofSourcePolicyNS::DofSource;
using vm_type = HyperelasticityNS::VectorsAndMatrices<TimeSchemeT>;
public:
//! Alias to unique_ptr.
using unique_ptr = std::unique_ptr<self>;
// \todo #501 Add a safety to ensure time step is constant! (model relies on this assumption).
using time_manager_type = TimeManager;//::HappyHeart::Private::TimeManager<TimeManagerNS::Policy::ConstantTimeStep>;
public:
/// \name Special members.
///@{
/*!
* \brief Constructor.
*
* \param[in] mpi Object that manages mpi informations.
* \param[in] time_manager Object in charge of time management.
* \param[in] god_of_dof God of dof into which the formulation works.
* \param[in] main_felt_space Finite element space upon which the variational formulation apply.
* \param[in] neumann_felt_space Finite element space upon which the Neumann condition apply.
* \param[in] solid_displacement Unknown considered; it should be a solid displacement.
* \param[in] dof_source_policy_args Supplementary arguments that might be required by the DofSourcePolicyT.
*
*/
explicit VariationalFormulationHyperElasticity(const Wrappers::Mpi& mpi,
const FEltSpace& main_felt_space,
const Unknown& solid_displacement,
const NumberingSubset& numbering_subset,
const time_manager_type& time_manager,
const GodOfDof& god_of_dof,
DirichletBoundaryCondition::vector_shared_ptr&& boundary_condition_list,
double factor,
const GlobalVector& dof_source);
//! Destructor.
~VariationalFormulationHyperElasticity() = default;
//! Copy constructor.
VariationalFormulationHyperElasticity(const VariationalFormulationHyperElasticity&) = delete;
//! Move constructor.
VariationalFormulationHyperElasticity(VariationalFormulationHyperElasticity&&) = delete;
//! Copy affectation.
VariationalFormulationHyperElasticity& operator=(const VariationalFormulationHyperElasticity&) = delete;
//! Move affectation.
VariationalFormulationHyperElasticity& operator=(VariationalFormulationHyperElasticity&&) = delete;
///@}
//! Prepare dynamic runs.
void PrepareDynamicRuns();
//! Update for next time step. (not called after each dynamic iteration).
void UpdateForNextTimeStep();
/*!
* \brief Get the only numbering subset relevant for this VariationalFormulation.
*
* There is a more generic accessor in the base class but is use is more unwieldy.
*/
const NumberingSubset& GetNumberingSubset() const;
private:
/// \name CRTP-required methods.
///@{
/*!
* \brief Specific initialisation for derived class attributes.
*
* \internal This method is called by base class method VariationalFormulation::Init().
*/
template<class InputParameterDataT>
void SupplInit(const InputParameterDataT& input_parameter_data);
/*!
* \brief Allocate the global matrices and vectors.
*/
void AllocateMatricesAndVectors();
//! Define the pointer function required to calculate the function required by the non-linear problem.
Wrappers::Petsc::Snes::SNESFunction ImplementSnesFunction() const;
//! Define the pointer function required to calculate the jacobian required by the non-linear problem.
Wrappers::Petsc::Snes::SNESJacobian ImplementSnesJacobian() const;
//! Define the pointer function required to view the results required by the non-linear problem.
Wrappers::Petsc::Snes::SNESViewer ImplementSnesViewer() const;
///@}
private:
/*!
* \brief Define the properties of all the global variational operators involved.
*
* \param[in] input_parameter_data Object which hold the values of all the parameters defined in
* the input file.
*/
template<class InputParameterDataT>
void DefineOperators(const InputParameterDataT& input_parameter_data);
/*!
* \brief Update the content of all the vectors and matrices relevant to the computation of the tangent
* and the residual.
*/
void UpdateVectorsAndMatrices();
//! An helper method of \a UpdateVectorsAndMatrices().
void AssembleVectorsAndMatrices();
/*!
* \brief Compute the matrix of the system.
*
* The index of the time iteration is used to know whether the static or dynamic system muste be computed.
*/
void ComputeTangent();
/*!
* \brief Compute the RHS of the system.
*
* The index of the time iteration is used to know whether the static or dynamic system muste be computed.
*/
void ComputeResidual();
void ComputeStaticTangent();
void ComputeStaticResidual();
void ComputeDynamicTangent();
void ComputeDynamicResidual();
//! Prepare vectors and matrices before assembling.
void PrepareMatricesAndVectorsBeforeAssembling();
//! Update current displacement. Already called in UpdateForNextTimeStep().
void UpdateDisplacementBetweenTimeStep();