/[escript]/trunk/finley/src/Quadrature.c
ViewVC logotype

Annotation of /trunk/finley/src/Quadrature.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1170 - (hide annotations)
Fri May 25 05:24:57 2007 UTC (13 years, 5 months ago) by btully
File MIME type: text/plain
File size: 19384 byte(s)
Added new methods for integration of Tri's and Tet's For Order 2 and 3.
NB: There is a conflict in the testing of the getNormal() function. On consultation with Lutz, we have concluded this is an issue with the test rather than the new methods.

1 jgs 150 /*
2 elspeth 616 ************************************************************
3     * Copyright 2006 by ACcESS MNRF *
4     * *
5     * http://www.access.edu.au *
6     * Primary Business: Queensland, Australia *
7     * Licensed under the Open Software License version 3.0 *
8     * http://www.opensource.org/licenses/osl-3.0.php *
9     * *
10     ************************************************************
11 jgs 150 */
12 jgs 82
13     /**************************************************************/
14    
15     /* Finley: */
16    
17     /**************************************************************/
18    
19     /* Author: gross@access.edu.au */
20     /* Version: $Id$ */
21    
22     /**************************************************************/
23    
24     #include "Quadrature.h"
25    
26     #define QUADNODES(_K_,_I_) quadNodes[INDEX2(_K_,_I_,DIM)]
27     #define QUADWEIGHTS(_I_) quadWeights[_I_]
28    
29     /**************************************************************/
30    
31     /* get a quadrature scheme with numQuadNodes quadrature nodes for the tri */
32     /* as a queezed scheme on a quad [0,1]^2 */
33    
34     void Finley_Quad_getNodesTri(int numQuadNodes,double* quadNodes,double* quadWeights) {
35     int i;
36     double Q1,Q2;
37     #define DIM 2
38    
39 btully 1170 /* the easy cases: */
40 jgs 82
41     if (numQuadNodes==1) {
42     QUADNODES(0,0)=1./3.;
43     QUADNODES(1,0)=1./3.;
44 btully 1170 QUADWEIGHTS(0)=1./2.;
45     } else if (numQuadNodes==3){
46     QUADNODES(0,0)=1./2.;
47     QUADNODES(1,0)=0.;
48     QUADWEIGHTS(0)=1./6.;
49     QUADNODES(0,1)=0.;
50     QUADNODES(1,1)=1./2.;
51     QUADWEIGHTS(1)=1./6.;
52     QUADNODES(0,2)=1./2.;
53     QUADNODES(1,2)=1./2.;
54     QUADWEIGHTS(2)=1./6.;
55     } else if (numQuadNodes==4){
56     QUADNODES(0,0)=1./3.;
57     QUADNODES(1,0)=1./3.;
58     QUADWEIGHTS(0)=-27./96.;
59     QUADNODES(0,1)=0.2;
60     QUADNODES(1,1)=0.2;
61     QUADWEIGHTS(1)=25./96.;
62     QUADNODES(0,2)=0.6;
63     QUADNODES(1,2)=0.2;
64     QUADWEIGHTS(2)=25./96.;
65     QUADNODES(0,3)=0.2;
66     QUADNODES(1,3)=0.6;
67     QUADWEIGHTS(3)=25./96.;
68 jgs 82 } else {
69    
70     /* get scheme on [0.1]^2 */
71     Finley_Quad_getNodesRec(numQuadNodes,quadNodes,quadWeights);
72 jgs 150 if (! Finley_noError()) return;
73 jgs 82
74     /* squeeze it: */
75    
76     for (i=0;i<numQuadNodes;i++) {
77     Q1=QUADNODES(0,i);
78     Q2=QUADNODES(1,i);
79     QUADWEIGHTS(i)=QUADWEIGHTS(i)*(1.-(1./2.)*(Q1+Q2));
80     QUADNODES(0,i)=Q1*(1.-(1./2.)*Q2);
81     QUADNODES(1,i)=Q2*(1.-(1./2.)*Q1);
82     }
83     }
84     #undef DIM
85     }
86    
87     /**************************************************************/
88    
89     /* get a quadrature scheme with numQuadNodes quadrature nodes for the tet */
90     /* as a queezed scheme on a hex [0,1]^3 */
91    
92     void Finley_Quad_getNodesTet(int numQuadNodes,double* quadNodes,double* quadWeights) {
93     int i;
94     double Q1,Q2,Q3,JA11,JA12,JA13,JA21,JA22,JA23,JA31,JA32,JA33,DET;
95     #define DIM 3
96    
97 btully 1170 /* the easy cases: */
98 jgs 82 if (numQuadNodes==1) {
99 btully 1170 QUADNODES(0,0)=0.25;
100     QUADNODES(1,0)=0.25;
101     QUADNODES(2,0)=0.25;
102 jgs 82 QUADWEIGHTS(0)=1./6.;
103 btully 1170 } else if (numQuadNodes==4){
104     double alpha=0.58541020;
105     double beta =0.13819660;
106     QUADNODES(0,0)=beta;
107     QUADNODES(1,0)=beta;
108     QUADNODES(2,0)=beta;
109     QUADWEIGHTS(0)=1./24.;
110     QUADNODES(0,1)=alpha;
111     QUADNODES(1,1)=beta;
112     QUADNODES(2,1)=beta;
113     QUADWEIGHTS(1)=1./24.;
114     QUADNODES(0,2)=beta;
115     QUADNODES(1,2)=alpha;
116     QUADNODES(2,2)=beta;
117     QUADWEIGHTS(2)=1./24.;
118     QUADNODES(0,3)=beta;
119     QUADNODES(1,3)=beta;
120     QUADNODES(2,3)=alpha;
121     QUADWEIGHTS(3)=1./24.;
122     } else if (numQuadNodes==5){
123     QUADNODES(0,0)=1./4.;
124     QUADNODES(1,0)=1./4.;
125     QUADNODES(2,0)=1./4.;
126     QUADWEIGHTS(0)=-2./15.;
127     QUADNODES(0,1)=1./6.;
128     QUADNODES(1,1)=1./6.;
129     QUADNODES(2,1)=1./6.;
130     QUADWEIGHTS(1)=3./40.;
131     QUADNODES(0,2)=1./2.;
132     QUADNODES(1,2)=1./6.;
133     QUADNODES(2,2)=1./6.;
134     QUADWEIGHTS(2)=3./40.;
135     QUADNODES(0,3)=1./6.;
136     QUADNODES(1,3)=1./2.;
137     QUADNODES(2,3)=1./6.;
138     QUADWEIGHTS(3)=3./40.;
139     QUADNODES(0,4)=1./6.;
140     QUADNODES(1,4)=1./6.;
141     QUADNODES(2,4)=1./2.;
142     QUADWEIGHTS(4)=3./40.;
143 jgs 82 } else {
144    
145     /* get scheme on [0.1]^3 */
146    
147     Finley_Quad_getNodesHex(numQuadNodes,quadNodes,quadWeights);
148 jgs 150 if (! Finley_noError()) return;
149 jgs 82
150     /* squeeze it: */
151    
152     for (i=0;i<numQuadNodes;i++) {
153     Q1=QUADNODES(0,i);
154     Q2=QUADNODES(1,i);
155     Q3=QUADNODES(2,i);
156    
157     JA11= (1./3.)*Q2*Q3-(1./2.)*(Q2+Q3) +1.;
158     JA12= (1./3.)*Q1*Q3-(1./2.)*Q1;
159     JA13= (1./3.)*Q1*Q2-(1./2.)*Q1;
160     JA21= (1./3.)*Q2*Q3-(1./2.)*Q2;
161     JA22= (1./3.)*Q1*Q3-(1./2.)*(Q1+Q3) +1.;
162     JA23= (1./3.)*Q1*Q2-(1./2.)*Q2;
163     JA31= (1./3.)*Q2*Q3-(1./2.)*Q3;
164     JA32= (1./3.)*Q1*Q3-(1./2.)*Q3;
165     JA33= (1./3.)*Q1*Q2-(1./2.)*(Q1+Q2) +1.;
166     DET=JA11*JA22*JA33+JA12*JA23*JA31+JA13*JA21*JA32-JA13*JA22*JA31-JA11*JA23*JA32-JA12*JA21*JA33;
167     quadWeights[i]=quadWeights[i]*ABS(DET);
168     QUADNODES(0,i)=Q1*((1./3.)*Q2*Q3-(1./2.)*(Q2+Q3)+1.);
169     QUADNODES(1,i)=Q2*((1./3.)*Q1*Q3-(1./2.)*(Q1+Q3)+1.);
170     QUADNODES(2,i)=Q3*((1./3.)*Q1*Q2-(1./2.)*(Q1+Q2)+1.);
171     }
172     }
173     #undef DIM
174     }
175    
176     /**************************************************************/
177    
178     /* get a quadrature scheme with numQuadNodes quadrature nodes for the quad [0.1]^2 */
179     /* as a X-product of a 1D scheme. */
180    
181     void Finley_Quad_getNodesRec(int numQuadNodes,double* quadNodes,double* quadWeights) {
182 jgs 150 char error_msg[LenErrorMsg_MAX];
183 jgs 82 int numQuadNodes1d,i,j,l;
184 gross 1028 double *quadNodes1d=NULL,*quadWeights1d=NULL;
185     bool_t set=FALSE;
186 jgs 82 #define DIM 2
187    
188 gross 1028 quadNodes1d=TMPMEMALLOC(numQuadNodes,double);
189     quadWeights1d=TMPMEMALLOC(numQuadNodes,double);
190     if (! ( Finley_checkPtr(quadNodes1d) || Finley_checkPtr(quadWeights1d) ) ) {
191     /* find numQuadNodes1d with numQuadNodes1d**2==numQuadNodes: */
192    
193     for (numQuadNodes1d=1;numQuadNodes1d<=MAX_numQuadNodesLine;numQuadNodes1d++) {
194     if (numQuadNodes1d*numQuadNodes1d==numQuadNodes) {
195 jgs 82
196 gross 1028 /* get 1D scheme: */
197    
198     Finley_Quad_getNodesLine(numQuadNodes1d,quadNodes1d,quadWeights1d);
199 jgs 82
200 gross 1028 /* make 2D scheme: */
201 jgs 82
202 gross 1028 if (Finley_noError()) {
203     l=0;
204     for (i=0;i<numQuadNodes1d;i++) {
205     for (j=0;j<numQuadNodes1d;j++) {
206     QUADNODES(0,l)=quadNodes1d[i];
207     QUADNODES(1,l)=quadNodes1d[j];
208     QUADWEIGHTS(l)=quadWeights1d[i]*quadWeights1d[j];
209     l++;
210     }
211     }
212     set=TRUE;
213     break;
214     }
215     }
216     }
217     if (!set) {
218     sprintf(error_msg,"Finley_Quad_getNodesRec: Illegal number of quadrature nodes %d on hexahedron.",numQuadNodes);
219     Finley_setError(VALUE_ERROR,error_msg);
220     }
221     TMPMEMFREE(quadNodes1d);
222     TMPMEMFREE(quadWeights1d);
223     }
224     #undef DIM
225 jgs 82 }
226    
227     /**************************************************************/
228    
229     /* get a quadrature scheme with numQuadNodes quadrature nodes for the hex [0.1]^3 */
230     /* as a X-product of a 1D scheme. */
231    
232     void Finley_Quad_getNodesHex(int numQuadNodes,double* quadNodes,double* quadWeights) {
233 jgs 150 char error_msg[LenErrorMsg_MAX];
234 jgs 82 int numQuadNodes1d,i,j,k,l;
235 gross 1028 double *quadNodes1d=NULL,*quadWeights1d=NULL;
236     bool_t set;
237 jgs 82 #define DIM 3
238    
239     /* find numQuadNodes1d with numQuadNodes1d**3==numQuadNodes: */
240    
241 gross 1028 quadNodes1d=TMPMEMALLOC(numQuadNodes,double);
242     quadWeights1d=TMPMEMALLOC(numQuadNodes,double);
243     if (! ( Finley_checkPtr(quadNodes1d) || Finley_checkPtr(quadWeights1d) ) ) {
244     for (numQuadNodes1d=1;numQuadNodes1d<=MAX_numQuadNodesLine;numQuadNodes1d++) {
245     if (numQuadNodes1d*numQuadNodes1d*numQuadNodes1d==numQuadNodes) {
246 jgs 82
247 gross 1028 /* get 1D scheme: */
248 jgs 82
249 gross 1028 Finley_Quad_getNodesLine(numQuadNodes1d,quadNodes1d,quadWeights1d);
250 jgs 82
251 gross 1028 /* make 3D scheme: */
252 jgs 82
253 gross 1028 if (Finley_noError()) {
254     l=0;
255     for (i=0;i<numQuadNodes1d;i++) {
256     for (j=0;j<numQuadNodes1d;j++) {
257     for (k=0;k<numQuadNodes1d;k++) {
258     QUADNODES(0,l)=quadNodes1d[i];
259     QUADNODES(1,l)=quadNodes1d[j];
260     QUADNODES(2,l)=quadNodes1d[k];
261     QUADWEIGHTS(l)=quadWeights1d[i]*quadWeights1d[j]*quadWeights1d[k];
262     l++;
263     }
264     }
265     }
266     set=TRUE;
267     break;
268     }
269     }
270     }
271     if (!set) {
272     sprintf(error_msg,"Finley_Quad_getNodesHex: Illegal number of quadrature nodes %d on hexahedron.",numQuadNodes);
273     Finley_setError(VALUE_ERROR,error_msg);
274     }
275     TMPMEMFREE(quadNodes1d);
276     TMPMEMFREE(quadWeights1d);
277 jgs 82 }
278     #undef DIM
279     }
280    
281     /**************************************************************/
282    
283     /* get a quadrature scheme with numQuadNodes quadrature nodes for a point. As there */
284     /* in no quadrature scheme for a point any value for numQuadNodes other than 0 throws */
285     /* an error. */
286    
287     void Finley_Quad_getNodesPoint(int numQuadNodes,double* quadNodes,double* quadWeights) {
288     if (numQuadNodes!=0) {
289 gross 964 Finley_setError(VALUE_ERROR,"Finley_Quad_getNodesPoint: There is no quadrature scheme on points.");
290 jgs 82 }
291     }
292    
293     /**************************************************************/
294    
295     /* get a quadrature scheme with numQuadNodes quadrature nodes on the line [0,1]: */
296     /* The nodes and weights are set from a table. */
297    
298     void Finley_Quad_getNodesLine(int numQuadNodes,double* quadNodes,double* quadWeights) {
299     switch(numQuadNodes) {
300     case 1:
301     quadNodes[0]=0.5;
302     quadWeights[0]=1.;
303     break;
304    
305     case 2:
306     quadNodes[0]=(1.-.577350269189626)/2.;
307     quadNodes[1]=(1.+.577350269189626)/2.;
308     quadWeights[0]=.5;
309     quadWeights[1]=.5;
310     break;
311    
312     case 3:
313     quadNodes[0]=(1.-.774596669241483)/2.;
314     quadNodes[1]=.5;
315     quadNodes[2]=(1.+.774596669241483)/2.;
316     quadWeights[0]=5./18.;
317     quadWeights[1]=4./ 9.;
318     quadWeights[2]=5./18.;
319     break;
320    
321     case 4:
322     quadNodes[0]=(1.-.861136311594053)/2.;
323     quadNodes[1]=(1.-.339981043584856)/2.;
324     quadNodes[2]=(1.+.339981043584856)/2.;
325     quadNodes[3]=(1.+.861136311594053)/2.;
326     quadWeights[0]=.347854845137454/2.;
327     quadWeights[1]=.652145154862546/2.;
328     quadWeights[2]=.652145154862546/2.;
329     quadWeights[3]=.347854845137454/2.;
330     break;
331    
332     case 5:
333     quadNodes[0]=(1.-.906179845938664)/2.;
334     quadNodes[1]=(1.-.538469310105683)/2.;
335     quadNodes[2]= .5;
336     quadNodes[3]=(1.+.538469310105683)/2.;
337     quadNodes[4]=(1.+.906179845938664)/2.;
338     quadWeights[0]=.236926885056189/2.;
339     quadWeights[1]=.478628670499366/2.;
340     quadWeights[2]=.568888888888889/2.;
341     quadWeights[3]=.478628670499366/2.;
342     quadWeights[4]=.236926885056189/2.;
343     break;
344    
345     case 6:
346     quadNodes[0]=(1.-.932469514203152)/2.;
347     quadNodes[1]=(1.-.661209386466265)/2.;
348     quadNodes[2]=(1.-.238619186083197)/2.;
349     quadNodes[3]=(1.+.238619186083197)/2.;
350     quadNodes[4]=(1.+.661209386466265)/2.;
351     quadNodes[5]=(1.+.932469514203152)/2.;
352     quadWeights[0]=.171324492379170/2.;
353     quadWeights[1]=.360761573048139/2.;
354     quadWeights[2]=.467913934572691/2.;
355     quadWeights[3]=.467913934572691/2.;
356     quadWeights[4]=.360761573048139/2.;
357     quadWeights[5]=.171324492379170/2.;
358     break;
359    
360     case 7:
361     quadNodes[0]=(1.-.949107912342759)/2.;
362     quadNodes[1]=(1.-.741531185599394)/2.;
363     quadNodes[2]=(1.-.405845151377397)/2.;
364     quadNodes[3]=0.5;
365     quadNodes[4]=(1.+.405845151377397)/2.;
366     quadNodes[5]=(1.+.741531185599394)/2.;
367     quadNodes[6]=(1.+.949107912342759)/2.;
368     quadWeights[0]= .129484966168870/2.;
369     quadWeights[1]= .279705391489277/2.;
370     quadWeights[2]= .381830050505119/2.;
371     quadWeights[3]= .417959183673469/2.;
372     quadWeights[4]= .381830050505119/2.;
373     quadWeights[5]= .279705391489277/2.;
374     quadWeights[6]= .129484966168870/2.;
375     break;
376    
377     case 8:
378     quadNodes[0]=(1.-.960289856497536)/2.;
379     quadNodes[1]=(1.-.796666477413627)/2.;
380     quadNodes[2]=(1.-.525532409916329)/2.;
381     quadNodes[3]=(1.-.183434642495650)/2.;
382     quadNodes[4]=(1.+.183434642495650)/2.;
383     quadNodes[5]=(1.+.525532409916329)/2.;
384     quadNodes[6]=(1.+.796666477413627)/2.;
385     quadNodes[7]=(1.+.960289856497536)/2.;
386     quadWeights[0]= .101228536290376/2.;
387     quadWeights[1]= .222381034453374/2.;
388     quadWeights[2]= .313706645877887/2.;
389     quadWeights[3]= .362683783378362/2.;
390     quadWeights[4]= .362683783378362/2.;
391     quadWeights[5]= .313706645877887/2.;
392     quadWeights[6]= .222381034453374/2.;
393     quadWeights[7]= .101228536290376/2.;
394     break;
395    
396     case 9:
397     quadNodes[0]=(1.-.968160239507626)/2.;
398     quadNodes[1]=(1.-.836031107326636)/2.;
399     quadNodes[2]=(1.-.613371432700590)/2.;
400     quadNodes[3]=(1.-.324253423403809)/2.;
401     quadNodes[4]= .5;
402     quadNodes[5]=(1.+.324253423403809)/2.;
403     quadNodes[6]=(1.+.613371432700590)/2.;
404     quadNodes[7]=(1.+.836031107326636)/2.;
405     quadNodes[8]=(1.+.968160239507626)/2.;
406     quadWeights[0]= .081274388361574/2.;
407     quadWeights[1]= .180648160694857/2.;
408     quadWeights[2]= .260610696402935/2.;
409     quadWeights[3]= .312347077040003/2.;
410     quadWeights[4]= .330239355001260/2.;
411     quadWeights[5]= .312347077040003/2.;
412     quadWeights[6]= .260610696402935/2.;
413     quadWeights[7]= .180648160694857/2.;
414     quadWeights[8]= .081274388361574/2.;
415     break;
416    
417     case 10:
418     quadNodes[0]=(1.-.973906528517172)/2.;
419     quadNodes[1]=(1.-.865063366688985)/2.;
420     quadNodes[2]=(1.-.679409568299024)/2.;
421     quadNodes[3]=(1.-.433395394129247)/2.;
422     quadNodes[4]=(1.-.148874338981631)/2.;
423     quadNodes[5]=(1.+.148874338981631)/2.;
424     quadNodes[6]=(1.+.433395394129247)/2.;
425     quadNodes[7]=(1.+.679409568299024)/2.;
426     quadNodes[8]=(1.+.865063366688985)/2.;
427     quadNodes[9]=(1.+.973906528517172)/2.;
428     quadWeights[0]= .066671344308688/2.;
429     quadWeights[1]= .149451349150581/2.;
430     quadWeights[2]= .219086362515982/2.;
431     quadWeights[3]= .269266719309996/2.;
432     quadWeights[4]= .295524224714753/2.;
433     quadWeights[5]= .295524224714753/2.;
434     quadWeights[6]= .269266719309996/2.;
435     quadWeights[7]= .219086362515982/2.;
436     quadWeights[8]= .149451349150581/2.;
437     quadWeights[9]= .066671344308688/2.;
438     break;
439    
440     default:
441 jgs 150 Finley_setError(VALUE_ERROR,"__FILE__: Negative intergration order.");
442 jgs 82 break;
443     }
444     }
445     /**************************************************************/
446     /* */
447     /* the following function are used define the meshes on the surface in the xy-plane */
448    
449     /* triangle surface on a tetrahedron */
450     void Finley_Quad_getNodesTriOnFace(int numQuadNodes,double* quadNodes,double* quadWeights) {
451     Finley_Quad_makeNodesOnFace(3,numQuadNodes,quadNodes,quadWeights,Finley_Quad_getNodesTri);
452     }
453     /* rectangular surface on a hexahedron */
454     void Finley_Quad_getNodesRecOnFace(int numQuadNodes,double* quadNodes,double* quadWeights) {
455     Finley_Quad_makeNodesOnFace(3,numQuadNodes,quadNodes,quadWeights,Finley_Quad_getNodesRec);
456     }
457     /* line surface on a triangle or rectangle */
458     void Finley_Quad_getNodesLineOnFace(int numQuadNodes,double* quadNodes,double* quadWeights) {
459     Finley_Quad_makeNodesOnFace(2,numQuadNodes,quadNodes,quadWeights,Finley_Quad_getNodesLine);
460     }
461     /* point surface on a line */
462     void Finley_Quad_getNodesPointOnFace(int numQuadNodes,double* quadNodes,double* quadWeights) {
463     Finley_Quad_makeNodesOnFace(1,numQuadNodes,quadNodes,quadWeights,Finley_Quad_getNodesPoint);
464     }
465    
466     void Finley_Quad_makeNodesOnFace(int dim, int numQuadNodes,double* quadNodes,double* quadWeights, Finley_Quad_getNodes getFaceNodes) {
467     int q,i;
468 gross 1028 double *quadNodesOnFace=NULL;
469     #define DIM dim
470     quadNodesOnFace=TMPMEMALLOC(numQuadNodes*(dim-1),double);
471 jgs 82
472 gross 1028 if (! Finley_checkPtr(quadNodesOnFace) ) {
473     getFaceNodes(numQuadNodes,quadNodesOnFace,quadWeights);
474     if (Finley_noError()) {
475     for (q=0;q<numQuadNodes;q++) {
476     for (i=0;i<dim-1;i++) QUADNODES(i,q)=quadNodesOnFace[INDEX2(i,q,dim-1)];
477     QUADNODES(dim-1,q)=0;
478     }
479     }
480     TMPMEMFREE(quadNodesOnFace);
481 jgs 82 }
482     #undef DIM
483     }
484    
485     /**************************************************************/
486    
487     /* The following functions Finley_Quad_getNumNodes* return the nmber of quadrature points needed to */
488     /* achieve a certain accuracy. Notice that for Tet and Tri the the order is increased */
489     /* to consider the accuracy reduction through the construction process. */
490    
491    
492     int Finley_Quad_getNumNodesPoint(int order) {
493     if (order <0 ) {
494 gross 964 Finley_setError(VALUE_ERROR,"Finley_Quad_getNumNodesPoint: Negative intergration order.");
495 jgs 82 return -1;
496     } else {
497     return 0;
498     }
499     }
500    
501     int Finley_Quad_getNumNodesLine(int order) {
502 jgs 150 char error_msg[LenErrorMsg_MAX];
503 jgs 82 if (order <0 ) {
504 gross 964 Finley_setError(VALUE_ERROR,"Finley_Quad_getNumNodesPoint: Negative intergration order.");
505 jgs 82 return -1;
506     } else {
507     if (order > 2*MAX_numQuadNodesLine-1) {
508 gross 964 sprintf(error_msg,"Finley_Quad_getNumNodesPoint: requested integration order %d on line is too large (>%d).",
509 jgs 82 order,2*MAX_numQuadNodesLine-1);
510 jgs 150 Finley_setError(VALUE_ERROR,error_msg);
511 jgs 82 return -1;
512     } else {
513 jgs 150 Finley_resetError();
514 jgs 82 return order/2+1;
515     }
516     }
517     }
518    
519     int Finley_Quad_getNumNodesTri(int order) {
520     int numQuadNodesLine;
521 gross 1072 if (order<=1) {
522 jgs 82 return 1;
523 btully 1170 } else if (order==2){
524     return 3;
525     } else if (order==3){
526     return 4;
527 jgs 82 } else {
528     numQuadNodesLine=Finley_Quad_getNumNodesLine(order+1);
529 jgs 150 if (Finley_noError()) {
530 jgs 82 return numQuadNodesLine*numQuadNodesLine;
531     } else {
532     return -1;
533     }
534     }
535     }
536    
537     int Finley_Quad_getNumNodesRec(int order) {
538     int numQuadNodesLine;
539     numQuadNodesLine=Finley_Quad_getNumNodesLine(order);
540 jgs 150 if (Finley_noError()) {
541 jgs 82 return numQuadNodesLine*numQuadNodesLine;
542     } else {
543     return -1;
544     }
545     }
546    
547     int Finley_Quad_getNumNodesTet(int order) {
548     int numQuadNodesLine;
549 gross 1072 if (order<=1) {
550     return 1;
551 btully 1170 } else if (order==2){
552     return 4;
553     } else if (order==3){
554     return 5;
555 jgs 82 } else {
556 gross 1072 numQuadNodesLine=Finley_Quad_getNumNodesLine(order+2);
557     if (Finley_noError()) {
558     return numQuadNodesLine*numQuadNodesLine*numQuadNodesLine;
559     } else {
560     return -1;
561     }
562 jgs 82 }
563     }
564    
565     int Finley_Quad_getNumNodesHex(int order) {
566     int numQuadNodesLine;
567     numQuadNodesLine=Finley_Quad_getNumNodesLine(order);
568 jgs 150 if (Finley_noError()) {
569 jgs 82 return numQuadNodesLine*numQuadNodesLine*numQuadNodesLine;
570     } else {
571     return -1;
572     }
573     }
574     /*
575     * $Log$
576 jgs 150 * Revision 1.2 2005/09/15 03:44:23 jgs
577     * Merge of development branch dev-02 back to main trunk on 2005-09-15
578 jgs 82 *
579 jgs 150 * Revision 1.1.1.1.6.1 2005/09/07 06:26:20 gross
580     * the solver from finley are put into the standalone package paso now
581     *
582     * Revision 1.1.1.1 2004/10/26 06:53:57 jgs
583     * initial import of project esys2
584     *
585 jgs 82 * Revision 1.2 2004/08/03 04:49:06 gross
586     * bug in Quadrature.c fixed
587     *
588     * Revision 1.1.1.1 2004/06/24 04:00:40 johng
589     * Initial version of eys using boost-python.
590     *
591     *
592     */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26