Reference documentation for deal.II version 9.0.0
notation.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2018 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_physics_notation_h
17 #define dealii_physics_notation_h
18 
19 #include <deal.II/base/exceptions.h>
20 #include <deal.II/base/numbers.h>
21 #include <deal.II/base/tensor.h>
23 #include <deal.II/lac/full_matrix.h>
24 #include <deal.II/lac/vector.h>
25 
26 #include <type_traits>
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 
31 namespace Physics
32 {
33  namespace Notation
34  {
260  namespace Kelvin
261  {
262 
267  << "The number of rows in the input matrix is " << arg1
268  << ", but needs to be either " << arg2
269  << " or " << arg3 << ".");
270 
271 
276  << "The number of rows in the input matrix is " << arg1
277  << ", but needs to be either " << arg2 << "," << arg3
278  << ", or " << arg4 << ".");
279 
280 
285  << "The number of columns in the input matrix is " << arg1
286  << ", but needs to be either " << arg2
287  << " or " << arg3 << ".");
288 
289 
294  << "The number of columns in the input matrix is " << arg1
295  << ", but needs to be either " << arg2 << "," << arg3
296  << ", or " << arg4 << ".");
297 
298 
303 
309  template<typename Number>
310  Vector<Number>
311  to_vector (const Number &s);
312 
313 
319  template<int dim, typename Number>
320  Vector<Number>
321  to_vector (const Tensor<0,dim,Number> &s);
322 
323 
329  template<int dim, typename Number>
330  Vector<Number>
331  to_vector (const Tensor<1,dim,Number> &v);
332 
333 
339  template<int dim, typename Number>
340  Vector<Number>
341  to_vector (const Tensor<2,dim,Number> &t);
342 
343 
350  template<int dim, typename Number>
351  Vector<Number>
353 
354 
360  template<typename Number>
362  to_matrix (const Number &s);
363 
364 
370  template<int dim, typename Number>
372  to_matrix (const Tensor<0,dim,Number> &s);
373 
374 
380  template<int dim, typename Number>
382  to_matrix (const Tensor<1,dim,Number> &v);
383 
384 
390  template<int dim, typename Number>
392  to_matrix (const Tensor<2,dim,Number> &t);
393 
394 
404  template<int dim, typename Number>
407 
434  template<int dim,
435  typename SubTensor1 = Tensor<2,dim>,
436  typename SubTensor2 = Tensor<1,dim>,
437  typename Number>
439  to_matrix (const Tensor<3,dim,Number> &t);
440 
441 
448  template<int dim, typename Number>
450  to_matrix (const Tensor<4,dim,Number> &t);
451 
452 
459  template<int dim, typename Number>
462 
464 
469 
473  template<typename Number>
474  void
475  to_tensor (const Vector<Number> &vec,
476  Number &s);
477 
478 
482  template<int dim, typename Number>
483  void
484  to_tensor (const Vector<Number> &vec,
486 
487 
491  template<int dim, typename Number>
492  void
493  to_tensor (const Vector<Number> &vec,
495 
496 
500  template<int dim, typename Number>
501  void
502  to_tensor (const Vector<Number> &vec,
504 
505 
509  template<int dim, typename Number>
510  void
511  to_tensor (const Vector<Number> &vec,
513 
514 
518  template<typename Number>
519  void
520  to_tensor (const FullMatrix<Number> &mtrx,
521  Number &s);
522 
523 
527  template<int dim, typename Number>
528  void
529  to_tensor (const FullMatrix<Number> &mtrx,
531 
532 
536  template<int dim, typename Number>
537  void
538  to_tensor (const FullMatrix<Number> &mtrx,
540 
541 
545  template<int dim, typename Number>
546  void
547  to_tensor (const FullMatrix<Number> &mtrx,
549 
550 
554  template<int dim, typename Number>
555  void
556  to_tensor (const FullMatrix<Number> &mtrx,
558 
559 
569  template<int dim, typename Number>
570  void
571  to_tensor (const FullMatrix<Number> &mtrx,
573 
574 
578  template<int dim, typename Number>
579  void
580  to_tensor (const FullMatrix<Number> &mtrx,
582 
583 
587  template<int dim, typename Number>
588  void
589  to_tensor (const FullMatrix<Number> &mtrx,
591 
592 
597  template<typename TensorType, typename Number>
598  TensorType
599  to_tensor (const Vector<Number> &vec);
600 
601 
606  template<typename TensorType, typename Number>
607  TensorType
608  to_tensor (const FullMatrix<Number> &vec);
610 
611  }
612 
613  }
614 }
615 
616 
617 #ifndef DOXYGEN
618 
619 
620 // ------------------------- inline functions ------------------------
621 
622 
623 namespace Physics
624 {
625  namespace Notation
626  {
627  namespace Kelvin
628  {
629 
630  namespace internal
631  {
632 
640  template<int dim>
641  std::pair<unsigned int, unsigned int>
642  indices_from_component (const unsigned int component_n,
643  const bool symmetric);
644 
645 
646  template<int dim>
647  std::pair<unsigned int, unsigned int>
648  indices_from_component (const unsigned int component_n,
649  const bool)
650  {
651  AssertThrow(false, ExcNotImplemented());
652  return std::pair<unsigned int, unsigned int> ();
653  }
654 
655 
656  template<>
657  inline
658  std::pair<unsigned int, unsigned int>
659  indices_from_component<1> (const unsigned int component_n,
660  const bool)
661  {
662  Assert(component_n < 1, ExcIndexRange(component_n,0,1));
663 
664  return std::make_pair(0u,0u);
665  }
666 
667 
668  template<>
669  inline
670  std::pair<unsigned int, unsigned int>
671  indices_from_component<2> (const unsigned int component_n,
672  const bool symmetric)
673  {
674  if (symmetric == true)
675  {
678  }
679  else
680  {
683  }
684 
685  static const unsigned int indices[4][2] =
686  {
687  {0,0}, {1,1},
688  {0,1}, {1,0}
689  };
690  return std::make_pair(indices[component_n][0], indices[component_n][1]);
691  }
692 
693 
694  template<>
695  inline
696  std::pair<unsigned int, unsigned int>
697  indices_from_component<3> (const unsigned int component_n,
698  const bool symmetric)
699  {
700  if (symmetric == true)
701  {
704  }
705  else
706  {
709  }
710 
711  static const unsigned int indices[9][2] =
712  {
713  {0,0}, {1,1}, {2,2},
714  {1,2}, {0,2}, {0,1},
715  {1,0}, {2,0}, {2,1}
716  };
717  return std::make_pair(indices[component_n][0], indices[component_n][1]);
718  }
719 
720 
725  template<int dim>
726  double
727  vector_component_factor (const unsigned int component_i,
728  const bool symmetric)
729  {
730  if (symmetric == false)
731  return 1.0;
732 
733  if (component_i < dim)
734  return 1.0;
735  else
736  return numbers::SQRT2;
737  }
738 
739 
744  template<int dim>
745  double
746  matrix_component_factor (const unsigned int component_i,
747  const unsigned int component_j,
748  const bool symmetric)
749  {
750  if (symmetric == false)
751  return 1.0;
752 
753  // This case check returns equivalent of this result:
754  // internal::vector_component_factor<dim>(component_i,symmetric)*internal::vector_component_factor<dim>(component_j,symmetric);
755  if (component_i < dim && component_j < dim)
756  return 1.0;
757  else if (component_i >= dim && component_j >= dim)
758  return 2.0;
759  else // ((component_i >= dim && component_j < dim) || (component_i < dim && component_j >= dim))
760  return numbers::SQRT2;
761  }
762 
763  }
764 
765 
766  template<typename Number>
767  Vector<Number>
768  to_vector (const Number &s)
769  {
770  Vector<Number> out (1);
771  const unsigned int n_rows = out.size();
772  for (unsigned int r=0; r<n_rows; ++r)
773  out(r) = s;
774  return out;
775  }
776 
777 
778  template<int dim, typename Number>
779  Vector<Number>
781  {
782  return to_vector(s.operator const Number &());
783  }
784 
785 
786  template<int dim, typename Number>
787  Vector<Number>
789  {
790  Vector<Number> out (v.n_independent_components);
791  const unsigned int n_rows = out.size();
792  for (unsigned int r=0; r<n_rows; ++r)
793  {
794  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,false);
795  Assert(indices.first < dim, ExcInternalError());
796  const unsigned int i = indices.first;
797  out(r) = v[i];
798  }
799  return out;
800  }
801 
802 
803  template<int dim, typename Number>
804  Vector<Number>
806  {
807  Vector<Number> out (t.n_independent_components);
808  const unsigned int n_rows = out.size();
809  for (unsigned int r=0; r<n_rows; ++r)
810  {
811  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,false);
812  Assert(indices.first < dim, ExcInternalError());
813  Assert(indices.second < dim, ExcInternalError());
814  const unsigned int i = indices.first;
815  const unsigned int j = indices.second;
816  out(r) = t[i][j];
817  }
818  return out;
819  }
820 
821 
822  template<int dim, typename Number>
823  Vector<Number>
825  {
826  Vector<Number> out (st.n_independent_components);
827  const unsigned int n_rows = out.size();
828  for (unsigned int r=0; r<n_rows; ++r)
829  {
830  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,true);
831  Assert(indices.first < dim, ExcInternalError());
832  Assert(indices.second < dim, ExcInternalError());
833  Assert(indices.second >= indices.first, ExcInternalError());
834  const unsigned int i = indices.first;
835  const unsigned int j = indices.second;
836 
837  const double factor = internal::vector_component_factor<dim>(r,true);
838 
839  out(r) = factor*st[i][j];
840  }
841  return out;
842  }
843 
844 
845  template<typename Number>
847  to_matrix (const Number &s)
848  {
849  FullMatrix<Number> out (1,1);
850  out(0,0) = s;
851  return out;
852  }
853 
854 
855  template<int dim, typename Number>
858  {
859  return to_matrix(s.operator const Number &());
860  }
861 
862 
863  template<int dim, typename Number>
866  {
868  const unsigned int n_rows = out.m();
869  const unsigned int n_cols = out.n();
870  for (unsigned int r=0; r<n_rows; ++r)
871  {
872  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,false);
873  Assert(indices.first < dim, ExcInternalError());
874  const unsigned int i = indices.first;
875 
876  for (unsigned int c=0; c<n_cols; ++c)
877  {
878  Assert(c < 1, ExcInternalError());
879  out(r,c) = v[i];
880  }
881  }
882  return out;
883  }
884 
885 
886  template<int dim, typename Number>
889  {
890  FullMatrix<Number> out (dim,dim);
891  const unsigned int n_rows = out.m();
892  const unsigned int n_cols = out.n();
893  for (unsigned int r=0; r<n_rows; ++r)
894  {
895  const std::pair<unsigned int,unsigned int> indices_i = internal::indices_from_component<dim>(r,false);
896  Assert(indices_i.first < dim, ExcInternalError());
897  Assert(indices_i.second < dim, ExcInternalError());
898  const unsigned int i = indices_i.first;
899 
900  for (unsigned int c=0; c<n_cols; ++c)
901  {
902  const std::pair<unsigned int,unsigned int> indices_j = internal::indices_from_component<dim>(c,false);
903  Assert(indices_j.first < dim, ExcInternalError());
904  Assert(indices_j.second < dim, ExcInternalError());
905  const unsigned int j = indices_j.second;
906 
907  out(r,c) = t[i][j];
908  }
909  }
910  return out;
911  }
912 
913 
914  template<int dim, typename Number>
917  {
918  return to_matrix(Tensor<2,dim,Number>(st));
919  }
920 
921 
922  namespace internal
923  {
924  namespace
925  {
926  template<typename TensorType>
927  struct is_rank_2_symmetric_tensor : std::false_type
928  {};
929 
930  template<int dim, typename Number>
931  struct is_rank_2_symmetric_tensor<SymmetricTensor<2,dim,Number> >
932  : std::true_type
933  {};
934  }
935  }
936 
937 
938  template<int dim, typename SubTensor1, typename SubTensor2, typename Number>
941  {
942  static_assert(
943  (SubTensor1::dimension == dim && SubTensor2::dimension == dim),
944  "Sub-tensor spatial dimension is different from those of the input tensor.");
945 
946  static_assert(
947  (SubTensor1::rank == 2 && SubTensor2::rank == 1) ||
948  (SubTensor1::rank == 1 && SubTensor2::rank == 2),
949  "Cannot build a rank 3 tensor from the given combination of sub-tensors.");
950 
951  FullMatrix<Number> out (SubTensor1::n_independent_components,SubTensor2::n_independent_components);
952  const unsigned int n_rows = out.m();
953  const unsigned int n_cols = out.n();
954 
955  if (SubTensor1::rank == 2 && SubTensor2::rank == 1)
956  {
957  const bool subtensor_is_rank_2_symmetric_tensor = internal::is_rank_2_symmetric_tensor<SubTensor1>::value;
958 
959  for (unsigned int r=0; r<n_rows; ++r)
960  {
961  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(r,subtensor_is_rank_2_symmetric_tensor);
962  Assert(indices_ij.first < dim, ExcInternalError());
963  Assert(indices_ij.second < dim, ExcInternalError());
964  if (subtensor_is_rank_2_symmetric_tensor)
965  {
966  Assert(indices_ij.second >= indices_ij.first, ExcInternalError());
967  }
968  const unsigned int i = indices_ij.first;
969  const unsigned int j = indices_ij.second;
970 
971  const double factor = internal::vector_component_factor<dim>(r,subtensor_is_rank_2_symmetric_tensor);
972 
973  for (unsigned int c=0; c<n_cols; ++c)
974  {
975  const std::pair<unsigned int,unsigned int> indices_k = internal::indices_from_component<dim>(c,false);
976  Assert(indices_k.first < dim, ExcInternalError());
977  const unsigned int k = indices_k.first;
978 
979  if (subtensor_is_rank_2_symmetric_tensor)
980  out(r,c) = factor*t[i][j][k];
981  else
982  out(r,c) = t[i][j][k];
983  }
984  }
985  }
986  else if (SubTensor1::rank == 1 && SubTensor2::rank == 2)
987  {
988  const bool subtensor_is_rank_2_symmetric_tensor = internal::is_rank_2_symmetric_tensor<SubTensor2>::value;
989 
990  for (unsigned int r=0; r<n_rows; ++r)
991  {
992  const std::pair<unsigned int,unsigned int> indices_k = internal::indices_from_component<dim>(r,false);
993  Assert(indices_k.first < dim, ExcInternalError());
994  const unsigned int k = indices_k.first;
995 
996  for (unsigned int c=0; c<n_cols; ++c)
997  {
998  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(c,subtensor_is_rank_2_symmetric_tensor);
999  Assert(indices_ij.first < dim, ExcInternalError());
1000  Assert(indices_ij.second < dim, ExcInternalError());
1001  if (subtensor_is_rank_2_symmetric_tensor)
1002  {
1003  Assert(indices_ij.second >= indices_ij.first, ExcInternalError());
1004  }
1005  const unsigned int i = indices_ij.first;
1006  const unsigned int j = indices_ij.second;
1007 
1008  if (subtensor_is_rank_2_symmetric_tensor)
1009  {
1010  const double factor = internal::vector_component_factor<dim>(c,subtensor_is_rank_2_symmetric_tensor);
1011  out(r,c) = factor*t[k][i][j];
1012  }
1013  else
1014  out(r,c) = t[k][i][j];
1015  }
1016  }
1017  }
1018  else
1019  {
1020  AssertThrow(false, ExcNotImplemented());
1021  }
1022 
1023  return out;
1024  }
1025 
1026 
1027  template<int dim, typename Number>
1029  to_matrix (const Tensor<4,dim,Number> &t)
1030  {
1033  const unsigned int n_rows = out.m();
1034  const unsigned int n_cols = out.n();
1035  for (unsigned int r=0; r<n_rows; ++r)
1036  {
1037  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(r,false);
1038  Assert(indices_ij.first < dim, ExcInternalError());
1039  Assert(indices_ij.second < dim, ExcInternalError());
1040  const unsigned int i = indices_ij.first;
1041  const unsigned int j = indices_ij.second;
1042 
1043  for (unsigned int c=0; c<n_cols; ++c)
1044  {
1045  const std::pair<unsigned int,unsigned int> indices_kl = internal::indices_from_component<dim>(c,false);
1046  Assert(indices_kl.first < dim, ExcInternalError());
1047  Assert(indices_kl.second < dim, ExcInternalError());
1048  const unsigned int k = indices_kl.first;
1049  const unsigned int l = indices_kl.second;
1050 
1051  out(r,c) = t[i][j][k][l];
1052  }
1053  }
1054  return out;
1055  }
1056 
1057 
1058  template<int dim, typename Number>
1061  {
1064  const unsigned int n_rows = out.m();
1065  const unsigned int n_cols = out.n();
1066  for (unsigned int r=0; r<n_rows; ++r)
1067  {
1068  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(r,true);
1069  Assert(indices_ij.first < dim, ExcInternalError());
1070  Assert(indices_ij.second < dim, ExcInternalError());
1071  Assert(indices_ij.second >= indices_ij.first, ExcInternalError());
1072  const unsigned int i = indices_ij.first;
1073  const unsigned int j = indices_ij.second;
1074 
1075  for (unsigned int c=0; c<n_cols; ++c)
1076  {
1077  const std::pair<unsigned int,unsigned int> indices_kl = internal::indices_from_component<dim>(c,true);
1078  Assert(indices_kl.first < dim, ExcInternalError());
1079  Assert(indices_kl.second < dim, ExcInternalError());
1080  Assert(indices_kl.second >= indices_kl.first, ExcInternalError());
1081  const unsigned int k = indices_kl.first;
1082  const unsigned int l = indices_kl.second;
1083 
1084  const double factor = internal::matrix_component_factor<dim>(r,c,true);
1085 
1086  out(r,c) = factor*st[i][j][k][l];
1087  }
1088  }
1089  return out;
1090  }
1091 
1092 
1093  template<typename Number>
1094  void
1095  to_tensor (const Vector<Number> &vec,
1096  Number &s)
1097  {
1098  Assert(vec.size() == 1, ExcDimensionMismatch(vec.size(), 1));
1099  s = vec(0);
1100  }
1101 
1102 
1103  template<int dim, typename Number>
1104  void
1105  to_tensor (const Vector<Number> &vec,
1107  {
1108  return to_tensor(vec, s.operator Number &());
1109  }
1110 
1111 
1112  template<int dim, typename Number>
1113  void
1114  to_tensor (const Vector<Number> &vec,
1116  {
1117  Assert(vec.size() == v.n_independent_components,
1119  const unsigned int n_rows = vec.size();
1120  for (unsigned int r=0; r<n_rows; ++r)
1121  {
1122  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,false);
1123  Assert(indices.first < dim, ExcInternalError());
1124  const unsigned int i = indices.first;
1125  v[i] = vec(r);
1126  }
1127  }
1128 
1129 
1130  template<int dim, typename Number>
1131  void
1132  to_tensor (const Vector<Number> &vec,
1134  {
1135  Assert(vec.size() == t.n_independent_components,
1137  const unsigned int n_rows = vec.size();
1138  for (unsigned int r=0; r<n_rows; ++r)
1139  {
1140  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,false);
1141  Assert(indices.first < dim, ExcInternalError());
1142  Assert(indices.second < dim, ExcInternalError());
1143  const unsigned int i = indices.first;
1144  const unsigned int j = indices.second;
1145  t[i][j] = vec(r);
1146  }
1147  }
1148 
1149 
1150  template<int dim, typename Number>
1151  void
1152  to_tensor (const Vector<Number> &vec,
1154  {
1155  Assert(vec.size() == st.n_independent_components,
1157  const unsigned int n_rows = vec.size();
1158  for (unsigned int r=0; r<n_rows; ++r)
1159  {
1160  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,true);
1161  Assert(indices.first < dim, ExcInternalError());
1162  Assert(indices.second < dim, ExcInternalError());
1163  Assert(indices.second >= indices.first, ExcInternalError());
1164  const unsigned int i = indices.first;
1165  const unsigned int j = indices.second;
1166 
1167  const double inv_factor = 1.0/internal::vector_component_factor<dim>(r,true);
1168 
1169  st[i][j] = inv_factor*vec(r);
1170  }
1171  }
1172 
1173 
1174  template<typename Number>
1175  void
1176  to_tensor (const FullMatrix<Number> &mtrx,
1177  Number &s)
1178  {
1179  Assert(mtrx.m() == 1, ExcDimensionMismatch(mtrx.m(), 1));
1180  Assert(mtrx.n() == 1, ExcDimensionMismatch(mtrx.n(), 1));
1181  Assert(mtrx.n_elements() == 1, ExcDimensionMismatch(mtrx.n_elements(), 1));
1182  s = mtrx(0,0);
1183  }
1184 
1185 
1186  template<int dim, typename Number>
1187  void
1188  to_tensor (const FullMatrix<Number> &mtrx,
1190  {
1191  return to_tensor(mtrx, s.operator Number &());
1192  }
1193 
1194 
1195  template<int dim, typename Number>
1196  void
1197  to_tensor (const FullMatrix<Number> &mtrx,
1199  {
1200  Assert(mtrx.m() == dim, ExcDimensionMismatch(mtrx.m(), dim));
1201  Assert(mtrx.n() == 1, ExcDimensionMismatch(mtrx.n(), 1));
1204 
1205  const unsigned int n_rows = mtrx.m();
1206  const unsigned int n_cols = mtrx.n();
1207  for (unsigned int r=0; r<n_rows; ++r)
1208  {
1209  const std::pair<unsigned int,unsigned int> indices = internal::indices_from_component<dim>(r,false);
1210  Assert(indices.first < dim, ExcInternalError());
1211  Assert(indices.second == 0, ExcInternalError());
1212  const unsigned int i = indices.first;
1213 
1214  for (unsigned int c=0; c<n_cols; ++c)
1215  {
1216  Assert(c < 1, ExcInternalError());
1217  v[i] = mtrx(r,c);
1218  }
1219  }
1220  }
1221 
1222 
1223  template<int dim, typename Number>
1224  void
1225  to_tensor (const FullMatrix<Number> &mtrx,
1227  {
1228  Assert(mtrx.m() == dim, ExcDimensionMismatch(mtrx.m(), dim));
1229  Assert(mtrx.n() == dim, ExcDimensionMismatch(mtrx.n(), dim));
1232 
1233  const unsigned int n_rows = mtrx.m();
1234  const unsigned int n_cols = mtrx.n();
1235  for (unsigned int r=0; r<n_rows; ++r)
1236  {
1237  const std::pair<unsigned int,unsigned int> indices_i = internal::indices_from_component<dim>(r,false);
1238  Assert(indices_i.first < dim, ExcInternalError());
1239  Assert(indices_i.second < dim, ExcInternalError());
1240  const unsigned int i = indices_i.first;
1241 
1242  for (unsigned int c=0; c<n_cols; ++c)
1243  {
1244  const std::pair<unsigned int,unsigned int> indices_j = internal::indices_from_component<dim>(c,false);
1245  Assert(indices_j.first < dim, ExcInternalError());
1246  Assert(indices_j.second < dim, ExcInternalError());
1247  const unsigned int j = indices_j.second;
1248 
1249  t[i][j] = mtrx(r,c);
1250  }
1251  }
1252 
1253  }
1254 
1255 
1256  template<int dim, typename Number>
1257  void
1258  to_tensor (const FullMatrix<Number> &mtrx,
1260  {
1261  // Its impossible to fit the (dim^2 + dim)/2 entries into a square matrix
1262  // We therefore assume that its been converted to a standard tensor format
1263  // using to_matrix (SymmetricTensor<2,dim,Number>) at some point...
1264  Assert(mtrx.m() == dim, ExcDimensionMismatch(mtrx.m(), dim));
1265  Assert(mtrx.n() == dim, ExcDimensionMismatch(mtrx.n(), dim));
1268 
1270  to_tensor(mtrx, tmp);
1271  st = symmetrize(tmp);
1272  Assert((Tensor<2,dim,Number>(st) - tmp).norm() < 1e-12,
1273  ExcMessage("The entries stored inside the matrix were not symmetric"));
1274  }
1275 
1276 
1277  template<int dim, typename Number>
1278  void
1279  to_tensor (const FullMatrix<Number> &mtrx,
1281  {
1286  mtrx.m(),
1294  mtrx.n(),
1298 
1299  const unsigned int n_rows = mtrx.m();
1300  const unsigned int n_cols = mtrx.n();
1302  {
1306  mtrx.m(),
1309  ));
1310 
1311  const bool subtensor_is_rank_2_symmetric_tensor
1313 
1314  for (unsigned int r=0; r<n_rows; ++r)
1315  {
1316  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(r,subtensor_is_rank_2_symmetric_tensor);
1317  Assert(indices_ij.first < dim, ExcInternalError());
1318  Assert(indices_ij.second < dim, ExcInternalError());
1319  if (subtensor_is_rank_2_symmetric_tensor)
1320  {
1321  Assert(indices_ij.second >= indices_ij.first, ExcInternalError());
1322  }
1323  const unsigned int i = indices_ij.first;
1324  const unsigned int j = indices_ij.second;
1325 
1326  const double inv_factor = 1.0/internal::vector_component_factor<dim>(r,subtensor_is_rank_2_symmetric_tensor);
1327 
1328  for (unsigned int c=0; c<n_cols; ++c)
1329  {
1330  const std::pair<unsigned int,unsigned int> indices_k = internal::indices_from_component<dim>(c,false);
1331  Assert(indices_k.first < dim, ExcInternalError());
1332  const unsigned int k = indices_k.first;
1333 
1334  if (subtensor_is_rank_2_symmetric_tensor)
1335  {
1336  t[i][j][k] = inv_factor*mtrx(r,c);
1337  t[j][i][k] = t[i][j][k];
1338  }
1339  else
1340  t[i][j][k] = mtrx(r,c);
1341  }
1342  }
1343  }
1344  else
1345  {
1351  mtrx.n(),
1354  ));
1355 
1356  const bool subtensor_is_rank_2_symmetric_tensor
1358 
1359  for (unsigned int r=0; r<n_rows; ++r)
1360  {
1361  const std::pair<unsigned int,unsigned int> indices_k = internal::indices_from_component<dim>(r,false);
1362  Assert(indices_k.first < dim, ExcInternalError());
1363  const unsigned int k = indices_k.first;
1364 
1365  for (unsigned int c=0; c<n_cols; ++c)
1366  {
1367  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(c,subtensor_is_rank_2_symmetric_tensor);
1368  Assert(indices_ij.first < dim, ExcInternalError());
1369  Assert(indices_ij.second < dim, ExcInternalError());
1370  if (subtensor_is_rank_2_symmetric_tensor)
1371  {
1372  Assert(indices_ij.second >= indices_ij.first, ExcInternalError());
1373  }
1374  const unsigned int i = indices_ij.first;
1375  const unsigned int j = indices_ij.second;
1376 
1377  if (subtensor_is_rank_2_symmetric_tensor)
1378  {
1379  const double inv_factor = 1.0/internal::vector_component_factor<dim>(c,subtensor_is_rank_2_symmetric_tensor);
1380  t[k][i][j] = inv_factor*mtrx(r,c);
1381  t[k][j][i] = t[k][i][j];
1382  }
1383  else
1384  t[k][i][j] = mtrx(r,c);
1385  }
1386  }
1387  }
1388  }
1389 
1390 
1391  template<int dim, typename Number>
1392  void
1393  to_tensor (const FullMatrix<Number> &mtrx,
1395  {
1402 
1403  const unsigned int n_rows = mtrx.m();
1404  const unsigned int n_cols = mtrx.n();
1405  for (unsigned int r=0; r<n_rows; ++r)
1406  {
1407  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(r,false);
1408  Assert(indices_ij.first < dim, ExcInternalError());
1409  Assert(indices_ij.second < dim, ExcInternalError());
1410  const unsigned int i = indices_ij.first;
1411  const unsigned int j = indices_ij.second;
1412 
1413  for (unsigned int c=0; c<n_cols; ++c)
1414  {
1415  const std::pair<unsigned int,unsigned int> indices_kl = internal::indices_from_component<dim>(c,false);
1416  Assert(indices_kl.first < dim, ExcInternalError());
1417  Assert(indices_kl.second < dim, ExcInternalError());
1418  const unsigned int k = indices_kl.first;
1419  const unsigned int l = indices_kl.second;
1420 
1421  t[i][j][k][l] = mtrx(r,c);
1422  }
1423  }
1424  }
1425 
1426 
1427  template<int dim, typename Number>
1428  void
1429  to_tensor (const FullMatrix<Number> &mtrx,
1431  {
1438 
1439  const unsigned int n_rows = mtrx.m();
1440  const unsigned int n_cols = mtrx.n();
1441  for (unsigned int r=0; r<n_rows; ++r)
1442  {
1443  const std::pair<unsigned int,unsigned int> indices_ij = internal::indices_from_component<dim>(r,false);
1444  Assert(indices_ij.first < dim, ExcInternalError());
1445  Assert(indices_ij.second < dim, ExcInternalError());
1446  const unsigned int i = indices_ij.first;
1447  const unsigned int j = indices_ij.second;
1448 
1449  for (unsigned int c=0; c<n_cols; ++c)
1450  {
1451  const std::pair<unsigned int,unsigned int> indices_kl = internal::indices_from_component<dim>(c,false);
1452  Assert(indices_kl.first < dim, ExcInternalError());
1453  Assert(indices_kl.second < dim, ExcInternalError());
1454  const unsigned int k = indices_kl.first;
1455  const unsigned int l = indices_kl.second;
1456 
1457  const double inv_factor = 1.0/internal::matrix_component_factor<dim>(r,c,true);
1458 
1459  st[i][j][k][l] = inv_factor*mtrx(r,c);
1460  }
1461  }
1462  }
1463 
1464 
1465  template<typename TensorType, typename Number>
1466  inline TensorType
1467  to_tensor (const Vector<Number> &vec)
1468  {
1469  TensorType out;
1470  to_tensor(vec, out);
1471  return out;
1472  }
1473 
1474 
1475  template<typename TensorType, typename Number>
1476  inline TensorType
1477  to_tensor (const FullMatrix<Number> &mtrx)
1478  {
1479  TensorType out;
1480  to_tensor(mtrx, out);
1481  return out;
1482  }
1483 
1484  }
1485  }
1486 }
1487 
1488 
1489 #endif // DOXYGEN
1490 
1491 
1492 DEAL_II_NAMESPACE_CLOSE
1493 
1494 #endif
Matrix is symmetric.
size_type m() const
static const double SQRT2
Definition: numbers.h:142
void to_tensor(const Vector< Number > &vec, Number &s)
static const unsigned int n_independent_components
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorColSize2(int arg1, int arg2, int arg3)
static const unsigned int n_independent_components
Definition: tensor.h:400
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorRowSize3(int arg1, int arg2, int arg3, int arg4)
Vector< Number > to_vector(const Number &s)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1221
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
size_type n_elements() const
static ::ExceptionBase & ExcMessage(std::string arg1)
size_type n() const
#define Assert(cond, exc)
Definition: exceptions.h:1142
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
FullMatrix< Number > to_matrix(const Number &s)
SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
double norm(const FEValuesBase< dim > &fe, const VectorSlice< const std::vector< std::vector< Tensor< 1, dim > > > > &Du)
Definition: divergence.h:534
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition: exceptions.h:382
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorRowSize2(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcNotImplemented()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:370
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotationExcFullMatrixToTensorColSize3(int arg1, int arg2, int arg3, int arg4)