Go to the documentation of this file.
35 namespace FiniteElement {
54 Eigen::MatrixXd
NB(
const Eigen::VectorXd&
N,
55 const int nDoFPerNode );
57 template <
int nDim,
int nNodes >
58 Eigen::Matrix< double, nDim, nDim * nNodes >
NB(
const Eigen::Matrix< double, 1, nNodes >&
N )
61 Eigen::Matrix< double, nDim, nDim* nNodes >
N_ = Eigen::Matrix< double, nDim, nDim * nNodes >::Zero();
63 for (
int j = 0;
j <
nDim;
j++ ) {
70 Eigen::MatrixXd
Jacobian(
const Eigen::MatrixXd& dN_dXi,
71 const Eigen::VectorXd& coordinates );
73 template <
int nDim,
int nNodes >
74 Eigen::Matrix< double, nDim, nDim >
Jacobian(
const Eigen::Matrix< double, nDim, nNodes >&
dNdXi,
75 const Eigen::Matrix< double, nDim * nNodes, 1 >& coordinates )
79 Eigen::Matrix< double, nDim, nDim >
J_ = Eigen::Matrix< double, nDim, nDim >::Zero();
80 for (
int i = 0;
i <
nDim;
i++ )
81 for (
int j = 0;
j <
nDim;
j++ )
93 using NSized = Eigen::Matrix< double, 1, nNodes >;
103 using NSized = Eigen::Matrix< double, 1, nNodes >;
111 namespace Spatial2D {
115 template <
int nNodes >
116 Eigen::Matrix< double, voigtSize, nNodes * nDim >
B(
const Eigen::Matrix< double, nDim, nNodes >& dNdX )
119 Eigen::Matrix< double, voigtSize, nNodes* nDim >
B_ = Eigen::Matrix< double, voigtSize, nNodes * nDim >::Zero();
121 B_( 0, 2 *
i ) = dNdX( 0,
i );
122 B_( 1, 2 *
i + 1 ) = dNdX( 1,
i );
123 B_( 2, 2 *
i ) = dNdX( 1,
i );
124 B_( 2, 2 *
i + 1 ) = dNdX( 0,
i );
129 template <
int nNodes >
130 Eigen::Matrix< double, voigtSize, nNodes * nDim >
BGreen(
const Eigen::Matrix< double, nDim, nNodes >& dNdX,
131 const Eigen::Matrix2d&
F )
135 Eigen::Matrix< double, voigtSize, nNodes* nDim >
B_ = Eigen::Matrix< double, voigtSize, nNodes * nDim >::Zero();
137 B_( 0, 2 *
i ) = dNdX( 0,
i ) *
F( 0, 0 );
138 B_( 0, 2 *
i + 1 ) = dNdX( 0,
i ) *
F( 1, 0 );
139 B_( 1, 2 *
i ) = dNdX( 1,
i ) *
F( 0, 1 );
140 B_( 1, 2 *
i + 1 ) = dNdX( 1,
i ) *
F( 1, 1 );
141 B_( 2, 2 *
i ) = dNdX( 0,
i ) *
F( 0, 1 ) + dNdX( 1,
i ) *
F( 0, 0 );
142 B_( 2, 2 *
i + 1 ) = dNdX( 0,
i ) *
F( 1, 1 ) + dNdX( 1,
i ) *
F( 1, 0 );
151 using NSized = Eigen::Matrix< double, 1, nNodes >;
164 using NSized = Eigen::Matrix< double, 1, nNodes >;
175 namespace Spatial3D {
179 template <
int nNodes >
180 Eigen::Matrix< double, voigtSize, nNodes * nDim >
B(
const Eigen::Matrix< double, nDim, nNodes >& dNdX )
191 Eigen::Matrix< double, voigtSize, nNodes* nDim >
B_ = Eigen::Matrix< double, voigtSize, nNodes * nDim >::Zero();
195 B_( 1,
nDim *
i + 1 ) = dNdX( 1,
i );
196 B_( 2,
nDim *
i + 2 ) = dNdX( 2,
i );
197 B_( 3,
nDim *
i + 0 ) = dNdX( 1,
i );
198 B_( 3,
nDim *
i + 1 ) = dNdX( 0,
i );
199 B_( 4,
nDim *
i + 0 ) = dNdX( 2,
i );
200 B_( 4,
nDim *
i + 2 ) = dNdX( 0,
i );
201 B_( 5,
nDim *
i + 1 ) = dNdX( 2,
i );
202 B_( 5,
nDim *
i + 2 ) = dNdX( 1,
i );
208 template <
int nNodes >
209 Eigen::Matrix< double, voigtSize, nNodes * nDim >
BGreen(
const Eigen::Matrix< double, nDim, nNodes >& dNdX,
215 Eigen::Matrix< double, voigtSize, nNodes* nDim >
B_ = Eigen::Matrix< double, voigtSize, nNodes * nDim >::Zero();
217 B_( 0,
nDim *
i ) = dNdX( 0,
i ) *
F( 0, 0 );
218 B_( 0,
nDim *
i + 1 ) = dNdX( 0,
i ) *
F( 1, 0 );
219 B_( 0,
nDim *
i + 2 ) = dNdX( 0,
i ) *
F( 2, 0 );
221 B_( 1,
nDim *
i ) = dNdX( 1,
i ) *
F( 0, 1 );
222 B_( 1,
nDim *
i + 1 ) = dNdX( 1,
i ) *
F( 1, 1 );
223 B_( 1,
nDim *
i + 2 ) = dNdX( 1,
i ) *
F( 2, 1 );
225 B_( 2,
nDim *
i ) = dNdX( 2,
i ) *
F( 0, 2 );
226 B_( 2,
nDim *
i + 1 ) = dNdX( 2,
i ) *
F( 1, 2 );
227 B_( 2,
nDim *
i + 2 ) = dNdX( 2,
i ) *
F( 2, 2 );
229 B_( 3,
nDim *
i ) = dNdX( 0,
i ) *
F( 0, 1 ) + dNdX( 1,
i ) *
F( 0, 0 );
230 B_( 3,
nDim *
i + 1 ) = dNdX( 0,
i ) *
F( 1, 1 ) + dNdX( 1,
i ) *
F( 1, 0 );
231 B_( 3,
nDim *
i + 2 ) = dNdX( 0,
i ) *
F( 2, 1 ) + dNdX( 1,
i ) *
F( 2, 0 );
233 B_( 4,
nDim *
i ) = dNdX( 0,
i ) *
F( 0, 2 ) + dNdX( 2,
i ) *
F( 0, 0 );
234 B_( 4,
nDim *
i + 1 ) = dNdX( 0,
i ) *
F( 1, 2 ) + dNdX( 2,
i ) *
F( 1, 0 );
235 B_( 4,
nDim *
i + 2 ) = dNdX( 0,
i ) *
F( 2, 2 ) + dNdX( 2,
i ) *
F( 2, 0 );
237 B_( 5,
nDim *
i ) = dNdX( 2,
i ) *
F( 0, 1 ) + dNdX( 1,
i ) *
F( 0, 2 );
238 B_( 5,
nDim *
i + 1 ) = dNdX( 2,
i ) *
F( 1, 1 ) + dNdX( 1,
i ) *
F( 1, 2 );
239 B_( 5,
nDim *
i + 2 ) = dNdX( 2,
i ) *
F( 2, 1 ) + dNdX( 1,
i ) *
F( 2, 2 );
248 using NSized = Eigen::Matrix< double, 1, nNodes >;
262 using NSized = Eigen::Matrix< double, 1, nNodes >;
275 using NSized = Eigen::Matrix< double, 1, nNodes >;
287 using NSized = Eigen::Matrix< double, 1, nNodes >;
327 int parentFaceNumber,
328 const Eigen::VectorXd& parentCoordinates );
343 Eigen::Ref< Eigen::VectorXd > ParentVector );
348 Eigen::Ref< Eigen::VectorXd > ParentVector );
350 Eigen::Ref< Eigen::MatrixXd > KParent );
354 namespace FiniteElement::Quadrature {
355 constexpr
double gp2 = 0.577350269189625764509;
356 constexpr
double gp3 = 0.774596669241483;
370 namespace Spatial1D {
375 { ( Eigen::VectorXd ( 1 ) << 0 ).finished(), 2.0 }
379 { ( Eigen::VectorXd ( 1 ) << -
gp2 ).finished(), 1.0 },
380 { ( Eigen::VectorXd ( 1 ) << +
gp2 ).finished(), 1.0 }
384 { ( Eigen::VectorXd ( 1 ) << -
gp3 ).finished(), 5./9 },
385 { ( Eigen::VectorXd ( 1 ) << 0. ).finished(), 8./9 },
386 { ( Eigen::VectorXd ( 1 ) << +
gp3 ).finished(), 5./9 }
392 namespace Spatial2D {
397 { Eigen::Vector2d::Zero(), 4. }
401 { ( Eigen::Vector2d () << +
gp2, +
gp2 ).finished(), 1.0 },
402 { ( Eigen::Vector2d () << -
gp2, +
gp2 ).finished(), 1.0 },
403 { ( Eigen::Vector2d () << -
gp2, -
gp2 ).finished(), 1.0 },
404 { ( Eigen::Vector2d () << +
gp2, -
gp2 ).finished(), 1.0 }
408 { ( Eigen::Vector2d () << 0, 0. ).finished(), 64./81},
409 { ( Eigen::Vector2d () << -
gp3, -
gp3 ).finished(), 25./81},
410 { ( Eigen::Vector2d () << +
gp3, -
gp3 ).finished(), 25./81},
411 { ( Eigen::Vector2d () << +
gp3, +
gp3 ).finished(), 25./81},
412 { ( Eigen::Vector2d () << -
gp3, +
gp3 ).finished(), 25./81},
413 { ( Eigen::Vector2d () << 0, -
gp3 ).finished(), 40./81},
414 { ( Eigen::Vector2d () <<
gp3, 0. ).finished(), 40./81},
415 { ( Eigen::Vector2d () << 0, +
gp3 ).finished(), 40./81},
416 { ( Eigen::Vector2d () << -
gp3, 0. ).finished(), 40./81},
423 namespace Spatial3D {
428 { Eigen::Vector3d::Zero(), 8.0 }
436 { (
Eigen::Vector3d() << (5-std::sqrt(5))/20, (5-std::sqrt(5))/20, (5-std::sqrt(5))/20 ).finished(), 1./24},
437 { (
Eigen::Vector3d() << (5-std::sqrt(5))/20, (5-std::sqrt(5))/20, (5+3*std::sqrt(5))/20 ).finished(), 1./24},
438 { (
Eigen::Vector3d() << (5-std::sqrt(5))/20, (5+3*std::sqrt(5))/20, (5-std::sqrt(5))/20 ).finished(), 1./24},
439 { (
Eigen::Vector3d() << (5+3*std::sqrt(5))/20, (5-std::sqrt(5))/20, (5-std::sqrt(5))/20 ).finished(), 1./24},
NSized N(const Eigen::Vector3d &xi)
Definition: MarmotFiniteElement.h:302
Fastor::Index< i_ > i
Definition: MarmotFastorTensorBasics.h:137
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:164
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:262
NSized N(double xi)
Definition: MarmotFiniteElement1D.cpp:41
constexpr int voigtSize
Definition: MarmotFiniteElement.h:113
const std::vector< QuadraturePointInfo > gaussPointList3x3
Definition: MarmotFiniteElement.h:407
Eigen::VectorXd computeScalarLoadVector()
Definition: MarmotFiniteElementBoundary.cpp:124
NSized N(const Eigen::Vector3d &xi)
Eigen::MatrixXd computeDSurfaceNormalVectorialLoadVector_dCoordinates()
Definition: MarmotFiniteElementBoundary.cpp:162
constexpr int nNodes
Definition: MarmotFiniteElement.h:246
BoundaryElement(ElementShapes parentShape, int nDim, int parentFaceNumber, const Eigen::VectorXd &parentCoordinates)
Definition: MarmotFiniteElementBoundary.cpp:10
@ FullIntegration
Definition: MarmotFiniteElement.h:358
Eigen::MatrixXd computeDScalarLoadVector_dCoordinates()
Definition: MarmotFiniteElementBoundary.cpp:139
Eigen::VectorXd condenseParentToBoundaryScalar(const Eigen::VectorXd &parentVector)
Definition: MarmotFiniteElementBoundary.cpp:229
void assembleIntoParentScalar(const Eigen::VectorXd &boundaryVector, Eigen::Ref< Eigen::VectorXd > ParentVector)
Definition: MarmotFiniteElementBoundary.cpp:242
Eigen::Matrix< double, voigtSize, nNodes *nDim > B(const Eigen::Matrix< double, nDim, nNodes > &dNdX)
Definition: MarmotFiniteElement.h:180
double weight
Definition: MarmotFiniteElement.h:303
Eigen::Matrix< double, nDim, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:288
const std::vector< QuadraturePointInfo > gaussPointList1x1x1
Definition: MarmotFiniteElement.h:427
@ N_
Definition: MarmotFastorTensorBasics.h:90
Eigen::VectorXd coordinates
Definition: MarmotFiniteElement.h:322
Eigen::Matrix< double, 1, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:94
std::vector< BoundaryElementQuadraturePoint > quadraturePoints
Definition: MarmotFiniteElement.h:318
@ ReducedIntegration
Definition: MarmotFiniteElement.h:358
Eigen::Matrix< double, nDim, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:152
Eigen::Matrix< double, nDim, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:165
@ Quad9
Definition: MarmotFiniteElement.h:41
constexpr double gp2
Definition: MarmotFiniteElement.h:355
Eigen::Matrix< double, nNodes *nDim, 1 > CoordinateSized
Definition: MarmotFiniteElement.h:247
@ Hexa20
Definition: MarmotFiniteElement.h:46
NSized N(const Eigen::Vector2d &xi)
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:275
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:103
@ B_
Definition: MarmotFastorTensorBasics.h:90
@ Bar3
Definition: MarmotFiniteElement.h:38
Eigen::VectorXd N
Definition: MarmotFiniteElement.h:306
void assembleIntoParentVectorial(const Eigen::VectorXd &boundaryVector, Eigen::Ref< Eigen::VectorXd > ParentVector)
Definition: MarmotFiniteElementBoundary.cpp:274
Eigen::VectorXd condenseParentToBoundaryVectorial(const Eigen::VectorXd &parentVector)
Definition: MarmotFiniteElementBoundary.cpp:261
constexpr int nNodes
Definition: MarmotFiniteElement.h:260
Definition: MarmotFiniteElement.h:297
Eigen::MatrixXd dNdXi
Definition: MarmotFiniteElement.h:307
@ Quad4
Definition: MarmotFiniteElement.h:39
constexpr int nDim
Definition: MarmotFiniteElement.h:112
constexpr int nNodes
Definition: MarmotFiniteElement.h:285
Eigen::Matrix< double, 3, 1 > Vector3d
Definition: MarmotTypedefs.h:42
int nParentCoordinates
Definition: MarmotFiniteElement.h:316
const std::vector< QuadraturePointInfo > gaussPointList1
Definition: MarmotFiniteElement.h:374
Eigen::Matrix< double, nNodes *nDim, 1 > CoordinateSized
Definition: MarmotFiniteElement.h:261
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:93
Eigen::Matrix< double, nNodes *nDim, 1 > CoordinateSized
Definition: MarmotFiniteElement.h:150
@ Tetra4
Definition: MarmotFiniteElement.h:43
Eigen::Matrix< double, voigtSize, nNodes *nDim > BGreen(const Eigen::Matrix< double, nDim, nNodes > &dNdX, const Eigen::Matrix3d &F)
Definition: MarmotFiniteElement.h:209
Eigen::VectorXd computeVectorialLoadVector(const Eigen::VectorXd &direction)
compute the element load vector for a unit vectorial load in a given direction.
Definition: MarmotFiniteElementBoundary.cpp:206
const std::vector< QuadraturePointInfo > gaussPointList2x2
Definition: MarmotFiniteElement.h:400
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:287
constexpr int nDim
Definition: MarmotFiniteElement.h:371
NSized N(const Eigen::Vector3d &xi)
constexpr int nNodes
Definition: MarmotFiniteElement.h:148
@ Hexa8
Definition: MarmotFiniteElement.h:45
constexpr int voigtSize
Definition: MarmotFiniteElement.h:177
dNdXiSized dNdXi(const Eigen::Vector2d &xi)
constexpr int nNodes
Definition: MarmotFiniteElement.h:161
const std::vector< QuadraturePointInfo > gaussPointList1x1
Definition: MarmotFiniteElement.h:396
Eigen::Matrix< double, nDim, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:249
constexpr int nDim
Definition: MarmotFiniteElement.h:176
Eigen::Matrix< double, 3, 3 > Matrix3d
Definition: MarmotTypedefs.h:40
Fastor::Index< j_ > j
Definition: MarmotFastorTensorBasics.h:182
@ Bar2
Definition: MarmotFiniteElement.h:37
Eigen::MatrixXd computeDVectorialLoadVector_dCoordinates(const Eigen::VectorXd &direction)
Definition: MarmotFiniteElementBoundary.cpp:221
This file includes functions needed for calculations with stress and strain tensors written in voigt ...
Definition: MarmotTesting.h:37
const std::vector< QuadraturePointInfo > gaussPointList3x3x3
Definition: MarmotFiniteElement.h:453
Eigen::Matrix< double, nNodes *nDim, 1 > CoordinateSized
Definition: MarmotFiniteElement.h:163
Eigen::VectorXd computeSurfaceNormalVectorialLoadVector()
compute the element load vector for a unit vectorial load normal to the surface.
Definition: MarmotFiniteElementBoundary.cpp:147
Eigen::Matrix< double, voigtSize, nNodes *nDim > B(const Eigen::Matrix< double, nDim, nNodes > &dNdX)
Definition: MarmotFiniteElement.h:116
Eigen::VectorXi expandNodeIndicesToCoordinateIndices(const Eigen::VectorXi &nodeIndices, int nDim)
@ Tetra10
Definition: MarmotFiniteElement.h:44
const std::vector< QuadraturePointInfo > gaussPointList2
Definition: MarmotFiniteElement.h:378
Eigen::Vector3i getBoundaryElementIndices(int faceID)
Definition: MarmotFiniteElement2D.cpp:131
double JxW
Definition: MarmotFiniteElement.h:304
NSized N(const Eigen::Vector3d &xi)
Eigen::VectorXd xi
Definition: MarmotFiniteElement.h:305
constexpr int nDim
Definition: MarmotFiniteElement.h:424
const std::vector< QuadraturePointInfo > gaussPointList3
Definition: MarmotFiniteElement.h:383
void assembleIntoParentStiffnessScalar(const Eigen::MatrixXd &KBoundary, Eigen::Ref< Eigen::MatrixXd > KParent)
Definition: MarmotFiniteElementBoundary.cpp:252
Eigen::MatrixXd NB(const Eigen::VectorXd &N, const int nDoFPerNode)
Eigen::VectorXd areaVector
Definition: MarmotFiniteElement.h:309
Eigen::VectorXi mapBoundaryToParentVectorial
Definition: MarmotFiniteElement.h:321
ElementShapes
Definition: MarmotFiniteElement.h:36
@ Quad16
Definition: MarmotFiniteElement.h:42
@ Hexa64
Definition: MarmotFiniteElement.h:48
Eigen::Vector3i getBoundaryElementIndices(int faceID)
Definition: MarmotFiniteElement3D.cpp:68
Eigen::Matrix< int, 8, 1 > Vector8i
Definition: MarmotTypedefs.h:47
dNdXiSized dNdXi(const Eigen::Vector2d &xi)
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:248
dNdXiSized dNdXi(const Eigen::Vector3d &xi)
constexpr int nNodes
Definition: MarmotFiniteElement.h:92
int getNumGaussPoints(Marmot::FiniteElement::ElementShapes shape, IntegrationTypes integrationType)
Definition: MarmotFiniteElementBasic.cpp:160
const std::vector< QuadraturePointInfo > & getGaussPointInfo(Marmot::FiniteElement::ElementShapes shape, IntegrationTypes integrationType)
Definition: MarmotFiniteElementBasic.cpp:100
Eigen::Matrix< double, nNodes *nDim, 1 > CoordinateSized
Definition: MarmotFiniteElement.h:286
Eigen::Vector2i getBoundaryElementIndices(int faceID)
Definition: MarmotFiniteElement2D.cpp:46
Eigen::Matrix< double, nDim, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:263
Eigen::Vector3i getBoundaryElementIndices(int faceID)
Definition: MarmotFiniteElement3D.cpp:177
ElementShapes getElementShapeByMetric(int nDim, int nNodes)
Definition: MarmotFiniteElementBasic.cpp:7
dNdXiSized dNdXi(double xi)
Definition: MarmotFiniteElement1D.cpp:22
const std::vector< QuadraturePointInfo > gaussPointListTetra10
Definition: MarmotFiniteElement.h:435
Marmot::Vector8i getBoundaryElementIndices(int faceID)
Definition: MarmotFiniteElement3D.cpp:419
@ Hexa27
Definition: MarmotFiniteElement.h:47
constexpr int nNodes
Definition: MarmotFiniteElement.h:102
ElementShapes boundaryShape
Definition: MarmotFiniteElement.h:314
dNdXiSized dNdXi(const Eigen::Vector3d &xi)
dNdXiSized dNdXi(const Eigen::Vector3d &xi)
Fastor::Index< k_ > k
Definition: MarmotFastorTensorBasics.h:195
int nNodes
Definition: MarmotFiniteElement.h:315
@ J_
Definition: MarmotFastorTensorBasics.h:90
Eigen::MatrixXd F(const Eigen::MatrixXd &J)
Eigen::Vector4i getBoundaryElementIndices(int faceID)
Definition: MarmotFiniteElement3D.cpp:261
IntegrationTypes
Definition: MarmotFiniteElement.h:358
constexpr double gp3
Definition: MarmotFiniteElement.h:356
Eigen::VectorXi mapBoundaryToParentScalar
Definition: MarmotFiniteElement.h:320
constexpr int nDim
Definition: MarmotFiniteElement.h:393
double weight
Definition: MarmotFiniteElement.h:362
Eigen::VectorXd xi
Definition: MarmotFiniteElement.h:361
dNdXiSized dNdXi(double xi)
Definition: MarmotFiniteElement1D.cpp:48
void assembleIntoParentStiffnessVectorial(const Eigen::MatrixXd &KBoundary, Eigen::Ref< Eigen::MatrixXd > KParent)
Definition: MarmotFiniteElementBoundary.cpp:284
Eigen::Matrix< double, 1, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:104
const std::vector< QuadraturePointInfo > gaussPointList2x2x2
Definition: MarmotFiniteElement.h:442
Eigen::MatrixXd dx_dXi
Definition: MarmotFiniteElement.h:308
dNdXiSized dNdXi(const Eigen::Vector3d &xi)
constexpr int nNodes
Definition: MarmotFiniteElement.h:273
@ Quad8
Definition: MarmotFiniteElement.h:40
Eigen::MatrixXd Jacobian(const Eigen::MatrixXd &dN_dXi, const Eigen::VectorXd &coordinates)
Eigen::Matrix< double, 1, nNodes > NSized
Definition: MarmotFiniteElement.h:151
NSized N(double xi)
Definition: MarmotFiniteElement1D.cpp:15
Definition: MarmotFiniteElement.h:360
const int nDim
Definition: MarmotFiniteElement.h:312
Eigen::Matrix< double, voigtSize, nNodes *nDim > BGreen(const Eigen::Matrix< double, nDim, nNodes > &dNdX, const Eigen::Matrix2d &F)
Definition: MarmotFiniteElement.h:130
NSized N(const Eigen::Vector2d &xi)
Eigen::Matrix< double, nDim, nNodes > dNdXiSized
Definition: MarmotFiniteElement.h:276
const std::vector< QuadraturePointInfo > gaussPointListTetra4
Definition: MarmotFiniteElement.h:431
void modifyCharElemLengthAbaqusLike(double &charElemLength, int intPoint)
Eigen::Matrix< double, nNodes *nDim, 1 > CoordinateSized
Definition: MarmotFiniteElement.h:274