ElasticityTensorTools.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 /****************************************************************/
7
8 #include "MooseTypes.h"
9 #include "PermutationTensor.h"
10 #include "RankFourTensor.h"
11
12 namespace ElasticityTensorTools
13 {
14
15 Real
16 elasticJacobian(const RankFourTensor & r4t,
17  unsigned int i,
18  unsigned int k,
21 {
22  // d(stress_ij*d(test)/dx_j)/du_k = d(C_ijmn*du_m/dx_n dtest/dx_j)/du_k (which is nonzero for m ==
23  // k)
24
25  const Real gt0 = grad_test(0);
26  const Real gt1 = grad_test(1);
27  const Real gt2 = grad_test(2);
28  const Real gp0 = grad_phi(0);
29  const Real gp1 = grad_phi(1);
30  const Real gp2 = grad_phi(2);
31
32  // clang-format off
33  // This is the algorithm that is unrolled below:
34  //
35  // Real sum = 0.0;
36  // for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
37  // for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
38  // sum += r4t(i, j, k, l) * grad_phi(l) * grad_test(j);
39  // return sum;
40
41  return
42  (
43  r4t(i,0,k,0) * gp0
44  + r4t(i,0,k,1) * gp1
45  + r4t(i,0,k,2) * gp2
46  ) * gt0
47  +
48  (
49  r4t(i,1,k,0) * gp0
50  + r4t(i,1,k,1) * gp1
51  + r4t(i,1,k,2) * gp2
52  ) * gt1
53  +
54  (
55  r4t(i,2,k,0) * gp0
56  + r4t(i,2,k,1) * gp1
57  + r4t(i,2,k,2) * gp2
58  ) * gt2;
59  // clang-format on
60 }
61
62 Real
63 elasticJacobianWC(const RankFourTensor & r4t,
64  unsigned int i,
65  unsigned int k,
67  Real phi)
68 {
69  // d(stress_ij*d(test)/dx_j)/dw_k = d(C_ijmn*eps_mnp*w_p*dtest/dx_j)/dw_k (only nonzero for p ==
70  // k)
71  Real sum = 0.0;
72  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
73  for (unsigned int m = 0; m < LIBMESH_DIM; ++m)
74  for (unsigned int n = 0; n < LIBMESH_DIM; ++n)
75  sum += r4t(i, j, m, n) * PermutationTensor::eps(m, n, k) * grad_test(j);
76  return sum * phi;
77 }
78
79 Real
80 momentJacobian(const RankFourTensor & r4t,
81  unsigned int i,
82  unsigned int k,
83  Real test,
85 {
86  // Jacobian entry: d(eps_ijm*stress_jm*test)/du_k = d(eps_ijm*C_jmln*du_l/dx_n*test)/du_k (only
87  // nonzero for l == k)
88  Real sum = 0.0;
89  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
90  for (unsigned int m = 0; m < LIBMESH_DIM; ++m)
91  for (unsigned int n = 0; n < LIBMESH_DIM; ++n)
92  sum += PermutationTensor::eps(i, j, m) * r4t(j, m, k, n) * grad_phi(n);
93  return test * sum;
94 }
95
96 Real
97 momentJacobianWC(const RankFourTensor & r4t, unsigned int i, unsigned int k, Real test, Real phi)
98 {
99  // Jacobian entry: d(eps_ijm*stress_jm*test)/dw_k = d(eps_ijm*C_jmln*eps_lnp*w_p*test)/dw_k (only
100  // nonzero for p ==k)
101  Real sum = 0.0;
102  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
103  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
104  for (unsigned int m = 0; m < LIBMESH_DIM; ++m)
105  for (unsigned int n = 0; n < LIBMESH_DIM; ++n)
106  sum +=
107  PermutationTensor::eps(i, j, m) * r4t(j, m, l, n) * PermutationTensor::eps(l, n, k);
108
109  return test * phi * sum;
110 }
111
112 Real
113 getIsotropicShearModulus(const RankFourTensor & elasticity_tensor)
114 {
115  const Real shear_modulus = elasticity_tensor(0, 1, 0, 1);
116  return shear_modulus;
117 }
118
119 Real
120 getIsotropicBulkModulus(const RankFourTensor & elasticity_tensor)
121 {
122  const Real shear_modulus = getIsotropicShearModulus(elasticity_tensor);
123  // dilatational modulus is defined as lambda plus two mu
124  const Real dilatational_modulus = elasticity_tensor(0, 0, 0, 0);
125  const Real lambda = dilatational_modulus - 2.0 * shear_modulus;
126  const Real bulk_modulus = lambda + 2.0 * shear_modulus / 3.0;
127  return bulk_modulus;
128 }
129
130 Real
131 getIsotropicYoungsModulus(const RankFourTensor & elasticity_tensor)
132 {
133  const Real shear_modulus = getIsotropicShearModulus(elasticity_tensor);
134  // dilatational modulus is defined as lambda plus two mu
135  const Real dilatational_modulus = elasticity_tensor(0, 0, 0, 0);
136  const Real lambda = dilatational_modulus - 2.0 * shear_modulus;
137  const Real youngs_modulus =
138  shear_modulus * (3.0 * lambda + 2.0 * shear_modulus) / (lambda + shear_modulus);
139  return youngs_modulus;
140 }
141
142 Real
143 getIsotropicPoissonsRatio(const RankFourTensor & elasticity_tensor)
144 {
145  const Real poissons_ratio = elasticity_tensor(1, 1, 0, 0) /
146  (elasticity_tensor(1, 1, 1, 1) + elasticity_tensor(1, 1, 0, 0));
147  return poissons_ratio;
148 }
149 }
Real elasticJacobianWC(const RankFourTensor &r4t, unsigned int i, unsigned int k, const RealGradient &grad_test, Real phi)
This is used for the standard kernel stress_ij*d(test)/dx_j, when varied wrt w_k (the cosserat rotati...
Real momentJacobianWC(const RankFourTensor &r4t, unsigned int i, unsigned int k, Real test, Real phi)
This is used for the moment-balancing kernel eps_ijk*stress_jk*test, when varied wrt w_k (the cossera...
Real getIsotropicBulkModulus(const RankFourTensor &elasticity_tensor)
Get the bulk modulus for an isotropic elasticity tensor param elasticity_tensor the tensor (must be i...