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