ESYS13  Revision_
BasisFunctions.h
Go to the documentation of this file.
00001 
00002 /*******************************************************
00003 *
00004 * Copyright (c) 2003-2012 by University of Queensland
00005 * Earth Systems Science Computational Center (ESSCC)
00006 * http://www.uq.edu.au/esscc
00007 *
00008 * Primary Business: Queensland, Australia
00009 * Licensed under the Open Software License version 3.0
00010 * http://www.opensource.org/licenses/osl-3.0.php
00011 *
00012 *******************************************************/
00013 
00014 
00015 /**************************************************************/
00016 
00017 /*   Finley: Reference elements */
00018 
00019 /**************************************************************/
00020 
00021 #ifndef INC_FINLEY_REFERENCEELEMENTS
00022 #define INC_FINLEY_REFERENCEELEMENTS
00023 
00024 
00025 /**************************************************************/
00026 
00027 #include "Finley.h"
00028 #include "ShapeFunctions.h"
00029 #include "Quadrature.h"
00030 
00031 /**************************************************************/
00032 
00033 /*     The ids of the allowed reference elements: */
00034 
00035 #define MAX_numNodes 64
00036 
00037 typedef enum {
00038   Point1,
00039   Line2,
00040   Line3,
00041   Line4,
00042   Tri3,
00043   Tri6,
00044   Tri9,
00045   Tri10,
00046   Rec4,
00047   Rec8,
00048   Rec9,
00049   Rec12,
00050   Rec16,
00051   Tet4,
00052   Tet10,
00053   Tet16,
00054   Hex8,
00055   Hex20,
00056   Hex27,
00057   Hex32,
00058   Line2Face,
00059   Line3Face,
00060   Line4Face,
00061   Tri3Face,
00062   Tri6Face,
00063   Tri9Face,
00064   Tri10Face,
00065   Rec4Face,
00066   Rec8Face,
00067   Rec9Face,
00068   Rec12Face,
00069   Rec16Face,
00070   Tet4Face, 
00071   Tet10Face, 
00072   Tet16Face,
00073   Hex8Face,
00074   Hex20Face, 
00075   Hex27Face, 
00076   Hex32Face, 
00077   Point1_Contact,
00078   Line2_Contact,
00079   Line3_Contact,
00080   Line4_Contact,
00081   Tri3_Contact,
00082   Tri6_Contact,
00083   Tri9_Contact,
00084   Tri10_Contact,
00085   Rec4_Contact,
00086   Rec8_Contact,
00087   Rec9_Contact,
00088   Rec12_Contact,
00089   Rec16_Contact,
00090   Line2Face_Contact,
00091   Line3Face_Contact,
00092   Line4Face_Contact,
00093   Tri3Face_Contact,
00094   Tri6Face_Contact,
00095   Tri9Face_Contact,
00096   Tri10Face_Contact,
00097   Rec4Face_Contact,
00098   Rec8Face_Contact,
00099   Rec9Face_Contact,
00100   Rec12Face_Contact,
00101   Rec16Face_Contact,
00102   Tet4Face_Contact, 
00103   Tet10Face_Contact, 
00104   Tet16Face_Contact,
00105   Hex8Face_Contact,
00106   Hex20Face_Contact,
00107   Hex27Face_Contact, 
00108   Hex32Face_Contact, 
00109   Line3Macro, 
00110   Tri6Macro, 
00111   Rec9Macro,
00112   Tet10Macro,
00113   Hex27Macro,
00114 
00115   NoType   /* marks end of list */
00116 } ElementTypeId;
00117 
00118 /**************************************************************/
00119 
00120 /*  this struct holds the definition of the reference element: */
00121 
00122 typedef struct Finley_ReferenceElementInfo {
00123   ElementTypeId TypeId;                      /* the id */
00124   char* Name;                                /* the name in text form e.g. Line1,Rec12,... */
00125   dim_t numLocalDim;                         /* local dimension of the element */
00126   dim_t numDim;                              /* dimension of the element */
00127   dim_t numNodes;                            /* number of nodes defining the element*/
00128   dim_t numShapes;                           /* number of shape functions, typically = numNodes*/
00129   dim_t numOrder;                            /* order of the shape functions */
00130   dim_t numVertices;                         /* number of vertices of the element */
00131   ElementTypeId LinearTypeId;                /* id of the linear version of the element */
00132   index_t linearNodes[MAX_numNodes];         /* gives the list of nodes defining the linear or macro element, typically it is linearNodes[i]=i */
00133   Finley_Shape_Function* getValues;          /* function to evaluate the shape functions at a set of points */
00134   Finley_Quad_getNodes* getQuadNodes;        /* function to set the quadrature points */
00135   Finley_Quad_getNumNodes* getNumQuadNodes;  /* function selects the number of quadrature nodes for a given accuracy order */
00136   
00137 /*********************************************************************************************************************************** */  
00138   dim_t numRelevantGeoNodes;                 /* number of nodes used to describe the geometry of the geometrically relevant part of the element
00139                                                 typically this is numNodes but for 'Face' elements where the quadrature points are defined on face of the element 
00140                         this is the number of nodes on the particular face. */
00141   index_t relevantGeoNodes[MAX_numNodes];    /* list to gather the geometrically relevant nodes (length used is numRelevantGeoNodes)
00142                                                 this list is used for the VTK interface */
00143   
00144   dim_t numNodesOnFace;                       /* if the element is allowed as a face element, numNodesOnFace defines the number of nodes defining the face */
00145                                               /* the following lists are only used for face elements defined by numNodesOnFace>0 */
00146   index_t faceNodes[MAX_numNodes];             /* list of the nodes defining the face */
00147   index_t shiftNodes[MAX_numNodes];           /* defines a permutation of the nodes which rotates the nodes on the face */
00148   index_t reverseNodes[MAX_numNodes];         /* reverses the order of the nodes on a face. The permutation has to keep 0 fixed. */
00149                                               /* shiftNodes={-1} or reverseNodes={-1} are ignored. */
00150 }  Finley_ReferenceElementInfo;
00151 
00152 /**************************************************************/
00153 
00154 /*  this struct holds the realization of a reference element */
00155 
00156 typedef struct Finley_ReferenceElement {
00157   Finley_ReferenceElementInfo* Type;     /* type of the reference element */
00158   int numQuadNodes;                /* number of quadrature points */
00159   double *QuadNodes;               /* coordinates of quadrature nodes */
00160   double *QuadWeights;             /* weights of the quadrature scheme */
00161   double *S;                       /* shape functions at quadrature nodes */
00162   double *dSdv;                    /* derivative of the shape functions at quadrature nodes */
00163 }  Finley_ReferenceElement;
00164 
00165 /**************************************************************/
00166 
00167 /*    interfaces: */
00168 
00169 Finley_ReferenceElement* Finley_ReferenceElement_alloc(ElementTypeId,int);
00170 void Finley_ReferenceElement_dealloc(Finley_ReferenceElement*);
00171 ElementTypeId Finley_ReferenceElement_getTypeId(char*);
00172 
00173 #endif /* #ifndef INC_FINLEY_REFERENCEELEMENTS */