Line data Source code
1 : // Copyright 2024, UChicago Argonne, LLC
2 : // All Rights Reserved
3 : // Software Name: NEML2 -- the New Engineering material Model Library, version 2
4 : // By: Argonne National Laboratory
5 : // OPEN SOURCE LICENSE (MIT)
6 : //
7 : // Permission is hereby granted, free of charge, to any person obtaining a copy
8 : // of this software and associated documentation files (the "Software"), to deal
9 : // in the Software without restriction, including without limitation the rights
10 : // to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
11 : // copies of the Software, and to permit persons to whom the Software is
12 : // furnished to do so, subject to the following conditions:
13 : //
14 : // The above copyright notice and this permission notice shall be included in
15 : // all copies or substantial portions of the Software.
16 : //
17 : // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
18 : // IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
19 : // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
20 : // AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
21 : // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
22 : // OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
23 : // THE SOFTWARE.
24 :
25 : #include "neml2/tensors/SSR4.h"
26 : #include "neml2/tensors/Scalar.h"
27 : #include "neml2/tensors/SR2.h"
28 : #include "neml2/tensors/R4.h"
29 : #include "neml2/tensors/R5.h"
30 : #include "neml2/tensors/SSFR5.h"
31 : #include "neml2/tensors/Rot.h"
32 : #include "neml2/tensors/SSSSR8.h"
33 : #include "neml2/tensors/R8.h"
34 : #include "neml2/tensors/assertions.h"
35 : #include "neml2/tensors/mandel_notation.h"
36 : #include "neml2/tensors/functions/bmv.h"
37 : #include "neml2/tensors/functions/linalg/inv.h"
38 :
39 : namespace neml2
40 : {
41 :
42 76 : SSR4::SSR4(const R4 & T)
43 152 : : SSR4(full_to_mandel(
44 152 : full_to_mandel((T + T.transpose_minor() + R4(T.transpose(0, 1)) + R4(T.transpose(2, 3))) /
45 : 4.0),
46 76 : 1))
47 : {
48 76 : }
49 :
50 : SSR4
51 90 : SSR4::identity(const TensorOptions & options)
52 : {
53 3960 : return SSR4::create({{1, 1, 1, 0, 0, 0},
54 : {1, 1, 1, 0, 0, 0},
55 : {1, 1, 1, 0, 0, 0},
56 : {0, 0, 0, 0, 0, 0},
57 : {0, 0, 0, 0, 0, 0},
58 : {0, 0, 0, 0, 0, 0}},
59 810 : options);
60 3420 : }
61 :
62 : SSR4
63 9 : SSR4::identity_C1(const TensorOptions & options)
64 : {
65 396 : return SSR4::create({{1, 0, 0, 0, 0, 0},
66 : {0, 1, 0, 0, 0, 0},
67 : {0, 0, 1, 0, 0, 0},
68 : {0, 0, 0, 0, 0, 0},
69 : {0, 0, 0, 0, 0, 0},
70 : {0, 0, 0, 0, 0, 0}},
71 81 : options);
72 342 : }
73 :
74 : SSR4
75 9 : SSR4::identity_C2(const TensorOptions & options)
76 : {
77 396 : return SSR4::create({{0, 1, 1, 0, 0, 0},
78 : {1, 0, 1, 0, 0, 0},
79 : {1, 1, 0, 0, 0, 0},
80 : {0, 0, 0, 0, 0, 0},
81 : {0, 0, 0, 0, 0, 0},
82 : {0, 0, 0, 0, 0, 0}},
83 81 : options);
84 342 : }
85 :
86 : SSR4
87 9 : SSR4::identity_C3(const TensorOptions & options)
88 : {
89 396 : return SSR4::create({{0, 0, 0, 0, 0, 0},
90 : {0, 0, 0, 0, 0, 0},
91 : {0, 0, 0, 0, 0, 0},
92 : {0, 0, 0, 1, 0, 0},
93 : {0, 0, 0, 0, 1, 0},
94 : {0, 0, 0, 0, 0, 1}},
95 81 : options);
96 342 : }
97 :
98 : SSR4
99 127 : SSR4::identity_sym(const TensorOptions & options)
100 : {
101 127 : return SSR4(at::eye(6, options), 0);
102 : }
103 :
104 : SSR4
105 15 : SSR4::identity_vol(const TensorOptions & options)
106 : {
107 30 : return SSR4::identity(options) / 3;
108 : }
109 :
110 : SSR4
111 43 : SSR4::identity_dev(const TensorOptions & options)
112 : {
113 43 : return SSR4::identity_sym(options) - SSR4::identity(options) / 3;
114 : }
115 :
116 : SSR4
117 5 : SSR4::isotropic_E_nu(const Scalar & E, const Scalar & nu)
118 : {
119 5 : neml_assert_broadcastable_dbg(E, nu);
120 :
121 5 : const auto zero = Scalar::zeros_like(E);
122 5 : const auto pf = E / ((1.0 + nu) * (1.0 - 2.0 * nu));
123 5 : const auto C1 = (1.0 - nu) * pf;
124 5 : const auto C2 = nu * pf;
125 5 : const auto C4 = (1.0 - 2.0 * nu) * pf;
126 :
127 10 : return SSR4::fill_C1_C2_C3(C1, C2, C4);
128 5 : }
129 :
130 : SSR4
131 3 : SSR4::isotropic_E_nu(const double & E, const double & nu, const TensorOptions & options)
132 : {
133 3 : return SSR4::isotropic_E_nu(Scalar(E, options), Scalar(nu, options));
134 : }
135 :
136 : SSR4
137 5 : SSR4::fill_C1_C2_C3(const Scalar & C1, const Scalar & C2, const Scalar & C3)
138 : {
139 5 : neml_assert_broadcastable_dbg(C1, C2, C3);
140 :
141 10 : return C1 * identity_C1(C1.options()) + C2 * identity_C2(C2.options()) +
142 15 : C3 * identity_C3(C3.options());
143 : }
144 :
145 : SSR4
146 0 : SSR4::fill_C1_C2_C3(const double & C1,
147 : const double & C2,
148 : const double & C3,
149 : const TensorOptions & options)
150 : {
151 0 : return SSR4::fill_C1_C2_C3(Scalar(C1, options), Scalar(C2, options), Scalar(C3, options));
152 : }
153 :
154 : SSSSR8
155 3 : SSR4::identity_map(const TensorOptions & options)
156 : {
157 3 : auto I = at::eye(6, options);
158 15 : return SSSSR8(at::einsum("ik,jl", {I, I}));
159 6 : }
160 :
161 : SSR4
162 55 : SSR4::rotate(const Rot & r) const
163 : {
164 55 : return R4(*this).rotate(r);
165 : }
166 :
167 : SSFR5
168 9 : SSR4::drotate(const Rot & r) const
169 : {
170 9 : auto dR = R4(*this).drotate(r);
171 18 : return full_to_mandel(full_to_mandel(dR), 1);
172 9 : }
173 :
174 : SSSSR8
175 6 : SSR4::drotate_self(const Rot & r) const
176 : {
177 6 : auto R = r.euler_rodrigues();
178 42 : auto Tsym = 0.25 * (at::einsum("...ma,...nb,...oc,...pd->...mnopabcd", {R, R, R, R}) +
179 54 : at::einsum("...mb,...na,...od,...pc->...mnopabcd", {R, R, R, R}) +
180 48 : at::einsum("...mb,...na,...oc,...pd->...mnopabcd", {R, R, R, R}) +
181 48 : at::einsum("...ma,...nb,...od,...pc->...mnopabcd", {R, R, R, R}));
182 12 : return SSSSR8(full_to_mandel(
183 24 : full_to_mandel(full_to_mandel(full_to_mandel(R8(Tsym, R.batch_dim()), 0), 1), 2), 3));
184 30 : }
185 :
186 : Scalar
187 324 : SSR4::operator()(Size i, Size j, Size k, Size l) const
188 : {
189 324 : const auto a = mandel_reverse_index[i][j];
190 324 : const auto b = mandel_reverse_index[k][l];
191 1296 : return base_index({a, b}) / (mandel_factor(a) * mandel_factor(b));
192 324 : }
193 :
194 : SSR4
195 47 : SSR4::inverse() const
196 : {
197 47 : return linalg::inv(*this);
198 : }
199 :
200 : SSSSR8
201 5 : SSR4::dinverse() const
202 : {
203 5 : auto SI = this->inverse();
204 25 : return SSSSR8(-at::einsum("...ik,...lj->...ijkl", {SI, SI}));
205 10 : }
206 :
207 : SSR4
208 0 : SSR4::transpose_minor() const
209 : {
210 0 : return *this;
211 : }
212 :
213 : SSR4
214 4 : SSR4::transpose_major() const
215 : {
216 4 : return TensorBase<SSR4>::base_transpose(0, 1);
217 : }
218 :
219 : SR2
220 27 : operator*(const SSR4 & a, const SR2 & b)
221 : {
222 27 : return SR2(bmv(a, b));
223 : }
224 :
225 : SR2
226 4 : operator*(const SR2 & a, const SSR4 & b)
227 : {
228 4 : return SR2(bmv(b.transpose_major(), a));
229 : }
230 :
231 : SSR4
232 18 : operator*(const SSR4 & a, const SSR4 & b)
233 : {
234 18 : return SSR4(at::matmul(a, b));
235 : }
236 : } // namespace neml2
|
