GMatElastoPlasticQPot 0.18.3
Loading...
Searching...
No Matches
Cartesian2d.hpp
Go to the documentation of this file.
1
9#ifndef GMATTENSOR_CARTESIAN2D_HPP
10#define GMATTENSOR_CARTESIAN2D_HPP
11
12#include "Cartesian2d.h"
13#include "detail.hpp"
14
15namespace GMatTensor {
16namespace Cartesian2d {
17
18namespace detail {
19using GMatTensor::detail::impl_A2;
20using GMatTensor::detail::impl_A4;
21} // namespace detail
22
24{
25 array_type::tensor<double, 2> ret = xt::random::randn<double>({2, 2});
26 return ret;
27}
28
30{
31 array_type::tensor<double, 4> ret = xt::random::randn<double>({2, 2, 2, 2});
32 return ret;
33}
34
36{
37 array_type::tensor<double, 2> ret = xt::zeros<double>({2, 2});
38 return ret;
39}
40
42{
43 array_type::tensor<double, 4> ret = xt::zeros<double>({2, 2, 2, 2});
44 return ret;
45}
46
48{
49 array_type::tensor<double, 2> ret = xt::empty<double>({2, 2});
50 pointer::I2(ret.data());
51 return ret;
52}
53
55{
56 array_type::tensor<double, 4> ret = xt::empty<double>({2, 2, 2, 2});
57 pointer::II(ret.data());
58 return ret;
59}
60
62{
63 array_type::tensor<double, 4> ret = xt::empty<double>({2, 2, 2, 2});
64 pointer::I4(ret.data());
65 return ret;
66}
67
69{
70 array_type::tensor<double, 4> ret = xt::empty<double>({2, 2, 2, 2});
71 pointer::I4rt(ret.data());
72 return ret;
73}
74
76{
77 array_type::tensor<double, 4> ret = xt::empty<double>({2, 2, 2, 2});
78 pointer::I4s(ret.data());
79 return ret;
80}
81
83{
84 array_type::tensor<double, 4> ret = xt::empty<double>({2, 2, 2, 2});
85 pointer::I4d(ret.data());
86 return ret;
87}
88
89template <class T, class R>
90inline void trace(const T& A, R& ret)
91{
92 return detail::impl_A2<T, 2>::ret0(A, ret, [](const auto& a) { return pointer::Trace(a); });
93}
94
95template <class T>
96inline auto Trace(const T& A)
97{
98 return detail::impl_A2<T, 2>::ret0(A, [](const auto& a) { return pointer::Trace(a); });
99}
100
101template <class T, class R>
102inline void hydrostatic(const T& A, R& ret)
103{
104 return detail::impl_A2<T, 2>::ret0(
105 A, ret, [](const auto& a) { return pointer::Hydrostatic(a); });
106}
107
108template <class T>
109inline auto Hydrostatic(const T& A)
110{
111 return detail::impl_A2<T, 2>::ret0(A, [](const auto& a) { return pointer::Hydrostatic(a); });
112}
113
114template <class T, class R>
115inline void A2_ddot_B2(const T& A, const T& B, R& ret)
116{
117 return detail::impl_A2<T, 2>::B2_ret0(
118 A, B, ret, [](const auto& a, const auto& b) { return pointer::A2_ddot_B2(a, b); });
119}
120
121template <class T>
122inline auto A2_ddot_B2(const T& A, const T& B)
123{
124 return detail::impl_A2<T, 2>::B2_ret0(
125 A, B, [](const auto& a, const auto& b) { return pointer::A2_ddot_B2(a, b); });
126}
127
128template <class T, class R>
129inline void A2s_ddot_B2s(const T& A, const T& B, R& ret)
130{
131 return detail::impl_A2<T, 2>::B2_ret0(
132 A, B, ret, [](const auto& a, const auto& b) { return pointer::A2s_ddot_B2s(a, b); });
133}
134
135template <class T>
136inline auto A2s_ddot_B2s(const T& A, const T& B)
137{
138 return detail::impl_A2<T, 2>::B2_ret0(
139 A, B, [](const auto& a, const auto& b) { return pointer::A2s_ddot_B2s(a, b); });
140}
141
142template <class T, class R>
143inline void norm_deviatoric(const T& A, R& ret)
144{
145 return detail::impl_A2<T, 2>::ret0(
146 A, ret, [](const auto& a) { return pointer::Norm_deviatoric(a); });
147}
148
149template <class T>
150inline auto Norm_deviatoric(const T& A)
151{
152 return detail::impl_A2<T, 2>::ret0(
153 A, [](const auto& a) { return pointer::Norm_deviatoric(a); });
154}
155
156template <class T, class R>
157inline void deviatoric(const T& A, R& ret)
158{
159 return detail::impl_A2<T, 2>::ret2(
160 A, ret, [](const auto& a, const auto& r) { return pointer::Hydrostatic_deviatoric(a, r); });
161}
162
163template <class T>
164inline auto Deviatoric(const T& A)
165{
166 return detail::impl_A2<T, 2>::ret2(
167 A, [](const auto& a, const auto& r) { return pointer::Hydrostatic_deviatoric(a, r); });
168}
169
170template <class T, class R>
171inline void sym(const T& A, R& ret)
172{
173 return detail::impl_A2<T, 2>::ret2(
174 A, ret, [](const auto& a, const auto& r) { return pointer::sym(a, r); });
175}
176
177template <class T>
178inline auto Sym(const T& A)
179{
180 return detail::impl_A2<T, 2>::ret2(
181 A, [](const auto& a, const auto& r) { return pointer::sym(a, r); });
182}
183
184template <class T, class R>
185inline void A2_dot_B2(const T& A, const T& B, R& ret)
186{
187 return detail::impl_A2<T, 2>::B2_ret2(
188 A, B, ret, [](const auto& a, const auto& b, const auto& r) {
189 return pointer::A2_dot_B2(a, b, r);
190 });
191}
192
193template <class T>
194inline auto A2_dot_B2(const T& A, const T& B)
195{
196 return detail::impl_A2<T, 2>::B2_ret2(A, B, [](const auto& a, const auto& b, const auto& r) {
197 return pointer::A2_dot_B2(a, b, r);
198 });
199}
200
201template <class T, class R>
202inline void A2_dyadic_B2(const T& A, const T& B, R& ret)
203{
204 return detail::impl_A2<T, 2>::B2_ret4(
205 A, B, ret, [](const auto& a, const auto& b, const auto& r) {
206 return pointer::A2_dyadic_B2(a, b, r);
207 });
208}
209
210template <class T>
211inline auto A2_dyadic_B2(const T& A, const T& B)
212{
213 return detail::impl_A2<T, 2>::B2_ret4(A, B, [](const auto& a, const auto& b, const auto& r) {
214 return pointer::A2_dyadic_B2(a, b, r);
215 });
216}
217
218template <class T, class U, class R>
219inline void A4_ddot_B2(const T& A, const U& B, R& ret)
220{
221 return detail::impl_A4<T, 2>::B2_ret2(
222 A, B, ret, [](const auto& a, const auto& b, const auto& r) {
223 return pointer::A4_ddot_B2(a, b, r);
224 });
225}
226
227template <class T, class U>
228inline auto A4_ddot_B2(const T& A, const U& B)
229{
230 return detail::impl_A4<T, 2>::B2_ret2(A, B, [](const auto& a, const auto& b, const auto& r) {
231 return pointer::A4_ddot_B2(a, b, r);
232 });
233}
234
235namespace pointer {
236
237template <class T>
238inline void O2(T* ret)
239{
240 std::fill(ret, ret + 4, T(0));
241}
242
243template <class T>
244inline void O4(T* ret)
245{
246 std::fill(ret, ret + 16, T(0));
247}
248
249template <class T>
250inline void I2(T* ret)
251{
252 ret[0] = 1.0;
253 ret[1] = 0.0;
254 ret[2] = 0.0;
255 ret[3] = 1.0;
256}
257
258template <class T>
259inline void II(T* ret)
260{
261 std::fill(ret, ret + 16, T(0));
262
263 for (size_t i = 0; i < 2; ++i) {
264 for (size_t j = 0; j < 2; ++j) {
265 for (size_t k = 0; k < 2; ++k) {
266 for (size_t l = 0; l < 2; ++l) {
267 if (i == j && k == l) {
268 ret[i * 8 + j * 4 + k * 2 + l] = 1.0;
269 }
270 }
271 }
272 }
273 }
274}
275
276template <class T>
277inline void I4(T* ret)
278{
279 std::fill(ret, ret + 16, T(0));
280
281 for (size_t i = 0; i < 2; ++i) {
282 for (size_t j = 0; j < 2; ++j) {
283 for (size_t k = 0; k < 2; ++k) {
284 for (size_t l = 0; l < 2; ++l) {
285 if (i == l && j == k) {
286 ret[i * 8 + j * 4 + k * 2 + l] = 1.0;
287 }
288 }
289 }
290 }
291 }
292}
293
294template <class T>
295inline void I4rt(T* ret)
296{
297 std::fill(ret, ret + 16, T(0));
298
299 for (size_t i = 0; i < 2; ++i) {
300 for (size_t j = 0; j < 2; ++j) {
301 for (size_t k = 0; k < 2; ++k) {
302 for (size_t l = 0; l < 2; ++l) {
303 if (i == k && j == l) {
304 ret[i * 8 + j * 4 + k * 2 + l] = 1.0;
305 }
306 }
307 }
308 }
309 }
310}
311
312template <class T>
313inline void I4s(T* ret)
314{
315 I4(ret);
316
317 std::array<double, 16> i4rt;
318 I4rt(&i4rt[0]);
319
320 std::transform(ret, ret + 16, &i4rt[0], ret, std::plus<T>());
321
322 std::transform(ret, ret + 16, ret, std::bind(std::multiplies<T>(), std::placeholders::_1, 0.5));
323}
324
325template <class T>
326inline void I4d(T* ret)
327{
328 I4s(ret);
329
330 std::array<double, 16> ii;
331 II(&ii[0]);
332
333 std::transform(
334 &ii[0], &ii[0] + 16, &ii[0], std::bind(std::multiplies<T>(), std::placeholders::_1, 0.5));
335
336 std::transform(ret, ret + 16, &ii[0], ret, std::minus<T>());
337}
338
339template <class T>
340inline T Trace(const T* A)
341{
342 return A[0] + A[3];
343}
344
345template <class T>
346inline T Hydrostatic(const T* A)
347{
348 return T(0.5) * Trace(A);
349}
350
351template <class T>
352inline void sym(const T* A, T* ret)
353{
354 ret[0] = A[0];
355 ret[1] = 0.5 * (A[1] + A[2]);
356 ret[2] = ret[1];
357 ret[3] = A[3];
358}
359
360template <class T>
361inline T Hydrostatic_deviatoric(const T* A, T* ret)
362{
363 T m = Hydrostatic(A);
364 ret[0] = A[0] - m;
365 ret[1] = A[1];
366 ret[2] = A[2];
367 ret[3] = A[3] - m;
368 return m;
369}
370
371template <class T>
372inline T Deviatoric_ddot_deviatoric(const T* A)
373{
374 T m = Hydrostatic(A);
375 return (A[0] - m) * (A[0] - m) + (A[3] - m) * (A[3] - m) + T(2) * A[1] * A[2];
376}
377
378template <class T>
379inline T Norm_deviatoric(const T* A)
380{
381 return std::sqrt(Deviatoric_ddot_deviatoric(A));
382}
383
384template <class T>
385inline T A2_ddot_B2(const T* A, const T* B)
386{
387 return A[0] * B[0] + A[3] * B[3] + A[1] * B[2] + A[2] * B[1];
388}
389
390template <class T>
391inline T A2s_ddot_B2s(const T* A, const T* B)
392{
393 return A[0] * B[0] + A[3] * B[3] + T(2) * A[1] * B[1];
394}
395
396template <class T>
397inline void A2_dyadic_B2(const T* A, const T* B, T* ret)
398{
399 for (size_t i = 0; i < 2; ++i) {
400 for (size_t j = 0; j < 2; ++j) {
401 for (size_t k = 0; k < 2; ++k) {
402 for (size_t l = 0; l < 2; ++l) {
403 ret[i * 8 + j * 4 + k * 2 + l] = A[i * 2 + j] * B[k * 2 + l];
404 }
405 }
406 }
407 }
408}
409
410template <class T>
411inline void A2_dot_B2(const T* A, const T* B, T* ret)
412{
413 ret[0] = A[1] * B[2] + A[0] * B[0];
414 ret[1] = A[0] * B[1] + A[1] * B[3];
415 ret[2] = A[2] * B[0] + A[3] * B[2];
416 ret[3] = A[2] * B[1] + A[3] * B[3];
417}
418
419template <class T>
420inline void A4_ddot_B2(const T* A, const T* B, T* ret)
421{
422 std::fill(ret, ret + 4, T(0));
423
424 for (size_t i = 0; i < 2; i++) {
425 for (size_t j = 0; j < 2; j++) {
426 for (size_t k = 0; k < 2; k++) {
427 for (size_t l = 0; l < 2; l++) {
428 ret[i * 2 + j] += A[i * 8 + j * 4 + k * 2 + l] * B[l * 2 + k];
429 }
430 }
431 }
432 }
433}
434
435} // namespace pointer
436
437} // namespace Cartesian2d
438} // namespace GMatTensor
439
440#endif
void I4d(T *ret)
See Cartesian2d::I4d()
void I4(T *ret)
See Cartesian2d::I4()
T Deviatoric_ddot_deviatoric(const T *A)
Double tensor contraction of the tensor's deviator.
void I2(T *ret)
See Cartesian2d::I2()
T A2_ddot_B2(const T *A, const T *B)
See Cartesian2d::A2_ddot_B2()
void I4rt(T *ret)
See Cartesian2d::I4rt()
T Trace(const T *A)
See Cartesian2d::Trace()
void A4_ddot_B2(const T *A, const T *B, T *ret)
See Cartesian2d::A4_ddot_B2()
void A2_dot_B2(const T *A, const T *B, T *ret)
See Cartesian2d::A2_dot_B2()
void A2_dyadic_B2(const T *A, const T *B, T *ret)
See Cartesian2d::A2_dyadic_B2()
void I4s(T *ret)
See Cartesian2d::I4s()
T A2s_ddot_B2s(const T *A, const T *B)
See Cartesian2d::A2s_ddot_B2s()
T Norm_deviatoric(const T *A)
See Cartesian2d::Norm_deviatoric()
void sym(const T *A, T *ret)
See Cartesian2d::Sym()
void II(T *ret)
See Cartesian2d::II()
T Hydrostatic_deviatoric(const T *A, T *ret)
Returns Cartesian2d::Hydrostatic() and computes Cartesian2d::Deviatoric()
T Hydrostatic(const T *A)
See Cartesian2d::Hydrostatic()
array_type::tensor< double, 4 > Random4()
Random 4th-order tensor (for example for use in testing).
Definition: Cartesian2d.hpp:29
array_type::tensor< double, 2 > O2()
2nd-order null tensor (all components equal to zero).
Definition: Cartesian2d.hpp:35
array_type::tensor< double, 2 > I2()
2nd-order identity tensor.
Definition: Cartesian2d.hpp:47
auto Deviatoric(const T &A)
Deviatoric part of a tensor:
auto Norm_deviatoric(const T &A)
Norm of the tensor's deviator:
auto Hydrostatic(const T &A)
Hydrostatic part of a tensor.
void hydrostatic(const T &A, R &ret)
Same as Hydrostatic() but writes to externally allocated output.
auto Sym(const T &A)
Symmetric part of a tensor:
array_type::tensor< double, 4 > II()
Result of the dyadic product of two 2nd-order identity tensors (see I2()).
Definition: Cartesian2d.hpp:54
void deviatoric(const T &A, R &ret)
Same as Deviatoric() but writes to externally allocated output.
array_type::tensor< double, 4 > I4()
Fourth order unit tensor.
Definition: Cartesian2d.hpp:61
auto A2_ddot_B2(const T &A, const T &B)
Double tensor contraction.
void sym(const T &A, R &ret)
Same as Sym() but writes to externally allocated output.
auto A4_ddot_B2(const T &A, const U &B)
Double tensor contraction.
void norm_deviatoric(const T &A, R &ret)
Same as Norm_deviatoric() but writes to externally allocated output.
void trace(const T &A, R &ret)
Same as Trace() but writes to externally allocated output.
Definition: Cartesian2d.hpp:90
array_type::tensor< double, 4 > I4s()
Fourth order symmetric projection.
Definition: Cartesian2d.hpp:75
auto Trace(const T &A)
Trace or 2nd-order tensor.
Definition: Cartesian2d.hpp:96
array_type::tensor< double, 2 > Random2()
Random 2nd-order tensor (for example for use in testing).
Definition: Cartesian2d.hpp:23
auto A2_dot_B2(const T &A, const T &B)
Dot-product (single tensor contraction)
array_type::tensor< double, 4 > O4()
4th-order null tensor (all components equal to zero).
Definition: Cartesian2d.hpp:41
array_type::tensor< double, 4 > I4rt()
Right-transposed fourth order unit tensor.
Definition: Cartesian2d.hpp:68
array_type::tensor< double, 4 > I4d()
Fourth order deviatoric projection.
Definition: Cartesian2d.hpp:82
auto A2s_ddot_B2s(const T &A, const T &B)
Same as A2_ddot_B2(const T& A, const T& B, R& ret) but for symmetric tensors.
auto A2_dyadic_B2(const T &A, const T &B)
Dyadic product.
xt::xtensor< T, N > tensor
Fixed (static) rank array.
Definition: config.h:86
Tensor products / operations.
Definition: Cartesian2d.h:20