www.mooseframework.org
SymmIsotropicElasticityTensor.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
8 
10  : SymmElasticityTensor(constant),
11  _lambda_set(false),
12  _mu_set(false),
13  _E_set(false),
14  _nu_set(false),
15  _k_set(false),
16  _lambda(0),
17  _mu(0),
18  _E(0),
19  _nu(0),
20  _k(0)
21 {
22 }
23 
24 void
26 {
27  _lambda = lambda;
28  _lambda_set = true;
29 }
30 
31 void
33 {
34  _mu = mu;
35  _mu_set = true;
36 }
37 
38 void
40 {
41  _E = E;
42  _E_set = true;
43 }
44 
45 void
47 {
48  _nu = nu;
49  _nu_set = true;
50 }
51 
52 void
54 {
55  _k = k;
56  _k_set = true;
57 }
58 
59 void
61 {
62  setMu(k);
63 }
64 
65 Real
67 {
68  return mu();
69 }
70 
71 Real
73 {
74  if (!_E_set)
75  mooseError("Youngs modulus not set");
76 
77  return _E;
78 }
79 
80 Real
82 {
83  if (!_mu_set)
84  {
85  mooseError("mu not set");
86  }
87  return _mu;
88 }
89 
90 void
92 {
93  if (_lambda_set && _mu_set) // First and second Lame
94  return;
95  else if (_lambda_set && _nu_set)
96  _mu = (_lambda * (1.0 - 2.0 * _nu)) / (2.0 * _nu);
97  else if (_lambda_set && _k_set)
98  _mu = (3.0 * (_k - _lambda)) / 2.0;
99  else if (_lambda_set && _E_set)
100  _mu = ((_E - 3.0 * _lambda) / 4.0) +
101  (std::sqrt((_E - 3.0 * _lambda) * (_E - 3.0 * _lambda) + 8.0 * _lambda * _E) / 4.0);
102  else if (_mu_set && _nu_set)
103  _lambda = (2.0 * _mu * _nu) / (1.0 - 2.0 * _nu);
104  else if (_mu_set && _k_set)
105  _lambda = (3.0 * _k - 2.0 * _mu) / 3.0;
106  else if (_mu_set && _E_set)
107  _lambda = ((2.0 * _mu - _E) * _mu) / (_E - 3.0 * _mu);
108  else if (_nu_set && _k_set)
109  {
110  _lambda = (3.0 * _k * _nu) / (1.0 + _nu);
111  _mu = (3.0 * _k * (1.0 - 2.0 * _nu)) / (2.0 * (1.0 + _nu));
112  }
113  else if (_E_set && _nu_set) // Young's Modulus and Poisson's Ratio
114  {
115  _lambda = (_nu * _E) / ((1.0 + _nu) * (1 - 2.0 * _nu));
116  _mu = _E / (2.0 * (1.0 + _nu));
117  }
118  else if (_E_set && _k_set)
119  {
120  _lambda = (3.0 * _k * (3.0 * _k - _E)) / (9.0 * _k - _E);
121  _mu = (3.0 * _E * _k) / (9.0 * _k - _E);
122  }
123  _lambda_set = true;
124  _mu_set = true;
125 }
126 
127 void
129 {
131 
132  const Real C12(_lambda);
133  const Real C44(_mu);
134  const Real C11(2 * C44 + C12);
135 
136  setEntries(C11, C12, C44);
137 }
138 
139 void
140 SymmIsotropicElasticityTensor::setEntries(Real C11, Real C12, Real C44)
141 {
142  _val[0] = _val[6] = _val[11] = C11;
143  _val[1] = _val[2] = _val[7] = C12;
144  _val[15] = _val[18] = _val[20] = C44;
145  _val[3] = _val[4] = _val[5] = 0;
146  _val[8] = _val[9] = _val[10] = 0;
147  _val[12] = _val[13] = _val[14] = 0;
148  _val[16] = _val[17] = 0;
149  _val[19] = 0;
150 }
151 
152 Real
154  const unsigned int j,
155  const RealGradient & test,
156  const RealGradient & phi) const
157 {
158  RealGradient b;
159  if (0 == i && 0 == j)
160  {
161  b(0) = _val[0] * phi(0);
162  b(1) = _val[15] * phi(1);
163  b(2) = _val[20] * phi(2);
164  }
165  else if (1 == i && 1 == j)
166  {
167  b(0) = _val[15] * phi(0);
168  b(1) = _val[6] * phi(1);
169  b(2) = _val[18] * phi(2);
170  }
171  else if (2 == i && 2 == j)
172  {
173  b(0) = _val[20] * phi(0);
174  b(1) = _val[18] * phi(1);
175  b(2) = _val[11] * phi(2);
176  }
177  else if (0 == i && 1 == j)
178  {
179  b(0) = _val[1] * phi(1);
180  b(1) = _val[15] * phi(0);
181  b(2) = 0;
182  }
183  else if (1 == i && 0 == j)
184  {
185  b(0) = _val[15] * phi(1);
186  b(1) = _val[1] * phi(0);
187  b(2) = 0;
188  }
189  else if (1 == i && 2 == j)
190  {
191  b(0) = 0;
192  b(1) = _val[7] * phi(2);
193  b(2) = _val[18] * phi(1);
194  }
195  else if (2 == i && 1 == j)
196  {
197  b(0) = 0;
198  b(1) = _val[18] * phi(2);
199  b(2) = _val[7] * phi(1);
200  }
201  else if (0 == i && 2 == j)
202  {
203  b(0) = _val[2] * phi(2);
204  b(1) = 0;
205  b(2) = _val[20] * phi(0);
206  }
207  else if (2 == i && 0 == j)
208  {
209  b(0) = _val[20] * phi(2);
210  b(1) = 0;
211  b(2) = _val[2] * phi(0);
212  }
213  else
214  {
215  mooseError("Wrong index in stiffness calculation");
216  }
217  return test * b;
218 }
219 
220 void
222 {
223  const Real xx = x.xx();
224  const Real yy = x.yy();
225  const Real zz = x.zz();
226  const Real xy = x.xy();
227  const Real yz = x.yz();
228  const Real zx = x.zx();
229 
230  b.xx() = _val[0] * xx + _val[1] * yy + _val[2] * zz;
231  b.yy() = _val[1] * xx + _val[6] * yy + _val[7] * zz;
232  b.zz() = _val[2] * xx + _val[7] * yy + _val[11] * zz;
233  b.xy() = 2 * _val[15] * xy;
234  b.yz() = 2 * _val[18] * yz;
235  b.zx() = 2 * _val[20] * zx;
236 
237  b.xx() += 2 * (_val[3] * xy + _val[4] * yz + _val[5] * zx);
238  b.yy() += 2 * (_val[8] * xy + _val[9] * yz + _val[10] * zx);
239  b.zz() += 2 * (_val[12] * xy + _val[13] * yz + _val[14] * zx);
240  b.xy() += _val[3] * xx + _val[8] * yy + _val[12] * zz;
241  b.yz() += _val[4] * xx + _val[9] * yy + _val[13] * zz;
242  b.zx() += _val[5] * xx + _val[10] * yy + _val[14] * zz;
243  b.yz() += 2 * _val[16] * xy;
244  b.zx() += 2 * _val[17] * xy + 2 * _val[19] * yz;
245 }
246 
247 void
248 SymmIsotropicElasticityTensor::adjustForCracking(const RealVectorValue & crack_flags)
249 {
250  const RealVectorValue & c(crack_flags);
251  const Real c0(c(0));
252  const Real c0_coupled(c0 < 1 ? 0 : 1);
253  const Real c1(c(1));
254  const Real c1_coupled(c1 < 1 ? 0 : 1);
255  const Real c2(c(2));
256  const Real c2_coupled(c2 < 1 ? 0 : 1);
257 
258  const Real c01(c0_coupled * c1_coupled);
259  const Real c02(c0_coupled * c2_coupled);
260  const Real c12(c1_coupled * c2_coupled);
261  const Real c012(c0_coupled * c12);
262 
263  const Real ym = _mu * (3 * _lambda + 2 * _mu) / (_lambda + _mu);
264 
265  // Assume Poisson's ratio goes to zero for the cracked direction.
266 
267  _val[0] = (c0 < 1 ? c0 * ym : _val[0]);
268  _val[1] *= c01;
269  _val[2] *= c02;
270  _val[3] *= c01;
271  _val[4] *= c012;
272  _val[5] *= c02;
273 
274  _val[6] = (c1 < 1 ? c1 * ym : _val[6]);
275  _val[7] *= c12;
276  _val[8] *= c01;
277  _val[9] *= c12;
278  _val[10] *= c012;
279 
280  _val[11] = (c2 < 1 ? c2 * ym : _val[11]);
281  _val[12] *= c012;
282  _val[13] *= c12;
283  _val[14] *= c02;
284 }
285 
286 void
288  const RealVectorValue & crack_flags)
289 {
290  const RealVectorValue & c = crack_flags;
291  const Real c0 = c(0);
292  const Real c0_coupled = (c0 < 1 ? 0 : 1);
293  const Real c1 = c(1);
294  const Real c1_coupled = (c1 < 1 ? 0 : 1);
295  const Real c2 = c(2);
296  const Real c2_coupled = (c2 < 1 ? 0 : 1);
297  const Real c01 = c0_coupled * c1_coupled;
298  const Real c02 = c0_coupled * c2_coupled;
299  const Real c12 = c1_coupled * c2_coupled;
300  const Real c012 = c0_coupled * c12;
301  adjustForCracking(crack_flags);
302  _val[15] *= c01;
303  _val[16] *= c012;
304  _val[17] *= c012;
305  _val[18] *= c12;
306  _val[19] *= c012;
307  _val[20] *= c02;
308 }
This class defines a basic set of capabilities any elasticity tensor should have. ...
Real youngsModulus() const
Return the youngs modulus.
void setBulkModulus(const Real k)
Set the Bulk Modulus.
Real yy() const
Definition: SymmTensor.h:130
virtual void multiply(const SymmTensor &x, SymmTensor &b) const
Real xx() const
Definition: SymmTensor.h:129
void setShearModulus(const Real k)
Set the shear modulus...
void setEntries(Real C11, Real C12, Real C44)
virtual void calculateEntries(unsigned int qp)
Fill in the matrix.
void setLambda(const Real lambda)
Set the first Lame Coefficient.
Real zz() const
Definition: SymmTensor.h:131
SymmIsotropicElasticityTensor(const bool constant=true)
void calculateLameCoefficients()
Calculates lambda and mu based on what has been set.
virtual Real stiffness(const unsigned i, const unsigned j, const RealGradient &test, const RealGradient &phi) const
Real xy() const
Definition: SymmTensor.h:132
Real yz() const
Definition: SymmTensor.h:133
Real zx() const
Definition: SymmTensor.h:134
Real shearModulus() const
Return the shear modulus...
void setYoungsModulus(const Real E)
Set the Young&#39;s Modulus.
virtual void adjustForCracking(const RealVectorValue &crack_flags)
virtual void adjustForCrackingWithShearRetention(const RealVectorValue &crack_flags)
void setMu(const Real mu)
Set the second Lame Coefficient.
void setPoissonsRatio(const Real nu)
Set Poissons Ratio.