Commit 14a38861 authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#820 Poromechanics/StVenantKirchhoff: law in progress.

parent a2953e61
......@@ -455,7 +455,12 @@ Result = {
-- Directory in which all the results will be written.
-- Expected format: "VALUE"
output_directory = '/Volumes/Data/${USER}/HappyHeart/Results/Poromechanics'
output_directory = '/Volumes/Data/${USER}/HappyHeart/Results/Poromechanics',
-- Enables to skip some printing in the console. Can be used to WriteSolution every n time.
-- Expected format: VALUE
display_value = 1
} -- Result
......@@ -729,6 +734,41 @@ Solid = {
} -- Solid
bulk_solid = {
-- How is given the volumic mass (as a constant, as a Lua function, per quadrature point, etc...)
-- Expected format: "VALUE"
-- Constraint: ops_in(v, {'constant', 'lua_function', 'at_quadrature_point', 'piecewise_constant_by_domain'})
nature = 'constant',
-- Value of the volumic mass in the case nature is 'constant'(and also initial value if nature is
-- 'at_quadrature_point'; irrelevant otherwise).
-- Expected format: VALUE
scalar_value = 1.e3,
-- Value of the volumic mass in the case nature is 'lua_function'(and also initial value if nature is
-- 'at_quadrature_point'; irrelevant otherwise).
-- Expected format: Function in Lua language, for instance:
-- -- function(arg1, arg2, arg3)
-- return arg1 + arg2 -
-- arg3
-- end
-- sin, cos and tan require a 'math.'
-- preffix.
-- If you do not wish to provide one, put anything you want (e.g. 'none'): the
-- content is not interpreted by Ops until an actual use of the underlying function.
lua_function = none,
-- Domain indices of the parameter in the case nature is 'piecewise_constant_by_domain'. The various
-- domains given here must not intersect.
-- Expected format: {VALUE1, VALUE2, ...}
piecewise_constant_domain_id = {},
-- Value of the parameter in the case nature is 'piecewise_constant_by_domain'.
-- Expected format: {VALUE1, VALUE2, ...}
piecewise_constant_domain_value = {}
} -- bulk_solid
-- Initial value of porosity at each dof.
......
MeshVersionFormatted 1
MeshVersionFormatted 2
Dimension
2
......
......@@ -679,6 +679,7 @@
BE794E841ABC6FCC003D9F2D /* DomainManager.hpp in Headers */ = {isa = PBXBuildFile; fileRef = BE794E811ABC6FCC003D9F2D /* DomainManager.hpp */; };
BE794E851ABC6FCC003D9F2D /* DomainManager.hxx in Headers */ = {isa = PBXBuildFile; fileRef = BE794E821ABC6FCC003D9F2D /* DomainManager.hxx */; };
BE83AF341C733EAD002C6FA3 /* StVenantKirchhoff.cpp in Sources */ = {isa = PBXBuildFile; fileRef = BE83AF311C733EAD002C6FA3 /* StVenantKirchhoff.cpp */; };
BE83AF3F1C73669D002C6FA3 /* BulkSolid.cpp in Sources */ = {isa = PBXBuildFile; fileRef = BE83AF3C1C73669D002C6FA3 /* BulkSolid.cpp */; };
BE8553A91BBD77BF00DB109E /* DofSource.cpp in Sources */ = {isa = PBXBuildFile; fileRef = BE8553A31BBD77BF00DB109E /* DofSource.cpp */; };
BE8553AA1BBD77BF00DB109E /* DofSource.hpp in Headers */ = {isa = PBXBuildFile; fileRef = BE8553A41BBD77BF00DB109E /* DofSource.hpp */; };
BE8553AB1BBD77BF00DB109E /* DofSource.hxx in Headers */ = {isa = PBXBuildFile; fileRef = BE8553A51BBD77BF00DB109E /* DofSource.hxx */; };
......@@ -4536,6 +4537,8 @@
BE83AF311C733EAD002C6FA3 /* StVenantKirchhoff.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; name = StVenantKirchhoff.cpp; path = HyperelasticLaw/StVenantKirchhoff.cpp; sourceTree = "<group>"; };
BE83AF321C733EAD002C6FA3 /* StVenantKirchhoff.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = StVenantKirchhoff.hpp; path = HyperelasticLaw/StVenantKirchhoff.hpp; sourceTree = "<group>"; };
BE83AF331C733EAD002C6FA3 /* StVenantKirchhoff.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = StVenantKirchhoff.hxx; path = HyperelasticLaw/StVenantKirchhoff.hxx; sourceTree = "<group>"; };
BE83AF3C1C73669D002C6FA3 /* BulkSolid.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; name = BulkSolid.cpp; path = InputParameter/BulkSolid.cpp; sourceTree = "<group>"; };
BE83AF3D1C73669D002C6FA3 /* BulkSolid.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = BulkSolid.hpp; path = InputParameter/BulkSolid.hpp; sourceTree = "<group>"; };
BE8553A31BBD77BF00DB109E /* DofSource.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = DofSource.cpp; sourceTree = "<group>"; };
BE8553A41BBD77BF00DB109E /* DofSource.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = DofSource.hpp; sourceTree = "<group>"; };
BE8553A51BBD77BF00DB109E /* DofSource.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = DofSource.hxx; sourceTree = "<group>"; };
......@@ -8064,6 +8067,8 @@
children = (
BEAEA8C41C71DFD30054FEF9 /* Porosity.cpp */,
BEAEA8C51C71DFD30054FEF9 /* Porosity.hpp */,
BE83AF3C1C73669D002C6FA3 /* BulkSolid.cpp */,
BE83AF3D1C73669D002C6FA3 /* BulkSolid.hpp */,
);
name = InputParameter;
sourceTree = "<group>";
......@@ -11568,6 +11573,7 @@
BE5842171C6B8C980040C694 /* Mass.cpp in Sources */,
BE0B0BE11C6A196C0038A9B9 /* main.cpp in Sources */,
BEF183281C6DE3A3008A6F1E /* Ale.cpp in Sources */,
BE83AF3F1C73669D002C6FA3 /* BulkSolid.cpp in Sources */,
);
runOnlyForDeploymentPostprocessing = 0;
};
......
......@@ -11,7 +11,6 @@
# include "Core/InputParameter/Crtp/Section.hpp"
# include "Core/InputParameter/Parameter/SpatialFunction.hpp"
# include "Core/InputParameter/Parameter/Private/ParameterUsualDescription.hpp"
# include "Core/InputParameter/Parameter/MaterialProperty/VolumicMass.hpp"
......@@ -36,22 +35,15 @@ namespace HappyHeart
//! Convenient alias.
using self = Fluid;
//! Friendship to section parent.
using parent = Crtp::Section<self, NoEnclosingSection>;
friend parent;
using VolumicMass = MaterialProperty::VolumicMass<self>;
struct Viscosity : public Crtp::Section<Viscosity, Fluid>
{
......
......@@ -134,7 +134,8 @@ namespace HappyHeart
const std::string& PiecewiseConstantByDomainId::DefaultValue()
{
return Utilities::EmptyString();
static std::string ret("{}");
return ret;
}
......@@ -159,12 +160,11 @@ namespace HappyHeart
const std::string& PiecewiseConstantByDomainValue::DefaultValue()
{
return Utilities::EmptyString();
static std::string ret("{}");
return ret;
}
} // namespace Private
......
......@@ -45,7 +45,6 @@ namespace HappyHeart
>
template<class InputParameterDataT>
void VariationalFormulation<DerivedT, SolverIndexT>
::Init(const InputParameterDataT& input_parameter_data)
{
auto& derived = static_cast<DerivedT&>(*this);
......@@ -446,7 +445,6 @@ namespace HappyHeart
unsigned int SolverIndexT
>
void VariationalFormulation<DerivedT, SolverIndexT>
::AllocateGlobalMatrix(const NumberingSubset& row_numbering_subset,
const NumberingSubset& col_numbering_subset)
{
......@@ -467,7 +465,6 @@ namespace HappyHeart
unsigned int SolverIndexT
>
void VariationalFormulation<DerivedT, SolverIndexT>
::AllocateGlobalVector(const NumberingSubset& numbering_subset)
{
auto& rhs_storage = GetNonCstRhsVectorStorage();
......@@ -625,10 +622,7 @@ namespace HappyHeart
// BOUNDARY CONDITIONS PRIVATE METHODS //
///////////////////////////////////////////////
template
<
class DerivedT,
......@@ -636,7 +630,6 @@ namespace HappyHeart
>
template<OperatorNS::Nature AppliedOnT>
void VariationalFormulation<DerivedT, SolverIndexT>
::ApplyEssentialBoundaryConditionPseudoElimination(const NumberingSubset& row_numbering_subset,
const NumberingSubset& col_numbering_subset)
{
......@@ -687,7 +680,6 @@ namespace HappyHeart
>
template<OperatorNS::Nature AppliedOnT>
[[noreturn]] void VariationalFormulation<DerivedT, SolverIndexT>
::ApplyEssentialBoundaryConditionPenalisation(const NumberingSubset& /*row_numbering_subset*/,
const NumberingSubset& /*col_numbering_subset*/)
{
......@@ -736,7 +728,6 @@ namespace HappyHeart
unsigned int SolverIndexT
>
const Wrappers::Petsc::Snes& VariationalFormulation<DerivedT, SolverIndexT>
::GetSnes() const
{
assert(!(!snes_));
......@@ -750,7 +741,6 @@ namespace HappyHeart
unsigned int SolverIndexT
>
Wrappers::Petsc::Snes& VariationalFormulation<DerivedT, SolverIndexT>
::GetNonCstSnes()
{
assert(!(!snes_));
......@@ -870,7 +860,6 @@ namespace HappyHeart
unsigned int SolverIndexT
>
const DirichletBoundaryCondition::vector_shared_ptr& VariationalFormulation<DerivedT, SolverIndexT>
::GetEssentialBoundaryConditionList() const noexcept
{
return boundary_condition_list_;
......
//
// BulkSolid.cpp
// HappyHeart
//
// Created by Sebastien Gilles on 16/02/16.
// Copyright © 2016 Inria. All rights reserved.
//
#include "Utilities/String/EmptyString.hpp"
#include "ModelInstances/UnderDevelopment/Poromechanics/InputParameter/BulkSolid.hpp"
namespace HappyHeart
{
namespace InputParameter
{
namespace PoromechanicsNS
{
const std::string& BulkSolid::GetName()
{
static std::string ret("bulk_solid");
return ret;
};
} // namespace PoromechanicsNS
} // namespace InputParameter
} // namespace HappyHeart
//
// BulkSolid.hpp
// HappyHeart
//
// Created by Sebastien Gilles on 16/02/16.
// Copyright © 2016 Inria. All rights reserved.
//
#ifndef HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_INPUT_PARAMETER_x_BULK_SOLID_HPP_
# define HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_INPUT_PARAMETER_x_BULK_SOLID_HPP_
# include "Core/InputParameter/Crtp/Section.hpp"
# include "Utilities/InputParameterList/Crtp/InputParameter.hpp"
# include "Core/InputParameter/Parameter/Private/ParameterUsualDescription.hpp"
namespace HappyHeart
{
namespace InputParameter
{
namespace PoromechanicsNS
{
struct BulkSolid : public Crtp::Section<BulkSolid, NoEnclosingSection>
{
//! Convenient alias.
using self = BulkSolid;
//! Friendship to section parent.
using parent = Crtp::Section<self, NoEnclosingSection>;
friend parent;
/*!
* \brief Return the name of the section in the input parameter.
*
*/
static const std::string& GetName();
/*!
* \brief Choose how is described the viscosity (through a scalar, a function, etc...)
*/
struct Nature : public Crtp::InputParameter<Nature, self, Private::Nature::storage_type>,
public Private::Nature
{ };
/*!
* \brief Scalar value. Irrelevant if nature is not scalar.
*/
struct Scalar : public Crtp::InputParameter<Scalar, self, Private::Scalar::storage_type>,
public Private::Scalar
{ };
/*!
* \brief Function that determines viscosity value. Irrelevant if nature is not lua_function.
*/
struct LuaFunction : public Crtp::InputParameter<LuaFunction, self, Private::LuaFunction::storage_type>,
public Private::LuaFunction
{ };
/*!
* \brief Piecewise Constant domain index.
*/
struct PiecewiseConstantByDomainId : public Crtp::InputParameter<PiecewiseConstantByDomainId, self, Private::PiecewiseConstantByDomainId::storage_type>,
public Private::PiecewiseConstantByDomainId
{ };
/*!
* \brief Piecewise Constant value by domain.
*/
struct PiecewiseConstantByDomainValue : public Crtp::InputParameter<PiecewiseConstantByDomainValue, self, Private::PiecewiseConstantByDomainValue::storage_type>,
public Private::PiecewiseConstantByDomainValue
{ };
//! Alias to the tuple of structs.
using section_content_type = std::tuple
<
Nature,
Scalar,
LuaFunction,
PiecewiseConstantByDomainId,
PiecewiseConstantByDomainValue
>;
private:
//! Content of the section.
section_content_type section_content;
}; // struct BulkSolid
} // namespace PoromechanicsNS
} // namespace InputParameter
} // namespace HappyHeart
#endif // HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_INPUT_PARAMETER_x_BULK_SOLID_HPP_
......@@ -21,6 +21,7 @@
# include "Core/InputParameter/Parameter/Solid/Solid.hpp"
# include "ModelInstances/UnderDevelopment/Poromechanics/InputParameter/Porosity.hpp"
# include "ModelInstances/UnderDevelopment/Poromechanics/InputParameter/BulkSolid.hpp"
namespace HappyHeart
......@@ -119,14 +120,16 @@ namespace HappyHeart
InputParameter::Petsc<EnumUnderlyingType(SolverIndex::main_solver)>,
// InputParameter::Petsc<EnumUnderlyingType(SolverIndex::dh_solver)>,
InputParameter::Result::OutputDirectory,
InputParameter::Result,
InputParameter::Solid,
InputParameter::Solid::PlaneStressStrain,
InputParameter::Fluid::Density,
InputParameter::PoromechanicsNS::Porosity
InputParameter::PoromechanicsNS::Porosity,
InputParameter::PoromechanicsNS::BulkSolid
>;
......
......@@ -36,14 +36,83 @@ namespace HappyHeart
const ScalarParameter& poisson_ratio,
const ScalarParameter& kappa1,
const ScalarParameter& kappa2,
const ScalarParameter& bulk)
const ScalarParameter& bulk,
// const Pa
const ScalarParameter& bulk_solid)
: without_poro_parent(young_modulus,
poisson_ratio,
kappa1,
kappa2,
bulk)
bulk),
bulk_solid_(bulk_solid)
{ }
double StVenantKirchhoff::W(const QuadraturePoint& quadrature_point, const GeometricElt& geom_elt) const
{
double ret = without_poro_parent::W(quadrature_point, geom_elt);
const double sqrt_I3 = std::sqrt(without_poro_parent::GetInvariant<3>());
const double bulk_solid = GetBulkSolid().GetValue(quadrature_point, geom_elt);
const double lame_lambda_value = without_poro_parent::GetLameLambda().GetValue(quadrature_point, geom_elt);
const double lame_mu_value = without_poro_parent::GetLameMu().GetValue(quadrature_point, geom_elt);
const double bulk = without_poro_parent::GetBulk().GetValue(quadrature_point, geom_elt);
// ret += bulk_solid
// * (sqrt_I3 -
//
// )
// ( bulkSolid*(1-isQuad)*( (I3^(1./2.)-fluidmass/rhoF)/(1.-phi0) - 1. - log( abs((I3^(1./2.)-fluidmass/rhoF)/(1.-phi0)) ) )
//
return ret;
}
double StVenantKirchhoff::FirstDerivativeWThirdInvariant(const QuadraturePoint& quadrature_point,
const GeometricElt& geom_elt) const
{
double ret = without_poro_parent::FirstDerivativeWThirdInvariant(quadrature_point, geom_elt);
return ret;
}
double StVenantKirchhoff::SecondDerivativeWThirdInvariant(const QuadraturePoint& quadrature_point,
const GeometricElt& geom_elt) const
{
double ret = without_poro_parent::SecondDerivativeWThirdInvariant(quadrature_point, geom_elt);
return ret;
}
double StVenantKirchhoff::SecondDerivativeWFirstAndThirdInvariant(const QuadraturePoint& quadrature_point,
const GeometricElt& geom_elt) const
{
double ret = without_poro_parent::SecondDerivativeWFirstAndThirdInvariant(quadrature_point, geom_elt);
return ret;
}
double StVenantKirchhoff::SecondDerivativeWSecondAndThirdInvariant(const QuadraturePoint& quadrature_point,
const GeometricElt& geom_elt) const
{
double ret = without_poro_parent::SecondDerivativeWSecondAndThirdInvariant(quadrature_point, geom_elt);
return ret;
}
......
......@@ -6,8 +6,8 @@
// Copyright © 2016 Inria. All rights reserved.
//
#ifndef __HappyHeart__StVenantKirchhoff__HPP
# define __HappyHeart__StVenantKirchhoff__HPP
#ifndef HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_LOCAL_VARIATIONAL_INSTANCES_x_HYPERELASTIC_LAW_x_ST_VENANT_KIRCHHOFF_HPP_
# define HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_LOCAL_VARIATIONAL_INSTANCES_x_HYPERELASTIC_LAW_x_ST_VENANT_KIRCHHOFF_HPP_
# include <memory>
# include <vector>
......@@ -69,8 +69,8 @@ namespace HappyHeart
const ScalarParameter& poisson_ratio,
const ScalarParameter& kappa1,
const ScalarParameter& kappa2,
const ScalarParameter& bulk);
const ScalarParameter& bulk,
const ScalarParameter& bulk_solid);
//! Destructor.
~StVenantKirchhoff() = default;
......@@ -126,6 +126,16 @@ namespace HappyHeart
double SecondDerivativeWSecondAndThirdInvariant(const QuadraturePoint& quadrature_point,
const GeometricElt& geom_elt) const;
private:
//! Bulk solid.
const ScalarParameter& GetBulkSolid() const noexcept;
private:
//! Bulk solid.
const ScalarParameter& bulk_solid_;
};
......@@ -146,4 +156,4 @@ namespace HappyHeart
# include "ModelInstances/UnderDevelopment/Poromechanics/LocalVariationalInstances/HyperelasticLaw/StVenantKirchhoff.hxx"
#endif /* defined(__HappyHeart__StVenantKirchhoff__HPP) */
#endif // HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_LOCAL_VARIATIONAL_INSTANCES_x_HYPERELASTIC_LAW_x_ST_VENANT_KIRCHHOFF_HPP_
......@@ -6,8 +6,8 @@
// Copyright © 2016 Inria. All rights reserved.
//
#ifndef __HappyHeart__StVenantKirchhoff__HXX
#define __HappyHeart__StVenantKirchhoff__HXX
#ifndef HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_LOCAL_VARIATIONAL_INSTANCES_x_HYPERELASTIC_LAW_x_ST_VENANT_KIRCHHOFF_HXX_
# define HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_LOCAL_VARIATIONAL_INSTANCES_x_HYPERELASTIC_LAW_x_ST_VENANT_KIRCHHOFF_HXX_
namespace HappyHeart
......@@ -26,6 +26,12 @@ namespace HappyHeart
{
inline const ScalarParameter& StVenantKirchhoff::GetBulkSolid() const noexcept
{
return bulk_solid_;
}
} // namespace HyperelasticLawNS
......@@ -39,4 +45,4 @@ namespace HappyHeart
} // namespace HappyHeart
#endif /* defined(__HappyHeart__StVenantKirchhoff__HXX) */
#endif // HAPPY_HEART_x_MODEL_INSTANCES_x_UNDER_DEVELOPMENT_x_POROMECHANICS_x_LOCAL_VARIATIONAL_INSTANCES_x_HYPERELASTIC_LAW_x_ST_VENANT_KIRCHHOFF_HXX_
......@@ -285,7 +285,7 @@ namespace HappyHeart
const auto size = porosity.GetSize(__FILE__, __LINE__);
const double fluid_density = GetFluidDensity().GetConstantValue();
assert(!Utilities::IsZero(fluid_density));
assert(!NumericNS::IsZero(fluid_density));
const double inv_fluid_density = 1. / fluid_density;
const double initial_porosity_value = GetInitialPorosityValue();
......