OGS
KelvinVector.h
Go to the documentation of this file.
1 
10 #pragma once
11 
12 #include <Eigen/Dense>
13 #include "BaseLib/Error.h"
14 
15 namespace MathLib
16 {
20 namespace KelvinVector
21 {
23 constexpr int kelvin_vector_dimensions(int const displacement_dim)
24 {
25  if (displacement_dim == 2)
26  {
27  return 4;
28  }
29  else if (displacement_dim == 3)
30  {
31  return 6;
32  }
33  OGS_FATAL(
34  "Cannot convert displacement dimension {} to kelvin vector dimension.",
35  displacement_dim);
36 }
37 
38 //
39 // Kelvin vector and matrix templates for given displacement dimension.
40 //
41 
45 template <int DisplacementDim>
47  Eigen::Matrix<double, kelvin_vector_dimensions(DisplacementDim), 1,
48  Eigen::ColMajor>;
49 
53 template <int DisplacementDim>
55  Eigen::Matrix<double, kelvin_vector_dimensions(DisplacementDim),
56  kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor>;
57 
61 template <int KelvinVectorSize>
62 struct Invariants final
63 {
64  static_assert(KelvinVectorSize == 4 || KelvinVectorSize == 6,
65  "KelvinVector invariants for vectors of size different than "
66  "4 or 6 is not allowed.");
69  static Eigen::Matrix<double, KelvinVectorSize, KelvinVectorSize> const
73  static Eigen::Matrix<double, KelvinVectorSize, KelvinVectorSize> const
76  static Eigen::Matrix<double, KelvinVectorSize, 1> const identity2;
77 
79  static double determinant(
80  Eigen::Matrix<double, KelvinVectorSize, 1> const& v);
81 
84  static double equivalentStress(
85  Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
86 
88  static double FrobeniusNorm(
89  Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
90 
93  static double J2(
94  Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
95 
98  static double J3(
99  Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
100 
102  static double trace(Eigen::Matrix<double, KelvinVectorSize, 1> const& v);
103 
106  static Eigen::Vector3d diagonal(
107  Eigen::Matrix<double, KelvinVectorSize, 1> const& v);
108 };
109 
110 //
111 // Inverses of a Kelvin vector.
112 //
113 
116 template <int KelvinVectorSize>
117 Eigen::Matrix<double, KelvinVectorSize, 1, Eigen::ColMajor, KelvinVectorSize, 1>
118 inverse(Eigen::Matrix<double,
119  KelvinVectorSize,
120  1,
121  Eigen::ColMajor,
122  KelvinVectorSize,
123  1> const& v);
124 
127 template <int KelvinVectorSize>
128 Eigen::Matrix<double, 3, 3> kelvinVectorToTensor(Eigen::Matrix<double,
129  KelvinVectorSize,
130  1,
131  Eigen::ColMajor,
132  KelvinVectorSize,
133  1> const& v);
134 
137 template <int DisplacementDim>
139  Eigen::Matrix<double, 3, 3> const& m);
140 
152 template <int KelvinVectorSize>
153 Eigen::Matrix<double, KelvinVectorSize, 1, Eigen::ColMajor, KelvinVectorSize, 1>
154 kelvinVectorToSymmetricTensor(Eigen::Matrix<double,
155  KelvinVectorSize,
156  1,
157  Eigen::ColMajor,
158  KelvinVectorSize,
159  1> const& v);
160 
168 template <typename Derived>
169 Eigen::Matrix<double, Eigen::MatrixBase<Derived>::RowsAtCompileTime, 1>
170 symmetricTensorToKelvinVector(Eigen::MatrixBase<Derived> const& v)
171 {
172  static_assert(
173  (Eigen::MatrixBase<Derived>::ColsAtCompileTime == 1) ||
174  (Eigen::MatrixBase<Derived>::ColsAtCompileTime == Eigen::Dynamic),
175  "KelvinVector must be a column vector");
176  if (v.cols() != 1)
177  {
178  OGS_FATAL(
179  "KelvinVector must be a column vector, but input has {:d} columns.",
180  v.cols());
181  }
182 
183  Eigen::Matrix<double, Eigen::MatrixBase<Derived>::RowsAtCompileTime, 1>
184  result;
185  if (v.rows() == 4)
186  {
187  result.resize(4, 1);
188  result << v[0], v[1], v[2], v[3] * std::sqrt(2.);
189  }
190  else if (v.rows() == 6)
191  {
192  result.resize(6, 1);
193  result << v[0], v[1], v[2], v[3] * std::sqrt(2.), v[4] * std::sqrt(2.),
194  v[5] * std::sqrt(2.);
195  }
196  else
197  {
198  OGS_FATAL(
199  "Symmetric tensor to Kelvin vector conversion expected an input "
200  "vector of size 4 or 6, but a vector of size {:d} was given.",
201  v.size());
202  }
203  return result;
204 }
205 
210 template <int DisplacementDim>
212  std::vector<double> const& values)
213 {
214  constexpr int kelvin_vector_size =
215  kelvin_vector_dimensions(DisplacementDim);
216 
217  if (values.size() != kelvin_vector_size)
218  {
219  OGS_FATAL(
220  "Symmetric tensor to Kelvin vector conversion expected an input "
221  "vector of size {:d}, but a vector of size {:d} was given.",
222  kelvin_vector_size, values.size());
223  }
224 
226  Eigen::Map<typename MathLib::KelvinVector::KelvinVectorType<
227  DisplacementDim> const>(
228  values.data(), kelvin_vector_dimensions(DisplacementDim), 1));
229 }
230 
235 template <int DisplacementDim>
237  Eigen::Matrix<double, DisplacementDim, DisplacementDim, Eigen::ColMajor,
238  DisplacementDim, DisplacementDim> const& transformation);
239 
240 } // namespace KelvinVector
241 } // namespace MathLib
242 
243 #include "KelvinVector-impl.h"
#define OGS_FATAL(...)
Definition: Error.h:26
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:48
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:23
Eigen::Matrix< double, 3, 3 > kelvinVectorToTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
Definition: KelvinVector.h:56
Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > inverse(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
KelvinVectorType< DisplacementDim > tensorToKelvin(Eigen::Matrix< double, 3, 3 > const &m)
KelvinMatrixType< DisplacementDim > fourthOrderRotationMatrix(Eigen::Matrix< double, DisplacementDim, DisplacementDim, Eigen::ColMajor, DisplacementDim, DisplacementDim > const &transformation)
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
Definition: KelvinVector.h:170
static double FrobeniusNorm(Eigen::Matrix< double, KelvinVectorSize, 1 > const &deviatoric_v)
Get the norm of the deviatoric stress.
static Eigen::Matrix< double, KelvinVectorSize, 1 > const identity2
Kelvin mapping of 2nd order identity tensor.
Definition: KelvinVector.h:76
static Eigen::Vector3d diagonal(Eigen::Matrix< double, KelvinVectorSize, 1 > const &v)
static double trace(Eigen::Matrix< double, KelvinVectorSize, 1 > const &v)
Trace of the corresponding tensor.
static Eigen::Matrix< double, KelvinVectorSize, KelvinVectorSize > const spherical_projection
Definition: KelvinVector.h:74
static double J2(Eigen::Matrix< double, KelvinVectorSize, 1 > const &deviatoric_v)
static Eigen::Matrix< double, KelvinVectorSize, KelvinVectorSize > const deviatoric_projection
Definition: KelvinVector.h:66
static double determinant(Eigen::Matrix< double, KelvinVectorSize, 1 > const &v)
Determinant of a matrix in Kelvin vector representation.
static double equivalentStress(Eigen::Matrix< double, KelvinVectorSize, 1 > const &deviatoric_v)
static double J3(Eigen::Matrix< double, KelvinVectorSize, 1 > const &deviatoric_v)