ESScript
Revision_4488
Main Page
Namespaces
Classes
Files
File List
File Members
finley
src
BasisFunctions.h
Go to the documentation of this file.
1
2
/*****************************************************************************
3
*
4
* Copyright (c) 2003-2013 by University of Queensland
5
* http://www.uq.edu.au
6
*
7
* Primary Business: Queensland, Australia
8
* Licensed under the Open Software License version 3.0
9
* http://www.opensource.org/licenses/osl-3.0.php
10
*
11
* Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12
* Development since 2012 by School of Earth Sciences
13
*
14
*****************************************************************************/
15
16
17
/************************************************************************************/
18
19
/* Finley: Reference elements */
20
21
/************************************************************************************/
22
23
#ifndef INC_FINLEY_REFERENCEELEMENTS
24
#define INC_FINLEY_REFERENCEELEMENTS
25
26
27
/************************************************************************************/
28
29
#include "
Finley.h
"
30
#include "
ShapeFunctions.h
"
31
#include "
Quadrature.h
"
32
33
/************************************************************************************/
34
35
/* The ids of the allowed reference elements: */
36
37
#define MAX_numNodes 64
38
39
typedef
enum
{
40
Point1
,
41
Line2
,
42
Line3
,
43
Line4
,
44
Tri3
,
45
Tri6
,
46
Tri9
,
47
Tri10
,
48
Rec4
,
49
Rec8
,
50
Rec9
,
51
Rec12
,
52
Rec16
,
53
Tet4
,
54
Tet10
,
55
Tet16
,
56
Hex8
,
57
Hex20
,
58
Hex27
,
59
Hex32
,
60
Line2Face
,
61
Line3Face
,
62
Line4Face
,
63
Tri3Face
,
64
Tri6Face
,
65
Tri9Face
,
66
Tri10Face
,
67
Rec4Face
,
68
Rec8Face
,
69
Rec9Face
,
70
Rec12Face
,
71
Rec16Face
,
72
Tet4Face
,
73
Tet10Face
,
74
Tet16Face
,
75
Hex8Face
,
76
Hex20Face
,
77
Hex27Face
,
78
Hex32Face
,
79
Point1_Contact
,
80
Line2_Contact
,
81
Line3_Contact
,
82
Line4_Contact
,
83
Tri3_Contact
,
84
Tri6_Contact
,
85
Tri9_Contact
,
86
Tri10_Contact
,
87
Rec4_Contact
,
88
Rec8_Contact
,
89
Rec9_Contact
,
90
Rec12_Contact
,
91
Rec16_Contact
,
92
Line2Face_Contact
,
93
Line3Face_Contact
,
94
Line4Face_Contact
,
95
Tri3Face_Contact
,
96
Tri6Face_Contact
,
97
Tri9Face_Contact
,
98
Tri10Face_Contact
,
99
Rec4Face_Contact
,
100
Rec8Face_Contact
,
101
Rec9Face_Contact
,
102
Rec12Face_Contact
,
103
Rec16Face_Contact
,
104
Tet4Face_Contact
,
105
Tet10Face_Contact
,
106
Tet16Face_Contact
,
107
Hex8Face_Contact
,
108
Hex20Face_Contact
,
109
Hex27Face_Contact
,
110
Hex32Face_Contact
,
111
Line3Macro
,
112
Tri6Macro
,
113
Rec9Macro
,
114
Tet10Macro
,
115
Hex27Macro
,
116
117
NoType
/* marks end of list */
118
}
ElementTypeId
;
119
120
/************************************************************************************/
121
122
/* this struct holds the definition of the reference element: */
123
124
typedef
struct
Finley_ReferenceElementInfo
{
125
ElementTypeId
TypeId
;
/* the id */
126
char
*
Name
;
/* the name in text form e.g. Line1,Rec12,... */
127
dim_t
numLocalDim
;
/* local dimension of the element */
128
dim_t
numDim
;
/* dimension of the element */
129
dim_t
numNodes
;
/* number of nodes defining the element*/
130
dim_t
numShapes
;
/* number of shape functions, typically = numNodes*/
131
dim_t
numOrder
;
/* order of the shape functions */
132
dim_t
numVertices
;
/* number of vertices of the element */
133
ElementTypeId
LinearTypeId
;
/* id of the linear version of the element */
134
index_t
linearNodes
[
MAX_numNodes
];
/* gives the list of nodes defining the linear or macro element, typically it is linearNodes[i]=i */
135
Finley_Shape_Function*
getValues
;
/* function to evaluate the shape functions at a set of points */
136
Finley_Quad_getNodes
*
getQuadNodes
;
/* function to set the quadrature points */
137
Finley_Quad_getNumNodes
*
getNumQuadNodes
;
/* function selects the number of quadrature nodes for a given accuracy order */
138
139
/********************************************************************************************************************************************************* */
140
dim_t
numRelevantGeoNodes
;
/* number of nodes used to describe the geometry of the geometrically relevant part of the element
141
typically this is numNodes but for 'Face' elements where the quadrature points are defined on face of the element
142
this is the number of nodes on the particular face. */
143
index_t
relevantGeoNodes
[
MAX_numNodes
];
/* list to gather the geometrically relevant nodes (length used is numRelevantGeoNodes)
144
this list is used for the VTK interface */
145
146
dim_t
numNodesOnFace
;
/* if the element is allowed as a face element, numNodesOnFace defines the number of nodes defining the face */
147
/* the following lists are only used for face elements defined by numNodesOnFace>0 */
148
index_t
faceNodes
[
MAX_numNodes
];
/* list of the nodes defining the face */
149
index_t
shiftNodes
[
MAX_numNodes
];
/* defines a permutation of the nodes which rotates the nodes on the face */
150
index_t
reverseNodes
[
MAX_numNodes
];
/* reverses the order of the nodes on a face. The permutation has to keep 0 fixed. */
151
/* shiftNodes={-1} or reverseNodes={-1} are ignored. */
152
}
Finley_ReferenceElementInfo
;
153
154
/************************************************************************************/
155
156
/* this struct holds the realization of a reference element */
157
158
typedef
struct
Finley_ReferenceElement
{
159
Finley_ReferenceElementInfo
*
Type
;
/* type of the reference element */
160
int
numQuadNodes
;
/* number of quadrature points */
161
double
*
QuadNodes
;
/* coordinates of quadrature nodes */
162
double
*
QuadWeights
;
/* weights of the quadrature scheme */
163
double
*
S
;
/* shape functions at quadrature nodes */
164
double
*
dSdv
;
/* derivative of the shape functions at quadrature nodes */
165
}
Finley_ReferenceElement
;
166
167
/************************************************************************************/
168
169
/* interfaces: */
170
171
Finley_ReferenceElement
*
Finley_ReferenceElement_alloc
(
ElementTypeId
,
int
);
172
void
Finley_ReferenceElement_dealloc
(
Finley_ReferenceElement
*);
173
ElementTypeId
Finley_ReferenceElement_getTypeId
(
char
*);
174
175
#endif
/* #ifndef INC_FINLEY_REFERENCEELEMENTS */
Generated on Fri Jun 28 2013 11:10:52 for ESScript by
1.8.1.2