/[escript]/branches/domexper/dudley/src/ReferenceElements.c
ViewVC logotype

Contents of /branches/domexper/dudley/src/ReferenceElements.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3130 - (show annotations)
Thu Sep 2 01:40:12 2010 UTC (8 years, 10 months ago) by jfenwick
File MIME type: text/plain
File size: 8427 byte(s)


1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2010 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
14
15 /**************************************************************/
16
17 /* Dudley: Reference elements */
18
19 /**************************************************************/
20
21 #include "ReferenceElements.h"
22
23 /**************************************************************
24
25 this list has been generated by generateReferenceElementList.py:.
26 */
27 Dudley_ReferenceElementInfo Dudley_ReferenceElement_InfoList[]={
28 { Point1, "Point1", 1, 1, { 0, 1 }, Point1,
29 { 0 }, PointQuad, Point1Shape, Point1Shape,
30 { 0 },
31 1, { 0 },
32 1, { 0 },
33 { 0 },
34 { -1 } },
35 { Line2, "Line2", 2, 1, { 0, 2 }, Line2,
36 { 0, 1 }, LineQuad, Line2Shape, Line2Shape,
37 { 0, 1 },
38 2, { 0, 1 },
39 2, { 0, 1 },
40 { 1, 0 },
41 { -1 } },
42 { Tri3, "Tri3", 3, 1, { 0, 3 }, Tri3,
43 { 0, 1, 2 }, TriQuad, Tri3Shape, Tri3Shape,
44 { 0, 1, 2 },
45 3, { 0, 1, 2 },
46 3, { 0, 1, 2 },
47 { 1, 2, 0 },
48 { 0, 2, 1 } },
49 { Tet4, "Tet4", 4, 1, { 0, 4 }, Tet4,
50 { 0, 1, 2, 3 }, TetQuad, Tet4Shape, Tet4Shape,
51 { 0, 1, 2, 3 },
52 4, { 0, 1, 2, 3 },
53 4, { 0, 1, 2, 3 },
54 { -1 },
55 { -1 } },
56 { Line2Face, "Line2Face", 2, 1, { 0, 2 }, Line2Face,
57 { 0, 1 }, PointQuad, Line2Shape, Line2Shape,
58 { 0, 1 },
59 1, { 0 },
60 1, { 0 },
61 { 0, 1, 2 },
62 { -1 } },
63 { Tri3Face, "Tri3Face", 3, 1, { 0, 3 }, Tri3Face,
64 { 0, 1, 2 }, LineQuad, Tri3Shape, Tri3Shape,
65 { 0, 1, 2 },
66 2, { 0, 1 },
67 2, { 0, 1 },
68 { 1, 0, 2 },
69 { -1 } },
70 { Tet4Face, "Tet4Face", 4, 1, { 0, 4 }, Tet4Face,
71 { 0, 1, 2, 3 }, TriQuad, Tet4Shape, Tet4Shape,
72 { 0, 1, 2, 3 },
73 3, { 0, 1, 2 },
74 4, { 0, 1, 2, 3 },
75 { 1, 2, 0, 3 },
76 { 0, 2, 1, 3 } },
77 { NoRef, "noElement", 0, 0, { 0 }, NoRef,
78 { 0 }, NoQuad, NoShape, NoShape,
79 { 0 },
80 -1, { 0 },
81 -1, { 0 },
82 { 0 },
83 { 0 } }
84
85 };
86
87 /**************************************************************
88
89 creates a ReferenceElement of type id and a given integration order
90
91 */
92
93
94 Dudley_ReferenceElement* Dudley_ReferenceElement_alloc(ElementTypeId id, int order)
95 {
96 dim_t nsub, numQuadNodes, numQuadNodes2;
97 double *quadWeights=NULL, *quadNodes=NULL, *quadWeights2=NULL, *quadNodes2=NULL;
98 Dudley_ReferenceElement *out=NULL;
99 Dudley_QuadInfo* quadscheme;
100 Dudley_ShapeFunctionInfo* parametrization, *basisfunction, *linearbasisfunction;
101 Dudley_ReferenceElementInfo *type, *linear_type;
102
103 type=Dudley_ReferenceElement_getInfo(id);
104 if (type == NULL) {
105 Dudley_setError(VALUE_ERROR,"Dudley_ReferenceElement_alloc: unable to identify element type.");
106 return NULL;
107 }
108 linear_type=Dudley_ReferenceElement_getInfo(type->LinearTypeId);
109 if (linear_type == NULL) {
110 Dudley_setError(VALUE_ERROR,"Dudley_ReferenceElement_alloc: unable to identify linear element type.");
111 return NULL;
112 }
113
114
115
116 out=MEMALLOC(1,Dudley_ReferenceElement);
117 if (Dudley_checkPtr(out)) return NULL;
118
119 out->reference_counter=0;
120 out->Parametrization=NULL;
121 out->BasisFunctions=NULL;
122 out->LinearBasisFunctions=NULL;
123 out->Type=type;
124 out->numNodes=out->Type->numNodes;
125 out->LinearType=linear_type;
126 out->numLinearNodes=out->LinearType->numNodes;
127 out->integrationOrder=-1;
128 out->DBasisFunctionDv=NULL;
129 out->DBasisFunctionDvShared=TRUE;
130
131 quadscheme=Dudley_QuadInfo_getInfo(out->Type->Quadrature);
132 parametrization=Dudley_ShapeFunction_getInfo(out->Type->Parametrization);
133 basisfunction=Dudley_ShapeFunction_getInfo(out->Type->BasisFunctions);
134 linearbasisfunction=Dudley_ShapeFunction_getInfo(Dudley_ReferenceElement_InfoList[out->Type->LinearTypeId].BasisFunctions);
135 nsub=out->Type->numSubElements;
136 out->numLocalDim=quadscheme->numDim;
137
138
139 /* set up the basic integration scheme
140 note that quadscheme->numDim is not necessarily the diemnsion of the element
141 */
142
143 if (order<0) order=MAX(2*basisfunction->numOrder,0);
144 out->integrationOrder=order;
145
146 numQuadNodes=quadscheme->getNumQuadNodes(order);
147
148 quadNodes=MEMALLOC(numQuadNodes*quadscheme->numDim*nsub,double);
149 quadWeights=MEMALLOC(numQuadNodes*nsub,double);
150 if ( !( Dudley_checkPtr(quadNodes) || Dudley_checkPtr(quadWeights) ) ) {
151
152 quadscheme->getQuadNodes(numQuadNodes, quadNodes, quadWeights);
153
154 /* set the basis functions on the quadrature points:
155 * note: Dudley_ShapeFunction_alloc will introduce 0. if the dimensions of the quadrature scheme and the dimension of the element don;t match.
156 */
157
158
159 /*
160 * before we can set the shape function for the parametrization the quadrature scheme needs to be replicated :
161 */
162
163 if (nsub>1) {
164 out->DBasisFunctionDv=MEMALLOC(numQuadNodes*nsub* (basisfunction->numShapes) * (basisfunction->numDim), double );
165 out->DBasisFunctionDvShared=FALSE;
166
167 out->BasisFunctions=Dudley_ShapeFunction_alloc(basisfunction->TypeId, quadscheme->numDim, numQuadNodes, quadNodes, quadWeights);
168 quadNodes2=MEMALLOC(numQuadNodes*quadscheme->numDim*nsub,double);
169 quadWeights2=MEMALLOC(numQuadNodes*nsub,double);
170 if ( !( Dudley_checkPtr(quadNodes2) || Dudley_checkPtr(quadWeights2) || Dudley_checkPtr(out->DBasisFunctionDv)) ) {
171
172 numQuadNodes2=nsub*out->BasisFunctions->numQuadNodes; // fake output of macro function
173 if (Dudley_noError()) {
174 out->Parametrization=Dudley_ShapeFunction_alloc(parametrization->TypeId, quadscheme->numDim, numQuadNodes2, quadNodes2, quadWeights2);
175 out->LinearBasisFunctions=Dudley_ShapeFunction_alloc(linearbasisfunction->TypeId, quadscheme->numDim, numQuadNodes2, quadNodes2, quadWeights2);
176 }
177
178 }
179 MEMFREE(quadNodes2);
180 MEMFREE(quadWeights2);
181
182 } else {
183 out->Parametrization=Dudley_ShapeFunction_alloc(parametrization->TypeId, quadscheme->numDim, numQuadNodes*nsub, quadNodes, quadWeights);
184 out->BasisFunctions=Dudley_ShapeFunction_alloc(basisfunction->TypeId, quadscheme->numDim, numQuadNodes, quadNodes, quadWeights);
185 out->LinearBasisFunctions=Dudley_ShapeFunction_alloc(linearbasisfunction->TypeId, quadscheme->numDim, numQuadNodes, quadNodes, quadWeights);
186 out->DBasisFunctionDv=out->BasisFunctions->dSdv;
187 out->DBasisFunctionDvShared=TRUE;
188 }
189 }
190 MEMFREE(quadNodes);
191 MEMFREE(quadWeights);
192 if (Dudley_noError()) {
193 return Dudley_ReferenceElement_reference(out);
194 } else {
195 Dudley_ReferenceElement_dealloc(out);
196 return NULL;
197 }
198 }
199
200 /**************************************************************/
201
202 void Dudley_ReferenceElement_dealloc(Dudley_ReferenceElement* in) {
203 if (in!=NULL) {
204 in->reference_counter--;
205 if (in->reference_counter<1) {
206 Dudley_ShapeFunction_dealloc(in->Parametrization);
207 Dudley_ShapeFunction_dealloc(in->BasisFunctions);
208 Dudley_ShapeFunction_dealloc(in->LinearBasisFunctions);
209 if (! in->DBasisFunctionDvShared) MEMFREE(in->DBasisFunctionDv);
210 MEMFREE(in);
211 }
212 }
213
214 }
215
216 /**************************************************************/
217
218 ElementTypeId Dudley_ReferenceElement_getTypeId(char* element_type) {
219 int ptr=0;
220 ElementTypeId out=NoRef;
221 while (Dudley_ReferenceElement_InfoList[ptr].TypeId!=NoRef && out==NoRef) {
222 if (strcmp(element_type,Dudley_ReferenceElement_InfoList[ptr].Name)==0) out=Dudley_ReferenceElement_InfoList[ptr].TypeId;
223 ptr++;
224 }
225 return out;
226 }
227
228 Dudley_ReferenceElement* Dudley_ReferenceElement_reference(Dudley_ReferenceElement* in) {
229 if (in!=NULL) ++(in->reference_counter);
230 return in;
231 }
232
233 Dudley_ReferenceElementInfo* Dudley_ReferenceElement_getInfo(ElementTypeId id)
234 {
235 int ptr=0;
236 Dudley_ReferenceElementInfo* out=NULL;
237 while (Dudley_ReferenceElement_InfoList[ptr].TypeId!=NoRef && out==NULL) {
238 if (Dudley_ReferenceElement_InfoList[ptr].TypeId==id) out=&(Dudley_ReferenceElement_InfoList[ptr]);
239 ptr++;
240 }
241 if (out==NULL) {
242 Dudley_setError(VALUE_ERROR,"Dudley_ReferenceElement_getInfo: cannot find requested reference element.");
243 }
244 return out;
245 }
246
247

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26