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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 4491 by jfenwick, Tue Apr 2 04:58:44 2013 UTC revision 4492 by caltinay, Tue Jul 2 01:44:11 2013 UTC
# Line 14  Line 14 
14  *****************************************************************************/  *****************************************************************************/
15    
16    
17  /************************************************************************************/  /****************************************************************************
18    
19  /*   Finley: Reference elements */    Finley: Reference elements
   
 /************************************************************************************/  
   
 #ifndef INC_FINLEY_REFERENCEELEMENTS  
 #define INC_FINLEY_REFERENCEELEMENTS  
20    
21    *****************************************************************************/
22    
23  /************************************************************************************/  #ifndef __FINLEY_REFERENCEELEMENTS_H__
24    #define __FINLEY_REFERENCEELEMENTS_H__
25    
26  #include "Finley.h"  #include "Finley.h"
27  #include "ShapeFunctions.h"  #include "ShapeFunctions.h"
28  #include "Quadrature.h"  #include "Quadrature.h"
29    
30  /************************************************************************************/  //  The ids of the allowed reference elements:
   
 /*     The ids of the allowed reference elements: */  
   
31  #define MAX_numNodes 64  #define MAX_numNodes 64
32  #define MAX_numSubElements 8  #define MAX_numSubElements 8
33  #define MAX_numSides 2  #define MAX_numSides 2
34    
35  typedef enum {  namespace finley {
   Finley_Point1,  
   Finley_Line2,  
   Finley_Line3,  
   Finley_Line4,  
   Finley_Tri3,  
   Finley_Tri6,  
   Finley_Tri9,  
   Finley_Tri10,  
   Finley_Rec4,  
   Finley_Rec8,  
   Finley_Rec9,  
   Finley_Rec12,  
   Finley_Rec16,  
   Finley_Tet4,  
   Finley_Tet10,  
   Finley_Tet16,  
   Finley_Hex8,  
   Finley_Hex20,  
   Finley_Hex27,  
   Finley_Hex32,  
   Finley_Line2Face,  
   Finley_Line3Face,  
   Finley_Line4Face,  
   Finley_Tri3Face,  
   Finley_Tri6Face,  
   Finley_Tri9Face,  
   Finley_Tri10Face,  
   Finley_Rec4Face,  
   Finley_Rec8Face,  
   Finley_Rec9Face,  
   Finley_Rec12Face,  
   Finley_Rec16Face,  
   Finley_Tet4Face,  
   Finley_Tet10Face,  
   Finley_Tet16Face,  
   Finley_Hex8Face,  
   Finley_Hex20Face,  
   Finley_Hex27Face,  
   Finley_Hex32Face,  
   Finley_Point1_Contact,  
   Finley_Line2_Contact,  
   Finley_Line3_Contact,  
   Finley_Line4_Contact,  
   Finley_Tri3_Contact,  
   Finley_Tri6_Contact,  
   Finley_Tri9_Contact,  
   Finley_Tri10_Contact,  
   Finley_Rec4_Contact,  
   Finley_Rec8_Contact,  
   Finley_Rec9_Contact,  
   Finley_Rec12_Contact,  
   Finley_Rec16_Contact,  
   Finley_Line2Face_Contact,  
   Finley_Line3Face_Contact,  
   Finley_Line4Face_Contact,  
   Finley_Tri3Face_Contact,  
   Finley_Tri6Face_Contact,  
   Finley_Tri9Face_Contact,  
   Finley_Tri10Face_Contact,  
   Finley_Rec4Face_Contact,  
   Finley_Rec8Face_Contact,  
   Finley_Rec9Face_Contact,  
   Finley_Rec12Face_Contact,  
   Finley_Rec16Face_Contact,  
   Finley_Tet4Face_Contact,  
   Finley_Tet10Face_Contact,  
   Finley_Tet16Face_Contact,  
   Finley_Hex8Face_Contact,  
   Finley_Hex20Face_Contact,  
   Finley_Hex27Face_Contact,  
   Finley_Hex32Face_Contact,  
   Finley_Line3Macro,  
   Finley_Tri6Macro,  
   Finley_Rec9Macro,  
   Finley_Tet10Macro,  
   Finley_Hex27Macro,  
   Finley_NoRef   /* marks end of list */  
 } Finley_ElementTypeId;  
   
 /************************************************************************************/  
   
 /*  this struct holds the definition of the reference element: */  
   
 typedef struct Finley_ReferenceElementInfo {  
   Finley_ElementTypeId TypeId;               /* the id */  
   const char* Name;                                /* the name in text form e.g. Line1,Rec12,... */  
   dim_t numNodes;                            /* number of nodes defining the element*/  
   dim_t numSubElements;                      /* number of subelements. >1 if macro elements are used. */  
   dim_t numSides;                            /* specifies the number of sides the element supports. This =2 if contact elements are used  
                                                 otherwise =1. */  
                                                   
   
   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*/                                  
   
     
   Finley_ElementTypeId LinearTypeId;         /* id of the linear version of the element */  
     
   index_t linearNodes[MAX_numNodes*MAX_numSides];  /* gives the list of nodes defining the linear or macro element */  
     
   Finley_QuadTypeId Quadrature;                /* quadrature scheme */  
   Finley_ShapeFunctionTypeId Parametrization;  /* shape function for parametrization of the element */  
   Finley_ShapeFunctionTypeId BasisFunctions;   /* shape function for the basis functions */  
   
   index_t subElementNodes[MAX_numNodes*MAX_numSides*MAX_numSubElements];         /* gives the list of nodes defining the subelements:  
                                                                         subElementNodes[INDEX2(i,s,BasisFunctions->numShape*numSides)] is the i-th node in the s-th subelement.*/  
 /********************************************************************************************************************************************************* */    
   dim_t numRelevantGeoNodes;                 /* number of nodes used to describe the geometry of the geometrically relevant part of the element  
                                                 typically this is numNodes but for 'Face' elements where the quadrature points are defined on face of the element  
                         this is the number of nodes on the particular face. */  
   index_t relevantGeoNodes[MAX_numNodes];    /* list to gather the geometrically relevant nodes (length used is numRelevantGeoNodes)  
                                                 this list is used for the VTK interface */  
     
   dim_t numNodesOnFace;                       /* if the element is allowed as a face element, numNodesOnFace defines the number of nodes defining the face */  
                                               /* the following lists are only used for face elements defined by numNodesOnFace>0 */  
   index_t faceNodes[MAX_numNodes];             /* list of the nodes defining the face */  
   index_t shiftNodes[MAX_numNodes];           /* defines a permutation of the nodes which rotates the nodes on the face */  
   index_t reverseNodes[MAX_numNodes];         /* reverses the order of the nodes on a face. The permutation has to keep 0 fixed. */  
                                               /* shiftNodes={-1} or reverseNodes={-1} are ignored. */  
 }  Finley_ReferenceElementInfo;  
   
   
 /************************************************************************************/  
   
 /*  this struct holds the realization of a reference element */  
   
 typedef struct Finley_ReferenceElement {  
     Finley_ReferenceElementInfo* Type;     /* type of the reference element */  
     Finley_ReferenceElementInfo* LinearType;     /* type of the linear reference element */  
     index_t reference_counter;         /* reference counter */  
         dim_t integrationOrder;                /* used integration order */  
     dim_t numNodes;  
         dim_t numLocalDim;  
     dim_t numLinearNodes;  
     Finley_ShapeFunction* Parametrization;  
     Finley_ShapeFunction* BasisFunctions;  
     Finley_ShapeFunction* LinearBasisFunctions;  
         double* DBasisFunctionDv;                              /* pointer to derivatives to basis function corresponding to the Parametrization quad points */  
         bool_t DBasisFunctionDvShared;                /* TRUE to indicate that DBasisFunctionDv is shared with another object which is managing it */  
   
 }  Finley_ReferenceElement;  
   
 /************************************************************************************/  
   
 /*    interfaces: */  
   
 Finley_ReferenceElement* Finley_ReferenceElement_alloc(Finley_ElementTypeId,int);  
 void Finley_ReferenceElement_dealloc(Finley_ReferenceElement*);  
 Finley_ElementTypeId Finley_ReferenceElement_getTypeId(char*);  
 Finley_ReferenceElement* Finley_ReferenceElement_reference(Finley_ReferenceElement* in);  
 Finley_ReferenceElementInfo* Finley_ReferenceElement_getInfo(Finley_ElementTypeId id);  
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      Hex27,
57      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      Hex27Face,
76      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      Hex20Face_Contact,
107      Hex27Face_Contact,
108      Hex32Face_Contact,
109      Line3Macro,
110      Tri6Macro,
111      Rec9Macro,
112      Tet10Macro,
113      Hex27Macro,
114      NoRef   /* marks end of list */
115    } ElementTypeId;
116    
117    
118    /// this struct holds the definition of the reference element
119    struct ReferenceElementInfo {
120        /// the type
121        ElementTypeId TypeId;
122        /// the name in text form e.g. "Line1", "Rec12", ...
123        const char* Name;
124        /// number of nodes defining the element
125        int numNodes;
126        /// number of subelements (>1 if macro elements are used)
127        int numSubElements;
128        /// the number of sides the element supports (=2 if contact elements are
129        /// used, otherwise =1)
130        int numSides;
131        /// offset to the side nodes: offsets[s]...offsets[s+1]-1 refer to the
132        /// nodes to be used for side s
133        int offsets[MAX_numSides+1];
134    
135        /// type id of the linear version of the element
136        ElementTypeId LinearTypeId;
137        /// stores the list of nodes defining the linear or macro element
138        int linearNodes[MAX_numNodes*MAX_numSides];
139        /// quadrature scheme
140        QuadTypeId Quadrature;
141        /// shape function for parametrization of the element
142        ShapeFunctionTypeId Parametrization;
143        /// shape function for the basis functions
144        ShapeFunctionTypeId BasisFunctions;
145    
146        /// the list of nodes defining the subelements, i.e.
147        /// subElementNodes[INDEX2(i,s,BasisFunctions->numShape*numSides)] is
148        /// the i-th node in the s-th subelement
149        int subElementNodes[MAX_numNodes*MAX_numSides*MAX_numSubElements];
150    
151        /// deprecated
152        int numRelevantGeoNodes;
153        int relevantGeoNodes[MAX_numNodes];
154    
155        /// if the element is allowed as a face element, numNodesOnFace defines
156        /// the number of nodes defining the face
157        int numNodesOnFace;
158    
159        // the following lists are only used for face elements defined by
160        // numNodesOnFace>0:
161    
162        /// list of the nodes defining the face
163        int faceNodes[MAX_numNodes];
164    
165        // shiftNodes={-1} or reverseNodes={-1} are ignored.
166        /// defines a permutation of the nodes which rotates the nodes on the face
167        int shiftNodes[MAX_numNodes];
168        /// reverses the order of the nodes on a face. The permutation has to keep
169        /// 0 fixed.
170        int reverseNodes[MAX_numNodes];
171    };
172    
173    
174    /// this struct holds the realization of a reference element
175    struct ReferenceElement {
176        /// type of the reference element
177        ReferenceElementInfo* Type;
178        /// type of the linear reference element
179        ReferenceElementInfo* LinearType;
180        /// reference counter
181        int reference_counter;
182        /// used integration order
183        int integrationOrder;
184        int numNodes;
185        int numLocalDim;
186        int numLinearNodes;
187        ShapeFunction* Parametrization;
188        ShapeFunction* BasisFunctions;
189        ShapeFunction* LinearBasisFunctions;
190        /// pointer to derivatives to basis function corresponding to the
191        /// Parametrization of quad points
192        double* DBasisFunctionDv;
193        /// TRUE to indicate that DBasisFunctionDv is shared with another object
194        /// which is managing it
195        bool_t DBasisFunctionDvShared;
196    };
197    
198    
199    /****** interfaces ******/
200    
201    ReferenceElement* ReferenceElement_alloc(ElementTypeId,int);
202    void ReferenceElement_dealloc(ReferenceElement*);
203    ElementTypeId ReferenceElement_getTypeId(char*);
204    ReferenceElement* ReferenceElement_reference(ReferenceElement* in);
205    ReferenceElementInfo* ReferenceElement_getInfo(ElementTypeId id);
206    
207    inline int ReferenceElement_getNumNodes(const ReferenceElement* in)
208    {
209        return in->Type->numNodes;
210    }
211    
212  #define Finley_ReferenceElement_getNumNodes(__in__) (__in__)->Type->numNodes  } // namespace finley
213    
214  #endif /* #ifndef INC_FINLEY_REFERENCEELEMENTS */  #endif // __FINLEY_REFERENCEELEMENTS_H__
215    

Legend:
Removed from v.4491  
changed lines
  Added in v.4492

  ViewVC Help
Powered by ViewVC 1.1.26