OGS
KelvinVector.cpp
Go to the documentation of this file.
1 
11 #include "KelvinVector.h"
12 
13 #include "BaseLib/Error.h"
14 
15 namespace MathLib
16 {
17 namespace KelvinVector
18 {
19 template <>
20 double Invariants<6>::determinant(Eigen::Matrix<double, 6, 1> const& v)
21 {
22  return v(0) * v(1) * v(2) + v(3) * v(4) * v(5) / std::sqrt(2.) -
23  v(3) * v(3) * v(2) / 2. - v(4) * v(4) * v(0) / 2. -
24  v(5) * v(5) * v(1) / 2.;
25 }
26 
27 template <>
28 double Invariants<4>::determinant(Eigen::Matrix<double, 4, 1> const& v)
29 {
30  return v(2) * (v(0) * v(1) - v(3) * v(3) / 2.);
31 }
32 
33 template <>
34 Eigen::Matrix<double, 4, 1, Eigen::ColMajor, 4, 1> inverse(
35  Eigen::Matrix<double, 4, 1, Eigen::ColMajor, 4, 1> const& v)
36 {
37  assert(Invariants<4>::determinant(v) != 0);
38 
39  Eigen::Matrix<double, 4, 1, Eigen::ColMajor, 4, 1> inv;
40  inv(0) = v(1) * v(2);
41  inv(1) = v(0) * v(2);
42  inv(2) = v(0) * v(1) - v(3) * v(3) / 2.;
43  inv(3) = -v(3) * v(2);
44  return inv / Invariants<4>::determinant(v);
45 }
46 
47 template <>
48 Eigen::Matrix<double, 6, 1, Eigen::ColMajor, 6, 1> inverse(
49  Eigen::Matrix<double, 6, 1, Eigen::ColMajor, 6, 1> const& v)
50 {
51  assert(Invariants<6>::determinant(v) != 0);
52 
53  Eigen::Matrix<double, 6, 1, Eigen::ColMajor, 6, 1> inv;
54  inv(0) = v(1) * v(2) - v(4) * v(4) / 2.;
55  inv(1) = v(0) * v(2) - v(5) * v(5) / 2.;
56  inv(2) = v(0) * v(1) - v(3) * v(3) / 2.;
57  inv(3) = v(4) * v(5) / std::sqrt(2.) - v(3) * v(2);
58  inv(4) = v(3) * v(5) / std::sqrt(2.) - v(4) * v(0);
59  inv(5) = v(4) * v(3) / std::sqrt(2.) - v(1) * v(5);
60  return inv / Invariants<6>::determinant(v);
61 }
62 
63 template <>
64 Eigen::Matrix<double, 3, 3> kelvinVectorToTensor(
65  Eigen::Matrix<double, 4, 1, Eigen::ColMajor, 4, 1> const& v)
66 {
67  Eigen::Matrix<double, 3, 3> m;
68  m << v[0], v[3] / std::sqrt(2.), 0, v[3] / std::sqrt(2.), v[1], 0, 0, 0,
69  v[2];
70  return m;
71 }
72 
73 template <>
74 Eigen::Matrix<double, 3, 3> kelvinVectorToTensor(
75  Eigen::Matrix<double, 6, 1, Eigen::ColMajor, 6, 1> const& v)
76 {
77  Eigen::Matrix<double, 3, 3> m;
78  m << v[0], v[3] / std::sqrt(2.), v[5] / std::sqrt(2.), v[3] / std::sqrt(2.),
79  v[1], v[4] / std::sqrt(2.), v[5] / std::sqrt(2.), v[4] / std::sqrt(2.),
80  v[2];
81  return m;
82 }
83 
84 template <>
85 Eigen::Matrix<double, 3, 3> kelvinVectorToTensor(Eigen::Matrix<double,
86  Eigen::Dynamic,
87  1,
88  Eigen::ColMajor,
89  Eigen::Dynamic,
90  1> const& v)
91 {
92  if (v.size() == 4)
93  {
94  Eigen::Matrix<double, 4, 1, Eigen::ColMajor, 4, 1> v4;
95  v4 << v[0], v[1], v[2], v[3];
96  return kelvinVectorToTensor(v4);
97  }
98  if (v.size() == 6)
99  {
100  Eigen::Matrix<double, 6, 1, Eigen::ColMajor, 6, 1> v6;
101  v6 << v[0], v[1], v[2], v[3], v[4], v[5];
102  return kelvinVectorToTensor(v6);
103  }
104  OGS_FATAL(
105  "Conversion of dynamic Kelvin vector of size {:d} to a tensor is not "
106  "possible. Kelvin vector must be of size 4 or 6.",
107  v.size());
108 }
109 
110 template <>
111 KelvinVectorType<2> tensorToKelvin<2>(Eigen::Matrix<double, 3, 3> const& m)
112 {
113  assert(std::abs(m(0, 1) - m(1, 0)) <
114  std::numeric_limits<double>::epsilon());
115  assert(m(0, 2) == 0);
116  assert(m(1, 2) == 0);
117  assert(m(2, 0) == 0);
118  assert(m(2, 1) == 0);
119 
120  KelvinVectorType<2> v;
121  v << m(0, 0), m(1, 1), m(2, 2), m(0, 1) * std::sqrt(2.);
122  return v;
123 }
124 
125 template <>
126 KelvinVectorType<3> tensorToKelvin<3>(Eigen::Matrix<double, 3, 3> const& m)
127 {
128  assert(std::abs(m(0, 1) - m(1, 0)) <
129  std::numeric_limits<double>::epsilon());
130  assert(std::abs(m(1, 2) - m(2, 1)) <
131  std::numeric_limits<double>::epsilon());
132  assert(std::abs(m(0, 2) - m(2, 0)) <
133  std::numeric_limits<double>::epsilon());
134 
135  KelvinVectorType<3> v;
136  v << m(0, 0), m(1, 1), m(2, 2), m(0, 1) * std::sqrt(2.),
137  m(1, 2) * std::sqrt(2.), m(0, 2) * std::sqrt(2.);
138  return v;
139 }
140 
141 template <>
142 Eigen::Matrix<double, 4, 1> kelvinVectorToSymmetricTensor(
143  Eigen::Matrix<double, 4, 1, Eigen::ColMajor, 4, 1> const& v)
144 {
145  Eigen::Matrix<double, 4, 1> m;
146  m << v[0], v[1], v[2], v[3] / std::sqrt(2.);
147  return m;
148 }
149 
150 template <>
151 Eigen::Matrix<double, 6, 1> kelvinVectorToSymmetricTensor(
152  Eigen::Matrix<double, 6, 1, Eigen::ColMajor, 6, 1> const& v)
153 {
154  Eigen::Matrix<double, 6, 1> m;
155  m << v[0], v[1], v[2], v[3] / std::sqrt(2.), v[4] / std::sqrt(2.),
156  v[5] / std::sqrt(2.);
157  return m;
158 }
159 
160 template <>
161 Eigen::Matrix<double, Eigen::Dynamic, 1, Eigen::ColMajor, Eigen::Dynamic, 1>
162 kelvinVectorToSymmetricTensor(Eigen::Matrix<double,
163  Eigen::Dynamic,
164  1,
165  Eigen::ColMajor,
166  Eigen::Dynamic,
167  1> const& v)
168 {
169  if (v.size() == 4)
170  {
171  return kelvinVectorToSymmetricTensor<4>(v);
172  }
173  if (v.size() == 6)
174  {
175  return kelvinVectorToSymmetricTensor<6>(v);
176  }
177  OGS_FATAL(
178  "Kelvin vector to tensor conversion expected an input vector of size 4 "
179  "or 6, but a vector of size {:d} was given.",
180  v.size());
181 }
182 
183 template <>
184 KelvinMatrixType<2> fourthOrderRotationMatrix<2>(
185  Eigen::Matrix<double, 2, 2, Eigen::ColMajor, 2, 2> const& transformation)
186 {
187  // 1-based index access for convenience.
188  auto Q = [&](int const i, int const j)
189  { return transformation(i - 1, j - 1); };
190 
191  MathLib::KelvinVector::KelvinMatrixType<2> R;
192  R << Q(1, 1) * Q(1, 1), Q(1, 2) * Q(1, 2), 0,
193  std::sqrt(2) * Q(1, 1) * Q(1, 2), Q(2, 1) * Q(2, 1), Q(2, 2) * Q(2, 2),
194  0, std::sqrt(2) * Q(2, 1) * Q(2, 2), 0, 0, 1, 0,
195  std::sqrt(2) * Q(1, 1) * Q(2, 1), std::sqrt(2) * Q(1, 2) * Q(2, 2), 0,
196  Q(1, 1) * Q(2, 2) + Q(1, 2) * Q(2, 1);
197  return R;
198 }
199 
200 template <>
201 KelvinMatrixType<3> fourthOrderRotationMatrix<3>(
202  Eigen::Matrix<double, 3, 3, Eigen::ColMajor, 3, 3> const& transformation)
203 {
204  // 1-based index access for convenience.
205  auto Q = [&](int const i, int const j)
206  { return transformation(i - 1, j - 1); };
207 
208  MathLib::KelvinVector::KelvinMatrixType<3> R;
209  R << Q(1, 1) * Q(1, 1), Q(1, 2) * Q(1, 2), Q(1, 3) * Q(1, 3),
210  std::sqrt(2) * Q(1, 1) * Q(1, 2), std::sqrt(2) * Q(1, 2) * Q(1, 3),
211  std::sqrt(2) * Q(1, 1) * Q(1, 3), Q(2, 1) * Q(2, 1), Q(2, 2) * Q(2, 2),
212  Q(2, 3) * Q(2, 3), std::sqrt(2) * Q(2, 1) * Q(2, 2),
213  std::sqrt(2) * Q(2, 2) * Q(2, 3), std::sqrt(2) * Q(2, 1) * Q(2, 3),
214  Q(3, 1) * Q(3, 1), Q(3, 2) * Q(3, 2), Q(3, 3) * Q(3, 3),
215  std::sqrt(2) * Q(3, 1) * Q(3, 2), std::sqrt(2) * Q(3, 2) * Q(3, 3),
216  std::sqrt(2) * Q(3, 1) * Q(3, 3), std::sqrt(2) * Q(1, 1) * Q(2, 1),
217  std::sqrt(2) * Q(1, 2) * Q(2, 2), std::sqrt(2) * Q(1, 3) * Q(2, 3),
218  Q(1, 1) * Q(2, 2) + Q(1, 2) * Q(2, 1),
219  Q(1, 2) * Q(2, 3) + Q(1, 3) * Q(2, 2),
220  Q(1, 1) * Q(2, 3) + Q(1, 3) * Q(2, 1), std::sqrt(2) * Q(2, 1) * Q(3, 1),
221  std::sqrt(2) * Q(2, 2) * Q(3, 2), std::sqrt(2) * Q(2, 3) * Q(3, 3),
222  Q(2, 1) * Q(3, 2) + Q(2, 2) * Q(3, 1),
223  Q(2, 2) * Q(3, 3) + Q(2, 3) * Q(3, 2),
224  Q(2, 1) * Q(3, 3) + Q(2, 3) * Q(3, 1), std::sqrt(2) * Q(1, 1) * Q(3, 1),
225  std::sqrt(2) * Q(1, 2) * Q(3, 2), std::sqrt(2) * Q(1, 3) * Q(3, 3),
226  Q(1, 1) * Q(3, 2) + Q(1, 2) * Q(3, 1),
227  Q(1, 2) * Q(3, 3) + Q(1, 3) * Q(3, 2),
228  Q(1, 1) * Q(3, 3) + Q(1, 3) * Q(3, 1);
229  return R;
230 }
231 } // namespace KelvinVector
232 } // namespace MathLib
OGS_FATAL
#define OGS_FATAL(...)
Definition: Error.h:26
MathLib::KelvinVector::kelvinVectorToSymmetricTensor
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
Definition: KelvinVector.cpp:142
MathLib::KelvinVector::tensorToKelvin< 2 >
KelvinVectorType< 2 > tensorToKelvin< 2 >(Eigen::Matrix< double, 3, 3 > const &m)
Definition: KelvinVector.cpp:111
MathLib::KelvinVector::tensorToKelvin< 3 >
KelvinVectorType< 3 > tensorToKelvin< 3 >(Eigen::Matrix< double, 3, 3 > const &m)
Definition: KelvinVector.cpp:126
MathLib::KelvinVector::fourthOrderRotationMatrix< 3 >
KelvinMatrixType< 3 > fourthOrderRotationMatrix< 3 >(Eigen::Matrix< double, 3, 3, Eigen::ColMajor, 3, 3 > const &transformation)
Definition: KelvinVector.cpp:201
MathLib::KelvinVector::kelvinVectorToTensor
Eigen::Matrix< double, 3, 3 > kelvinVectorToTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
Definition: KelvinVector.cpp:64
MathLib::KelvinVector::KelvinMatrixType
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
Definition: KelvinVector.h:56
MathLib::KelvinVector::inverse
Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > inverse(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
Definition: KelvinVector.cpp:34
MathLib::KelvinVector::fourthOrderRotationMatrix< 2 >
KelvinMatrixType< 2 > fourthOrderRotationMatrix< 2 >(Eigen::Matrix< double, 2, 2, Eigen::ColMajor, 2, 2 > const &transformation)
Definition: KelvinVector.cpp:184
MathLib::KelvinVector::Invariants::determinant
static double determinant(Eigen::Matrix< double, KelvinVectorSize, 1 > const &v)
Determinant of a matrix in Kelvin vector representation.