Line data Source code
1 : //
2 : // The isotropic tensor is defined to be
3 : //
4 : // .. math::
5 : //
6 : // C_{ijkl} = \mu(d_{il}d_{jk} + d_{ik}d_{jl}) + \lambda d_{ij} d_{kl}
7 : //
8 : // The inverse ("compliance") tensor is
9 : //
10 : // .. math::
11 : //
12 : // S_{ijkl} = ((1+\nu)/2E)(d_{il}d_{jl} + d_{ik}d_{jl}) - (\nu/E)d_{ij}d_{kl}
13 : //
14 : // Replacing E, nu with Lame constants gives:
15 : //
16 : // .. math::
17 : //
18 : // S_{ijkl} = (1/(4 \mu))(d_{il}d_{jl} + d_{ik}d_{jl}) + (\lambda/(2*\mu*(3\lambda+2\mu))) * d_{ij} * d_{kl}
19 : //
20 : // For reference: http://solidmechanics.org/text/Chapter3_2/Chapter3_2.htm
21 : // https://en.wikipedia.org/wiki/Lam%C3%A9_parameters
22 : //
23 :
24 : #ifndef SET_MATRIX4_ISOTROPIC_H
25 : #define SET_MATRIX4_ISOTROPIC_H
26 :
27 : #include "Util/Util.H"
28 : #include "Base.H"
29 :
30 : namespace Set
31 : {
32 : template<>
33 : class Matrix4<AMREX_SPACEDIM,Sym::Isotropic>
34 : {
35 : Set::Scalar lambda=NAN, mu=NAN;
36 : public:
37 52221969 : AMREX_GPU_HOST_DEVICE Matrix4() {};
38 88305 : AMREX_GPU_HOST_DEVICE Matrix4(Set::Scalar a_lambda, Set::Scalar a_mu) : lambda(a_lambda), mu(a_mu) {};
39 :
40 : /// Note: for the Isotropic Matrix4 this routine works for **retrieval only**!
41 : /// If you try to assign a value using this with, say
42 : ///
43 : /// isotropicmatrix4(i,j,k,l) = 8.0
44 : ///
45 : /// you will get a `lvalue required as left operand of assignment` compile error.
46 : /// You should probably consider using a lower symmetry operator.
47 : AMREX_FORCE_INLINE
48 : Scalar operator () (const int i, const int j, const int k, const int l) const
49 : {
50 3880 : Set::Scalar ret = 0.0;
51 3880 : if (i==k && j==l) ret += mu;
52 3880 : if (i==l && j==k) ret += mu;
53 3880 : if (i==j && k==l) ret += lambda;
54 3880 : return ret;
55 : }
56 : void Randomize()
57 : {
58 : lambda = Util::Random();
59 : mu = Util::Random();
60 : }
61 0 : void Print (std::ostream& os )
62 : {
63 0 : os << "lambda = " << lambda << " mu = " << mu;
64 0 : }
65 450711 : Set::Scalar Lambda () const
66 : {
67 450711 : return lambda;
68 : }
69 450711 : Set::Scalar Mu () const
70 : {
71 450711 : return mu;
72 : }
73 : Set::Scalar Youngs() const
74 : {
75 : return mu * ( 3.0*lambda + 2.0*mu ) / (lambda + mu);
76 : }
77 : Set::Scalar Nu() const
78 : {
79 : return 0.5 * lambda / (lambda + mu);
80 : }
81 562 : static Matrix4<AMREX_SPACEDIM,Sym::Isotropic> Zero()
82 : {
83 562 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> zero;
84 562 : zero.lambda = 0.0;
85 562 : zero.mu = 0.0;
86 562 : return zero;
87 : }
88 :
89 :
90 79 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> Inverse() const
91 : {
92 79 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> inv;
93 79 : inv.mu = 1./4./mu;
94 79 : inv.lambda = lambda / (2*mu*(3*lambda + 2*mu));
95 79 : return inv;
96 : }
97 : friend Set::Matrix operator * (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Matrix &b);
98 : friend Set::Vector operator * (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Matrix3 &b);
99 : friend Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator - (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &b);
100 : friend Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator * (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Scalar &b);
101 : friend Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator * (const Set::Scalar &b, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a);
102 : friend Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator / (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Scalar &b);
103 : friend Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator * (const Set::Scalar &b, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a);
104 : friend bool operator == (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a,const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &b);
105 : friend Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator + (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a,const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &b);
106 :
107 : //AMREX_GPU_HOST_DEVICE void operator = (Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a) {lambda = a.lambda; mu = a.mu;}
108 126627 : AMREX_GPU_HOST_DEVICE void operator += (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a) {lambda += a.lambda; mu += a.mu;}
109 : AMREX_GPU_HOST_DEVICE void operator -= (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a) {lambda -= a.lambda; mu -= a.mu;}
110 : AMREX_GPU_HOST_DEVICE void operator *= (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a) {lambda *= a.lambda; mu *= a.mu;}
111 : AMREX_GPU_HOST_DEVICE void operator /= (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a) {lambda /= a.lambda; mu /= a.mu;}
112 : AMREX_GPU_HOST_DEVICE void operator *= (const Set::Scalar &alpha) {lambda *= alpha; mu *= alpha;}
113 : AMREX_GPU_HOST_DEVICE void operator /= (const Set::Scalar &alpha) {lambda /= alpha; mu /= alpha;}
114 :
115 : Set::Scalar Norm()
116 : {
117 : return std::sqrt(lambda*lambda + mu*mu);
118 : }
119 :
120 0 : bool contains_nan() const
121 : {
122 0 : if (std::isnan(lambda)) return true;
123 0 : if (std::isnan(mu)) return true;
124 0 : return false;
125 : }
126 : };
127 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
128 : Set::Matrix operator * (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Matrix &b)
129 : {
130 50021357 : Set::Matrix ret;
131 :
132 : #if AMREX_SPACEDIM == 2
133 35046273 : ret(0,0) = (a.lambda + 2.*a.mu) * b(0,0) + a.lambda *b(1,1);
134 35046273 : ret(1,1) = a.lambda * b(0,0) + (a.lambda + 2.*a.mu)*b(1,1);
135 35046273 : ret(0,1) = a.mu*(b(0,1) + b(1,0)); ret(1,0) = ret(0,1);
136 :
137 : #elif AMREX_SPACEDIM == 3
138 14975084 : ret(0,0) = (a.lambda + 2.*a.mu) * b(0,0) + a.lambda *b(1,1) + a.lambda *b(2,2);
139 14975084 : ret(1,1) = a.lambda * b(0,0) + (a.lambda + 2.*a.mu)*b(1,1) + a.lambda *b(2,2);
140 14975084 : ret(2,2) = a.lambda * b(0,0) + a.lambda *b(1,1) + (a.lambda + 2.*a.mu)*b(2,2);
141 14975084 : ret(1,2) = a.mu*(b(1,2) + b(2,1)); ret(2,1) = ret(1,2);
142 14975084 : ret(2,0) = a.mu*(b(2,0) + b(0,2)); ret(0,2) = ret(2,0);
143 14975084 : ret(0,1) = a.mu*(b(0,1) + b(1,0)); ret(1,0) = ret(0,1);
144 :
145 : #endif
146 50021357 : return ret;
147 : }
148 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
149 : Set::Matrix operator * (const Set::Matrix &b, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a)
150 : {return a*b;}
151 :
152 :
153 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
154 : Set::Vector operator * (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Matrix3 &b)
155 : {
156 11898313 : Set::Vector ret = Set::Vector::Zero();
157 :
158 37311639 : for (int i = 0; i < AMREX_SPACEDIM; i++)
159 81090078 : for (int j=0; j < AMREX_SPACEDIM; j++)
160 222707008 : ret(i) += a.mu*(b(i,j,j) + b(j,i,j)) + a.lambda*b(j,j,i);
161 :
162 11898313 : return ret;
163 : }
164 :
165 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
166 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator * (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Scalar &b)
167 : {
168 25223672 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> ret;// = Set::Vector::Zero();
169 25223672 : ret.mu = a.mu * b;
170 25223672 : ret.lambda = a.lambda * b;
171 25223672 : return ret;
172 : }
173 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
174 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator / (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Set::Scalar &b)
175 : {
176 25041668 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> ret;// = Set::Vector::Zero();
177 25041668 : ret.mu = a.mu / b;
178 25041668 : ret.lambda = a.lambda / b;
179 25041668 : return ret;
180 : }
181 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
182 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator * (const Set::Scalar &b, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a)
183 : {
184 0 : return a*b;
185 : }
186 :
187 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
188 : bool operator == (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &b)
189 : {
190 20 : if (a.mu != b.mu) return false;
191 20 : if (a.lambda != b.lambda) return false;
192 20 : return true;
193 : }
194 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
195 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator + (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a, const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &b)
196 : {
197 70959 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> ret;// = Set::Vector::Zero();
198 70959 : ret.mu = a.mu + b.mu;
199 70959 : ret.lambda = a.lambda + b.lambda;
200 70959 : return ret;
201 : }
202 :
203 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
204 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> operator - (const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &a,
205 : const Matrix4<AMREX_SPACEDIM,Sym::Isotropic> &b)
206 : {
207 25042830 : Matrix4<AMREX_SPACEDIM,Sym::Isotropic> ret = a;
208 25042830 : ret.mu -= b.mu;
209 25042830 : ret.lambda -= b.lambda;
210 25042830 : return ret;
211 : }
212 :
213 :
214 : }
215 : #endif
|
