OGS
KelvinVector.h
Go to the documentation of this file.
1
10#pragma once
11
12#include <Eigen/Core>
13
14#include "BaseLib/Error.h"
15
16namespace MathLib
17{
21namespace KelvinVector
22{
24constexpr int kelvin_vector_dimensions(int const displacement_dim)
25{
26 if (displacement_dim == 2)
27 {
28 return 4;
29 }
30 else if (displacement_dim == 3)
31 {
32 return 6;
33 }
35 "Cannot convert displacement dimension {} to kelvin vector dimension.",
36 displacement_dim);
37}
38
39//
40// Kelvin vector and matrix templates for given displacement dimension.
41//
42
46template <int DisplacementDim>
48 Eigen::Matrix<double, kelvin_vector_dimensions(DisplacementDim), 1,
49 Eigen::ColMajor>;
50
54template <int DisplacementDim>
56 Eigen::Matrix<double, kelvin_vector_dimensions(DisplacementDim),
57 kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor>;
58
62template <int KelvinVectorSize>
63struct Invariants final
64{
65 static_assert(KelvinVectorSize == 4 || KelvinVectorSize == 6,
66 "KelvinVector invariants for vectors of size different than "
67 "4 or 6 is not allowed.");
70 static Eigen::Matrix<double, KelvinVectorSize, KelvinVectorSize> const
74 static Eigen::Matrix<double, KelvinVectorSize, KelvinVectorSize> const
77 static Eigen::Matrix<double, KelvinVectorSize, 1> const identity2;
78
80 static double determinant(
81 Eigen::Matrix<double, KelvinVectorSize, 1> const& v);
82
85 static double equivalentStress(
86 Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
87
89 static double FrobeniusNorm(
90 Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
91
94 static double J2(
95 Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
96
99 static double J3(
100 Eigen::Matrix<double, KelvinVectorSize, 1> const& deviatoric_v);
101
103 static double trace(Eigen::Matrix<double, KelvinVectorSize, 1> const& v);
104
107 static Eigen::Vector3d diagonal(
108 Eigen::Matrix<double, KelvinVectorSize, 1> const& v);
109};
110
111//
112// Inverses of a Kelvin vector.
113//
114
117template <int KelvinVectorSize>
118Eigen::Matrix<double, KelvinVectorSize, 1, Eigen::ColMajor, KelvinVectorSize, 1>
119inverse(Eigen::Matrix<double,
120 KelvinVectorSize,
121 1,
122 Eigen::ColMajor,
123 KelvinVectorSize,
124 1> const& v);
125
128template <int KelvinVectorSize>
129Eigen::Matrix<double, 3, 3> kelvinVectorToTensor(Eigen::Matrix<double,
130 KelvinVectorSize,
131 1,
132 Eigen::ColMajor,
133 KelvinVectorSize,
134 1> const& v);
135
138template <int DisplacementDim>
140 Eigen::Matrix<double, 3, 3> const& m);
141
153template <int KelvinVectorSize>
154Eigen::Matrix<double, KelvinVectorSize, 1, Eigen::ColMajor, KelvinVectorSize, 1>
155kelvinVectorToSymmetricTensor(Eigen::Matrix<double,
156 KelvinVectorSize,
157 1,
158 Eigen::ColMajor,
159 KelvinVectorSize,
160 1> const& v);
161
169template <typename Derived>
170Eigen::Matrix<double, Eigen::MatrixBase<Derived>::RowsAtCompileTime, 1>
171symmetricTensorToKelvinVector(Eigen::MatrixBase<Derived> const& v)
172{
173 static_assert(
174 (Eigen::MatrixBase<Derived>::ColsAtCompileTime == 1) ||
175 (Eigen::MatrixBase<Derived>::ColsAtCompileTime == Eigen::Dynamic),
176 "KelvinVector must be a column vector");
177 if (v.cols() != 1)
178 {
179 OGS_FATAL(
180 "KelvinVector must be a column vector, but input has {:d} columns.",
181 v.cols());
182 }
183
184 Eigen::Matrix<double, Eigen::MatrixBase<Derived>::RowsAtCompileTime, 1>
185 result;
186 if (v.rows() == 4)
187 {
188 result.resize(4, 1);
189 result << v[0], v[1], v[2], v[3] * std::sqrt(2.);
190 }
191 else if (v.rows() == 6)
192 {
193 result.resize(6, 1);
194 result << v[0], v[1], v[2], v[3] * std::sqrt(2.), v[4] * std::sqrt(2.),
195 v[5] * std::sqrt(2.);
196 }
197 else
198 {
199 OGS_FATAL(
200 "Symmetric tensor to Kelvin vector conversion expected an input "
201 "vector of size 4 or 6, but a vector of size {:d} was given.",
202 v.size());
203 }
204 return result;
205}
206
211template <int DisplacementDim>
213 std::vector<double> const& values)
214{
215 constexpr int kelvin_vector_size =
216 kelvin_vector_dimensions(DisplacementDim);
217
218 if (values.size() != kelvin_vector_size)
219 {
220 OGS_FATAL(
221 "Symmetric tensor to Kelvin vector conversion expected an input "
222 "vector of size {:d}, but a vector of size {:d} was given.",
223 kelvin_vector_size, values.size());
224 }
225
227 Eigen::Map<typename MathLib::KelvinVector::KelvinVectorType<
228 DisplacementDim> const>(
229 values.data(), kelvin_vector_dimensions(DisplacementDim), 1));
230}
231
235template <int DisplacementDim>
236Eigen::Matrix<double, DisplacementDim,
237 kelvin_vector_dimensions(DisplacementDim)>
238liftVectorToKelvin(Eigen::Matrix<double, DisplacementDim, 1> const& v);
239
242template <int DisplacementDim>
243Eigen::Matrix<double, DisplacementDim, 1> reduceKelvinToVector(
244 Eigen::Matrix<double, DisplacementDim,
245 kelvin_vector_dimensions(DisplacementDim)> const& m);
246
251template <int DisplacementDim>
253 Eigen::Matrix<double, DisplacementDim, DisplacementDim, Eigen::ColMajor,
254 DisplacementDim, DisplacementDim> const& transformation);
255
256} // namespace KelvinVector
257} // namespace MathLib
258
259#include "KelvinVector-impl.h"
#define OGS_FATAL(...)
Definition: Error.h:26
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:49
Eigen::Matrix< double, DisplacementDim, 1 > reduceKelvinToVector(Eigen::Matrix< double, DisplacementDim, kelvin_vector_dimensions(DisplacementDim)> const &m)
KelvinMatrixType< DisplacementDim > fourthOrderRotationMatrix(Eigen::Matrix< double, DisplacementDim, DisplacementDim, Eigen::ColMajor, DisplacementDim, DisplacementDim > const &transformation)
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:24
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
Definition: KelvinVector.h:171
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:57
Eigen::Matrix< double, DisplacementDim, kelvin_vector_dimensions(DisplacementDim)> liftVectorToKelvin(Eigen::Matrix< double, DisplacementDim, 1 > const &v)
KelvinVectorType< DisplacementDim > tensorToKelvin(Eigen::Matrix< double, 3, 3 > const &m)
Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > inverse(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
static const double v
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > KelvinVector
Definition: Base.h:22
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:77
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:75
static double J2(Eigen::Matrix< double, KelvinVectorSize, 1 > const &deviatoric_v)
static Eigen::Matrix< double, KelvinVectorSize, KelvinVectorSize > const deviatoric_projection
Definition: KelvinVector.h:71
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)