Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "MooseTypes.h"
11 : #include "PermutationTensor.h"
12 : #include "RankFourTensor.h"
13 : #include "ElasticityTensorTools.h"
14 :
15 : #include "libmesh/vector_value.h"
16 :
17 : namespace ElasticityTensorTools
18 : {
19 :
20 : Real
21 2302056944 : elasticJacobian(const RankFourTensor & r4t,
22 : unsigned int i,
23 : unsigned int k,
24 : const RealGradient & grad_test,
25 : const RealGradient & grad_phi)
26 : {
27 : // 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 ==
28 : // k)
29 :
30 2302056944 : const Real gt0 = grad_test(0);
31 2302056944 : const Real gt1 = grad_test(1);
32 2302056944 : const Real gt2 = grad_test(2);
33 2302056944 : const Real gp0 = grad_phi(0);
34 2302056944 : const Real gp1 = grad_phi(1);
35 2302056944 : const Real gp2 = grad_phi(2);
36 :
37 : // clang-format off
38 : // This is the algorithm that is unrolled below:
39 : //
40 : // Real sum = 0.0;
41 : // for (const auto j: make_range(Moose::dim))
42 : // for (const auto l: make_range(Moose::dim))
43 : // sum += r4t(i, j, k, l) * grad_phi(l) * grad_test(j);
44 : // return sum;
45 :
46 : return
47 : (
48 2302056944 : r4t(i,0,k,0) * gp0
49 2302056944 : + r4t(i,0,k,1) * gp1
50 2302056944 : + r4t(i,0,k,2) * gp2
51 2302056944 : ) * gt0
52 2302056944 : +
53 : (
54 2302056944 : r4t(i,1,k,0) * gp0
55 2302056944 : + r4t(i,1,k,1) * gp1
56 2302056944 : + r4t(i,1,k,2) * gp2
57 2302056944 : ) * gt1
58 : +
59 : (
60 2302056944 : r4t(i,2,k,0) * gp0
61 2302056944 : + r4t(i,2,k,1) * gp1
62 2302056944 : + r4t(i,2,k,2) * gp2
63 2302056944 : ) * gt2;
64 : // clang-format on
65 : }
66 :
67 : Real
68 13320192 : elasticJacobianWC(const RankFourTensor & r4t,
69 : unsigned int i,
70 : unsigned int k,
71 : const RealGradient & grad_test,
72 : Real phi)
73 : {
74 : // 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 ==
75 : // k)
76 : Real sum = 0.0;
77 53280768 : for (const auto j : make_range(Moose::dim))
78 159842304 : for (const auto m : make_range(Moose::dim))
79 479526912 : for (const auto n : make_range(Moose::dim))
80 359645184 : sum += r4t(i, j, m, n) * PermutationTensor::eps(m, n, k) * grad_test(j);
81 13320192 : return sum * phi;
82 : }
83 :
84 : Real
85 13320192 : momentJacobian(const RankFourTensor & r4t,
86 : unsigned int i,
87 : unsigned int k,
88 : Real test,
89 : const RealGradient & grad_phi)
90 : {
91 : // Jacobian entry: d(eps_ijm*stress_jm*test)/du_k = d(eps_ijm*C_jmln*du_l/dx_n*test)/du_k (only
92 : // nonzero for l == k)
93 : Real sum = 0.0;
94 53280768 : for (const auto j : make_range(Moose::dim))
95 159842304 : for (const auto m : make_range(Moose::dim))
96 479526912 : for (const auto n : make_range(Moose::dim))
97 359645184 : sum += PermutationTensor::eps(i, j, m) * r4t(j, m, k, n) * grad_phi(n);
98 13320192 : return test * sum;
99 : }
100 :
101 : Real
102 9453568 : momentJacobianWC(const RankFourTensor & r4t, unsigned int i, unsigned int k, Real test, Real phi)
103 : {
104 : // Jacobian entry: d(eps_ijm*stress_jm*test)/dw_k = d(eps_ijm*C_jmln*eps_lnp*w_p*test)/dw_k (only
105 : // nonzero for p ==k)
106 : Real sum = 0.0;
107 37814272 : for (const auto j : make_range(Moose::dim))
108 113442816 : for (const auto l : make_range(Moose::dim))
109 340328448 : for (const auto m : make_range(Moose::dim))
110 1020985344 : for (const auto n : make_range(Moose::dim))
111 765739008 : sum +=
112 765739008 : PermutationTensor::eps(i, j, m) * r4t(j, m, l, n) * PermutationTensor::eps(l, n, k);
113 :
114 9453568 : return test * phi * sum;
115 : }
116 :
117 : void
118 7817472 : toVoigtNotationIndexConversion(int k, int & a, int & b)
119 : {
120 7817472 : if (k < 3 && k >= 0)
121 3908736 : a = b = k;
122 3908736 : else if (k == 3)
123 : {
124 1302912 : a = 1;
125 1302912 : b = 2;
126 : }
127 2605824 : else if (k == 4)
128 : {
129 1302912 : a = 0;
130 1302912 : b = 2;
131 : }
132 1302912 : else if (k == 5)
133 : {
134 1302912 : a = 0;
135 1302912 : b = 1;
136 : }
137 : else
138 0 : mooseError("\nIndex out of bound while converting from tensor to voigt notation in "
139 : "toVoigtNotationIndexConversion");
140 7817472 : }
141 :
142 : // MOOSE uses the stress tensor with components ordered as (11, 22, 33, 12, 23, 13) in voigt form.
143 : // This is in conflict with the voigt notation used in literature which has components ordered as
144 : // (11, 22, 33, 23, 13, 12). Whenever an operation involving the voigt forms of stress and
145 : // elasticity tensor has to be performed the voigt form of the elasticity tensor should be built
146 : // taking into account this difference in component-ordering for stress. The
147 : // toMooseVoigtNotationIndexConversion function facilitates that.
148 :
149 : void
150 529128864 : toMooseVoigtNotationIndexConversion(int k, int & a, int & b)
151 : {
152 529128864 : if (k < 3 && k >= 0)
153 264564432 : a = b = k;
154 264564432 : else if (k == 3)
155 : {
156 88188144 : a = 0;
157 88188144 : b = 1;
158 : }
159 176376288 : else if (k == 4)
160 : {
161 88188144 : a = 1;
162 88188144 : b = 2;
163 : }
164 88188144 : else if (k == 5)
165 : {
166 88188144 : a = 0;
167 88188144 : b = 2;
168 : }
169 : else
170 0 : mooseError("\nIndex out of bound while converting from tensor to MOOSE voigt notation in "
171 : "toMooseVoigtNotationIndexConversion");
172 529128864 : }
173 : }
|