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

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

Parent Directory Parent Directory | Revision Log Revision Log


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