Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F70997149
FE_Quadrature.h
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Mon, Jul 8, 16:53
Size
31 KB
Mime Type
text/x-c
Expires
Wed, Jul 10, 16:53 (1 d, 23 h)
Engine
blob
Format
Raw Data
Handle
18805373
Attached To
rLAMMPS lammps
FE_Quadrature.h
View Options
#ifndef FE_QUADRATURE_H
#define FE_QUADRATURE_H
#include <vector>
#include "MatrixLibrary.h"
#include "ATC_TypeDefs.h"
namespace
ATC
{
/**
* @class FE_Quadrature
* @brief Stores quadrature schemes
*/
struct
FE_Quadrature
{
int
numIPs
;
int
numFaceIPs
;
DENS_MAT
ipCoords
;
std
::
vector
<
DENS_MAT
>
ipFaceCoords
;
DENS_MAT
ipFace2DCoords
;
DENS_VEC
ipWeights
;
DENS_VEC
ipFaceWeights
;
FE_Quadrature
(
FeEltGeometry
geo
,
FeIntQuadrature
quad
)
{
if
(
geo
==
HEXA
)
{
switch
(
quad
)
{
// Degenerate 1 point quadrature
case
GAUSS1
:
{
// Set number of IPs and face IPs
numIPs
=
1
;
numFaceIPs
=
1
;
// Size matrices and vectors thereof accordingly
ipCoords
.
resize
(
3
,
numIPs
);
ipFaceCoords
.
assign
(
6
,
DENS_MAT
(
3
,
numFaceIPs
));
ipFace2DCoords
.
resize
(
2
,
numFaceIPs
);
ipWeights
.
reset
(
numIPs
);
ipFaceWeights
.
reset
(
numFaceIPs
);
// "Matrix" of integration point location
ipCoords
(
0
,
0
)
=
0.0
;
ipCoords
(
1
,
0
)
=
0.0
;
ipCoords
(
2
,
0
)
=
0.0
;
// Integration point for each face
ipFaceCoords
[
0
](
0
,
0
)
=
-
1.0
;
ipFaceCoords
[
3
](
0
,
0
)
=
0.0
;
ipFaceCoords
[
0
](
1
,
0
)
=
0.0
;
ipFaceCoords
[
3
](
1
,
0
)
=
1.0
;
ipFaceCoords
[
0
](
2
,
0
)
=
0.0
;
ipFaceCoords
[
3
](
2
,
0
)
=
0.0
;
//
ipFaceCoords
[
1
](
0
,
0
)
=
1.0
;
ipFaceCoords
[
4
](
0
,
0
)
=
0.0
;
ipFaceCoords
[
1
](
1
,
0
)
=
0.0
;
ipFaceCoords
[
4
](
1
,
0
)
=
0.0
;
ipFaceCoords
[
1
](
2
,
0
)
=
0.0
;
ipFaceCoords
[
4
](
2
,
0
)
=
-
1.0
;
//
ipFaceCoords
[
2
](
0
,
0
)
=
0.0
;
ipFaceCoords
[
5
](
0
,
0
)
=
0.0
;
ipFaceCoords
[
2
](
1
,
0
)
=
-
1.0
;
ipFaceCoords
[
5
](
1
,
0
)
=
0.0
;
ipFaceCoords
[
2
](
2
,
0
)
=
0.0
;
ipFaceCoords
[
5
](
2
,
0
)
=
1.0
;
// 2D integration scheme for the faces
ipFace2DCoords
(
0
,
0
)
=
0.0
;
ipFace2DCoords
(
1
,
0
)
=
0.0
;
// Integration point weights
ipWeights
=
8.0
;
// Face integration point weights
ipFaceWeights
=
4.0
;
break
;
}
// 8 point quadratures
case
NODAL
:
case
GAUSS2
:
{
// Set number of IPs and face IPs
numIPs
=
8
;
numFaceIPs
=
4
;
// Size matrices and vectors thereof accordingly
ipCoords
.
resize
(
3
,
numIPs
);
ipFaceCoords
.
assign
(
6
,
DENS_MAT
(
3
,
numFaceIPs
));
ipFace2DCoords
.
resize
(
2
,
numFaceIPs
);
ipWeights
.
reset
(
numIPs
);
ipFaceWeights
.
reset
(
numFaceIPs
);
// Dictates difference in node locations for nodal/GAUSS2
double
a
=
1.0
/
sqrt
(
3.0
);
if
(
quad
==
NODAL
)
a
=
1.0
;
// Matrix of integration point locations & follows local
// conn
ipCoords
(
0
,
0
)
=
-
a
;
ipCoords
(
0
,
4
)
=
-
a
;
ipCoords
(
1
,
0
)
=
-
a
;
ipCoords
(
1
,
4
)
=
-
a
;
ipCoords
(
2
,
0
)
=
-
a
;
ipCoords
(
2
,
4
)
=
a
;
//
ipCoords
(
0
,
1
)
=
a
;
ipCoords
(
0
,
5
)
=
a
;
ipCoords
(
1
,
1
)
=
-
a
;
ipCoords
(
1
,
5
)
=
-
a
;
ipCoords
(
2
,
1
)
=
-
a
;
ipCoords
(
2
,
5
)
=
a
;
//
ipCoords
(
0
,
2
)
=
a
;
ipCoords
(
0
,
6
)
=
a
;
ipCoords
(
1
,
2
)
=
a
;
ipCoords
(
1
,
6
)
=
a
;
ipCoords
(
2
,
2
)
=
-
a
;
ipCoords
(
2
,
6
)
=
a
;
//
ipCoords
(
0
,
3
)
=
-
a
;
ipCoords
(
0
,
7
)
=
-
a
;
ipCoords
(
1
,
3
)
=
a
;
ipCoords
(
1
,
7
)
=
a
;
ipCoords
(
2
,
3
)
=
-
a
;
ipCoords
(
2
,
7
)
=
a
;
// Integration points by face
ipFaceCoords
[
0
](
0
,
0
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
0
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
3
](
1
,
0
)
=
1
;
ipFaceCoords
[
0
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
3
](
2
,
0
)
=
-
a
;
//
ipFaceCoords
[
0
](
0
,
1
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
1
)
=
a
;
ipFaceCoords
[
0
](
1
,
1
)
=
a
;
ipFaceCoords
[
3
](
1
,
1
)
=
1
;
ipFaceCoords
[
0
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
3
](
2
,
1
)
=
-
a
;
//
ipFaceCoords
[
0
](
0
,
2
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
2
)
=
a
;
ipFaceCoords
[
0
](
1
,
2
)
=
a
;
ipFaceCoords
[
3
](
1
,
2
)
=
1
;
ipFaceCoords
[
0
](
2
,
2
)
=
a
;
ipFaceCoords
[
3
](
2
,
2
)
=
a
;
//
ipFaceCoords
[
0
](
0
,
3
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
0
](
1
,
3
)
=
-
a
;
ipFaceCoords
[
3
](
1
,
3
)
=
1
;
ipFaceCoords
[
0
](
2
,
3
)
=
a
;
ipFaceCoords
[
3
](
2
,
3
)
=
a
;
ipFaceCoords
[
1
](
0
,
0
)
=
1
;
ipFaceCoords
[
4
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
1
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
4
](
2
,
0
)
=
-
1
;
//
ipFaceCoords
[
1
](
0
,
1
)
=
1
;
ipFaceCoords
[
4
](
0
,
1
)
=
a
;
ipFaceCoords
[
1
](
1
,
1
)
=
a
;
ipFaceCoords
[
4
](
1
,
1
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
4
](
2
,
1
)
=
-
1
;
//
ipFaceCoords
[
1
](
0
,
2
)
=
1
;
ipFaceCoords
[
4
](
0
,
2
)
=
a
;
ipFaceCoords
[
1
](
1
,
2
)
=
a
;
ipFaceCoords
[
4
](
1
,
2
)
=
a
;
ipFaceCoords
[
1
](
2
,
2
)
=
a
;
ipFaceCoords
[
4
](
2
,
2
)
=
-
1
;
//
ipFaceCoords
[
1
](
0
,
3
)
=
1
;
ipFaceCoords
[
4
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
1
](
1
,
3
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
3
)
=
a
;
ipFaceCoords
[
1
](
2
,
3
)
=
a
;
ipFaceCoords
[
4
](
2
,
3
)
=
-
1
;
ipFaceCoords
[
2
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
5
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
2
](
1
,
0
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
2
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
5
](
2
,
0
)
=
1
;
//
ipFaceCoords
[
2
](
0
,
1
)
=
a
;
ipFaceCoords
[
5
](
0
,
1
)
=
a
;
ipFaceCoords
[
2
](
1
,
1
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
1
)
=
-
a
;
ipFaceCoords
[
2
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
5
](
2
,
1
)
=
1
;
//
ipFaceCoords
[
2
](
0
,
2
)
=
a
;
ipFaceCoords
[
5
](
0
,
2
)
=
a
;
ipFaceCoords
[
2
](
1
,
2
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
2
)
=
a
;
ipFaceCoords
[
2
](
2
,
2
)
=
a
;
ipFaceCoords
[
5
](
2
,
2
)
=
1
;
//
ipFaceCoords
[
2
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
5
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
2
](
1
,
3
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
3
)
=
a
;
ipFaceCoords
[
2
](
2
,
3
)
=
a
;
ipFaceCoords
[
5
](
2
,
3
)
=
1
;
// Integration points for all faces ignoring the
// redundant dim
ipFace2DCoords
(
0
,
0
)
=
-
a
;
ipFace2DCoords
(
0
,
2
)
=
a
;
ipFace2DCoords
(
1
,
0
)
=
-
a
;
ipFace2DCoords
(
1
,
2
)
=
a
;
//
ipFace2DCoords
(
0
,
1
)
=
a
;
ipFace2DCoords
(
0
,
3
)
=
-
a
;
ipFace2DCoords
(
1
,
1
)
=
-
a
;
ipFace2DCoords
(
1
,
3
)
=
a
;
// Integration point weights
ipWeights
=
1.0
;
// Face integration point weights
ipFaceWeights
=
1.0
;
break
;
}
// 6 point "face" quadrature
case
FACE
:
{
printf
(
"using face quad!
\n
"
);
// Set number of IPs and face IPs
numIPs
=
6
;
numFaceIPs
=
4
;
// Size matrices and vectors thereof accordingly
ipCoords
.
resize
(
3
,
numIPs
);
ipFaceCoords
.
assign
(
6
,
DENS_MAT
(
3
,
numFaceIPs
));
ipFace2DCoords
.
resize
(
2
,
numFaceIPs
);
ipWeights
.
reset
(
numIPs
);
ipFaceWeights
.
reset
(
numFaceIPs
);
// Use GAUSS2 for faces for now...
double
a
=
1.0
/
sqrt
(
3.0
);
// Matrix of integration point locations
ipCoords
(
0
,
0
)
=
1.0
;
ipCoords
(
0
,
3
)
=
0.0
;
ipCoords
(
1
,
0
)
=
0.0
;
ipCoords
(
1
,
3
)
=
-
1.0
;
ipCoords
(
2
,
0
)
=
0.0
;
ipCoords
(
2
,
3
)
=
0.0
;
//
ipCoords
(
0
,
1
)
=
-
1.0
;
ipCoords
(
0
,
4
)
=
0.0
;
ipCoords
(
1
,
1
)
=
0.0
;
ipCoords
(
1
,
4
)
=
0.0
;
ipCoords
(
2
,
1
)
=
0.0
;
ipCoords
(
2
,
4
)
=
1.0
;
//
ipCoords
(
0
,
2
)
=
0.0
;
ipCoords
(
0
,
5
)
=
0.0
;
ipCoords
(
1
,
2
)
=
1.0
;
ipCoords
(
1
,
5
)
=
0.0
;
ipCoords
(
2
,
2
)
=
0.0
;
ipCoords
(
2
,
5
)
=
-
1.0
;
// Integration points by face
ipFaceCoords
[
0
](
0
,
0
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
0
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
3
](
1
,
0
)
=
1
;
ipFaceCoords
[
0
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
3
](
2
,
0
)
=
-
a
;
//
ipFaceCoords
[
0
](
0
,
1
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
1
)
=
-
a
;
ipFaceCoords
[
0
](
1
,
1
)
=
a
;
ipFaceCoords
[
3
](
1
,
1
)
=
1
;
ipFaceCoords
[
0
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
3
](
2
,
1
)
=
a
;
//
ipFaceCoords
[
0
](
0
,
2
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
2
)
=
a
;
ipFaceCoords
[
0
](
1
,
2
)
=
a
;
ipFaceCoords
[
3
](
1
,
2
)
=
1
;
ipFaceCoords
[
0
](
2
,
2
)
=
a
;
ipFaceCoords
[
3
](
2
,
2
)
=
a
;
//
ipFaceCoords
[
0
](
0
,
3
)
=
-
1
;
ipFaceCoords
[
3
](
0
,
3
)
=
a
;
ipFaceCoords
[
0
](
1
,
3
)
=
-
a
;
ipFaceCoords
[
3
](
1
,
3
)
=
1
;
ipFaceCoords
[
0
](
2
,
3
)
=
a
;
ipFaceCoords
[
3
](
2
,
3
)
=
-
a
;
ipFaceCoords
[
1
](
0
,
0
)
=
1
;
ipFaceCoords
[
4
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
1
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
4
](
2
,
0
)
=
-
1
;
//
ipFaceCoords
[
1
](
0
,
1
)
=
1
;
ipFaceCoords
[
4
](
0
,
1
)
=
a
;
ipFaceCoords
[
1
](
1
,
1
)
=
a
;
ipFaceCoords
[
4
](
1
,
1
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
4
](
2
,
1
)
=
-
1
;
//
ipFaceCoords
[
1
](
0
,
2
)
=
1
;
ipFaceCoords
[
4
](
0
,
2
)
=
a
;
ipFaceCoords
[
1
](
1
,
2
)
=
a
;
ipFaceCoords
[
4
](
1
,
2
)
=
a
;
ipFaceCoords
[
1
](
2
,
2
)
=
a
;
ipFaceCoords
[
4
](
2
,
2
)
=
-
1
;
//
ipFaceCoords
[
1
](
0
,
3
)
=
1
;
ipFaceCoords
[
4
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
1
](
1
,
3
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
3
)
=
a
;
ipFaceCoords
[
1
](
2
,
3
)
=
a
;
ipFaceCoords
[
4
](
2
,
3
)
=
-
1
;
ipFaceCoords
[
2
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
5
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
2
](
1
,
0
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
2
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
5
](
2
,
0
)
=
1
;
//
ipFaceCoords
[
2
](
0
,
1
)
=
-
a
;
ipFaceCoords
[
5
](
0
,
1
)
=
a
;
ipFaceCoords
[
2
](
1
,
1
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
1
)
=
-
a
;
ipFaceCoords
[
2
](
2
,
1
)
=
a
;
ipFaceCoords
[
5
](
2
,
1
)
=
1
;
//
ipFaceCoords
[
2
](
0
,
2
)
=
a
;
ipFaceCoords
[
5
](
0
,
2
)
=
a
;
ipFaceCoords
[
2
](
1
,
2
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
2
)
=
a
;
ipFaceCoords
[
2
](
2
,
2
)
=
a
;
ipFaceCoords
[
5
](
2
,
2
)
=
1
;
//
ipFaceCoords
[
2
](
0
,
3
)
=
a
;
ipFaceCoords
[
5
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
2
](
1
,
3
)
=
-
1
;
ipFaceCoords
[
5
](
1
,
3
)
=
a
;
ipFaceCoords
[
2
](
2
,
3
)
=
-
a
;
ipFaceCoords
[
5
](
2
,
3
)
=
1
;
// Integration points for all faces ignoring the
// redundant dim
ipFace2DCoords
(
0
,
0
)
=
-
a
;
ipFace2DCoords
(
0
,
2
)
=
a
;
ipFace2DCoords
(
1
,
0
)
=
-
a
;
ipFace2DCoords
(
1
,
2
)
=
a
;
//
ipFace2DCoords
(
0
,
1
)
=
a
;
ipFace2DCoords
(
0
,
3
)
=
-
a
;
ipFace2DCoords
(
1
,
1
)
=
-
a
;
ipFace2DCoords
(
1
,
3
)
=
a
;
// Integration point weights
ipWeights
=
4.0
/
3.0
;
// Face integration point weights
ipFaceWeights
=
1.0
;
break
;
}
// 27 point quadrature
case
GAUSS3
:
{
// Set number of IPs and face IPs
numIPs
=
27
;
numFaceIPs
=
9
;
// Size matrices and vectors thereof accordingly
ipCoords
.
resize
(
3
,
numIPs
);
ipFaceCoords
.
assign
(
6
,
DENS_MAT
(
3
,
numFaceIPs
));
ipFace2DCoords
.
resize
(
2
,
numFaceIPs
);
ipWeights
.
reset
(
numIPs
);
ipFaceWeights
.
reset
(
numFaceIPs
);
double
a
=
sqrt
(
3.0
/
5.0
);
// Matrix of integration point locations & follows local
// conn
ipCoords
(
0
,
0
)
=
-
a
;
ipCoords
(
0
,
16
)
=
-
a
;
ipCoords
(
1
,
0
)
=
-
a
;
ipCoords
(
1
,
16
)
=
-
a
;
ipCoords
(
2
,
0
)
=
-
a
;
ipCoords
(
2
,
16
)
=
0
;
//
ipCoords
(
0
,
1
)
=
a
;
ipCoords
(
0
,
17
)
=
-
a
;
ipCoords
(
1
,
1
)
=
-
a
;
ipCoords
(
1
,
17
)
=
a
;
ipCoords
(
2
,
1
)
=
-
a
;
ipCoords
(
2
,
17
)
=
0
;
//
ipCoords
(
0
,
2
)
=
a
;
ipCoords
(
0
,
18
)
=
a
;
ipCoords
(
1
,
2
)
=
a
;
ipCoords
(
1
,
18
)
=
-
a
;
ipCoords
(
2
,
2
)
=
-
a
;
ipCoords
(
2
,
18
)
=
0
;
//
ipCoords
(
0
,
3
)
=
-
a
;
ipCoords
(
0
,
19
)
=
a
;
ipCoords
(
1
,
3
)
=
a
;
ipCoords
(
1
,
19
)
=
a
;
ipCoords
(
2
,
3
)
=
-
a
;
ipCoords
(
2
,
19
)
=
0
;
ipCoords
(
0
,
4
)
=
-
a
;
ipCoords
(
0
,
20
)
=
0
;
ipCoords
(
1
,
4
)
=
-
a
;
ipCoords
(
1
,
20
)
=
0
;
ipCoords
(
2
,
4
)
=
a
;
ipCoords
(
2
,
20
)
=
-
a
;
ipCoords
(
0
,
5
)
=
a
;
ipCoords
(
0
,
21
)
=
0
;
ipCoords
(
1
,
5
)
=
-
a
;
ipCoords
(
1
,
21
)
=
0
;
ipCoords
(
2
,
5
)
=
a
;
ipCoords
(
2
,
21
)
=
a
;
ipCoords
(
0
,
6
)
=
a
;
ipCoords
(
0
,
22
)
=
0
;
ipCoords
(
1
,
6
)
=
a
;
ipCoords
(
1
,
22
)
=
-
a
;
ipCoords
(
2
,
6
)
=
a
;
ipCoords
(
2
,
22
)
=
0
;
ipCoords
(
0
,
7
)
=
-
a
;
ipCoords
(
0
,
23
)
=
0
;
ipCoords
(
1
,
7
)
=
a
;
ipCoords
(
1
,
23
)
=
a
;
ipCoords
(
2
,
7
)
=
a
;
ipCoords
(
2
,
23
)
=
0
;
ipCoords
(
0
,
8
)
=
0
;
ipCoords
(
0
,
24
)
=
-
a
;
ipCoords
(
1
,
8
)
=
-
a
;
ipCoords
(
1
,
24
)
=
0
;
ipCoords
(
2
,
8
)
=
-
a
;
ipCoords
(
2
,
24
)
=
0
;
ipCoords
(
0
,
9
)
=
0
;
ipCoords
(
0
,
25
)
=
a
;
ipCoords
(
1
,
9
)
=
-
a
;
ipCoords
(
1
,
25
)
=
0
;
ipCoords
(
2
,
9
)
=
a
;
ipCoords
(
2
,
25
)
=
0
;
ipCoords
(
0
,
10
)
=
0
;
ipCoords
(
0
,
26
)
=
0
;
ipCoords
(
1
,
10
)
=
a
;
ipCoords
(
1
,
26
)
=
0
;
ipCoords
(
2
,
10
)
=
-
a
;
ipCoords
(
2
,
26
)
=
0
;
ipCoords
(
0
,
11
)
=
0
;
ipCoords
(
1
,
11
)
=
a
;
ipCoords
(
2
,
11
)
=
a
;
ipCoords
(
0
,
12
)
=
-
a
;
ipCoords
(
1
,
12
)
=
0
;
ipCoords
(
2
,
12
)
=
-
a
;
ipCoords
(
0
,
13
)
=
-
a
;
ipCoords
(
1
,
13
)
=
0
;
ipCoords
(
2
,
13
)
=
a
;
ipCoords
(
0
,
14
)
=
a
;
ipCoords
(
1
,
14
)
=
0
;
ipCoords
(
2
,
14
)
=
-
a
;
ipCoords
(
0
,
15
)
=
a
;
ipCoords
(
1
,
15
)
=
0
;
ipCoords
(
2
,
15
)
=
a
;
// Integration points by face
ipFaceCoords
[
0
](
0
,
0
)
=
-
1
;
ipFaceCoords
[
0
](
0
,
5
)
=
-
1
;
ipFaceCoords
[
0
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
0
](
1
,
5
)
=
0
;
ipFaceCoords
[
0
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
0
](
2
,
5
)
=
a
;
//
ipFaceCoords
[
0
](
0
,
1
)
=
-
1
;
ipFaceCoords
[
0
](
0
,
6
)
=
-
1
;
ipFaceCoords
[
0
](
1
,
1
)
=
a
;
ipFaceCoords
[
0
](
1
,
6
)
=
-
a
;
ipFaceCoords
[
0
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
0
](
2
,
6
)
=
0
;
//
ipFaceCoords
[
0
](
0
,
2
)
=
-
1
;
ipFaceCoords
[
0
](
0
,
7
)
=
-
1
;
ipFaceCoords
[
0
](
1
,
2
)
=
a
;
ipFaceCoords
[
0
](
1
,
7
)
=
a
;
ipFaceCoords
[
0
](
2
,
2
)
=
a
;
ipFaceCoords
[
0
](
2
,
7
)
=
0
;
//
ipFaceCoords
[
0
](
0
,
3
)
=
-
1
;
ipFaceCoords
[
0
](
0
,
8
)
=
-
1
;
ipFaceCoords
[
0
](
1
,
3
)
=
-
a
;
ipFaceCoords
[
0
](
1
,
8
)
=
0
;
ipFaceCoords
[
0
](
2
,
3
)
=
a
;
ipFaceCoords
[
0
](
2
,
8
)
=
0
;
//
ipFaceCoords
[
0
](
0
,
4
)
=
-
1
;
ipFaceCoords
[
0
](
1
,
4
)
=
0
;
ipFaceCoords
[
0
](
2
,
4
)
=
-
a
;
ipFaceCoords
[
1
](
0
,
0
)
=
1
;
ipFaceCoords
[
1
](
0
,
5
)
=
1
;
ipFaceCoords
[
1
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
1
](
1
,
5
)
=
0
;
ipFaceCoords
[
1
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
5
)
=
a
;
//
ipFaceCoords
[
1
](
0
,
1
)
=
1
;
ipFaceCoords
[
1
](
0
,
6
)
=
1
;
ipFaceCoords
[
1
](
1
,
1
)
=
a
;
ipFaceCoords
[
1
](
1
,
6
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
1
)
=
-
a
;
ipFaceCoords
[
1
](
2
,
6
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
2
)
=
1
;
ipFaceCoords
[
1
](
0
,
7
)
=
1
;
ipFaceCoords
[
1
](
1
,
2
)
=
a
;
ipFaceCoords
[
1
](
1
,
7
)
=
a
;
ipFaceCoords
[
1
](
2
,
2
)
=
a
;
ipFaceCoords
[
1
](
2
,
7
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
3
)
=
1
;
ipFaceCoords
[
1
](
0
,
8
)
=
1
;
ipFaceCoords
[
1
](
1
,
3
)
=
-
a
;
ipFaceCoords
[
1
](
1
,
8
)
=
0
;
ipFaceCoords
[
1
](
2
,
3
)
=
a
;
ipFaceCoords
[
1
](
2
,
8
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
4
)
=
1
;
ipFaceCoords
[
1
](
1
,
4
)
=
0
;
ipFaceCoords
[
1
](
2
,
4
)
=
-
a
;
ipFaceCoords
[
2
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
2
](
0
,
5
)
=
0
;
ipFaceCoords
[
2
](
1
,
0
)
=
-
1
;
ipFaceCoords
[
2
](
1
,
5
)
=
-
1
;
ipFaceCoords
[
2
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
2
](
2
,
5
)
=
a
;
//
ipFaceCoords
[
2
](
0
,
1
)
=
-
a
;
ipFaceCoords
[
2
](
0
,
6
)
=
-
a
;
ipFaceCoords
[
2
](
1
,
1
)
=
-
1
;
ipFaceCoords
[
2
](
1
,
6
)
=
-
1
;
ipFaceCoords
[
2
](
2
,
1
)
=
a
;
ipFaceCoords
[
2
](
2
,
6
)
=
0
;
//
ipFaceCoords
[
2
](
0
,
2
)
=
a
;
ipFaceCoords
[
2
](
0
,
7
)
=
a
;
ipFaceCoords
[
2
](
1
,
2
)
=
-
1
;
ipFaceCoords
[
2
](
1
,
7
)
=
-
1
;
ipFaceCoords
[
2
](
2
,
2
)
=
a
;
ipFaceCoords
[
2
](
2
,
7
)
=
0
;
//
ipFaceCoords
[
2
](
0
,
3
)
=
a
;
ipFaceCoords
[
2
](
0
,
8
)
=
0
;
ipFaceCoords
[
2
](
1
,
3
)
=
-
1
;
ipFaceCoords
[
2
](
1
,
8
)
=
-
1
;
ipFaceCoords
[
2
](
2
,
3
)
=
-
a
;
ipFaceCoords
[
2
](
2
,
8
)
=
0
;
//
ipFaceCoords
[
2
](
0
,
4
)
=
0
;
ipFaceCoords
[
2
](
1
,
4
)
=
-
1
;
ipFaceCoords
[
2
](
2
,
4
)
=
-
a
;
ipFaceCoords
[
3
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
3
](
0
,
5
)
=
0
;
ipFaceCoords
[
3
](
1
,
0
)
=
1
;
ipFaceCoords
[
3
](
1
,
5
)
=
1
;
ipFaceCoords
[
3
](
2
,
0
)
=
-
a
;
ipFaceCoords
[
3
](
2
,
5
)
=
a
;
//
ipFaceCoords
[
3
](
0
,
1
)
=
-
a
;
ipFaceCoords
[
3
](
0
,
6
)
=
-
a
;
ipFaceCoords
[
3
](
1
,
1
)
=
1
;
ipFaceCoords
[
3
](
1
,
6
)
=
1
;
ipFaceCoords
[
3
](
2
,
1
)
=
a
;
ipFaceCoords
[
3
](
2
,
6
)
=
0
;
//
ipFaceCoords
[
3
](
0
,
2
)
=
a
;
ipFaceCoords
[
3
](
0
,
7
)
=
a
;
ipFaceCoords
[
3
](
1
,
2
)
=
1
;
ipFaceCoords
[
3
](
1
,
7
)
=
1
;
ipFaceCoords
[
3
](
2
,
2
)
=
a
;
ipFaceCoords
[
3
](
2
,
7
)
=
0
;
//
ipFaceCoords
[
3
](
0
,
3
)
=
a
;
ipFaceCoords
[
3
](
0
,
8
)
=
0
;
ipFaceCoords
[
3
](
1
,
3
)
=
1
;
ipFaceCoords
[
3
](
1
,
8
)
=
1
;
ipFaceCoords
[
3
](
2
,
3
)
=
-
a
;
ipFaceCoords
[
3
](
2
,
8
)
=
0
;
//
ipFaceCoords
[
3
](
0
,
4
)
=
0
;
ipFaceCoords
[
3
](
1
,
4
)
=
1
;
ipFaceCoords
[
3
](
2
,
4
)
=
-
a
;
ipFaceCoords
[
4
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
4
](
0
,
5
)
=
0
;
ipFaceCoords
[
4
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
5
)
=
a
;
ipFaceCoords
[
4
](
2
,
0
)
=
-
1
;
ipFaceCoords
[
4
](
2
,
5
)
=
-
1
;
//
ipFaceCoords
[
4
](
0
,
1
)
=
a
;
ipFaceCoords
[
4
](
0
,
6
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
1
)
=
-
a
;
ipFaceCoords
[
4
](
1
,
6
)
=
0
;
ipFaceCoords
[
4
](
2
,
1
)
=
-
1
;
ipFaceCoords
[
4
](
2
,
6
)
=
-
1
;
//
ipFaceCoords
[
4
](
0
,
2
)
=
a
;
ipFaceCoords
[
4
](
0
,
7
)
=
a
;
ipFaceCoords
[
4
](
1
,
2
)
=
a
;
ipFaceCoords
[
4
](
1
,
7
)
=
0
;
ipFaceCoords
[
4
](
2
,
2
)
=
-
1
;
ipFaceCoords
[
4
](
2
,
7
)
=
-
1
;
//
ipFaceCoords
[
4
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
4
](
0
,
8
)
=
0
;
ipFaceCoords
[
4
](
1
,
3
)
=
a
;
ipFaceCoords
[
4
](
1
,
8
)
=
0
;
ipFaceCoords
[
4
](
2
,
3
)
=
-
1
;
ipFaceCoords
[
4
](
2
,
8
)
=
-
1
;
//
ipFaceCoords
[
4
](
0
,
4
)
=
0
;
ipFaceCoords
[
4
](
1
,
4
)
=
-
a
;
ipFaceCoords
[
4
](
2
,
4
)
=
-
1
;
ipFaceCoords
[
5
](
0
,
0
)
=
-
a
;
ipFaceCoords
[
5
](
0
,
5
)
=
0
;
ipFaceCoords
[
5
](
1
,
0
)
=
-
a
;
ipFaceCoords
[
5
](
1
,
5
)
=
a
;
ipFaceCoords
[
5
](
2
,
0
)
=
1
;
ipFaceCoords
[
5
](
2
,
5
)
=
1
;
//
ipFaceCoords
[
5
](
0
,
1
)
=
a
;
ipFaceCoords
[
5
](
0
,
6
)
=
-
a
;
ipFaceCoords
[
5
](
1
,
1
)
=
-
a
;
ipFaceCoords
[
5
](
1
,
6
)
=
0
;
ipFaceCoords
[
5
](
2
,
1
)
=
1
;
ipFaceCoords
[
5
](
2
,
6
)
=
1
;
//
ipFaceCoords
[
5
](
0
,
2
)
=
a
;
ipFaceCoords
[
5
](
0
,
7
)
=
a
;
ipFaceCoords
[
5
](
1
,
2
)
=
a
;
ipFaceCoords
[
5
](
1
,
7
)
=
0
;
ipFaceCoords
[
5
](
2
,
2
)
=
1
;
ipFaceCoords
[
5
](
2
,
7
)
=
1
;
//
ipFaceCoords
[
5
](
0
,
3
)
=
-
a
;
ipFaceCoords
[
5
](
0
,
8
)
=
0
;
ipFaceCoords
[
5
](
1
,
3
)
=
a
;
ipFaceCoords
[
5
](
1
,
8
)
=
0
;
ipFaceCoords
[
5
](
2
,
3
)
=
1
;
ipFaceCoords
[
5
](
2
,
8
)
=
1
;
//
ipFaceCoords
[
5
](
0
,
4
)
=
0
;
ipFaceCoords
[
5
](
1
,
4
)
=
-
a
;
ipFaceCoords
[
5
](
2
,
4
)
=
1
;
// Integration points for all faces ignoring the
// redundant dim
ipFace2DCoords
(
0
,
0
)
=
-
a
;
ipFace2DCoords
(
0
,
5
)
=
0
;
ipFace2DCoords
(
1
,
0
)
=
-
a
;
ipFace2DCoords
(
1
,
5
)
=
a
;
//
ipFace2DCoords
(
0
,
1
)
=
a
;
ipFace2DCoords
(
0
,
6
)
=
-
a
;
ipFace2DCoords
(
1
,
1
)
=
-
a
;
ipFace2DCoords
(
1
,
6
)
=
0
;
//
ipFace2DCoords
(
0
,
2
)
=
a
;
ipFace2DCoords
(
0
,
7
)
=
a
;
ipFace2DCoords
(
1
,
2
)
=
a
;
ipFace2DCoords
(
1
,
7
)
=
0
;
//
ipFace2DCoords
(
0
,
3
)
=
-
a
;
ipFace2DCoords
(
0
,
8
)
=
0
;
ipFace2DCoords
(
1
,
3
)
=
a
;
ipFace2DCoords
(
1
,
8
)
=
0
;
//
ipFace2DCoords
(
0
,
4
)
=
0
;
ipFace2DCoords
(
1
,
4
)
=
-
a
;
// Integration point weights
for
(
int
i
=
0
;
i
<
numIPs
;
++
i
)
{
if
(
i
<
8
)
ipWeights
[
i
]
=
125.0
/
729.0
;
else
if
(
i
<
20
)
ipWeights
[
i
]
=
200.0
/
729.0
;
else
if
(
i
<
26
)
ipWeights
[
i
]
=
320.0
/
729.0
;
else
ipWeights
[
i
]
=
512.0
/
729.0
;
}
// Face integration point weights
for
(
int
i
=
0
;
i
<
numFaceIPs
;
++
i
)
{
if
(
i
<
4
)
ipFaceWeights
[
i
]
=
25.0
/
81.0
;
else
if
(
i
<
8
)
ipFaceWeights
[
i
]
=
40.0
/
81.0
;
else
ipFaceWeights
[
i
]
=
64.0
/
81.0
;
}
break
;
}
// Error
default:
{
throw
ATC_Error
(
"Unrecognized quadrature type "
"for element type HEXA."
);
}
}
}
else
if
(
geo
==
TETRA
)
{
switch
(
quad
)
{
// 4 point quadratures
case
NODAL
:
case
GAUSS2
:
{
// Set number of IPs and face IPs
numIPs
=
4
;
numFaceIPs
=
3
;
// Size matrices and vectors thereof accordingly
ipCoords
.
resize
(
4
,
numIPs
);
ipFaceCoords
.
assign
(
4
,
DENS_MAT
(
3
,
numFaceIPs
));
ipFace2DCoords
.
resize
(
3
,
numFaceIPs
);
ipWeights
.
reset
(
numIPs
);
ipFaceWeights
.
reset
(
numFaceIPs
);
double
v1
,
v2
,
a1
,
a2
;
v1
=
0.585410196624969
;
v2
=
0.138196601125011
;
a1
=
2.0
/
3.0
;
a2
=
1.0
/
6.0
;
if
(
quad
==
NODAL
)
{
// Nodal quadrature
v1
=
1
;
v2
=
0
;
a1
=
1
;
a2
=
0
;
}
// Integration point coordinates
ipCoords
.
resize
(
4
,
numIPs
);
ipCoords
(
0
,
0
)
=
v2
;
ipCoords
(
0
,
2
)
=
v2
;
ipCoords
(
1
,
0
)
=
v2
;
ipCoords
(
1
,
2
)
=
v1
;
ipCoords
(
2
,
0
)
=
v2
;
ipCoords
(
2
,
2
)
=
v2
;
ipCoords
(
3
,
0
)
=
v1
;
ipCoords
(
3
,
2
)
=
v2
;
//
ipCoords
(
0
,
1
)
=
v1
;
ipCoords
(
0
,
3
)
=
v2
;
ipCoords
(
1
,
1
)
=
v2
;
ipCoords
(
1
,
3
)
=
v2
;
ipCoords
(
2
,
1
)
=
v2
;
ipCoords
(
2
,
3
)
=
v1
;
ipCoords
(
3
,
1
)
=
v2
;
ipCoords
(
3
,
3
)
=
v2
;
// Integration points by face
// ...face 0 ...face 2
ipFaceCoords
[
0
](
0
,
0
)
=
a1
;
ipFaceCoords
[
2
](
0
,
0
)
=
a2
;
ipFaceCoords
[
0
](
1
,
0
)
=
a2
;
ipFaceCoords
[
2
](
1
,
0
)
=
0
;
ipFaceCoords
[
0
](
2
,
0
)
=
a2
;
ipFaceCoords
[
2
](
2
,
0
)
=
a2
;
//
ipFaceCoords
[
0
](
0
,
1
)
=
a2
;
ipFaceCoords
[
2
](
0
,
1
)
=
a1
;
ipFaceCoords
[
0
](
1
,
1
)
=
a1
;
ipFaceCoords
[
2
](
1
,
1
)
=
0
;
ipFaceCoords
[
0
](
2
,
1
)
=
a2
;
ipFaceCoords
[
2
](
2
,
1
)
=
a2
;
//
ipFaceCoords
[
0
](
0
,
2
)
=
a2
;
ipFaceCoords
[
2
](
0
,
2
)
=
a2
;
ipFaceCoords
[
0
](
1
,
2
)
=
a2
;
ipFaceCoords
[
2
](
1
,
2
)
=
0
;
ipFaceCoords
[
0
](
2
,
2
)
=
a1
;
ipFaceCoords
[
2
](
2
,
2
)
=
a1
;
// ...face 1 ...face 3
ipFaceCoords
[
1
](
0
,
0
)
=
0
;
ipFaceCoords
[
3
](
0
,
0
)
=
a2
;
ipFaceCoords
[
1
](
1
,
0
)
=
a1
;
ipFaceCoords
[
3
](
1
,
0
)
=
a2
;
ipFaceCoords
[
1
](
2
,
0
)
=
a2
;
ipFaceCoords
[
3
](
2
,
0
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
1
)
=
0
;
ipFaceCoords
[
3
](
0
,
1
)
=
a2
;
ipFaceCoords
[
1
](
1
,
1
)
=
a2
;
ipFaceCoords
[
3
](
1
,
1
)
=
a1
;
ipFaceCoords
[
1
](
2
,
1
)
=
a2
;
ipFaceCoords
[
3
](
2
,
1
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
2
)
=
0
;
ipFaceCoords
[
3
](
0
,
2
)
=
a1
;
ipFaceCoords
[
1
](
1
,
2
)
=
a2
;
ipFaceCoords
[
3
](
1
,
2
)
=
a2
;
ipFaceCoords
[
1
](
2
,
2
)
=
a1
;
ipFaceCoords
[
3
](
2
,
2
)
=
0
;
// 2D integration points for faces
ipFace2DCoords
(
0
,
0
)
=
a2
;
ipFace2DCoords
(
0
,
2
)
=
a2
;
ipFace2DCoords
(
1
,
0
)
=
a2
;
ipFace2DCoords
(
1
,
2
)
=
a1
;
ipFace2DCoords
(
2
,
0
)
=
a1
;
ipFace2DCoords
(
2
,
2
)
=
a2
;
//
ipFace2DCoords
(
0
,
1
)
=
a1
;
ipFace2DCoords
(
1
,
1
)
=
a2
;
ipFace2DCoords
(
2
,
1
)
=
a2
;
// Integration point weights
ipWeights
=
(
1.0
/
6.0
)
/
numIPs
;
// Integration point face weights
ipFaceWeights
=
(
1.0
/
2.0
)
/
numFaceIPs
;
break
;
}
case
GAUSS3
:
{
// Set number of IPs and face IPs
numIPs
=
5
;
numFaceIPs
=
4
;
// Size matrices and vectors thereof accordingly
ipCoords
.
resize
(
4
,
numIPs
);
ipFaceCoords
.
assign
(
4
,
DENS_MAT
(
3
,
numFaceIPs
));
ipFace2DCoords
.
resize
(
3
,
numFaceIPs
);
ipWeights
.
reset
(
numIPs
);
ipFaceWeights
.
reset
(
numFaceIPs
);
double
v1
,
v2
,
v3
,
a1
,
a2
,
a3
;
/* These weights for calculating Gaussian Quadrature
* points are taken from the paper "Integration Points
* for Triangles and Tetrahedra Obtained from the
* Gaussian Quadrature Points for a Line" by K. S. Sunder
* and R. A. Cookson, Computers and Structures, Vol 21,
* No. 5, 1985. */
v1
=
1.0
/
4.0
;
v2
=
1.0
/
2.0
;
v3
=
1.0
/
6.0
;
a1
=
1.0
/
3.0
;
a2
=
3.0
/
5.0
;
a3
=
1.0
/
5.0
;
// Integration point coordinates
ipCoords
.
resize
(
4
,
numIPs
);
ipCoords
(
0
,
0
)
=
v1
;
ipCoords
(
1
,
0
)
=
v1
;
ipCoords
(
2
,
0
)
=
v1
;
ipCoords
(
3
,
0
)
=
v1
;
//
ipCoords
(
0
,
1
)
=
v3
;
ipCoords
(
0
,
3
)
=
v3
;
ipCoords
(
1
,
1
)
=
v3
;
ipCoords
(
1
,
3
)
=
v2
;
ipCoords
(
2
,
1
)
=
v3
;
ipCoords
(
2
,
3
)
=
v3
;
ipCoords
(
3
,
1
)
=
v2
;
ipCoords
(
3
,
3
)
=
v3
;
//
ipCoords
(
0
,
2
)
=
v2
;
ipCoords
(
0
,
4
)
=
v3
;
ipCoords
(
1
,
2
)
=
v3
;
ipCoords
(
1
,
4
)
=
v3
;
ipCoords
(
2
,
2
)
=
v3
;
ipCoords
(
2
,
4
)
=
v2
;
ipCoords
(
3
,
2
)
=
v3
;
ipCoords
(
3
,
4
)
=
v3
;
// Integration points by face
// ...face 0 ...face 2
ipFaceCoords
[
0
](
0
,
0
)
=
a1
;
ipFaceCoords
[
2
](
0
,
0
)
=
a1
;
ipFaceCoords
[
0
](
1
,
0
)
=
a1
;
ipFaceCoords
[
2
](
1
,
0
)
=
0
;
ipFaceCoords
[
0
](
2
,
0
)
=
a1
;
ipFaceCoords
[
2
](
2
,
0
)
=
a1
;
//
ipFaceCoords
[
0
](
0
,
1
)
=
a2
;
ipFaceCoords
[
2
](
0
,
1
)
=
a3
;
ipFaceCoords
[
0
](
1
,
1
)
=
a3
;
ipFaceCoords
[
2
](
1
,
1
)
=
0
;
ipFaceCoords
[
0
](
2
,
1
)
=
a3
;
ipFaceCoords
[
2
](
2
,
1
)
=
a3
;
//
ipFaceCoords
[
0
](
0
,
2
)
=
a3
;
ipFaceCoords
[
2
](
0
,
2
)
=
a2
;
ipFaceCoords
[
0
](
1
,
2
)
=
a2
;
ipFaceCoords
[
2
](
1
,
2
)
=
0
;
ipFaceCoords
[
0
](
2
,
2
)
=
a3
;
ipFaceCoords
[
2
](
2
,
2
)
=
a3
;
//
ipFaceCoords
[
0
](
0
,
3
)
=
a3
;
ipFaceCoords
[
2
](
0
,
3
)
=
a3
;
ipFaceCoords
[
0
](
1
,
3
)
=
a3
;
ipFaceCoords
[
2
](
1
,
3
)
=
0
;
ipFaceCoords
[
0
](
2
,
3
)
=
a2
;
ipFaceCoords
[
2
](
2
,
3
)
=
a2
;
// ...face 1 ...face 3
ipFaceCoords
[
1
](
0
,
0
)
=
0
;
ipFaceCoords
[
3
](
0
,
0
)
=
a1
;
ipFaceCoords
[
1
](
1
,
0
)
=
a1
;
ipFaceCoords
[
3
](
1
,
0
)
=
a1
;
ipFaceCoords
[
1
](
2
,
0
)
=
a1
;
ipFaceCoords
[
3
](
2
,
0
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
1
)
=
0
;
ipFaceCoords
[
3
](
0
,
1
)
=
a3
;
ipFaceCoords
[
1
](
1
,
1
)
=
a2
;
ipFaceCoords
[
3
](
1
,
1
)
=
a3
;
ipFaceCoords
[
1
](
2
,
1
)
=
a3
;
ipFaceCoords
[
3
](
2
,
1
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
2
)
=
0
;
ipFaceCoords
[
3
](
0
,
2
)
=
a3
;
ipFaceCoords
[
1
](
1
,
2
)
=
a3
;
ipFaceCoords
[
3
](
1
,
2
)
=
a2
;
ipFaceCoords
[
1
](
2
,
2
)
=
a3
;
ipFaceCoords
[
3
](
2
,
2
)
=
0
;
//
ipFaceCoords
[
1
](
0
,
3
)
=
0
;
ipFaceCoords
[
3
](
0
,
3
)
=
a2
;
ipFaceCoords
[
1
](
1
,
3
)
=
a3
;
ipFaceCoords
[
3
](
1
,
3
)
=
a3
;
ipFaceCoords
[
1
](
2
,
3
)
=
a2
;
ipFaceCoords
[
3
](
2
,
3
)
=
0
;
// 2D integration points for faces
//
ipFace2DCoords
(
0
,
0
)
=
a1
;
ipFace2DCoords
(
0
,
2
)
=
a2
;
ipFace2DCoords
(
1
,
0
)
=
a1
;
ipFace2DCoords
(
1
,
2
)
=
a3
;
ipFace2DCoords
(
2
,
0
)
=
a1
;
ipFace2DCoords
(
2
,
2
)
=
a3
;
//
ipFace2DCoords
(
0
,
1
)
=
a3
;
ipFace2DCoords
(
0
,
3
)
=
a3
;
ipFace2DCoords
(
1
,
1
)
=
a3
;
ipFace2DCoords
(
1
,
3
)
=
a2
;
ipFace2DCoords
(
2
,
1
)
=
a2
;
ipFace2DCoords
(
2
,
3
)
=
a3
;
for
(
int
i
=
0
;
i
<
numIPs
;
++
i
)
{
if
(
i
<
1
)
ipWeights
[
i
]
=
-
4.0
/
5.0
;
else
ipWeights
[
i
]
=
9.0
/
20.0
;
}
// Face integration point weights
for
(
int
i
=
0
;
i
<
numFaceIPs
;
++
i
)
{
if
(
i
<
1
)
ipFaceWeights
[
i
]
=
-
9.0
/
16.0
;
else
ipFaceWeights
[
i
]
=
25.0
/
48.0
;
}
break
;
}
// Error
default:
{
throw
ATC_Error
(
"Unrecognized quadrature type "
"for element type TETRA."
);
}
}
}
}
};
};
// namespace ATC
#endif
// FE_QUADRATURE_H
Event Timeline
Log In to Comment