deal.II version GIT relicensing-2173-gae8fc9d14b 2024-11-24 06:40:00+00:00
\(\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\}}\)
Loading...
Searching...
No Matches
notation.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2017 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_physics_notation_h
16#define dealii_physics_notation_h
17
18#include <deal.II/base/config.h>
19
23#include <deal.II/base/tensor.h>
24
26#include <deal.II/lac/vector.h>
27
28#include <type_traits>
29
31
32
33namespace Physics
34{
35 namespace Notation
36 {
258 namespace Kelvin
259 {
264 int,
265 int,
266 int,
267 << "The number of rows in the input matrix is " << arg1
268 << ", but needs to be either " << arg2 << " or " << arg3
269 << '.');
270
271
276 int,
277 int,
278 int,
279 int,
280 << "The number of rows in the input matrix is " << arg1
281 << ", but needs to be either " << arg2 << ',' << arg3
282 << ", or " << arg4 << '.');
283
284
289 int,
290 int,
291 int,
292 << "The number of columns in the input matrix is " << arg1
293 << ", but needs to be either " << arg2 << " or " << arg3
294 << '.');
295
296
301 int,
302 int,
303 int,
304 int,
305 << "The number of columns in the input matrix is " << arg1
306 << ", but needs to be either " << arg2 << ',' << arg3
307 << ", or " << arg4 << '.');
308
309
320 template <typename Number>
322 to_vector(const Number &s);
323
324
330 template <int dim, typename Number>
333
334
340 template <int dim, typename Number>
343
344
350 template <int dim, typename Number>
353
354
361 template <int dim, typename Number>
364
365
371 template <typename Number>
373 to_matrix(const Number &s);
374
375
381 template <int dim, typename Number>
384
385
391 template <int dim, typename Number>
394
395
401 template <int dim, typename Number>
404
405
415 template <int dim, typename Number>
418
451 template <int dim,
452 typename SubTensor1 = Tensor<2, dim>,
453 typename SubTensor2 = Tensor<1, dim>,
454 typename Number>
457
458
465 template <int dim, typename Number>
468
469
477 template <int dim, typename Number>
480
491 template <typename Number>
492 void
493 to_tensor(const Vector<Number> &vec, Number &s);
494
495
499 template <int dim, typename Number>
500 void
502
503
507 template <int dim, typename Number>
508 void
510
511
515 template <int dim, typename Number>
516 void
518
519
523 template <int dim, typename Number>
524 void
526
527
531 template <typename Number>
532 void
533 to_tensor(const FullMatrix<Number> &mtrx, Number &s);
534
535
539 template <int dim, typename Number>
540 void
542
543
547 template <int dim, typename Number>
548 void
550
551
555 template <int dim, typename Number>
556 void
558
559
563 template <int dim, typename Number>
564 void
567
568
578 template <int dim, typename Number>
579 void
581
582
586 template <int dim, typename Number>
587 void
589
590
594 template <int dim, typename Number>
595 void
598
599
604 template <typename TensorType, typename Number>
605 TensorType
607
608
613 template <typename TensorType, typename Number>
614 TensorType
618 } // namespace Kelvin
619
620 } // namespace Notation
621} // namespace Physics
622
623
624#ifndef DOXYGEN
625
626
627// ------------------------- inline functions ------------------------
628
629
630namespace Physics
631{
632 namespace Notation
633 {
634 namespace Kelvin
635 {
636 namespace internal
637 {
645 template <int dim>
646 std::pair<unsigned int, unsigned int>
647 indices_from_component(const unsigned int component_n,
648 const bool symmetric);
649
650
651 template <int dim>
652 std::pair<unsigned int, unsigned int>
653 indices_from_component(const unsigned int /*component_n*/, const bool)
654 {
656 return std::make_pair(0u, 0u);
657 }
658
659
660 template <>
661 inline std::pair<unsigned int, unsigned int>
662 indices_from_component<1>(const unsigned int component_n, const bool)
663 {
664 AssertIndexRange(component_n, 1);
665
666 return std::make_pair(0u, 0u);
667 }
668
669
670 template <>
671 inline std::pair<unsigned int, unsigned int>
672 indices_from_component<2>(const unsigned int component_n,
673 const bool symmetric)
674 {
675 if (symmetric == true)
676 {
677 Assert(
679 ExcIndexRange(component_n,
680 0,
682 }
683 else
684 {
686 ExcIndexRange(component_n,
687 0,
689 }
690
691 static const unsigned int indices[4][2] = {{0, 0},
692 {1, 1},
693 {0, 1},
694 {1, 0}};
695 return std::make_pair(indices[component_n][0],
696 indices[component_n][1]);
697 }
698
699
700 template <>
701 inline std::pair<unsigned int, unsigned int>
702 indices_from_component<3>(const unsigned int component_n,
703 const bool symmetric)
704 {
705 if (symmetric == true)
706 {
707 Assert(
709 ExcIndexRange(component_n,
710 0,
712 }
713 else
714 {
716 ExcIndexRange(component_n,
717 0,
719 }
720
721 static const unsigned int indices[9][2] = {{0, 0},
722 {1, 1},
723 {2, 2},
724 {1, 2},
725 {0, 2},
726 {0, 1},
727 {1, 0},
728 {2, 0},
729 {2, 1}};
730 return std::make_pair(indices[component_n][0],
731 indices[component_n][1]);
732 }
733
734
739 template <int dim>
740 double
741 vector_component_factor(const unsigned int component_i,
742 const bool symmetric)
743 {
744 if (symmetric == false)
745 return 1.0;
746
747 if (component_i < dim)
748 return 1.0;
749 else
750 return numbers::SQRT2;
751 }
752
753
758 template <int dim>
759 double
760 matrix_component_factor(const unsigned int component_i,
761 const unsigned int component_j,
762 const bool symmetric)
763 {
764 if (symmetric == false)
765 return 1.0;
766
767 // This case check returns equivalent of this result:
768 // internal::vector_component_factor<dim>(component_i,symmetric)*internal::vector_component_factor<dim>(component_j,symmetric);
769 if (component_i < dim && component_j < dim)
770 return 1.0;
771 else if (component_i >= dim && component_j >= dim)
772 return 2.0;
773 else // ((component_i >= dim && component_j < dim) || (component_i <
774 // dim && component_j >= dim))
775 return numbers::SQRT2;
776 }
777
778 } // namespace internal
779
780
781 template <typename Number>
783 to_vector(const Number &s)
784 {
785 Vector<Number> out(1);
786 const unsigned int n_rows = out.size();
787 for (unsigned int r = 0; r < n_rows; ++r)
788 out(r) = s;
789 return out;
790 }
791
792
793 template <int dim, typename Number>
796 {
797 return to_vector(s.operator const Number &());
798 }
799
800
801 template <int dim, typename Number>
804 {
806 const unsigned int n_rows = out.size();
807 for (unsigned int r = 0; r < n_rows; ++r)
808 {
809 const std::pair<unsigned int, unsigned int> indices =
810 internal::indices_from_component<dim>(r, false);
811 Assert(indices.first < dim, ExcInternalError());
812 const unsigned int i = indices.first;
813 out(r) = v[i];
814 }
815 return out;
816 }
817
818
819 template <int dim, typename Number>
822 {
824 const unsigned int n_rows = out.size();
825 for (unsigned int r = 0; r < n_rows; ++r)
826 {
827 const std::pair<unsigned int, unsigned int> indices =
828 internal::indices_from_component<dim>(r, false);
829 Assert(indices.first < dim, ExcInternalError());
830 Assert(indices.second < dim, ExcInternalError());
831 const unsigned int i = indices.first;
832 const unsigned int j = indices.second;
833 out(r) = t[i][j];
834 }
835 return out;
836 }
837
838
839 template <int dim, typename Number>
842 {
844 const unsigned int n_rows = out.size();
845 for (unsigned int r = 0; r < n_rows; ++r)
846 {
847 const std::pair<unsigned int, unsigned int> indices =
848 internal::indices_from_component<dim>(r, true);
849 Assert(indices.first < dim, ExcInternalError());
850 Assert(indices.second < dim, ExcInternalError());
851 Assert(indices.second >= indices.first, ExcInternalError());
852 const unsigned int i = indices.first;
853 const unsigned int j = indices.second;
854
855 const double factor =
856 internal::vector_component_factor<dim>(r, true);
857
858 out(r) = factor * st[i][j];
859 }
860 return out;
861 }
862
863
864 template <typename Number>
866 to_matrix(const Number &s)
867 {
868 FullMatrix<Number> out(1, 1);
869 out(0, 0) = s;
870 return out;
871 }
872
873
874 template <int dim, typename Number>
877 {
878 return to_matrix(s.operator const Number &());
879 }
880
881
882 template <int dim, typename Number>
885 {
887 const unsigned int n_rows = out.m();
888 const unsigned int n_cols = out.n();
889 for (unsigned int r = 0; r < n_rows; ++r)
890 {
891 const std::pair<unsigned int, unsigned int> indices =
892 internal::indices_from_component<dim>(r, false);
893 Assert(indices.first < dim, ExcInternalError());
894 const unsigned int i = indices.first;
895
896 for (unsigned int c = 0; c < n_cols; ++c)
897 {
898 Assert(c < 1, ExcInternalError());
899 out(r, c) = v[i];
900 }
901 }
902 return out;
903 }
904
905
906 template <int dim, typename Number>
909 {
910 FullMatrix<Number> out(dim, dim);
911 const unsigned int n_rows = out.m();
912 const unsigned int n_cols = out.n();
913 for (unsigned int r = 0; r < n_rows; ++r)
914 {
915 const std::pair<unsigned int, unsigned int> indices_i =
916 internal::indices_from_component<dim>(r, false);
917 Assert(indices_i.first < dim, ExcInternalError());
918 Assert(indices_i.second < dim, ExcInternalError());
919 const unsigned int i = indices_i.first;
920
921 for (unsigned int c = 0; c < n_cols; ++c)
922 {
923 const std::pair<unsigned int, unsigned int> indices_j =
924 internal::indices_from_component<dim>(c, false);
925 Assert(indices_j.first < dim, ExcInternalError());
926 Assert(indices_j.second < dim, ExcInternalError());
927 const unsigned int j = indices_j.second;
928
929 out(r, c) = t[i][j];
930 }
931 }
932 return out;
933 }
934
935
936 template <int dim, typename Number>
939 {
941 }
942
943
944 namespace internal
945 {
946 template <typename TensorType>
947 struct is_rank_2_symmetric_tensor : std::false_type
948 {};
949
950 template <int dim, typename Number>
951 struct is_rank_2_symmetric_tensor<SymmetricTensor<2, dim, Number>>
952 : std::true_type
953 {};
954 } // namespace internal
955
956
957 template <int dim,
958 typename SubTensor1,
959 typename SubTensor2,
960 typename Number>
963 {
964 static_assert(
965 (SubTensor1::dimension == dim && SubTensor2::dimension == dim),
966 "Sub-tensor spatial dimension is different from those of the input tensor.");
967
968 static_assert(
969 (SubTensor1::rank == 2 && SubTensor2::rank == 1) ||
970 (SubTensor1::rank == 1 && SubTensor2::rank == 2),
971 "Cannot build a rank 3 tensor from the given combination of sub-tensors.");
972
973 FullMatrix<Number> out(SubTensor1::n_independent_components,
974 SubTensor2::n_independent_components);
975 const unsigned int n_rows = out.m();
976 const unsigned int n_cols = out.n();
977
978 if (SubTensor1::rank == 2 && SubTensor2::rank == 1)
979 {
980 const bool subtensor_is_rank_2_symmetric_tensor =
981 internal::is_rank_2_symmetric_tensor<SubTensor1>::value;
982
983 for (unsigned int r = 0; r < n_rows; ++r)
984 {
985 const std::pair<unsigned int, unsigned int> indices_ij =
986 internal::indices_from_component<dim>(
987 r, subtensor_is_rank_2_symmetric_tensor);
988 Assert(indices_ij.first < dim, ExcInternalError());
989 Assert(indices_ij.second < dim, ExcInternalError());
990 if (subtensor_is_rank_2_symmetric_tensor)
991 {
992 Assert(indices_ij.second >= indices_ij.first,
994 }
995 const unsigned int i = indices_ij.first;
996 const unsigned int j = indices_ij.second;
997
998 const double factor = internal::vector_component_factor<dim>(
999 r, subtensor_is_rank_2_symmetric_tensor);
1000
1001 for (unsigned int c = 0; c < n_cols; ++c)
1002 {
1003 const std::pair<unsigned int, unsigned int> indices_k =
1004 internal::indices_from_component<dim>(c, false);
1005 Assert(indices_k.first < dim, ExcInternalError());
1006 const unsigned int k = indices_k.first;
1007
1008 if (subtensor_is_rank_2_symmetric_tensor)
1009 out(r, c) = factor * t[i][j][k];
1010 else
1011 out(r, c) = t[i][j][k];
1012 }
1013 }
1014 }
1015 else if (SubTensor1::rank == 1 && SubTensor2::rank == 2)
1016 {
1017 const bool subtensor_is_rank_2_symmetric_tensor =
1018 internal::is_rank_2_symmetric_tensor<SubTensor2>::value;
1019
1020 for (unsigned int r = 0; r < n_rows; ++r)
1021 {
1022 const std::pair<unsigned int, unsigned int> indices_k =
1023 internal::indices_from_component<dim>(r, false);
1024 Assert(indices_k.first < dim, ExcInternalError());
1025 const unsigned int k = indices_k.first;
1026
1027 for (unsigned int c = 0; c < n_cols; ++c)
1028 {
1029 const std::pair<unsigned int, unsigned int> indices_ij =
1030 internal::indices_from_component<dim>(
1031 c, subtensor_is_rank_2_symmetric_tensor);
1032 Assert(indices_ij.first < dim, ExcInternalError());
1033 Assert(indices_ij.second < dim, ExcInternalError());
1034 if (subtensor_is_rank_2_symmetric_tensor)
1035 {
1036 Assert(indices_ij.second >= indices_ij.first,
1038 }
1039 const unsigned int i = indices_ij.first;
1040 const unsigned int j = indices_ij.second;
1041
1042 if (subtensor_is_rank_2_symmetric_tensor)
1043 {
1044 const double factor =
1045 internal::vector_component_factor<dim>(
1046 c, subtensor_is_rank_2_symmetric_tensor);
1047 out(r, c) = factor * t[k][i][j];
1048 }
1049 else
1050 out(r, c) = t[k][i][j];
1051 }
1052 }
1053 }
1054 else
1055 {
1057 }
1058
1059 return out;
1060 }
1061
1062
1063 template <int dim, typename Number>
1066 {
1070 const unsigned int n_rows = out.m();
1071 const unsigned int n_cols = out.n();
1072 for (unsigned int r = 0; r < n_rows; ++r)
1073 {
1074 const std::pair<unsigned int, unsigned int> indices_ij =
1075 internal::indices_from_component<dim>(r, false);
1076 Assert(indices_ij.first < dim, ExcInternalError());
1077 Assert(indices_ij.second < dim, ExcInternalError());
1078 const unsigned int i = indices_ij.first;
1079 const unsigned int j = indices_ij.second;
1080
1081 for (unsigned int c = 0; c < n_cols; ++c)
1082 {
1083 const std::pair<unsigned int, unsigned int> indices_kl =
1084 internal::indices_from_component<dim>(c, false);
1085 Assert(indices_kl.first < dim, ExcInternalError());
1086 Assert(indices_kl.second < dim, ExcInternalError());
1087 const unsigned int k = indices_kl.first;
1088 const unsigned int l = indices_kl.second;
1089
1090 out(r, c) = t[i][j][k][l];
1091 }
1092 }
1093 return out;
1094 }
1095
1096
1097 template <int dim, typename Number>
1100 {
1104 const unsigned int n_rows = out.m();
1105 const unsigned int n_cols = out.n();
1106 for (unsigned int r = 0; r < n_rows; ++r)
1107 {
1108 const std::pair<unsigned int, unsigned int> indices_ij =
1109 internal::indices_from_component<dim>(r, true);
1110 Assert(indices_ij.first < dim, ExcInternalError());
1111 Assert(indices_ij.second < dim, ExcInternalError());
1112 Assert(indices_ij.second >= indices_ij.first, ExcInternalError());
1113 const unsigned int i = indices_ij.first;
1114 const unsigned int j = indices_ij.second;
1115
1116 for (unsigned int c = 0; c < n_cols; ++c)
1117 {
1118 const std::pair<unsigned int, unsigned int> indices_kl =
1119 internal::indices_from_component<dim>(c, true);
1120 Assert(indices_kl.first < dim, ExcInternalError());
1121 Assert(indices_kl.second < dim, ExcInternalError());
1122 Assert(indices_kl.second >= indices_kl.first,
1124 const unsigned int k = indices_kl.first;
1125 const unsigned int l = indices_kl.second;
1126
1127 const double factor =
1128 internal::matrix_component_factor<dim>(r, c, true);
1129
1130 out(r, c) = factor * st[i][j][k][l];
1131 }
1132 }
1133 return out;
1134 }
1135
1136
1137 template <typename Number>
1138 void
1139 to_tensor(const Vector<Number> &vec, Number &s)
1140 {
1141 Assert(vec.size() == 1, ExcDimensionMismatch(vec.size(), 1));
1142 s = vec(0);
1143 }
1144
1145
1146 template <int dim, typename Number>
1147 void
1149 {
1150 return to_tensor(vec, s.operator Number &());
1151 }
1152
1153
1154 template <int dim, typename Number>
1155 void
1157 {
1160 const unsigned int n_rows = vec.size();
1161 for (unsigned int r = 0; r < n_rows; ++r)
1162 {
1163 const std::pair<unsigned int, unsigned int> indices =
1164 internal::indices_from_component<dim>(r, false);
1165 Assert(indices.first < dim, ExcInternalError());
1166 const unsigned int i = indices.first;
1167 v[i] = vec(r);
1168 }
1169 }
1170
1171
1172 template <int dim, typename Number>
1173 void
1175 {
1178 const unsigned int n_rows = vec.size();
1179 for (unsigned int r = 0; r < n_rows; ++r)
1180 {
1181 const std::pair<unsigned int, unsigned int> indices =
1182 internal::indices_from_component<dim>(r, false);
1183 Assert(indices.first < dim, ExcInternalError());
1184 Assert(indices.second < dim, ExcInternalError());
1185 const unsigned int i = indices.first;
1186 const unsigned int j = indices.second;
1187 t[i][j] = vec(r);
1188 }
1189 }
1190
1191
1192 template <int dim, typename Number>
1193 void
1195 {
1198 const unsigned int n_rows = vec.size();
1199 for (unsigned int r = 0; r < n_rows; ++r)
1200 {
1201 const std::pair<unsigned int, unsigned int> indices =
1202 internal::indices_from_component<dim>(r, true);
1203 Assert(indices.first < dim, ExcInternalError());
1204 Assert(indices.second < dim, ExcInternalError());
1205 Assert(indices.second >= indices.first, ExcInternalError());
1206 const unsigned int i = indices.first;
1207 const unsigned int j = indices.second;
1208
1209 const double inv_factor =
1210 1.0 / internal::vector_component_factor<dim>(r, true);
1211
1212 st[i][j] = inv_factor * vec(r);
1213 }
1214 }
1215
1216
1217 template <typename Number>
1218 void
1219 to_tensor(const FullMatrix<Number> &mtrx, Number &s)
1220 {
1221 Assert(mtrx.m() == 1, ExcDimensionMismatch(mtrx.m(), 1));
1222 Assert(mtrx.n() == 1, ExcDimensionMismatch(mtrx.n(), 1));
1223 Assert(mtrx.n_elements() == 1,
1224 ExcDimensionMismatch(mtrx.n_elements(), 1));
1225 s = mtrx(0, 0);
1226 }
1227
1228
1229 template <int dim, typename Number>
1230 void
1232 {
1233 return to_tensor(mtrx, s.operator Number &());
1234 }
1235
1236
1237 template <int dim, typename Number>
1238 void
1240 {
1241 Assert(mtrx.m() == dim, ExcDimensionMismatch(mtrx.m(), dim));
1242 Assert(mtrx.n() == 1, ExcDimensionMismatch(mtrx.n(), 1));
1243 Assert(mtrx.n_elements() == v.n_independent_components,
1244 ExcDimensionMismatch(mtrx.n_elements(),
1246
1247 const unsigned int n_rows = mtrx.m();
1248 const unsigned int n_cols = mtrx.n();
1249 for (unsigned int r = 0; r < n_rows; ++r)
1250 {
1251 const std::pair<unsigned int, unsigned int> indices =
1252 internal::indices_from_component<dim>(r, false);
1253 Assert(indices.first < dim, ExcInternalError());
1254 Assert(indices.second == 0, ExcInternalError());
1255 const unsigned int i = indices.first;
1256
1257 for (unsigned int c = 0; c < n_cols; ++c)
1258 {
1259 Assert(c < 1, ExcInternalError());
1260 v[i] = mtrx(r, c);
1261 }
1262 }
1263 }
1264
1265
1266 template <int dim, typename Number>
1267 void
1269 {
1270 Assert(mtrx.m() == dim, ExcDimensionMismatch(mtrx.m(), dim));
1271 Assert(mtrx.n() == dim, ExcDimensionMismatch(mtrx.n(), dim));
1272 Assert(mtrx.n_elements() == t.n_independent_components,
1273 ExcDimensionMismatch(mtrx.n_elements(),
1275
1276 const unsigned int n_rows = mtrx.m();
1277 const unsigned int n_cols = mtrx.n();
1278 for (unsigned int r = 0; r < n_rows; ++r)
1279 {
1280 const std::pair<unsigned int, unsigned int> indices_i =
1281 internal::indices_from_component<dim>(r, false);
1282 Assert(indices_i.first < dim, ExcInternalError());
1283 Assert(indices_i.second < dim, ExcInternalError());
1284 const unsigned int i = indices_i.first;
1285
1286 for (unsigned int c = 0; c < n_cols; ++c)
1287 {
1288 const std::pair<unsigned int, unsigned int> indices_j =
1289 internal::indices_from_component<dim>(c, false);
1290 Assert(indices_j.first < dim, ExcInternalError());
1291 Assert(indices_j.second < dim, ExcInternalError());
1292 const unsigned int j = indices_j.second;
1293
1294 t[i][j] = mtrx(r, c);
1295 }
1296 }
1297 }
1298
1299
1300 template <int dim, typename Number>
1301 void
1302 to_tensor(const FullMatrix<Number> &mtrx,
1304 {
1305 // Its impossible to fit the (dim^2 + dim)/2 entries into a square
1306 // matrix We therefore assume that its been converted to a standard
1307 // tensor format using to_matrix (SymmetricTensor<2,dim,Number>) at some
1308 // point...
1309 Assert(mtrx.m() == dim, ExcDimensionMismatch(mtrx.m(), dim));
1310 Assert(mtrx.n() == dim, ExcDimensionMismatch(mtrx.n(), dim));
1311 Assert((mtrx.n_elements() ==
1314 mtrx.n_elements(),
1316
1318 to_tensor(mtrx, tmp);
1319 st = symmetrize(tmp);
1320 Assert((Tensor<2, dim, Number>(st) - tmp).norm() < 1e-12,
1321 ExcMessage(
1322 "The entries stored inside the matrix were not symmetric"));
1323 }
1324
1325
1326 template <int dim, typename Number>
1327 void
1329 {
1331 (mtrx.m() ==
1333 (mtrx.m() ==
1336 mtrx.m(),
1341 (mtrx.n() ==
1343 (mtrx.n() ==
1346 mtrx.n(),
1350
1351 const unsigned int n_rows = mtrx.m();
1352 const unsigned int n_cols = mtrx.n();
1354 {
1355 Assert(
1357 (mtrx.m() ==
1360 mtrx.m(),
1363
1364 const bool subtensor_is_rank_2_symmetric_tensor =
1365 (mtrx.m() ==
1367
1368 for (unsigned int r = 0; r < n_rows; ++r)
1369 {
1370 const std::pair<unsigned int, unsigned int> indices_ij =
1371 internal::indices_from_component<dim>(
1372 r, subtensor_is_rank_2_symmetric_tensor);
1373 Assert(indices_ij.first < dim, ExcInternalError());
1374 Assert(indices_ij.second < dim, ExcInternalError());
1375 if (subtensor_is_rank_2_symmetric_tensor)
1376 {
1377 Assert(indices_ij.second >= indices_ij.first,
1379 }
1380 const unsigned int i = indices_ij.first;
1381 const unsigned int j = indices_ij.second;
1382
1383 const double inv_factor =
1384 1.0 / internal::vector_component_factor<dim>(
1385 r, subtensor_is_rank_2_symmetric_tensor);
1386
1387 for (unsigned int c = 0; c < n_cols; ++c)
1388 {
1389 const std::pair<unsigned int, unsigned int> indices_k =
1390 internal::indices_from_component<dim>(c, false);
1391 Assert(indices_k.first < dim, ExcInternalError());
1392 const unsigned int k = indices_k.first;
1393
1394 if (subtensor_is_rank_2_symmetric_tensor)
1395 {
1396 t[i][j][k] = inv_factor * mtrx(r, c);
1397 t[j][i][k] = t[i][j][k];
1398 }
1399 else
1400 t[i][j][k] = mtrx(r, c);
1401 }
1402 }
1403 }
1404 else
1405 {
1406 Assert(
1410 Assert(
1412 (mtrx.n() ==
1415 mtrx.n(),
1418
1419 const bool subtensor_is_rank_2_symmetric_tensor =
1420 (mtrx.n() ==
1422
1423 for (unsigned int r = 0; r < n_rows; ++r)
1424 {
1425 const std::pair<unsigned int, unsigned int> indices_k =
1426 internal::indices_from_component<dim>(r, false);
1427 Assert(indices_k.first < dim, ExcInternalError());
1428 const unsigned int k = indices_k.first;
1429
1430 for (unsigned int c = 0; c < n_cols; ++c)
1431 {
1432 const std::pair<unsigned int, unsigned int> indices_ij =
1433 internal::indices_from_component<dim>(
1434 c, subtensor_is_rank_2_symmetric_tensor);
1435 Assert(indices_ij.first < dim, ExcInternalError());
1436 Assert(indices_ij.second < dim, ExcInternalError());
1437 if (subtensor_is_rank_2_symmetric_tensor)
1438 {
1439 Assert(indices_ij.second >= indices_ij.first,
1441 }
1442 const unsigned int i = indices_ij.first;
1443 const unsigned int j = indices_ij.second;
1444
1445 if (subtensor_is_rank_2_symmetric_tensor)
1446 {
1447 const double inv_factor =
1448 1.0 / internal::vector_component_factor<dim>(
1449 c, subtensor_is_rank_2_symmetric_tensor);
1450 t[k][i][j] = inv_factor * mtrx(r, c);
1451 t[k][j][i] = t[k][i][j];
1452 }
1453 else
1454 t[k][i][j] = mtrx(r, c);
1455 }
1456 }
1457 }
1458 }
1459
1460
1461 template <int dim, typename Number>
1462 void
1464 {
1471 Assert(mtrx.n_elements() == t.n_independent_components,
1472 ExcDimensionMismatch(mtrx.n_elements(),
1474
1475 const unsigned int n_rows = mtrx.m();
1476 const unsigned int n_cols = mtrx.n();
1477 for (unsigned int r = 0; r < n_rows; ++r)
1478 {
1479 const std::pair<unsigned int, unsigned int> indices_ij =
1480 internal::indices_from_component<dim>(r, false);
1481 Assert(indices_ij.first < dim, ExcInternalError());
1482 Assert(indices_ij.second < dim, ExcInternalError());
1483 const unsigned int i = indices_ij.first;
1484 const unsigned int j = indices_ij.second;
1485
1486 for (unsigned int c = 0; c < n_cols; ++c)
1487 {
1488 const std::pair<unsigned int, unsigned int> indices_kl =
1489 internal::indices_from_component<dim>(c, false);
1490 Assert(indices_kl.first < dim, ExcInternalError());
1491 Assert(indices_kl.second < dim, ExcInternalError());
1492 const unsigned int k = indices_kl.first;
1493 const unsigned int l = indices_kl.second;
1494
1495 t[i][j][k][l] = mtrx(r, c);
1496 }
1497 }
1498 }
1499
1500
1501 template <int dim, typename Number>
1502 void
1503 to_tensor(const FullMatrix<Number> &mtrx,
1505 {
1506 Assert((mtrx.m() ==
1509 mtrx.m(),
1511 Assert((mtrx.n() ==
1514 mtrx.n(),
1516 Assert(mtrx.n_elements() == st.n_independent_components,
1517 ExcDimensionMismatch(mtrx.n_elements(),
1519
1520 const unsigned int n_rows = mtrx.m();
1521 const unsigned int n_cols = mtrx.n();
1522 for (unsigned int r = 0; r < n_rows; ++r)
1523 {
1524 const std::pair<unsigned int, unsigned int> indices_ij =
1525 internal::indices_from_component<dim>(r, false);
1526 Assert(indices_ij.first < dim, ExcInternalError());
1527 Assert(indices_ij.second < dim, ExcInternalError());
1528 const unsigned int i = indices_ij.first;
1529 const unsigned int j = indices_ij.second;
1530
1531 for (unsigned int c = 0; c < n_cols; ++c)
1532 {
1533 const std::pair<unsigned int, unsigned int> indices_kl =
1534 internal::indices_from_component<dim>(c, false);
1535 Assert(indices_kl.first < dim, ExcInternalError());
1536 Assert(indices_kl.second < dim, ExcInternalError());
1537 const unsigned int k = indices_kl.first;
1538 const unsigned int l = indices_kl.second;
1539
1540 const double inv_factor =
1541 1.0 / internal::matrix_component_factor<dim>(r, c, true);
1542
1543 st[i][j][k][l] = inv_factor * mtrx(r, c);
1544 }
1545 }
1546 }
1547
1548
1549 template <typename TensorType, typename Number>
1550 inline TensorType
1551 to_tensor(const Vector<Number> &vec)
1552 {
1553 TensorType out;
1554 to_tensor(vec, out);
1555 return out;
1556 }
1557
1558
1559 template <typename TensorType, typename Number>
1560 inline TensorType
1561 to_tensor(const FullMatrix<Number> &mtrx)
1562 {
1563 TensorType out;
1564 to_tensor(mtrx, out);
1565 return out;
1566 }
1567
1568 } // namespace Kelvin
1569 } // namespace Notation
1570} // namespace Physics
1571
1572
1573#endif // DOXYGEN
1574
1575
1577
1578#endif
size_type n() const
size_type m() const
static constexpr unsigned int n_independent_components
static constexpr unsigned int n_independent_components
Definition tensor.h:497
virtual size_type size() const override
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorRowSize2(int arg1, int arg2, int arg3)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition exceptions.h:580
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorRowSize3(int arg1, int arg2, int arg3, int arg4)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:557
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorColSize3(int arg1, int arg2, int arg3, int arg4)
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorColSize2(int arg1, int arg2, int arg3)
#define AssertThrow(cond, exc)
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &Du)
Definition divergence.h:471
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
void to_tensor(const Vector< Number > &vec, Number &s)
FullMatrix< Number > to_matrix(const Number &s)
Vector< Number > to_vector(const Number &s)
static constexpr double SQRT2
Definition numbers.h:269
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)