ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/outil/src/ot_quadrature_gauss.cpp
Revision: 885
Committed: Tue May 23 21:07:21 2017 UTC (8 years, 1 month ago) by francois
Original Path: magic/lib/outil/src/ot_quadrature_gauss.cpp
File size: 19513 byte(s)
Log Message: 
points de gauss pour les pentaèdres . Cas des 6 points de gauss

File Contents

# User Rev Content
1 francois 283 //---------------------------------------------------------------------------
2    
3     #include <math.h>
4     #include <stdio.h>
5     #define EPS 3.0e-11
6    
7     #pragma hdrstop
8    
9     #include "gestionversion.h"
10     #include "ot_quadrature_gauss.h"
11     //---------------------------------------------------------------------------
12     #pragma package(smart_init)
13    
14     void OT_QUADRATURE_GAUSS::gauss_legendre_points(double x1, double x2, double x[], double w[], int n)
15     {
16     int m,j,i;
17     double z1,z,xm,xl,pp,p3,p2,p1;
18     m=(n+1)/2;
19     xm=0.5*(x2+x1);
20     xl=0.5*(x2-x1);
21     for (i=1;i<=m;i++) {
22     z=cos(M_PI*(i-0.25)/(n+0.5));
23     do {
24     p1=1;
25     p2=0;
26     for (j=1;j<=n;j++) {
27     p3=p2;
28     p2=p1;
29     p1=((2*j-1)*z*p2-(j-1)*p3)/j;
30     }
31     pp=n*(z*p1-p2)/(z*z-1);
32     z1=z;
33     z=z1-p1/pp;
34     } while (fabs(z-z1) > EPS);
35     x[i]=xm-xl*z;
36     x[n+1-i]=xm+xl*z;
37     w[i]=2*xl/((1-z*z)*pp*pp);
38     w[n+1-i]=w[i];
39     }
40     }
41 francois 754
42    
43    
44    
45    
46     int OT_POINTS_GAUSS::get_nb_point_seg(int degre)
47     {
48     if (degre<2) return 1;
49     if (degre<4) return 2;
50     if (degre<6) return 3;
51     if (degre<8) return 4;
52     if (degre<10) return 5;
53     if (degre<12) return 6;
54     if (degre<14) return 7;
55     return 0;
56     }
57    
58 francois 757
59    
60     int OT_POINTS_GAUSS::get_degre_gauss_seg(int num)
61     {
62     if (num<2) return 1;
63     if (num<3) return 3;
64     if (num<4) return 7;
65     if (num<5) return 9;
66     if (num<6) return 11;
67     if (num<7) return 13;
68     if (num<8) return 15;
69     return 0;
70     }
71    
72    
73    
74 francois 754 int OT_POINTS_GAUSS::get_nb_point_tri(int degre)
75     {
76     if (degre<2) return 1;
77     if (degre<3) return 3;
78     if (degre<4) return 4;
79     if (degre<5) return 6;
80     if (degre<6) return 7;
81 francois 757 return 0;
82     }
83 francois 754
84 francois 757 int OT_POINTS_GAUSS::get_degre_gauss_tri(int num)
85     {
86     if (num<2) return 1;
87     if (num<4) return 2;
88     if (num<6) return 3;
89     if (num<7) return 4;
90     if (num<8) return 5;
91     return 0;
92     }
93    
94 francois 754 int OT_POINTS_GAUSS::get_nb_point_quad(int degre)
95     {
96     if (degre<3) return 3;
97     if (degre<4) return 4;
98     if (degre<6) return 7;
99     return 0;
100     }
101 francois 757
102     int OT_POINTS_GAUSS::get_degre_gauss_quad(int num)
103     {
104     if (num<4) return 2;
105     if (num<5) return 3;
106     if (num<8) return 5;
107     return 0;
108     }
109    
110    
111 francois 754 int OT_POINTS_GAUSS::get_nb_point_quad_prod(int degre)
112     {
113     if (degre<2) return 1;
114     if (degre<4) return 4;
115     if (degre<6) return 9;
116     return 0;
117     }
118 francois 757
119     int OT_POINTS_GAUSS::get_degre_gauss_quad_prod(int num)
120     {
121     if (num<2) return 1;
122     if (num<5) return 3;
123     if (num<10) return 5;
124     return 0;
125     }
126    
127 francois 754 int OT_POINTS_GAUSS::get_nb_point_tetra(int degre)
128     {
129     if (degre<2) return 1;
130     if (degre<3) return 4;
131     if (degre<4) return 5;
132     if (degre<6) return 15;
133 francois 757 return 0;
134     }
135 francois 754
136 francois 757 int OT_POINTS_GAUSS::get_degre_gauss_tetra(int num)
137     {
138     if (num<2) return 1;
139     if (num<5) return 2;
140     if (num<6) return 3;
141     if (num<16) return 5;
142     return 0;
143     }
144    
145 francois 754 int OT_POINTS_GAUSS::get_nb_point_hexa(int degre)
146     {
147     if (degre<3) return 4;
148     if (degre<4) return 6;
149     if (degre<6) return 14;
150     return 0;
151     }
152    
153 francois 757 int OT_POINTS_GAUSS::get_degre_gauss_hexa(int num)
154     {
155     if (num<5) return 2;
156     if (num<7) return 3;
157     if (num<15) return 5;
158     return 0;
159     }
160    
161 francois 754 int OT_POINTS_GAUSS::get_nb_point_hexa_prod(int degre)
162     {
163     if (degre<2) return 1;
164     if (degre<4) return 8;
165     if (degre<6) return 27;
166     return 0;
167     }
168    
169 francois 757 int OT_POINTS_GAUSS::get_degre_gauss_hexa_prod(int num)
170     {
171     if (num<2) return 1;
172     if (num<9) return 3;
173     if (num<28) return 5;
174     return 0;
175     }
176    
177 francois 876 int OT_POINTS_GAUSS::get_nb_point_penta(int degre)
178     {
179     if (degre<3) return 6;
180     if (degre<4) return 8;
181     if (degre<6) return 21;
182     return 0;
183     }
184 francois 757
185 francois 876 int OT_POINTS_GAUSS::get_degre_gauss_penta(int num)
186     {
187     if (num<7) return 2;
188     if (num<9) return 3;
189     if (num<22) return 5;
190     return 0;
191     }
192 francois 757
193    
194 francois 876
195 francois 754 void OT_POINTS_GAUSS::get_pt_gauss_seg(int degre,int num,double &w,double *u)
196     {
197     if (degre<2)
198     {
199     if (num==0) {w=2;u[0]=0;}
200     return;
201     }
202     if (degre<4)
203     {
204     if (num==0) {w=1;u[0]=0.577350269189626;}
205     if (num==1) {w=1;u[0]=-0.577350269189626;}
206     return;
207     }
208     if (degre<6)
209     {
210     if (num==0) {w=0.555555555555556;u[0]=-0.774596669241483;}
211     if (num==1) {w=0.888888888888889;u[0]=0.;}
212     if (num==2) {w=0.555555555555556;u[0]=0.774596669241483;}
213     return;
214     }
215     if (degre<8)
216     {
217     if (num==0) {w=0.652145154862546;u[0]=0.339981043584856;}
218     if (num==1) {w=0.652145154862546;u[0]=-0.339981043584856;}
219     if (num==2) {w=0.347854845137454;u[0]=0.861136311594053;}
220     if (num==3) {w=0.347854845137454;u[0]=-0.861136311594053;}
221     return;
222     }
223     // limite conformite code aster
224     if (degre<10)
225     {
226     if (num==0) {w=0.568888888888889;u[0]=0.;}
227     if (num==1) {w=0.478628670499366;u[0]=-0.538469310105683;}
228     if (num==2) {w=0.478628670499366;u[0]=0.538469310105683;}
229     if (num==3) {w=0.236926885056189;u[0]=-0.906179845938664;}
230     if (num==4) {w=0.236926885056189;u[0]=0.906179845938664;}
231     return;
232     }
233     if (degre<12)
234     {
235     if (num==0) {w=0.467913934572691;u[0]=-0.238619186083197;}
236     if (num==1) {w=0.467913934572691;u[0]=0.238619186083197;}
237     if (num==2) {w=0.360761573048139;u[0]=-0.661209386466265;}
238     if (num==3) {w=0.360761573048139;u[0]=0.661209386466265;}
239     if (num==4) {w=0.171324492379170;u[0]=-0.932469514203152;}
240     if (num==5) {w=0.171324492379170;u[0]=0.932469514203152;}
241     return;
242     }
243     if (degre<14)
244     {
245     if (num==0) {w=0.417959183673469;u[0]=0;}
246     if (num==1) {w=0.381830050505119;u[0]=-0.405845151377397;}
247     if (num==2) {w=0.381830050505119;u[0]=0.405845151377397;}
248     if (num==3) {w=0.279705391489277;u[0]=-0.741531185599394;}
249     if (num==4) {w=0.279705391489277;u[0]=0.741531185599394;}
250     if (num==5) {w=0.129484966168870;u[0]=-0.949107912342759;}
251     if (num==6) {w=0.129484966168870;u[0]=0.949107912342759;}
252     return;
253     }
254     }
255    
256    
257     void OT_POINTS_GAUSS::get_pt_gauss_tri(int degre, int num, double& w, double* uv)
258     {
259     if (degre<2)
260     {
261     if (num==0) {w=0.5;uv[0]=0.333333333333333;uv[1]=0.333333333333333;}
262     return;
263     }
264     if (degre<3)
265     {
266     if (num==0) {w=0.166666666666667;uv[0]=0.166666666666667;uv[1]=0.166666666666667;}
267     if (num==1) {w=0.166666666666667;uv[0]=0.666666666666667;uv[1]=0.166666666666667;}
268     if (num==2) {w=0.166666666666667;uv[0]=0.166666666666667;uv[1]=0.666666666666667;}
269     return;
270     }
271     if (degre<4)
272     {
273     if (num==0) {w=0.260416666666667;uv[0]=0.2;uv[1]=0.2;}
274     if (num==1) {w=0.260416666666667;uv[0]=0.6;uv[1]=0.2;}
275     if (num==2) {w=0.260416666666667;uv[0]=0.2;uv[1]=0.6;}
276     if (num==3) {w=-0.28125;uv[0]=0.333333333333333;uv[1]=0.333333333333333;}
277     return;
278     }
279     if (degre<5)
280     {
281     if (num==0) {w=0.054975871827661;uv[0]=0.091576213509771;uv[1]=0.091576213509771;}
282     if (num==1) {w=0.054975871827661;uv[0]=0.816847572980458;uv[1]=0.091576213509771;}
283     if (num==2) {w=0.054975871827661;uv[0]=0.091576213509771;uv[1]=0.816847572980458;}
284     if (num==3) {w=0.111690794839005;uv[0]=0.445948490915965;uv[1]=0.445948490915965;}
285     if (num==4) {w=0.111690794839005;uv[0]=0.445948490915965;uv[1]=0.108103018168070;}
286     if (num==5) {w=0.111690794839005;uv[0]=0.108103018168070;uv[1]=0.445948490915965;}
287     return;
288     }
289     if (degre<6)
290     {
291     if (num==0) {w=0.1125;uv[0]=0.333333333333333;uv[1]=0.333333333333333;}
292     if (num==1) {w=0.066197076394253;uv[0]=0.470142064105115;uv[1]=0.470142064105115;}
293     if (num==2) {w=0.066197076394253;uv[0]=0.059715871789770;uv[1]=0.470142064105115;}
294     if (num==3) {w=0.066197076394253;uv[0]=0.470142064105115;uv[1]=0.059715871789770;}
295     if (num==4) {w=0.062969590272413;uv[0]=0.101286507323456;uv[1]=0.101286507323456;}
296     if (num==5) {w=0.062969590272413;uv[0]=0.797426985353088;uv[1]=0.101286507323456;}
297     if (num==6) {w=0.062969590272413;uv[0]=0.101286507323456;uv[1]=0.797426985353088;}
298     return;
299     }
300     }
301    
302 francois 757
303    
304    
305    
306 francois 754 void OT_POINTS_GAUSS::get_pt_gauss_tet(int degre, int num, double& w, double* uvw)
307     {
308     if (degre<2)
309     {
310     if (num==0) {w=0.166666666666667;uvw[0]=0.25;uvw[1]=0.25;uvw[2]=0.25;}
311     return;
312     }
313     if (degre<3)
314     {
315     if (num==0) {w=0.041666666666667;uvw[0]=0.138196601125011;uvw[1]=0.138196601125011;uvw[2]=0.138196601125011;}
316     if (num==1) {w=0.041666666666667;uvw[0]=0.138196601125011;uvw[1]=0.138196601125011;uvw[2]=0.585410196624968;}
317     if (num==2) {w=0.041666666666667;uvw[0]=0.138196601125011;uvw[1]=0.585410196624968;uvw[2]=0.138196601125011;}
318     if (num==3) {w=0.041666666666667;uvw[0]=0.585410196624968;uvw[1]=0.138196601125011;uvw[2]=0.138196601125011;}
319     return;
320     }
321     if (degre<4)
322     {
323     if (num==0) {w=-0.133333333333333;uvw[0]=0.25;uvw[1]=0.25;uvw[2]=0.25;}
324     if (num==1) {w=0.075;uvw[0]=0.166666666666667;uvw[1]=0.166666666666667;uvw[2]=0.166666666666667;}
325     if (num==2) {w=0.075;uvw[0]=0.166666666666667;uvw[1]=0.166666666666667;uvw[2]=0.5;}
326     if (num==3) {w=0.075;uvw[0]=0.166666666666667;uvw[1]=0.5;uvw[2]=0.166666666666667;}
327     if (num==4) {w=0.075;uvw[0]=0.5;uvw[1]=0.166666666666667;uvw[2]=0.166666666666667;}
328     return;
329     }
330     if (degre<6)
331     {
332     double a=0.25;
333     double b1=0.3197936278;
334     double b2=0.0919710781;
335     double c2=0.7240867658;
336     double c1=0.0406191165;
337     double d=0.0563508327;
338     double e=0.4436491673;
339     double w1=0.0197530864;
340     double w3=0.011989514;
341     double w2=0.0115113679;
342     double w4=0.0088183422;
343     if (num==0) {w=w1;uvw[0]=a;uvw[1]=a;uvw[2]=a;}
344     if (num==1) {w=w2;uvw[0]=b1;uvw[1]=b1;uvw[2]=b1;}
345     if (num==2) {w=w2;uvw[0]=b1;uvw[1]=b1;uvw[2]=c1;}
346     if (num==3) {w=w2;uvw[0]=b1;uvw[1]=c1;uvw[2]=b1;}
347     if (num==4) {w=w2;uvw[0]=c1;uvw[1]=b1;uvw[2]=b1;}
348     if (num==5) {w=w3;uvw[0]=b2;uvw[1]=b2;uvw[2]=b2;}
349     if (num==6) {w=w3;uvw[0]=b2;uvw[1]=b2;uvw[2]=c2;}
350     if (num==7) {w=w3;uvw[0]=b2;uvw[1]=c2;uvw[2]=b2;}
351     if (num==8) {w=w3;uvw[0]=c2;uvw[1]=b2;uvw[2]=b2;}
352     if (num==9) {w=w4;uvw[0]=d;uvw[1]=d;uvw[2]=e;}
353     if (num==10) {w=w4;uvw[0]=d;uvw[1]=e;uvw[2]=d;}
354     if (num==11) {w=w4;uvw[0]=e;uvw[1]=d;uvw[2]=d;}
355     if (num==12) {w=w4;uvw[0]=d;uvw[1]=e;uvw[2]=e;}
356     if (num==13) {w=w4;uvw[0]=e;uvw[1]=d;uvw[2]=e;}
357     if (num==14) {w=w4;uvw[0]=e;uvw[1]=e;uvw[2]=d;}
358     return;
359     }
360     }
361    
362     void OT_POINTS_GAUSS::get_pt_gauss_qua(int degre, int num, double& w, double* uv)
363     {
364     if (degre<3)
365     {
366     if (num==0) {w=1.333333333333333;uv[0]=0.816496580927726;uv[1]=0.;}
367     if (num==1) {w=1.333333333333333;uv[0]=-0.408248290463863;uv[1]=-0.707106781186548;}
368     if (num==2) {w=1.333333333333333;uv[0]=-0.408248290463863;uv[1]=0.707106781186548;}
369     return;
370     }
371     if (degre<4)
372     {
373     if (num==0) {w=1.;uv[0]=0.816496580927726;uv[1]=0.;}
374     if (num==1) {w=1.;uv[0]=-0.816496580927726;uv[1]=0.;}
375     if (num==2) {w=1.;uv[0]=0.;uv[1]=0.816496580927726;}
376     if (num==3) {w=1.;uv[0]=0.;uv[1]=-0.816496580927726;}
377     return;
378     }
379     if (degre<6)
380     {
381     if (num==0) {w=1.142857142857140;uv[0]=0.;uv[1]=0.;}
382     if (num==1) {w=0.317460317460317;uv[0]=0.;uv[1]=0.966091783079296;}
383     if (num==2) {w=0.317460317460317;uv[0]=0.;uv[1]=-0.966091783079296;}
384     if (num==3) {w=0.555555555555556;uv[0]=0.774596669241483;uv[1]=-0.577350269189626;}
385     if (num==4) {w=0.555555555555556;uv[0]=-0.774596669241483;uv[1]=-0.577350269189626;}
386     if (num==5) {w=0.555555555555556;uv[0]=-0.774596669241483;uv[1]=0.577350269189626;}
387     if (num==6) {w=0.555555555555556;uv[0]=0.774596669241483;uv[1]=0.577350269189626;}
388     return;
389     }
390     }
391    
392     void OT_POINTS_GAUSS::get_pt_gauss_qua_prod(int degre, int num, double& w, double* uv)
393     {
394     if (degre<2)
395     {
396     if (num==0) {w=4.0;uv[0]=0.;uv[1]=0.;}
397     return;
398     }
399     if (degre<4)
400     {
401     if (num==0) {w=1.;uv[0]=-0.577350269189626;uv[1]=-0.577350269189626;}
402 francois 876 if (num==1) {w=1.;uv[0]=0.577350269189626;uv[1]=-0.577350269189626;}
403     if (num==2) {w=1.;uv[0]=0.577350269189626;uv[1]=0.577350269189626;}
404     if (num==3) {w=1.;uv[0]=-0.577350269189626;uv[1]=0.577350269189626;}
405 francois 754 return;
406     }
407     if (degre<6)
408     {
409     double a=0.774596669241483;
410     if (num==0) {w=25./81.;uv[0]=-a;uv[1]=-a;}
411     if (num==1) {w=25./81.;uv[0]=a;uv[1]=-a;}
412     if (num==2) {w=25./81.;uv[0]=a;uv[1]=a;}
413     if (num==3) {w=25./81.;uv[0]=-a;uv[1]=a;}
414     if (num==4) {w=40./81.;uv[0]=0.;uv[1]=-a;}
415     if (num==5) {w=40./81.;uv[0]=a;uv[1]=0.;}
416     if (num==6) {w=40./81.;uv[0]=0.;uv[1]=a;}
417     if (num==7) {w=40./81.;uv[0]=-a;uv[1]=0.;}
418     if (num==8) {w=64./81.;uv[0]=0.;uv[1]=0.;}
419     return;
420     }
421     }
422    
423    
424     void OT_POINTS_GAUSS::get_pt_gauss_hex_prod(int degre, int num, double& w, double* uvw)
425     {
426     if (degre<2)
427     {
428     if (num==0) {w=8.0;uvw[0]=0.;uvw[1]=0.;uvw[2]=0.;}
429     return;
430     }
431     if (degre<4)
432     {
433     double a=0.577350269189626;
434     if (num==0) {w=1.0;uvw[0]=-a;uvw[1]=-a;uvw[2]=-a;}
435     if (num==1) {w=1.0;uvw[0]=-a;uvw[1]=-a;uvw[2]=a;}
436     if (num==2) {w=1.0;uvw[0]=-a;uvw[1]=a;uvw[2]=-a;}
437     if (num==3) {w=1.0;uvw[0]=-a;uvw[1]=a;uvw[2]=a;}
438     if (num==4) {w=1.0;uvw[0]=a;uvw[1]=-a;uvw[2]=-a;}
439     if (num==5) {w=1.0;uvw[0]=a;uvw[1]=-a;uvw[2]=a;}
440     if (num==6) {w=1.0;uvw[0]=a;uvw[1]=a;uvw[2]=-a;}
441     if (num==7) {w=1.0;uvw[0]=a;uvw[1]=a;uvw[2]=a;}
442     return;
443     }
444     if (degre<6)
445     {
446     double a=0.774596669241483;
447     double c1=0.555555555555556;
448     double c2=0.888888888888889;
449     if (num==0) {w=c1*c1*c1;uvw[0]=-a;uvw[1]=-a;uvw[2]=-a;}
450     if (num==1) {w=c1*c1*c2;uvw[0]=-a;uvw[1]=-a;uvw[2]=0.;}
451     if (num==2) {w=c1*c1*c1;uvw[0]=-a;uvw[1]=-a;uvw[2]=a;}
452     if (num==3) {w=c1*c1*c2;uvw[0]=-a;uvw[1]=0.;uvw[2]=-a;}
453     if (num==4) {w=c1*c2*c2;uvw[0]=-a;uvw[1]=0.;uvw[2]=0.;}
454     if (num==5) {w=c1*c1*c2;uvw[0]=-a;uvw[1]=0.;uvw[2]=a;}
455     if (num==6) {w=c1*c1*c1;uvw[0]=-a;uvw[1]=a;uvw[2]=-a;}
456     if (num==7) {w=c1*c1*c2;uvw[0]=-a;uvw[1]=a;uvw[2]=0.;}
457     if (num==8) {w=c1*c1*c1;uvw[0]=-a;uvw[1]=a;uvw[2]=a;}
458     if (num==9) {w=c1*c1*c2;uvw[0]=0.;uvw[1]=-a;uvw[2]=-a;}
459     if (num==10) {w=c1*c2*c2;uvw[0]=0.;uvw[1]=-a;uvw[2]=0.;}
460     if (num==11) {w=c1*c1*c2;uvw[0]=0.;uvw[1]=-a;uvw[2]=a;}
461     if (num==12) {w=c1*c2*c2;uvw[0]=0.;uvw[1]=0.;uvw[2]=-a;}
462     if (num==13) {w=c2*c2*c2;uvw[0]=0.;uvw[1]=0.;uvw[2]=0.;}
463     if (num==14) {w=c1*c2*c2;uvw[0]=0.;uvw[1]=0.;uvw[2]=a;}
464     if (num==15) {w=c1*c1*c2;uvw[0]=0.;uvw[1]=a;uvw[2]=-a;}
465     if (num==16) {w=c1*c2*c2;uvw[0]=0.;uvw[1]=a;uvw[2]=0.;}
466     if (num==17) {w=c1*c1*c2;uvw[0]=0.;uvw[1]=a;uvw[2]=a;}
467     if (num==18) {w=c1*c1*c1;uvw[0]=a;uvw[1]=-a;uvw[2]=-a;}
468     if (num==19) {w=c1*c1*c2;uvw[0]=a;uvw[1]=-a;uvw[2]=0.;}
469     if (num==20) {w=c1*c1*c1;uvw[0]=a;uvw[1]=-a;uvw[2]=a;}
470     if (num==21) {w=c1*c1*c2;uvw[0]=a;uvw[1]=0.;uvw[2]=-a;}
471     if (num==22) {w=c1*c2*c2;uvw[0]=a;uvw[1]=0.;uvw[2]=0.;}
472     if (num==23) {w=c1*c1*c2;uvw[0]=a;uvw[1]=0.;uvw[2]=a;}
473     if (num==24) {w=c1*c1*c1;uvw[0]=a;uvw[1]=a;uvw[2]=-a;}
474     if (num==25) {w=c1*c1*c2;uvw[0]=a;uvw[1]=a;uvw[2]=0.;}
475     if (num==26) {w=c1*c1*c1;uvw[0]=a;uvw[1]=a;uvw[2]=a;}
476     return;
477     }
478     }
479    
480    
481    
482     void OT_POINTS_GAUSS::get_pt_gauss_hex(int degre, int num, double& w, double* uvw)
483     {
484     if (degre<3)
485     {
486     if (num==0) {w=2.;uvw[0]=0.;uvw[1]=0.816496580927726;uvw[2]=-0,577350269189626;}
487     if (num==1) {w=2.;uvw[0]=0.;uvw[1]=-0.816496580927726;uvw[2]=-0,577350269189626;}
488     if (num==2) {w=2.;uvw[0]=-0.816496580927726;uvw[1]=0.;uvw[2]=0,577350269189626;}
489 francois 876 if (num==3) {w=2.;uvw[0]=0.816496580927726;uvw[1]=0.;uvw[2]=0,577350269189626;}
490 francois 754 return;
491     }
492     if (degre<4)
493     {
494     if (num==0) {w=1.333333333333333;uvw[0]=1.;uvw[1]=0.;uvw[2]=0.;}
495     if (num==1) {w=1.333333333333333;uvw[0]=-1.;uvw[1]=0.;uvw[2]=0.;}
496     if (num==2) {w=1.333333333333333;uvw[0]=0.;uvw[1]=1.;uvw[2]=0.;}
497     if (num==3) {w=1.333333333333333;uvw[0]=0.;uvw[1]=-1.;uvw[2]=0.;}
498     if (num==4) {w=1.333333333333333;uvw[0]=0.;uvw[1]=0.;uvw[2]=1.;}
499     if (num==5) {w=1.333333333333333;uvw[0]=0.;uvw[1]=0.;uvw[2]=-1.;}
500     return;
501     }
502     if (degre<6)
503     {
504     double a=0.795822425754221;
505     double b=0.758786910639328;
506     double w1=0.886426592797784;
507     double w2=0.335180055401662;
508     if (num==0) {w=w1;uvw[0]=a;uvw[1]=0.;uvw[2]=0.;}
509     if (num==1) {w=w1;uvw[0]=-a;uvw[1]=0.;uvw[2]=0.;}
510     if (num==2) {w=w1;uvw[0]=0.;uvw[1]=a;uvw[2]=0.;}
511     if (num==3) {w=w1;uvw[0]=0.;uvw[1]=-a;uvw[2]=0.;}
512     if (num==4) {w=w1;uvw[0]=0.;uvw[1]=0.;uvw[2]=a;}
513     if (num==5) {w=w1;uvw[0]=0.;uvw[1]=0.;uvw[2]=-a;}
514     if (num==6) {w=w2;uvw[0]=b;uvw[1]=b;uvw[2]=b;}
515     if (num==7) {w=w2;uvw[0]=b;uvw[1]=-b;uvw[2]=-b;}
516     if (num==8) {w=w2;uvw[0]=b;uvw[1]=b;uvw[2]=-b;}
517     if (num==9) {w=w2;uvw[0]=b;uvw[1]=-b;uvw[2]=b;}
518     if (num==10) {w=w2;uvw[0]=-b;uvw[1]=b;uvw[2]=b;}
519     if (num==11) {w=w2;uvw[0]=-b;uvw[1]=-b;uvw[2]=-b;}
520     if (num==12) {w=w2;uvw[0]=-b;uvw[1]=b;uvw[2]=-b;}
521     if (num==13) {w=w2;uvw[0]=-b;uvw[1]=-b;uvw[2]=b;}
522     return;
523     }
524    
525 francois 876 }
526    
527     void OT_POINTS_GAUSS::get_pt_gauss_penta(int degre, int num, double& w, double* uvw)
528     {
529     if (degre<3)
530     {
531 francois 885 double a=1./sqrt(3.);
532     if (num==0) {w=1./6.;uvw[0]=-a;uvw[1]=0.5;uvw[2]=0.5;}
533     if (num==1) {w=1./6.;uvw[0]=-a;uvw[1]=0.;uvw[2]=0.5;}
534     if (num==2) {w=1./6.;uvw[0]=-a;uvw[1]=0.5;uvw[2]=0.;}
535     if (num==3) {w=1./6.;uvw[0]=a;uvw[1]=0.5;uvw[2]=0.5;}
536     if (num==4) {w=1./6.;uvw[0]=a;uvw[1]=0.;uvw[2]=0.5;}
537     if (num==5) {w=1./6.;uvw[0]=a;uvw[1]=0.5;uvw[2]=0.;}
538 francois 876 return;
539     }
540     if (degre<4)
541     {
542     double a=1./sqrt(3);
543     if (num==0) {w=-27./96.;uvw[0]=-a;uvw[1]=1./3.;uvw[2]=1./3.;}
544     if (num==1) {w=25./96.;uvw[0]=-a;uvw[1]=0.6;uvw[2]=0.2;}
545     if (num==2) {w=25./96.;uvw[0]=-a;uvw[1]=0.2;uvw[2]=0.6;}
546     if (num==3) {w=25./96.;uvw[0]=-a;uvw[1]=0.2;uvw[2]=0.2;}
547     if (num==4) {w=-27./96.;uvw[0]=a;uvw[1]=1./3.;uvw[2]=1./3.;}
548     if (num==5) {w=25./96.;uvw[0]=a;uvw[1]=0.6;uvw[2]=0.2;}
549     if (num==6) {w=25./96.;uvw[0]=a;uvw[1]=0.2;uvw[2]=0.6;}
550     if (num==7) {w=25./96.;uvw[0]=a;uvw[1]=0.2;uvw[2]=0.2;}
551     return;
552     }
553     if (degre<6)
554     {
555     double alpha=sqrt(3./5.);
556     double c1=5./9.;
557     double c2=8./9.;
558     double a=(6+sqrt(15))/21;
559     double b=(6-sqrt(15))/21;
560     if (num==0) {w=c1*9./80.;uvw[0]=-alpha;uvw[1]=1./3.;uvw[2]=1./3.;}
561     if (num==1) {w=c1*(155.+sqrt(15))/2400.;uvw[0]=-alpha;uvw[1]=a;uvw[2]=a;}
562     if (num==2) {w=c1*(155.+sqrt(15))/2400.;uvw[0]=-alpha;uvw[1]=1.-2.*a;uvw[2]=a;}
563     if (num==3) {w=c1*(155.+sqrt(15))/2400.;uvw[0]=-alpha;uvw[1]=a;uvw[2]=1.-2.*a;}
564     if (num==4) {w=c1*(155.-sqrt(15))/2400.;uvw[0]=-alpha;uvw[1]=b;uvw[2]=b;}
565     if (num==5) {w=c1*(155.-sqrt(15))/2400.;uvw[0]=-alpha;uvw[1]=1.-2.*b;uvw[2]=b;}
566     if (num==6) {w=c1*(155.-sqrt(15))/2400.;uvw[0]=-alpha;uvw[1]=b;uvw[2]=1.-2.*b;}
567     if (num==7) {w=c2*9./80.;uvw[0]=0.;uvw[1]=1./3.;uvw[2]=1./3.;}
568     if (num==8) {w=c2*(155.+sqrt(15))/2400.;uvw[0]=0.;uvw[1]=a;uvw[2]=a;}
569     if (num==9) {w=c2*(155.+sqrt(15))/2400.;uvw[0]=0.;uvw[1]=1.-2.*a;uvw[2]=a;}
570     if (num==10) {w=c2*(155.+sqrt(15))/2400.;uvw[0]=0.;uvw[1]=a;uvw[2]=1.-2.*a;}
571     if (num==11) {w=c2*(155.-sqrt(15))/2400.;uvw[0]=0.;uvw[1]=b;uvw[2]=b;}
572     if (num==12) {w=c2*(155.-sqrt(15))/2400.;uvw[0]=0.;uvw[1]=1.-2.*b;uvw[2]=b;}
573     if (num==13) {w=c2*(155.-sqrt(15))/2400.;uvw[0]=0.;uvw[1]=b;uvw[2]=1.-2.*b;}
574     if (num==14) {w=c1*9./80.;uvw[0]=alpha;uvw[1]=1./3.;uvw[2]=1./3.;}
575     if (num==15) {w=c1*(155.+sqrt(15))/2400.;uvw[0]=alpha;uvw[1]=a;uvw[2]=a;}
576     if (num==16) {w=c1*(155.+sqrt(15))/2400.;uvw[0]=alpha;uvw[1]=1.-2.*a;uvw[2]=a;}
577     if (num==17) {w=c1*(155.+sqrt(15))/2400.;uvw[0]=alpha;uvw[1]=a;uvw[2]=1.-2.*a;}
578     if (num==18) {w=c1*(155.-sqrt(15))/2400.;uvw[0]=alpha;uvw[1]=b;uvw[2]=b;}
579     if (num==19) {w=c1*(155.-sqrt(15))/2400.;uvw[0]=alpha;uvw[1]=1.-2.*b;uvw[2]=b;}
580     if (num==20) {w=c1*(155.-sqrt(15))/2400.;uvw[0]=alpha;uvw[1]=b;uvw[2]=1.-2.*b;}
581    
582    
583     return;
584     }
585    
586 francois 754 }