escript
Revision_
|
00001 00002 /***************************************************************************** 00003 * 00004 * Copyright (c) 2003-2014 by University of Queensland 00005 * http://www.uq.edu.au 00006 * 00007 * Primary Business: Queensland, Australia 00008 * Licensed under the Open Software License version 3.0 00009 * http://www.opensource.org/licenses/osl-3.0.php 00010 * 00011 * Development until 2012 by Earth Systems Science Computational Center (ESSCC) 00012 * Development 2012-2013 by School of Earth Sciences 00013 * Development from 2014 by Centre for Geoscience Computing (GeoComp) 00014 * 00015 *****************************************************************************/ 00016 00017 00018 /**************************************************************************** 00019 00020 Finley: Reference elements 00021 00022 *****************************************************************************/ 00023 00024 #ifndef __FINLEY_REFERENCEELEMENTS_H__ 00025 #define __FINLEY_REFERENCEELEMENTS_H__ 00026 00027 #include "Finley.h" 00028 #include "ShapeFunctions.h" 00029 #include "Quadrature.h" 00030 00031 // The ids of the allowed reference elements: 00032 #define MAX_numNodes 64 00033 #define MAX_numSubElements 8 00034 #define MAX_numSides 2 00035 00036 namespace finley { 00037 00038 typedef enum { 00039 Point1, 00040 Line2, 00041 Line3, 00042 Line4, 00043 Tri3, 00044 Tri6, 00045 Tri9, 00046 Tri10, 00047 Rec4, 00048 Rec8, 00049 Rec9, 00050 Rec12, 00051 Rec16, 00052 Tet4, 00053 Tet10, 00054 Tet16, 00055 Hex8, 00056 Hex20, 00057 Hex27, 00058 Hex32, 00059 Line2Face, 00060 Line3Face, 00061 Line4Face, 00062 Tri3Face, 00063 Tri6Face, 00064 Tri9Face, 00065 Tri10Face, 00066 Rec4Face, 00067 Rec8Face, 00068 Rec9Face, 00069 Rec12Face, 00070 Rec16Face, 00071 Tet4Face, 00072 Tet10Face, 00073 Tet16Face, 00074 Hex8Face, 00075 Hex20Face, 00076 Hex27Face, 00077 Hex32Face, 00078 Point1_Contact, 00079 Line2_Contact, 00080 Line3_Contact, 00081 Line4_Contact, 00082 Tri3_Contact, 00083 Tri6_Contact, 00084 Tri9_Contact, 00085 Tri10_Contact, 00086 Rec4_Contact, 00087 Rec8_Contact, 00088 Rec9_Contact, 00089 Rec12_Contact, 00090 Rec16_Contact, 00091 Line2Face_Contact, 00092 Line3Face_Contact, 00093 Line4Face_Contact, 00094 Tri3Face_Contact, 00095 Tri6Face_Contact, 00096 Tri9Face_Contact, 00097 Tri10Face_Contact, 00098 Rec4Face_Contact, 00099 Rec8Face_Contact, 00100 Rec9Face_Contact, 00101 Rec12Face_Contact, 00102 Rec16Face_Contact, 00103 Tet4Face_Contact, 00104 Tet10Face_Contact, 00105 Tet16Face_Contact, 00106 Hex8Face_Contact, 00107 Hex20Face_Contact, 00108 Hex27Face_Contact, 00109 Hex32Face_Contact, 00110 Line3Macro, 00111 Tri6Macro, 00112 Rec9Macro, 00113 Tet10Macro, 00114 Hex27Macro, 00115 NoRef // marks end of list 00116 } ElementTypeId; 00117 00118 00120 struct ReferenceElementInfo { 00122 ElementTypeId TypeId; 00124 const char* Name; 00126 int numNodes; 00128 int numSubElements; 00131 int numSides; 00134 int offsets[MAX_numSides+1]; 00135 00137 ElementTypeId LinearTypeId; 00139 int linearNodes[MAX_numNodes*MAX_numSides]; 00141 QuadTypeId Quadrature; 00143 ShapeFunctionTypeId Parametrization; 00145 ShapeFunctionTypeId BasisFunctions; 00146 00150 int subElementNodes[MAX_numNodes*MAX_numSides*MAX_numSubElements]; 00151 00153 int numRelevantGeoNodes; 00154 int relevantGeoNodes[MAX_numNodes]; 00155 00158 int numNodesOnFace; 00159 00160 // the following lists are only used for face elements defined by 00161 // numNodesOnFace>0: 00162 00164 int faceNodes[MAX_numNodes]; 00165 00166 // shiftNodes={-1} or reverseNodes={-1} are ignored. 00168 int shiftNodes[MAX_numNodes]; 00171 int reverseNodes[MAX_numNodes]; 00172 }; 00173 00174 00176 struct ReferenceElement { 00178 ReferenceElement(ElementTypeId id, int order); 00179 00181 ~ReferenceElement(); 00182 00184 static const ReferenceElementInfo* getInfo(ElementTypeId id); 00185 00187 static ElementTypeId getTypeId(const char*); 00188 00190 int getNumNodes() const { return Type->numNodes; } 00191 00193 const ReferenceElementInfo* Type; 00195 const ReferenceElementInfo* LinearType; 00197 int integrationOrder; 00198 int numNodes; 00199 int numLocalDim; 00200 int numLinearNodes; 00201 const_ShapeFunction_ptr Parametrization; 00202 const_ShapeFunction_ptr BasisFunctions; 00203 const_ShapeFunction_ptr LinearBasisFunctions; 00206 double* DBasisFunctionDv; 00209 bool DBasisFunctionDvShared; 00210 }; 00211 00212 typedef boost::shared_ptr<ReferenceElement> ReferenceElement_ptr; 00213 typedef boost::shared_ptr<const ReferenceElement> const_ReferenceElement_ptr; 00214 00215 } // namespace finley 00216 00217 #endif // __FINLEY_REFERENCEELEMENTS_H__ 00218