Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.4.1
\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}} \newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=} \newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]} \newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
transformations.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2016 - 2021 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_transformations_h
17#define dealii_transformations_h
18
19#include <deal.II/base/config.h>
20
22#include <deal.II/base/tensor.h>
23
25
26
27namespace Physics
28{
29 namespace Transformations
30 {
35 namespace Rotations
36 {
41
57 template <typename Number>
59 rotation_matrix_2d(const Number &angle);
60
61
90 template <typename Number>
92 rotation_matrix_3d(const Tensor<1, 3, Number> &axis, const Number &angle);
93
99 template <typename Number>
100 DEAL_II_DEPRECATED_EARLY Tensor<2, 3, Number>
101 rotation_matrix_3d(const Point<3, Number> &axis, const Number &angle);
102
104
105 } // namespace Rotations
106
123 namespace Contravariant
124 {
129
143 template <int dim, typename Number>
146 const Tensor<2, dim, Number> &F);
147
162 template <int dim, typename Number>
165 const Tensor<2, dim, Number> &F);
166
182 template <int dim, typename Number>
185 const Tensor<2, dim, Number> & F);
186
201 template <int dim, typename Number>
204 const Tensor<2, dim, Number> &F);
205
221 template <int dim, typename Number>
224 const Tensor<2, dim, Number> & F);
225
227
232
246 template <int dim, typename Number>
249 const Tensor<2, dim, Number> &F);
250
265 template <int dim, typename Number>
268 const Tensor<2, dim, Number> &F);
269
284 template <int dim, typename Number>
287 const Tensor<2, dim, Number> & F);
288
303 template <int dim, typename Number>
306 const Tensor<2, dim, Number> &F);
307
322 template <int dim, typename Number>
325 const Tensor<2, dim, Number> & F);
326
328 } // namespace Contravariant
329
348 namespace Covariant
349 {
354
368 template <int dim, typename Number>
371 const Tensor<2, dim, Number> &F);
372
387 template <int dim, typename Number>
390 const Tensor<2, dim, Number> &F);
391
407 template <int dim, typename Number>
410 const Tensor<2, dim, Number> & F);
411
426 template <int dim, typename Number>
429 const Tensor<2, dim, Number> &F);
430
446 template <int dim, typename Number>
449 const Tensor<2, dim, Number> & F);
450
452
457
471 template <int dim, typename Number>
474 const Tensor<2, dim, Number> &F);
475
490 template <int dim, typename Number>
493 const Tensor<2, dim, Number> &F);
494
509 template <int dim, typename Number>
512 const Tensor<2, dim, Number> & F);
513
528 template <int dim, typename Number>
531 const Tensor<2, dim, Number> &F);
532
547 template <int dim, typename Number>
550 const Tensor<2, dim, Number> & F);
551
553 } // namespace Covariant
554
560 namespace Piola
561 {
566
582 template <int dim, typename Number>
585 const Tensor<2, dim, Number> &F);
586
602 template <int dim, typename Number>
605 const Tensor<2, dim, Number> &F);
606
623 template <int dim, typename Number>
626 const Tensor<2, dim, Number> & F);
627
644 template <int dim, typename Number>
647 const Tensor<2, dim, Number> &F);
648
666 template <int dim, typename Number>
669 const Tensor<2, dim, Number> & F);
670
672
677
693 template <int dim, typename Number>
696 const Tensor<2, dim, Number> &F);
697
713 template <int dim, typename Number>
716 const Tensor<2, dim, Number> &F);
717
733 template <int dim, typename Number>
736 const Tensor<2, dim, Number> & F);
737
754 template <int dim, typename Number>
757 const Tensor<2, dim, Number> &F);
758
775 template <int dim, typename Number>
778 const Tensor<2, dim, Number> & F);
779
781 } // namespace Piola
782
787
810 template <int dim, typename Number>
813 const Tensor<2, dim, Number> &F);
814
816
821
832 template <int dim, typename Number>
835 const Tensor<2, dim, Number> &B);
836
848 template <int dim, typename Number>
851 const Tensor<2, dim, Number> &B);
852
864 template <int dim, typename Number>
867 const Tensor<2, dim, Number> & B);
868
879 template <int dim, typename Number>
882 const Tensor<2, dim, Number> &B);
883
895 template <int dim, typename Number>
898 const Tensor<2, dim, Number> & B);
899
901
902 } // namespace Transformations
903} // namespace Physics
904
905
906
907#ifndef DOXYGEN
908
909
910
911template <typename Number>
914{
915 const Number rotation[2][2] = {{std::cos(angle), -std::sin(angle)},
916 {std::sin(angle), std::cos(angle)}};
917 return Tensor<2, 2>(rotation);
918}
919
920
921
922template <typename Number>
925 const Tensor<1, 3, Number> &axis,
926 const Number & angle)
927{
928 Assert(std::abs(axis.norm() - 1.0) < 1e-9,
929 ExcMessage("The supplied axial vector is not a unit vector."));
930 const Number c = std::cos(angle);
931 const Number s = std::sin(angle);
932 const Number t = 1. - c;
933 const Number rotation[3][3] = {{t * axis[0] * axis[0] + c,
934 t * axis[0] * axis[1] - s * axis[2],
935 t * axis[0] * axis[2] + s * axis[1]},
936 {t * axis[0] * axis[1] + s * axis[2],
937 t * axis[1] * axis[1] + c,
938 t * axis[1] * axis[2] - s * axis[0]},
939 {t * axis[0] * axis[2] - s * axis[1],
940 t * axis[1] * axis[2] + s * axis[0],
941 t * axis[2] * axis[2] + c}};
942 return Tensor<2, 3, Number>(rotation);
943}
944
945
946
947template <typename Number>
950 const Point<3, Number> &axis,
951 const Number & angle)
952{
953 return rotation_matrix_3d(static_cast<Tensor<1, 3, Number>>(axis), angle);
954}
955
956
957
958template <int dim, typename Number>
961 const Tensor<1, dim, Number> &V,
962 const Tensor<2, dim, Number> &F)
963{
965}
966
967
968
969template <int dim, typename Number>
972 const Tensor<2, dim, Number> &T,
973 const Tensor<2, dim, Number> &F)
974{
976}
977
978
979
980template <int dim, typename Number>
984 const Tensor<2, dim, Number> & F)
985{
987}
988
989
990
991template <int dim, typename Number>
994 const Tensor<4, dim, Number> &H,
995 const Tensor<2, dim, Number> &F)
996{
998}
999
1000
1001
1002template <int dim, typename Number>
1006 const Tensor<2, dim, Number> & F)
1007{
1009}
1010
1011
1012
1013template <int dim, typename Number>
1016 const Tensor<1, dim, Number> &v,
1017 const Tensor<2, dim, Number> &F)
1018{
1020}
1021
1022
1023
1024template <int dim, typename Number>
1027 const Tensor<2, dim, Number> &t,
1028 const Tensor<2, dim, Number> &F)
1029{
1031}
1032
1033
1034
1035template <int dim, typename Number>
1039 const Tensor<2, dim, Number> & F)
1040{
1042}
1043
1044
1045
1046template <int dim, typename Number>
1049 const Tensor<4, dim, Number> &h,
1050 const Tensor<2, dim, Number> &F)
1051{
1053}
1054
1055
1056
1057template <int dim, typename Number>
1061 const Tensor<2, dim, Number> & F)
1062{
1064}
1065
1066
1067
1068template <int dim, typename Number>
1071 const Tensor<1, dim, Number> &V,
1072 const Tensor<2, dim, Number> &F)
1073{
1075 transpose(invert(F)));
1076}
1077
1078
1079
1080template <int dim, typename Number>
1083 const Tensor<2, dim, Number> &T,
1084 const Tensor<2, dim, Number> &F)
1085{
1087 transpose(invert(F)));
1088}
1089
1090
1091
1092template <int dim, typename Number>
1096 const Tensor<2, dim, Number> & F)
1097{
1099 transpose(invert(F)));
1100}
1101
1102
1103
1104template <int dim, typename Number>
1107 const Tensor<4, dim, Number> &H,
1108 const Tensor<2, dim, Number> &F)
1109{
1111 transpose(invert(F)));
1112}
1113
1114
1115
1116template <int dim, typename Number>
1120 const Tensor<2, dim, Number> & F)
1121{
1123 transpose(invert(F)));
1124}
1125
1126
1127
1128template <int dim, typename Number>
1131 const Tensor<2, dim, Number> &F)
1132{
1134}
1135
1136
1137
1138template <int dim, typename Number>
1141 const Tensor<2, dim, Number> &F)
1142{
1144}
1145
1146
1147
1148template <int dim, typename Number>
1152 const Tensor<2, dim, Number> & F)
1153{
1155}
1156
1157
1158
1159template <int dim, typename Number>
1162 const Tensor<2, dim, Number> &F)
1163{
1165}
1166
1167
1168
1169template <int dim, typename Number>
1173 const Tensor<2, dim, Number> & F)
1174{
1176}
1177
1178
1179
1180template <int dim, typename Number>
1183 const Tensor<2, dim, Number> &F)
1184{
1185 return Number(1.0 / determinant(F)) * Contravariant::push_forward(V, F);
1186}
1187
1188
1189
1190template <int dim, typename Number>
1193 const Tensor<2, dim, Number> &F)
1194{
1195 return Number(1.0 / determinant(F)) * Contravariant::push_forward(T, F);
1196}
1197
1198
1199
1200template <int dim, typename Number>
1204 const Tensor<2, dim, Number> & F)
1205{
1206 return Number(1.0 / determinant(F)) * Contravariant::push_forward(T, F);
1207}
1208
1209
1210
1211template <int dim, typename Number>
1214 const Tensor<2, dim, Number> &F)
1215{
1216 return Number(1.0 / determinant(F)) * Contravariant::push_forward(H, F);
1217}
1218
1219
1220
1221template <int dim, typename Number>
1225 const Tensor<2, dim, Number> & F)
1226{
1227 return Number(1.0 / determinant(F)) * Contravariant::push_forward(H, F);
1228}
1229
1230
1231
1232template <int dim, typename Number>
1235 const Tensor<2, dim, Number> &F)
1236{
1237 return Number(determinant(F)) * Contravariant::pull_back(v, F);
1238}
1239
1240
1241
1242template <int dim, typename Number>
1245 const Tensor<2, dim, Number> &F)
1246{
1247 return Number(determinant(F)) * Contravariant::pull_back(t, F);
1248}
1249
1250
1251
1252template <int dim, typename Number>
1256 const Tensor<2, dim, Number> & F)
1257{
1258 return Number(determinant(F)) * Contravariant::pull_back(t, F);
1259}
1260
1261
1262
1263template <int dim, typename Number>
1266 const Tensor<2, dim, Number> &F)
1267{
1268 return Number(determinant(F)) * Contravariant::pull_back(h, F);
1269}
1270
1271
1272
1273template <int dim, typename Number>
1277 const Tensor<2, dim, Number> & F)
1278{
1279 return Number(determinant(F)) * Contravariant::pull_back(h, F);
1280}
1281
1282
1283
1284template <int dim, typename Number>
1287 const Tensor<2, dim, Number> &F)
1288{
1289 return cofactor(F) * N;
1290}
1291
1292
1293template <int dim, typename Number>
1296 const Tensor<2, dim, Number> &B)
1297{
1298 return contract<1, 0>(B, V);
1299}
1300
1301
1302
1303template <int dim, typename Number>
1306 const Tensor<2, dim, Number> &B)
1307{
1308 return contract<1, 0>(B, contract<1, 1>(T, B));
1309}
1310
1311
1312
1313template <int dim, typename Number>
1317 const Tensor<2, dim, Number> & B)
1318{
1320 for (unsigned int i = 0; i < dim; ++i)
1321 for (unsigned int J = 0; J < dim; ++J)
1322 // Loop over I but complex.h defines a macro I, so use I_ instead
1323 for (unsigned int I_ = 0; I_ < dim; ++I_)
1324 tmp_1[i][J] += B[i][I_] * T[I_][J];
1325
1327 for (unsigned int i = 0; i < dim; ++i)
1328 for (unsigned int j = i; j < dim; ++j)
1329 for (unsigned int J = 0; J < dim; ++J)
1330 out[i][j] += B[j][J] * tmp_1[i][J];
1331
1332 return out;
1333}
1334
1335
1336
1337template <int dim, typename Number>
1340 const Tensor<2, dim, Number> &B)
1341{
1342 // This contraction order and indexing might look a bit dubious, so a
1343 // quick explanation as to what's going on is probably in order:
1344 //
1345 // When the contract() function operates on the inner indices, the
1346 // result has the inner index and outer index transposed, i.e.
1347 // contract<2,1>(H,F) implies
1348 // T_{IJLk} = (H_{IJMN} F_{mM}) \delta_{mL} \delta_{Nk}
1349 // rather than T_{IJkL} (the desired result).
1350 // So, in effect, contraction of the 3rd (inner) index with F as the
1351 // second argument results in its transposition with respect to its
1352 // adjacent neighbor. This is due to the position of the argument F,
1353 // leading to the free index being on the right hand side of the result.
1354 // However, given that we can do two transformations from the LHS of H
1355 // and two from the right we can undo the otherwise erroneous
1356 // swapping of the outer indices upon application of the second
1357 // sets of contractions.
1358 //
1359 // Note: Its significantly quicker (in 3d) to push forward
1360 // each index individually
1361 return contract<1, 1>(
1362 B, contract<1, 1>(B, contract<2, 1>(contract<2, 1>(H, B), B)));
1363}
1364
1365
1366
1367template <int dim, typename Number>
1371 const Tensor<2, dim, Number> & B)
1372{
1373 // The first and last transformation operations respectively
1374 // break and recover the symmetry properties of the tensors.
1375 // We also want to perform a minimal number of operations here
1376 // and avoid some complications related to the transposition of
1377 // tensor indices when contracting inner indices using the contract()
1378 // function. (For an explanation of the contraction operations,
1379 // please see the note in the equivalent function for standard
1380 // Tensors.) So what we'll do here is manually perform the first
1381 // and last contractions that break/recover the tensor symmetries
1382 // on the inner indices, and use the contract() function only on
1383 // the outer indices.
1384 //
1385 // Note: Its significantly quicker (in 3d) to push forward
1386 // each index individually
1387
1388 // Push forward (inner) index 1
1390 // Loop over I but complex.h defines a macro I, so use I_ instead
1391 for (unsigned int I_ = 0; I_ < dim; ++I_)
1392 for (unsigned int j = 0; j < dim; ++j)
1393 for (unsigned int K = 0; K < dim; ++K)
1394 for (unsigned int L = 0; L < dim; ++L)
1395 for (unsigned int J = 0; J < dim; ++J)
1396 tmp[I_][j][K][L] += B[j][J] * H[I_][J][K][L];
1397
1398 // Push forward (outer) indices 0 and 3
1399 tmp = contract<1, 0>(B, contract<3, 1>(tmp, B));
1400
1401 // Push forward (inner) index 2
1403 for (unsigned int i = 0; i < dim; ++i)
1404 for (unsigned int j = i; j < dim; ++j)
1405 for (unsigned int k = 0; k < dim; ++k)
1406 for (unsigned int l = k; l < dim; ++l)
1407 for (unsigned int K = 0; K < dim; ++K)
1408 out[i][j][k][l] += B[k][K] * tmp[i][j][K][l];
1409
1410 return out;
1411}
1412
1413#endif // DOXYGEN
1414
1416
1417#endif
Definition: point.h:111
Definition: tensor.h:503
numbers::NumberTraits< Number >::real_type norm() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcMessage(std::string arg1)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Tensor< 1, dim, Number > pull_back(const Tensor< 1, dim, Number > &v, const Tensor< 2, dim, Number > &F)
Tensor< 1, dim, Number > push_forward(const Tensor< 1, dim, Number > &V, const Tensor< 2, dim, Number > &F)
Tensor< 1, dim, Number > push_forward(const Tensor< 1, dim, Number > &V, const Tensor< 2, dim, Number > &F)
Tensor< 1, dim, Number > pull_back(const Tensor< 1, dim, Number > &v, const Tensor< 2, dim, Number > &F)
Tensor< 1, dim, Number > push_forward(const Tensor< 1, dim, Number > &V, const Tensor< 2, dim, Number > &F)
Tensor< 1, dim, Number > pull_back(const Tensor< 1, dim, Number > &v, const Tensor< 2, dim, Number > &F)
Tensor< 2, 3, Number > rotation_matrix_3d(const Tensor< 1, 3, Number > &axis, const Number &angle)
Tensor< 2, 2, Number > rotation_matrix_2d(const Number &angle)
Tensor< 1, dim, Number > nansons_formula(const Tensor< 1, dim, Number > &N, const Tensor< 2, dim, Number > &F)
Tensor< 1, dim, Number > basis_transformation(const Tensor< 1, dim, Number > &V, const Tensor< 2, dim, Number > &B)
::VectorizedArray< Number, width > cos(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sin(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
constexpr Number determinant(const SymmetricTensor< 2, dim, Number > &)
constexpr SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)
constexpr Tensor< 2, dim, Number > cofactor(const Tensor< 2, dim, Number > &t)
Definition: tensor.h:2926