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 1780662926 : 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 1780662926 : const Real gt0 = grad_test(0);
31 1780662926 : const Real gt1 = grad_test(1);
32 1780662926 : const Real gt2 = grad_test(2);
33 1780662926 : const Real gp0 = grad_phi(0);
34 1780662926 : const Real gp1 = grad_phi(1);
35 1780662926 : 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 1780662926 : r4t(i,0,k,0) * gp0
49 1780662926 : + r4t(i,0,k,1) * gp1
50 1780662926 : + r4t(i,0,k,2) * gp2
51 1780662926 : ) * gt0
52 1780662926 : +
53 : (
54 1780662926 : r4t(i,1,k,0) * gp0
55 1780662926 : + r4t(i,1,k,1) * gp1
56 1780662926 : + r4t(i,1,k,2) * gp2
57 1780662926 : ) * gt1
58 : +
59 : (
60 1780662926 : r4t(i,2,k,0) * gp0
61 1780662926 : + r4t(i,2,k,1) * gp1
62 1780662926 : + r4t(i,2,k,2) * gp2
63 1780662926 : ) * gt2;
64 : // clang-format on
65 : }
66 :
67 : Real
68 8368128 : 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 33472512 : for (const auto j : make_range(Moose::dim))
78 100417536 : for (const auto m : make_range(Moose::dim))
79 301252608 : for (const auto n : make_range(Moose::dim))
80 225939456 : sum += r4t(i, j, m, n) * PermutationTensor::eps(m, n, k) * grad_test(j);
81 8368128 : return sum * phi;
82 : }
83 :
84 : Real
85 8368128 : 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 33472512 : for (const auto j : make_range(Moose::dim))
95 100417536 : for (const auto m : make_range(Moose::dim))
96 301252608 : for (const auto n : make_range(Moose::dim))
97 225939456 : sum += PermutationTensor::eps(i, j, m) * r4t(j, m, k, n) * grad_phi(n);
98 8368128 : return test * sum;
99 : }
100 :
101 : Real
102 6932480 : 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 27729920 : for (const auto j : make_range(Moose::dim))
108 83189760 : for (const auto l : make_range(Moose::dim))
109 249569280 : for (const auto m : make_range(Moose::dim))
110 748707840 : for (const auto n : make_range(Moose::dim))
111 561530880 : sum +=
112 561530880 : PermutationTensor::eps(i, j, m) * r4t(j, m, l, n) * PermutationTensor::eps(l, n, k);
113 :
114 6932480 : return test * phi * sum;
115 : }
116 :
117 : void
118 6252480 : toVoigtNotationIndexConversion(int k, int & a, int & b)
119 : {
120 6252480 : if (k < 3 && k >= 0)
121 3126240 : a = b = k;
122 3126240 : else if (k == 3)
123 : {
124 1042080 : a = 1;
125 1042080 : b = 2;
126 : }
127 2084160 : else if (k == 4)
128 : {
129 1042080 : a = 0;
130 1042080 : b = 2;
131 : }
132 1042080 : else if (k == 5)
133 : {
134 1042080 : a = 0;
135 1042080 : b = 1;
136 : }
137 : else
138 0 : mooseError("\nIndex out of bound while converting from tensor to voigt notation in "
139 : "toVoigtNotationIndexConversion");
140 6252480 : }
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 400194144 : toMooseVoigtNotationIndexConversion(int k, int & a, int & b)
151 : {
152 400194144 : if (k < 3 && k >= 0)
153 200097072 : a = b = k;
154 200097072 : else if (k == 3)
155 : {
156 66699024 : a = 0;
157 66699024 : b = 1;
158 : }
159 133398048 : else if (k == 4)
160 : {
161 66699024 : a = 1;
162 66699024 : b = 2;
163 : }
164 66699024 : else if (k == 5)
165 : {
166 66699024 : a = 0;
167 66699024 : b = 2;
168 : }
169 : else
170 0 : mooseError("\nIndex out of bound while converting from tensor to MOOSE voigt notation in "
171 : "toMooseVoigtNotationIndexConversion");
172 400194144 : }
173 : }
|
