Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <limits>
8 : #include <memory>
9 : #include <pup.h>
10 :
11 : #include "DataStructures/Tensor/TypeAliases.hpp"
12 : #include "Options/String.hpp"
13 : #include "PointwiseFunctions/Elasticity/ConstitutiveRelations/ConstitutiveRelation.hpp"
14 : #include "Utilities/Gsl.hpp"
15 : #include "Utilities/Serialization/CharmPupable.hpp"
16 : #include "Utilities/TMPL.hpp"
17 :
18 : /// \cond
19 : class DataVector;
20 : /// \endcond
21 :
22 : namespace Elasticity {
23 : namespace ConstitutiveRelations {
24 :
25 : /*!
26 : * \brief An isotropic and homogeneous material
27 : *
28 : * \details For an isotropic and homogeneous material the linear constitutive
29 : * relation \f$T^{ij}=-Y^{ijkl}S_{kl}\f$ reduces to
30 : *
31 : * \f[
32 : * Y^{ijkl} = \lambda \delta^{ij}\delta^{kl} + \mu
33 : * \left(\delta^{ik}\delta^{jl} + \delta^{il}\delta^{jk}\right) \\
34 : * \implies \quad T^{ij} = -\lambda \mathrm{Tr}(S) \delta^{ij} - 2\mu S^{ij}
35 : * \f]
36 : *
37 : * with the _Lamé parameter_ \f$\lambda\f$ and the _shear modulus_ (or
38 : * _rigidity_) \f$\mu\f$. In the parametrization chosen in this implementation
39 : * we use the _bulk modulus_ (or _incompressibility_)
40 : *
41 : * \f[
42 : * K=\lambda + \frac{2}{3}\mu
43 : * \f]
44 : *
45 : * instead of the Lamé parameter. In this parametrization the
46 : * stress-strain relation
47 : *
48 : * \f[
49 : * T^{ij} = -K \mathrm{Tr}(S) \delta^{ij} - 2\mu\left(S^{ij} -
50 : * \frac{1}{3}\mathrm{Tr}(S)\delta^{ij}\right)
51 : * \f]
52 : *
53 : * decomposes into a scalar and a traceless part (Eq. 11.18 in
54 : * \cite ThorneBlandford2017). Parameters also often used in this context are
55 : * the _Young's modulus_
56 : *
57 : * \f[
58 : * E=\frac{9K\mu}{3K+\mu}=\frac{\mu(3\lambda+2\mu)}{\lambda+\mu}
59 : * \f]
60 : *
61 : * and the _Poisson ratio_
62 : *
63 : * \f[
64 : * \nu=\frac{3K-2\mu}{2(3K+\mu)}=\frac{\lambda}{2(\lambda+\mu)}\text{.}
65 : * \f]
66 : *
67 : * Inversely, these relations read:
68 : *
69 : * \f[
70 : * K =\frac{E}{3(1-2\nu)}, \quad
71 : * \lambda =\frac{E\nu}{(1+\nu)(1-2\nu)}, \quad
72 : * \mu =\frac{E}{2(1+\nu)}
73 : * \f]
74 : *
75 : * **In two dimensions** this implementation reduces to the plane-stress
76 : * approximation. We assume that all stresses apply in the plane of the
77 : * computational domain, which corresponds to scenarios of in-plane stretching
78 : * and shearing of thin slabs of material. Since orthogonal stresses vanish as
79 : * \f$T^{i3}=0=T^{3i}\f$ we find \f$\mathrm{Tr}(S)=\frac{2\mu}{\lambda +
80 : * 2\mu}\mathrm{Tr}^{(2)}(S)\f$, where \f$\mathrm{Tr}^{(2)}\f$ denotes that the
81 : * trace only applies to the two dimensions within the plane. The constitutive
82 : * relation thus reduces to
83 : *
84 : * \f{align}
85 : * T^{ij}=&-\frac{2\lambda\mu}{\lambda + 2\mu}\mathrm{Tr}^{(2)}(S)\delta^{ij} -
86 : * 2\mu S^{ij} \\
87 : * =&-\frac{E\nu}{1-\nu^2}\mathrm{Tr}^{(2)}(S)\delta^{ij} -
88 : * \frac{E}{1+\nu}S^{ij}
89 : * \f}
90 : *
91 : * which is non-zero only in the directions of the plane. Since the stresses
92 : * are also assumed to be constant along the thickness of the plane
93 : * \f$\partial_3T^{ij}=0\f$ the elasticity problem \f$-\partial_i T^{ij}=F^j\f$
94 : * reduces to two dimensions.
95 : */
96 : template <size_t Dim>
97 1 : class IsotropicHomogeneous : public ConstitutiveRelation<Dim> {
98 : public:
99 0 : static constexpr size_t volume_dim = Dim;
100 :
101 0 : struct BulkModulus {
102 0 : using type = double;
103 0 : static constexpr Options::String help = {
104 : "The incompressibility of the material"};
105 0 : static type lower_bound() { return 0.0; }
106 : };
107 :
108 0 : struct ShearModulus {
109 0 : using type = double;
110 0 : static constexpr Options::String help = {"The rigidity of the material"};
111 0 : static type lower_bound() { return 0.0; }
112 : };
113 :
114 0 : using options = tmpl::list<BulkModulus, ShearModulus>;
115 :
116 0 : static constexpr Options::String help = {
117 : "A constitutive relation that describes an isotropic, homogeneous "
118 : "material in terms of two elastic moduli. These bulk and shear moduli "
119 : "indicate the material's resistance to volume and shape changes, "
120 : "respectively. Both are measured in units of stress, typically Pascals."};
121 :
122 0 : IsotropicHomogeneous() = default;
123 0 : IsotropicHomogeneous(const IsotropicHomogeneous&) = default;
124 0 : IsotropicHomogeneous& operator=(const IsotropicHomogeneous&) = default;
125 0 : IsotropicHomogeneous(IsotropicHomogeneous&&) = default;
126 0 : IsotropicHomogeneous& operator=(IsotropicHomogeneous&&) = default;
127 0 : ~IsotropicHomogeneous() override = default;
128 :
129 0 : IsotropicHomogeneous(double bulk_modulus, double shear_modulus);
130 :
131 1 : std::unique_ptr<ConstitutiveRelation<Dim>> get_clone() const override;
132 :
133 : /// The constitutive relation that characterizes the elastic properties of a
134 : /// material
135 1 : void stress(gsl::not_null<tnsr::II<DataVector, Dim>*> stress,
136 : const tnsr::ii<DataVector, Dim>& strain,
137 : const tnsr::I<DataVector, Dim>& x) const override;
138 :
139 : /// The bulk modulus (or incompressibility) \f$K\f$
140 1 : double bulk_modulus() const;
141 : /// The shear modulus (or rigidity) \f$\mu\f$
142 1 : double shear_modulus() const;
143 : /// The Lamé parameter \f$\lambda\f$
144 1 : double lame_parameter() const;
145 : /// The Young's modulus \f$E\f$
146 1 : double youngs_modulus() const;
147 : /// The Poisson ratio \f$\nu\f$
148 1 : double poisson_ratio() const;
149 :
150 : // clang-tidy: no runtime references
151 0 : void pup(PUP::er& /*p*/) override; // NOLINT
152 :
153 0 : explicit IsotropicHomogeneous(CkMigrateMessage* /*unused*/) {}
154 :
155 0 : WRAPPED_PUPable_decl_base_template( // NOLINT
156 : SINGLE_ARG(ConstitutiveRelation<Dim>), IsotropicHomogeneous);
157 :
158 : private:
159 0 : double bulk_modulus_ = std::numeric_limits<double>::signaling_NaN();
160 0 : double shear_modulus_ = std::numeric_limits<double>::signaling_NaN();
161 : };
162 :
163 : template <size_t Dim>
164 0 : bool operator==(const IsotropicHomogeneous<Dim>& lhs,
165 : const IsotropicHomogeneous<Dim>& rhs);
166 : template <size_t Dim>
167 0 : bool operator!=(const IsotropicHomogeneous<Dim>& lhs,
168 : const IsotropicHomogeneous<Dim>& rhs);
169 :
170 : /// \cond
171 : template <size_t Dim>
172 : PUP::able::PUP_ID IsotropicHomogeneous<Dim>::my_PUP_ID = 0;
173 : /// \endcond
174 :
175 : } // namespace ConstitutiveRelations
176 : } // namespace Elasticity
|