version 0.4.1
material.hh
Go to the documentation of this file.
1// SPDX-FileCopyrightText: 2021-2025 The Ikarus Developers mueller@ibb.uni-stuttgart.de
2// SPDX-License-Identifier: LGPL-3.0-or-later
3
9#pragma once
10
11#include <dune/common/classname.hh>
12#include <dune/python/common/typeregistry.hh>
13#include <dune/python/pybind11/eigen.h>
14#include <dune/python/pybind11/pybind11.h>
15#include <dune/python/pybind11/stl.h>
16
19
20#define MAKE_MaterialFunction(clsName, materialName, functionname, vecSize) \
21 clsName.def( \
22 #functionname, \
23 [](materialName& self, StrainTags straintag, Eigen::Ref<const Eigen::Vector<double, vecSize>> eVoigt_) { \
24 if constexpr (not Concepts::IsMaterial<Materials::LinearElasticityT, materialName>) { \
25 Eigen::Vector<double, vecSize> eVoigt = eVoigt_; \
26 if (straintag == StrainTags::rightCauchyGreenTensor) \
27 return self.template functionname<StrainTags::rightCauchyGreenTensor>(eVoigt); \
28 else if (straintag == StrainTags::greenLagrangian) \
29 return self.template functionname<StrainTags::greenLagrangian>(eVoigt); \
30 else if (straintag == StrainTags::linear) \
31 DUNE_THROW(Dune::MathError, "Passing linear strain to " + std::string(#materialName) + \
32 " does not makes sense use LinearElastic class"); \
33 else if (straintag == StrainTags::displacementGradient) \
34 DUNE_THROW(Dune::MathError, \
35 "Passing displacementGradient strain in 6d Voigt notation does not make any sense!"); \
36 else if (straintag == StrainTags::deformationGradient) \
37 DUNE_THROW(Dune::MathError, \
38 "Passing deformationGradient strain in 6d Voigt notation does not make any sense!"); \
39 else \
40 DUNE_THROW(Dune::MathError, toString(straintag) + "is not a valid strain tag."); \
41 } else { \
42 Eigen::Vector<double, vecSize> eVoigt = eVoigt_; /* Linear elastic path */ \
43 if (straintag == StrainTags::linear) \
44 return self.template functionname<StrainTags::linear>(eVoigt); \
45 else \
46 DUNE_THROW(Dune::MathError, "Linear elastic material only accepts linear strains!"); \
47 } \
48 __builtin_unreachable(); \
49 }, \
50 "StrainName"_a, "strainVector"_a);
51
52namespace Ikarus::Python {
53
54namespace Impl {
55 // Function to extract and convert parameters using a conversion strategy
56 template <typename T>
57 LamesFirstParameterAndShearModulus convertMaterialParameters(const pybind11::kwargs& kwargs,
58 const std::string& param1, const std::string& param2) {
59 auto converter =
60 convertLameConstants(T{kwargs[param1.c_str()].cast<double>(), kwargs[param2.c_str()].cast<double>()});
61
62 // This is necessary as the converter can't convert to a parameter already present due to compile-time constraints
63 double lamesFirst = [&]() {
64 if constexpr (requires { converter.toLamesFirstParameter(); })
65 return converter.toLamesFirstParameter();
66 else
67 return kwargs["Lambda"].cast<double>();
68 }();
69 double shearModulus = [&]() {
70 if constexpr (requires { converter.toShearModulus(); })
71 return converter.toShearModulus();
72 else
73 return kwargs["mu"].cast<double>();
74 }();
75 return {lamesFirst, shearModulus};
76 }
77
78 Ikarus::LamesFirstParameterAndShearModulus extractMaterialParameters(const pybind11::kwargs& kwargs) {
79 // clang-format off
80 static const std::map<std::array<std::string, 2>, std::function<LamesFirstParameterAndShearModulus(const pybind11::kwargs&)>> conversionMap = {
81 {{"E", "nu"}, [](const auto& kw){ return convertMaterialParameters<YoungsModulusAndPoissonsRatio>(kw, "E", "nu"); }},
82 {{"E", "mu"}, [](const auto& kw){ return convertMaterialParameters<YoungsModulusAndShearModulus>(kw, "E", "mu"); }},
83 {{"E", "K"}, [](const auto& kw){ return convertMaterialParameters<YoungsModulusAndBulkModulus>(kw, "E", "K"); }},
84 {{"E", "Lambda"}, [](const auto& kw){ return convertMaterialParameters<YoungsModulusAndLamesFirstParameter>(kw, "E", "Lambda"); }},
85 {{"K", "Lambda"}, [](const auto& kw){ return convertMaterialParameters<BulkModulusAndLamesFirstParameter>(kw, "K", "Lambda"); }},
86 {{"Lambda", "mu"}, [](const auto& kw){ return LamesFirstParameterAndShearModulus{kw["Lambda"].template cast<double>(), kw["mu"].template cast<double>()}; }}
87 };
88 // clang-format on
89
90 if (kwargs.size() != 2)
91 DUNE_THROW(Dune::IOError, "The number of material parameters passed to the material should be 2");
92
93 for (const auto& [materialParameters, parameterConverter] : conversionMap) {
94 const auto [firstPar, secondPar] = materialParameters;
95 if (kwargs.contains(firstPar) && kwargs.contains(secondPar)) {
96 return parameterConverter(kwargs);
97 }
98 }
99
100 DUNE_THROW(Dune::IOError,
101 "No suitable combination of material parameters found, valid combinations are: (E, nu), (E, mu), (E, "
102 "K), (E, Lambda), (K, Lambda), (Lambda, nu)");
103 }
104} // namespace Impl
105
106template <class Material, size_t vecSize, class... options>
107void registerMaterial(pybind11::handle scope, pybind11::class_<Material, options...> cls) {
108 using pybind11::operator""_a;
109 namespace py = pybind11;
110
111 cls.def(pybind11::init([](const py::kwargs& kwargs) {
112 auto matParameter = Impl::extractMaterialParameters(kwargs);
113 return new Material(matParameter);
114 }));
115
116 std::string materialname = Material::name();
117
118 MAKE_MaterialFunction(cls, Material, storedEnergy, vecSize);
119 MAKE_MaterialFunction(cls, Material, stresses, vecSize);
120 MAKE_MaterialFunction(cls, Material, tangentModuli, vecSize);
121
122 using PlaneStressClass = decltype(Materials::planeStress(std::declval<Material>()));
123 auto includes = Dune::Python::IncludeFiles{"ikarus/finiteelements/mechanics/materials.hh"};
124 auto pS = Dune::Python::insertClass<PlaneStressClass>(
125 scope, std::string("PlaneStress_") + materialname,
126 Dune::Python::GenerateTypeName(
127 "Ikarus::Materials::VanishingStress<std::array<Ikarus::Materials::MatrixIndexPair, 3ul >"
128 "{"
129 "{Ikarus::Materials::MatrixIndexPair{2ul, 1ul}, Ikarus::Materials::MatrixIndexPair{2ul, 0ul},"
130 "Ikarus::Materials::MatrixIndexPair{2ul, 2ul}}}," +
131 Dune::className<Material>() + ">"),
132 includes)
133 .first;
134 MAKE_MaterialFunction(pS, PlaneStressClass, storedEnergy, 3);
135 MAKE_MaterialFunction(pS, PlaneStressClass, stresses, 3);
136 MAKE_MaterialFunction(pS, PlaneStressClass, tangentModuli, 3);
137 MAKE_MaterialFunction(pS, PlaneStressClass, storedEnergy, 6);
138 MAKE_MaterialFunction(pS, PlaneStressClass, stresses, 6);
139 MAKE_MaterialFunction(pS, PlaneStressClass, tangentModuli, 6);
140
141 cls.def("asPlaneStress", [](Material& self) {
142 return Materials::planeStress(self);
143 }); /* no keep_alive since planeStress copies the material */
144
145 using PlaneStrainClass = decltype(Materials::planeStrain(std::declval<Material>()));
146 auto pStrain = Dune::Python::insertClass<PlaneStrainClass>(
147 scope, std::string("PlaneStrain_") + materialname,
148 Dune::Python::GenerateTypeName(
149 "Ikarus::Materials::VanishingStrain<std::array<Ikarus::Materials::MatrixIndexPair, "
150 "3ul>{{Ikarus::Materials::MatrixIndexPair{2ul, 1ul},"
151 "Ikarus::Materials::MatrixIndexPair{2ul,0ul}, Ikarus::Materials::MatrixIndexPair{"
152 "2ul, 2ul}}}," +
153 Dune::className<Material>() + ">"),
154 includes)
155 .first;
156 MAKE_MaterialFunction(pStrain, PlaneStrainClass, storedEnergy, 3);
157 MAKE_MaterialFunction(pStrain, PlaneStrainClass, stresses, 3);
158 MAKE_MaterialFunction(pStrain, PlaneStrainClass, tangentModuli, 3);
159 MAKE_MaterialFunction(pStrain, PlaneStrainClass, storedEnergy, 6);
160 MAKE_MaterialFunction(pStrain, PlaneStrainClass, stresses, 6);
161 MAKE_MaterialFunction(pStrain, PlaneStrainClass, tangentModuli, 6);
162
163 cls.def("asPlaneStrain", [](Material& self) {
164 return Materials::planeStrain(self);
165 }); /* no keep_alive since planeStrain copies the material */
166 using ShellMaterialClass = decltype(Materials::shellMaterial(std::declval<Material>()));
167 auto shellmaterial = Dune::Python::insertClass<ShellMaterialClass>(
168 scope, std::string("Shell_") + materialname,
169 Dune::Python::GenerateTypeName(
170 "Ikarus::Materials::VanishingStress<std::array<Ikarus::Materials::MatrixIndexPair,"
171 "1ul>{{Ikarus::Materials::MatrixIndexPair{2ul, 2ul}}}," +
172 Dune::className<Material>() + ">"),
173 includes)
174 .first;
175 MAKE_MaterialFunction(shellmaterial, ShellMaterialClass, storedEnergy, 5);
176 MAKE_MaterialFunction(shellmaterial, ShellMaterialClass, stresses, 5);
177 MAKE_MaterialFunction(shellmaterial, ShellMaterialClass, tangentModuli, 5);
178 MAKE_MaterialFunction(shellmaterial, ShellMaterialClass, storedEnergy, 6);
179 MAKE_MaterialFunction(shellmaterial, ShellMaterialClass, stresses, 6);
180 MAKE_MaterialFunction(shellmaterial, ShellMaterialClass, tangentModuli, 6);
181
182 cls.def("asShellMaterial", [](Material& self) {
183 return Materials::shellMaterial(self);
184 }); /* no keep_alive since shellMaterial copies the material */
185 using BeamMaterialClass = decltype(Materials::beamMaterial(std::declval<Material>()));
186 auto beammaterial = Dune::Python::insertClass<BeamMaterialClass>(
187 scope, std::string("Beam_") + materialname,
188 Dune::Python::GenerateTypeName(
189 "Ikarus::Materials::VanishingStress<std::array<Ikarus::Materials::MatrixIndexPair, "
190 "2ul>{{Materials::MatrixIndexPair{1, 1},Ikarus::Materials::MatrixIndexPair{2ul, 2ul}}}," +
191 Dune::className<Material>() + ">"),
192 includes)
193 .first;
194 MAKE_MaterialFunction(beammaterial, BeamMaterialClass, storedEnergy, 4);
195 MAKE_MaterialFunction(beammaterial, BeamMaterialClass, stresses, 4);
196 MAKE_MaterialFunction(beammaterial, BeamMaterialClass, tangentModuli, 4);
197 MAKE_MaterialFunction(beammaterial, BeamMaterialClass, storedEnergy, 6);
198 MAKE_MaterialFunction(beammaterial, BeamMaterialClass, stresses, 6);
199 MAKE_MaterialFunction(beammaterial, BeamMaterialClass, tangentModuli, 6);
200 cls.def("asBeamMaterial", [](Material& self) {
201 return Materials::beamMaterial(self);
202 }); /* no keep_alive since beamMaterial copies the material */
203}
204
205#define MAKE_MATERIAL_REGISTRY_FUNCTION(name, vecSize) \
206 template <class Material, class... options> \
207 void register##name(pybind11::handle scope, pybind11::class_<Material, options...> cls) { \
208 Ikarus::Python::registerMaterial<Material, vecSize>(scope, cls); \
209 }
210
214} // namespace Ikarus::Python
#define MAKE_MaterialFunction(clsName, materialName, functionname, vecSize)
Definition: material.hh:20
ConvertLameConstants< YoungsModulusAndPoissonsRatio > convertLameConstants(const YoungsModulusAndPoissonsRatio &valuePair)
Definition: physicshelper.hh:250
void init(int argc, char **argv, bool enableFileLogger=true)
Initializes the Ikarus framework.
Definition: init.hh:82
NeoHookeT< double > NeoHooke
Alias for NeoHookeT with double as the default scalar type.
Definition: neohooke.hh:160
auto planeStrain(const MaterialImpl &mat)
Factory function to create a VanishingStrain material for plane strain conditions.
Definition: vanishingstrain.hh:197
StVenantKirchhoffT< double > StVenantKirchhoff
Alias for StVenantKirchhoffT with double as the default scalar type.
Definition: svk.hh:182
auto shellMaterial(const MaterialImpl &mat, typename MaterialImpl::ScalarType tol=1e-8)
Factory function to create a VanishingStress material for a shell material with zero normal stress co...
Definition: vanishingstress.hh:255
LinearElasticityT< double > LinearElasticity
Convenience typedef for LinearElasticity with double as ScalarType.
Definition: linearelasticity.hh:130
auto planeStress(const MaterialImpl &mat, typename MaterialImpl::ScalarType tol=1e-8)
Factory function to create a VanishingStress material for plane stress conditions.
Definition: vanishingstress.hh:242
auto beamMaterial(const MaterialImpl &mat, typename MaterialImpl::ScalarType tol=1e-8)
Factory function to create a VanishingStress material for a beam material with two zero normal stress...
Definition: vanishingstress.hh:268
Definition: flatassembler.hh:21
void registerMaterial(pybind11::handle scope, pybind11::class_< Material, options... > cls)
Definition: material.hh:107
MAKE_MATERIAL_REGISTRY_FUNCTION(LinearElasticity, 6)
Definition: physicshelper.hh:54
Header file for material models in Ikarus finite element mechanics.
Several concepts.