/[escript]/trunk/finley/src/ReferenceElements.h
ViewVC logotype

Annotation of /trunk/finley/src/ReferenceElements.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1811 - (hide annotations)
Thu Sep 25 23:11:13 2008 UTC (11 years, 5 months ago) by ksteube
File MIME type: text/plain
File size: 5621 byte(s)
Copyright updated in all files

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

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26