ESScript  Revision_4488
ReferenceElements.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 #define MAX_numSubElements 8
39 #define MAX_numSides 2
40 
41 typedef enum {
118  Finley_NoRef /* marks end of list */
120 
121 /************************************************************************************/
122 
123 /* this struct holds the definition of the reference element: */
124 
125 typedef struct Finley_ReferenceElementInfo {
127  const char* Name; /* the name in text form e.g. Line1,Rec12,... */
128  dim_t numNodes; /* number of nodes defining the element*/
129  dim_t numSubElements; /* number of subelements. >1 if macro elements are used. */
130  dim_t numSides; /* specifies the number of sides the element supports. This =2 if contact elements are used
131  otherwise =1. */
132 
133 
134  index_t offsets[MAX_numSides+1]; /* offset to the side nodes: offsets[s]...offset[s+1]-1 refers to the nodes to be used for side s*/
135 
136 
137  Finley_ElementTypeId LinearTypeId; /* id of the linear version of the element */
138 
139  index_t linearNodes[MAX_numNodes*MAX_numSides]; /* gives the list of nodes defining the linear or macro element */
140 
141  Finley_QuadTypeId Quadrature; /* quadrature scheme */
142  Finley_ShapeFunctionTypeId Parametrization; /* shape function for parametrization of the element */
143  Finley_ShapeFunctionTypeId BasisFunctions; /* shape function for the basis functions */
144 
145  index_t subElementNodes[MAX_numNodes*MAX_numSides*MAX_numSubElements]; /* gives the list of nodes defining the subelements:
146  subElementNodes[INDEX2(i,s,BasisFunctions->numShape*numSides)] is the i-th node in the s-th subelement.*/
147 /********************************************************************************************************************************************************* */
148  dim_t numRelevantGeoNodes; /* number of nodes used to describe the geometry of the geometrically relevant part of the element
149  typically this is numNodes but for 'Face' elements where the quadrature points are defined on face of the element
150  this is the number of nodes on the particular face. */
151  index_t relevantGeoNodes[MAX_numNodes]; /* list to gather the geometrically relevant nodes (length used is numRelevantGeoNodes)
152  this list is used for the VTK interface */
153 
154  dim_t numNodesOnFace; /* if the element is allowed as a face element, numNodesOnFace defines the number of nodes defining the face */
155  /* the following lists are only used for face elements defined by numNodesOnFace>0 */
156  index_t faceNodes[MAX_numNodes]; /* list of the nodes defining the face */
157  index_t shiftNodes[MAX_numNodes]; /* defines a permutation of the nodes which rotates the nodes on the face */
158  index_t reverseNodes[MAX_numNodes]; /* reverses the order of the nodes on a face. The permutation has to keep 0 fixed. */
159  /* shiftNodes={-1} or reverseNodes={-1} are ignored. */
161 
162 
163 /************************************************************************************/
164 
165 /* this struct holds the realization of a reference element */
166 
167 typedef struct Finley_ReferenceElement {
168  Finley_ReferenceElementInfo* Type; /* type of the reference element */
169  Finley_ReferenceElementInfo* LinearType; /* type of the linear reference element */
170  index_t reference_counter; /* reference counter */
171  dim_t integrationOrder; /* used integration order */
178  double* DBasisFunctionDv; /* pointer to derivatives to basis function corresponding to the Parametrization quad points */
179  bool_t DBasisFunctionDvShared; /* TRUE to indicate that DBasisFunctionDv is shared with another object which is managing it */
180 
182 
183 /************************************************************************************/
184 
185 /* interfaces: */
186 
192 
193 
194 #define Finley_ReferenceElement_getNumNodes(__in__) (__in__)->Type->numNodes
195 
196 #endif /* #ifndef INC_FINLEY_REFERENCEELEMENTS */
197