OGS
MaterialLib::Solids::MFront Namespace Reference

Namespaces

namespace  detail
 

Classes

struct  DeformationGradient
 Meta data for deformation gradient. More...
 
struct  ForcesGradsCombinations
 
struct  GreenLagrangeStrain
 Meta data for Green-Lagrange strain. More...
 
struct  LiquidPressure
 
struct  MaterialStateVariablesMFront
 
class  MFront
 
struct  MFrontConfig
 
class  MFrontGeneric
 
class  OGSMFrontTangentOperatorBlocksView
 
struct  OGSMFrontTangentOperatorData
 Used for disambiguation with MFront's thermodynamic forces data. More...
 
struct  OGSMFrontThermodynamicForcesData
 Used for disambiguation with MFront's tangent operator blocks data. More...
 
class  OGSMFrontThermodynamicForcesView
 
struct  Saturation
 
struct  SecondPiolaKirchhoffStress
 
struct  Strain
 Meta data for strain. More...
 
struct  Stress
 
struct  Temperature
 
struct  Variable
 

Functions

template<int DisplacementDim>
std::unique_ptr< MechanicsBase< DisplacementDim > > createMFront (std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, std::optional< ParameterLib::CoordinateSystem > const &local_coordinate_system, BaseLib::ConfigTree const &config)
 
template std::unique_ptr< MechanicsBase< 2 > > createMFront< 2 > (std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, std::optional< ParameterLib::CoordinateSystem > const &local_coordinate_system, BaseLib::ConfigTree const &config)
 
template std::unique_ptr< MechanicsBase< 3 > > createMFront< 3 > (std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, std::optional< ParameterLib::CoordinateSystem > const &local_coordinate_system, BaseLib::ConfigTree const &config)
 
MFrontConfig createMFrontConfig (int const displacement_dim, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config, bool const library_path_is_relative_to_prj_file)
 
template<int DisplacementDim, typename Gradients , typename TDynForces , typename ExtStateVars >
std::unique_ptr< MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars > > createMFrontGeneric (std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, std::optional< ParameterLib::CoordinateSystem > const &local_coordinate_system, BaseLib::ConfigTree const &config, bool const library_path_is_relative_to_prj_file)
 
const char * varTypeToString (int v)
 
int getEquivalentPlasticStrainOffset (mgis::behaviour::Behaviour const &b)
 
template<int DisplacementDim>
OGSMFrontTangentOperatorData tangentOperatorDataMFrontToOGS (std::vector< double > const &mfront_data, std::optional< MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > const &Q, mgis::behaviour::Behaviour const &behaviour)
 
template OGSMFrontTangentOperatorData tangentOperatorDataMFrontToOGS< 2 > (std::vector< double > const &mfront_data, std::optional< MathLib::KelvinVector::KelvinMatrixType< 2 > > const &Q, mgis::behaviour::Behaviour const &behaviour)
 
template OGSMFrontTangentOperatorData tangentOperatorDataMFrontToOGS< 3 > (std::vector< double > const &mfront_data, std::optional< MathLib::KelvinVector::KelvinMatrixType< 3 > > const &Q, mgis::behaviour::Behaviour const &behaviour)
 
template<typename Derived >
constexpr auto eigenSwap45View (Eigen::MatrixBase< Derived > const &matrix)
 
template<int DisplacementDim, typename Derived >
constexpr auto ogsTensorToMFrontTensor (Eigen::MatrixBase< Derived > const &matrix)
 

Variables

static constexpr Strain strain
 Instance that can be used for overload resolution/template type deduction.
 
static constexpr GreenLagrangeStrain green_lagrange_strain
 Instance that can be used for overload resolution/template type deduction.
 
static constexpr DeformationGradient deformation_gradient
 Instance that can be used for overload resolution/template type deduction.
 
static constexpr LiquidPressure liquid_pressure
 
static constexpr Stress stress
 
static constexpr SecondPiolaKirchhoffStress second_piola_kirchhoff_stress
 
static constexpr Saturation saturation
 
static constexpr Temperature temperature
 

Function Documentation

◆ createMFront()

template<int DisplacementDim>
std::unique_ptr< MechanicsBase< DisplacementDim > > MaterialLib::Solids::MFront::createMFront ( std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
std::optional< ParameterLib::CoordinateSystem > const & local_coordinate_system,
BaseLib::ConfigTree const & config )

Definition at line 26 of file CreateMFront.cpp.

31{
32 bool const library_path_is_relative_to_prj_file = true;
33 auto conf = createMFrontConfig(DisplacementDim, parameters, config,
34 library_path_is_relative_to_prj_file);
35
36 return std::make_unique<MFront<DisplacementDim>>(
37 std::move(conf.behaviour), std::move(conf.material_properties),
38 std::move(conf.state_variables_initial_properties),
39 local_coordinate_system);
40}
MFrontConfig createMFrontConfig(int const displacement_dim, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config, bool const library_path_is_relative_to_prj_file)

References createMFrontConfig().

◆ createMFront< 2 >()

template std::unique_ptr< MechanicsBase< 2 > > MaterialLib::Solids::MFront::createMFront< 2 > ( std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
std::optional< ParameterLib::CoordinateSystem > const & local_coordinate_system,
BaseLib::ConfigTree const & config )

◆ createMFront< 3 >()

template std::unique_ptr< MechanicsBase< 3 > > MaterialLib::Solids::MFront::createMFront< 3 > ( std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
std::optional< ParameterLib::CoordinateSystem > const & local_coordinate_system,
BaseLib::ConfigTree const & config )

◆ createMFrontConfig()

MFrontConfig MaterialLib::Solids::MFront::createMFrontConfig ( int const displacement_dim,
std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
BaseLib::ConfigTree const & config,
bool const library_path_is_relative_to_prj_file )
Input File Parameter
material__solid__constitutive_relation__type
Input File Parameter
material__solid__constitutive_relation__MFront__library
Input File Parameter
material__solid__constitutive_relation__MFront__behaviour

Definition at line 324 of file CreateMFrontGeneric.cpp.

329{
330 INFO("### MFRONT ########################################################");
331
333 config.checkConfigParameter("type", "MFront");
334
335 auto const library_name =
337 config.getConfigParameterOptional<std::string>("library");
338 auto const lib_path =
339 library_name ? (library_path_is_relative_to_prj_file
341 *library_name)
342 : *library_name)
343 : "libOgsMFrontBehaviour";
344
345 mgis::behaviour::Hypothesis hypothesis;
346 if (displacement_dim == 2)
347 {
348 // TODO support the axial symmetry modelling hypothesis.
349 WARN(
350 "The model is defined in 2D. On the material level currently a "
351 "plane strain setting is used. In particular it is not checked if "
352 "axial symmetry or plane stress are assumed. Special material "
353 "behaviour for these settings is currently not supported.");
354 hypothesis = mgis::behaviour::Hypothesis::PLANESTRAIN;
355 }
356 else if (displacement_dim == 3)
357 {
358 hypothesis = mgis::behaviour::Hypothesis::TRIDIMENSIONAL;
359 }
360 else
361 {
362 OGS_FATAL("Displacement dim {} is not supported.", displacement_dim);
363 }
364
365 auto behaviour = loadBehaviour(
366 lib_path,
368 config.getConfigParameter<std::string>("behaviour"),
369 hypothesis);
370
371 INFO("Behaviour: `{:s}'.", behaviour.behaviour);
372 INFO("Hypothesis: `{:s}'.", mgis::behaviour::toString(hypothesis));
373 INFO("Source: `{:s}'.", behaviour.source);
374 INFO("TFEL version: `{:s}'.", behaviour.tfel_version);
375 INFO("Behaviour type: `{:s}'.", btypeToString(behaviour.btype));
376 INFO("Kinematic: `{:s}'.", toString(behaviour.kinematic));
377 INFO("Symmetry: `{:s}'.", toString(behaviour.symmetry));
378
379 varInfo("Mat. props.", behaviour.mps, hypothesis);
380 varInfo("Gradients", behaviour.gradients, hypothesis);
381 varInfo("Thdyn. forces", behaviour.thermodynamic_forces, hypothesis);
382 varInfo("Int. StVars.", behaviour.isvs, hypothesis);
383 varInfo("Ext. StVars.", behaviour.esvs, hypothesis);
384
385 // TODO read parameters from prj file, not yet (2018-11-05) supported by
386 // MGIS library.
387 varInfo("Real-valued parameters", behaviour.params);
388 varInfo("Integer parameters", behaviour.iparams);
389 varInfo("Unsigned parameters", behaviour.usparams);
390
391 INFO("#Tangent operator blocks: {}.", behaviour.to_blocks.size());
392 for (auto const& [var1, var2] : behaviour.to_blocks)
393 {
394 INFO(" --> ({}, {}).", var1.name, var2.name);
395 }
396
397 std::vector<ParameterLib::Parameter<double> const*> material_properties =
398 readMaterialProperties(behaviour, hypothesis, parameters, config);
399
400 std::map<std::string, ParameterLib::Parameter<double> const*>
401 state_variables_initial_properties =
402 readStateVariablesInitialValueProperties(behaviour, hypothesis,
403 parameters, config);
404
405 INFO("### MFRONT END ####################################################");
406
407 return {std::move(behaviour), std::move(material_properties),
408 std::move(state_variables_initial_properties)};
409}
#define OGS_FATAL(...)
Definition Error.h:26
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:40
std::string const & getProjectDirectory()
Returns the directory where the prj file resides.
std::string joinPaths(std::string const &pathA, std::string const &pathB)
std::vector< ParameterLib::Parameter< double > const * > readMaterialProperties(mgis::behaviour::Behaviour const &behaviour, mgis::behaviour::Hypothesis const &hypothesis, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config)
void varInfo(std::string const &msg, std::vector< mgis::behaviour::Variable > const &vars, mgis::behaviour::Hypothesis hypothesis)
Prints info about MFront variables.
const char * toString(mgis::behaviour::Behaviour::Kinematic kin)
std::map< std::string, ParameterLib::Parameter< double > const * > readStateVariablesInitialValueProperties(mgis::behaviour::Behaviour const &behaviour, mgis::behaviour::Hypothesis const &hypothesis, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, BaseLib::ConfigTree const &config)
mgis::behaviour::Behaviour loadBehaviour(std::string const &lib_path, std::string behaviour_name, mgis::behaviour::Hypothesis const hypothesis)

References BaseLib::ConfigTree::checkConfigParameter(), BaseLib::ConfigTree::getConfigParameter(), BaseLib::ConfigTree::getConfigParameterOptional(), BaseLib::getProjectDirectory(), INFO(), BaseLib::joinPaths(), OGS_FATAL, and WARN().

Referenced by createMFront(), and createMFrontGeneric().

◆ createMFrontGeneric()

template<int DisplacementDim, typename Gradients , typename TDynForces , typename ExtStateVars >
std::unique_ptr< MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars > > MaterialLib::Solids::MFront::createMFrontGeneric ( std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
std::optional< ParameterLib::CoordinateSystem > const & local_coordinate_system,
BaseLib::ConfigTree const & config,
bool const library_path_is_relative_to_prj_file )

Definition at line 36 of file CreateMFrontGeneric.h.

42{
43 auto conf = createMFrontConfig(DisplacementDim, parameters, config,
44 library_path_is_relative_to_prj_file);
45
46 return std::make_unique<
48 std::move(conf.behaviour), std::move(conf.material_properties),
49 std::move(conf.state_variables_initial_properties),
50 local_coordinate_system);
51}

References createMFrontConfig().

Referenced by ProcessLib::LargeDeformation::ConstitutiveRelations::createMFrontGeneric(), and ProcessLib::ThermoRichardsMechanics::ConstitutiveStressSaturation_StrainPressureTemperature::createMFrontGeneric().

◆ eigenSwap45View()

template<typename Derived >
constexpr auto MaterialLib::Solids::MFront::eigenSwap45View ( Eigen::MatrixBase< Derived > const & matrix)
constexpr

Converts between OGS' and MFront's Kelvin vectors and matrices and vice versa. Numbering of the Kelvin vectors and matrices in the two cases is: MFront: 11 22 33 12 13 23 OGS: 11 22 33 12 23 13 The function swaps the 4th and 5th row and column of a matrix.

Definition at line 36 of file MFrontGeneric.h.

37{
38 using Matrix =
39 Eigen::Matrix<typename Derived::Scalar, Derived::RowsAtCompileTime,
40 Derived::ColsAtCompileTime>;
41
42 return Matrix::NullaryExpr(
43 matrix.rows(), matrix.cols(),
44 [&m = matrix.derived()](Eigen::Index const row, Eigen::Index const col)
45 {
46 constexpr std::ptrdiff_t result[6] = {0, 1, 2, 3, 5, 4};
47 return m(result[row], result[col]);
48 });
49}

Referenced by MaterialLib::Solids::MFront::detail::SetGradient< DisplacementDim >::operator()(), and tangentOperatorDataMFrontToOGS().

◆ getEquivalentPlasticStrainOffset()

int MaterialLib::Solids::MFront::getEquivalentPlasticStrainOffset ( mgis::behaviour::Behaviour const & b)

Definition at line 29 of file MFrontGeneric.cpp.

30{
31 return mgis::behaviour::contains(b.isvs, "EquivalentPlasticStrain")
32 ? mgis::behaviour::getVariableOffset(
33 b.isvs, "EquivalentPlasticStrain", b.hypothesis)
34 : -1;
35}

Referenced by MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::MFrontGeneric().

◆ ogsTensorToMFrontTensor()

template<int DisplacementDim, typename Derived >
constexpr auto MaterialLib::Solids::MFront::ogsTensorToMFrontTensor ( Eigen::MatrixBase< Derived > const & matrix)
constexpr

Converts between OGS' and MFront's tensors, which are represented as vectors. An OGS tensor 11 12 13 21 22 23 31 32 33 0 1 2 3 4 5 6 7 8 is converted to MFront tensor 11 22 33 12 21 13 31 23 32 0 4 8 1 3 2 6 5 7.

Definition at line 59 of file MFrontGeneric.h.

60{
61 using Matrix =
62 Eigen::Matrix<typename Derived::Scalar, Derived::RowsAtCompileTime,
63 Derived::ColsAtCompileTime>;
64
65 if constexpr (DisplacementDim == 2)
66 {
67 return Matrix::NullaryExpr(
68 matrix.rows(), matrix.cols(),
69 [&m = matrix.derived()](Eigen::Index const row,
70 Eigen::Index const col)
71 {
72 constexpr std::ptrdiff_t result[5] = {0, 3, 4, 1, 2};
73 return m(result[row], result[col]);
74 });
75 }
76 if constexpr (DisplacementDim == 3)
77 {
78 return Matrix::NullaryExpr(
79 matrix.rows(), matrix.cols(),
80 [&m = matrix.derived()](Eigen::Index const row,
81 Eigen::Index const col)
82 {
83 constexpr std::ptrdiff_t result[9] = {0, 4, 8, 1, 3,
84 2, 6, 5, 7};
85 return m(result[row], result[col]);
86 });
87 }
88}

References ogsTensorToMFrontTensor().

Referenced by ogsTensorToMFrontTensor().

◆ tangentOperatorDataMFrontToOGS()

template<int DisplacementDim>
OGSMFrontTangentOperatorData MaterialLib::Solids::MFront::tangentOperatorDataMFrontToOGS ( std::vector< double > const & mfront_data,
std::optional< MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > const & Q,
mgis::behaviour::Behaviour const & behaviour )

Transforms MFront's to OGS's tangent operator data.

Essentially swaps some off-diagonal components of symmetric tensors and applies the given rotation of the frame of reference Q.

Definition at line 38 of file MFrontGeneric.cpp.

43{
44 using VT = mgis::behaviour::Variable::Type;
47
48 std::vector<double> ogs_data(mfront_data.size());
49
50 std::size_t offset = 0;
51 constexpr auto kv_size =
53
54 for (auto const& [var1, var2] : behaviour.to_blocks)
55 {
56 std::size_t size;
57
58 auto const* const d_in = mfront_data.data() + offset;
59 auto* const d_out = ogs_data.data() + offset;
60
61 if (var1.type == VT::SCALAR && var2.type == VT::SCALAR)
62 {
63 size = 1;
64 *d_out = *d_in;
65 }
66 else if (var1.type == VT::SCALAR && var2.type == VT::STENSOR)
67 {
68 assert(getVariableSize(var2, behaviour.hypothesis) == kv_size);
69 size = kv_size;
70
71 if (Q)
72 {
74 "Coordinate frame rotation not yet implemented for "
75 "dScalar/dSTensor.");
76 }
77
78 Eigen::Map<KV>{d_out} = eigenSwap45View(Eigen::Map<const KV>{d_in});
79 }
80 else if (var1.type == VT::STENSOR && var2.type == VT::SCALAR)
81 {
82 assert(getVariableSize(var1, behaviour.hypothesis) == kv_size);
83 size = kv_size;
84
85 if (Q)
86 {
88 "Coordinate frame rotation not yet implemented for "
89 "dSTensor/dScalar.");
90 }
91
92 Eigen::Map<KV>{d_out} = eigenSwap45View(Eigen::Map<const KV>{d_in});
93 }
94 else if (var1.type == VT::STENSOR && var2.type == VT::STENSOR)
95 {
96 assert(getVariableSize(var1, behaviour.hypothesis) == kv_size);
97 assert(getVariableSize(var2, behaviour.hypothesis) == kv_size);
98 size = kv_size * kv_size;
99
100 if (Q)
101 {
102 Eigen::Map<KM>{d_out} =
103 *Q * eigenSwap45View(Eigen::Map<const KM>{d_in}) *
104 Q->transpose();
105 }
106 else
107 {
108 Eigen::Map<KM>{d_out} =
109 eigenSwap45View(Eigen::Map<const KM>{d_in});
110 }
111 }
112 else
113 {
114 OGS_FATAL("unsupported variable type combination");
115 }
116
117 offset += size;
118 }
119
120 return {std::move(ogs_data)};
121}
constexpr auto eigenSwap45View(Eigen::MatrixBase< Derived > const &matrix)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
constexpr int size(int const displacement_dim)
Vectorized tensor size for given displacement dimension.

References eigenSwap45View(), MathLib::KelvinVector::kelvin_vector_dimensions(), and OGS_FATAL.

◆ tangentOperatorDataMFrontToOGS< 2 >()

template OGSMFrontTangentOperatorData MaterialLib::Solids::MFront::tangentOperatorDataMFrontToOGS< 2 > ( std::vector< double > const & mfront_data,
std::optional< MathLib::KelvinVector::KelvinMatrixType< 2 > > const & Q,
mgis::behaviour::Behaviour const & behaviour )

◆ tangentOperatorDataMFrontToOGS< 3 >()

template OGSMFrontTangentOperatorData MaterialLib::Solids::MFront::tangentOperatorDataMFrontToOGS< 3 > ( std::vector< double > const & mfront_data,
std::optional< MathLib::KelvinVector::KelvinMatrixType< 3 > > const & Q,
mgis::behaviour::Behaviour const & behaviour )

◆ varTypeToString()

const char * MaterialLib::Solids::MFront::varTypeToString ( int v)

Definition at line 14 of file MFrontGeneric.cpp.

15{
16 using V = mgis::behaviour::Variable;
17 if (v == V::SCALAR)
18 return "SCALAR";
19 if (v == V::VECTOR)
20 return "VECTOR";
21 if (v == V::STENSOR)
22 return "STENSOR";
23 if (v == V::TENSOR)
24 return "TENSOR";
25
26 OGS_FATAL("Unknown variable type {}.", v);
27}

References OGS_FATAL.

Referenced by MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::MFrontGeneric(), MaterialLib::Solids::MFront::detail::SetGradient< DisplacementDim >::operator()(), and anonymous_namespace{CreateMFrontGeneric.cpp}::varInfo().

Variable Documentation

◆ deformation_gradient

constexpr DeformationGradient MaterialLib::Solids::MFront::deformation_gradient
staticconstexpr

Instance that can be used for overload resolution/template type deduction.

Definition at line 130 of file Variable.h.

◆ green_lagrange_strain

constexpr GreenLagrangeStrain MaterialLib::Solids::MFront::green_lagrange_strain
staticconstexpr

Instance that can be used for overload resolution/template type deduction.

Definition at line 112 of file Variable.h.

Referenced by ProcessLib::LargeDeformation::ConstitutiveRelations::SolidMechanicsModel< DisplacementDim >::eval().

◆ liquid_pressure

constexpr LiquidPressure MaterialLib::Solids::MFront::liquid_pressure
staticconstexpr

◆ saturation

constexpr Saturation MaterialLib::Solids::MFront::saturation
staticconstexpr

◆ second_piola_kirchhoff_stress

constexpr SecondPiolaKirchhoffStress MaterialLib::Solids::MFront::second_piola_kirchhoff_stress
staticconstexpr

◆ strain

constexpr Strain MaterialLib::Solids::MFront::strain
staticconstexpr

◆ stress

◆ temperature

constexpr Temperature MaterialLib::Solids::MFront::temperature
staticconstexpr