xref: /phasta/phSolver/common/symtet.c (revision 7acde132a6def0fe2daaec0d1a712dff0e5c6636)
1 #include <FCMangle.h>
2 #define symtet FortranCInterface_GLOBAL_(symtet, SYMTET)
3 
4 typedef double DARR4[4];
5 
6 int tetIntPnt(int,DARR4**,double **);
7 
8 
symtet(int * n1,double pt[][4],double wt[],int * err)9 void symtet(int *n1, double pt[][4], double wt[], int *err)
10 {
11   double *lwt;
12   DARR4 *lpt;
13   int i,j;
14   *err = tetIntPnt(*n1, &lpt, &lwt);
15   for(i=0; i < *n1; i++) {
16     wt[i] = lwt[i];
17     for(j=0; j <4; j++)
18       pt[i][j]=lpt[i][j];
19   }
20 }
21 
22 /*$Id$*/
23 #include <stdio.h>
24 
25 /* constants for 4-point rule */
26 #define  a4  0.5854101966249685
27 #define  b4  0.1381966011250150
28 
29 /* constants for 5-point rule */
30 #define  a5  0.5000000000000000
31 #define  b5  0.1666666666666667
32 
33 /* constants for 16-point rule */
34 #define  a16  0.0503737941001228
35 #define  b16  0.0665420686332923
36 #define  c16  0.7716429020672371
37 #define  d16  0.0761190326442543
38 #define  e16  0.1197005277978019
39 #define  f16  0.0718316452676693
40 #define  g16  0.4042339134672644
41 
42 /* constants for 17-point rule */
43 #define  a17  0.1884185567365411
44 #define  b17  0.0670385837260428
45 #define  c17  0.0452855923632739
46 #define  p17  0.7316369079576180
47 #define  q17  0.0894543640141273
48 #define  e17  0.1325810999384657
49 #define  f17  0.0245400397290300
50 #define  g17  0.4214394310662522
51 
52 
53 /* constants for 29-point rule */
54 #define  a29  0.0904012904601475
55 #define  b29  0.0191198342789912
56 #define  c29  0.0436149384066657
57 #define  d29  0.0258116759619916
58 #define  p29  0.8277192480479295
59 #define  q29  0.0574269173173568
60 #define  e29  0.0513518841255634
61 #define  f29  0.4860510285706072
62 #define  g29  0.2312985436519147
63 #define  h29  0.2967538129690260
64 #define  i29  0.6081079894015281
65 #define  j29  0.0475690988147229
66 
67 #define a25   0.2500000000000000
68 #define a26  -0.8000000000000000
69 #define a27   0.4500000000000000
70 
71 /* typedef double DARR4[4] ; */
72 
73   static double rstw1[][4] = {{a25,a25,a25,a25}};
74   static double twt1[] = {1.000000000000000000000};
75 
76   static double rstw4[][4] = {{a4,b4,b4,b4},{b4,a4,b4,b4},{b4,b4,a4,b4},
77                               {b4,b4,b4,a4}};
78   static double twt4[] = {a25,a25,a25,a25};
79 
80   static double rstw5[][4] = {{a25,a25,a25,a25},
81                              {a5,b5,b5,b5},
82                              {b5,a5,b5,b5},
83                              {b5,b5,a5,b5},
84                              {b5,b5,b5,a5}};
85   static double twt5[] = {a26,a27,a27,a27,a27};
86 
87   static double rstw16[][4] = {{c16,d16,d16,d16},
88                               {d16,c16,d16,d16},
89                               {d16,d16,c16,d16},
90                               {d16,d16,d16,c16},
91 
92                               {e16,f16,g16,g16},
93                               {f16,e16,g16,g16},
94                               {e16,g16,g16,f16},
95                               {f16,g16,g16,e16},
96                               {g16,g16,e16,f16},
97                               {g16,g16,f16,e16},
98                               {g16,e16,f16,g16},
99                               {g16,f16,e16,g16},
100                               {e16,g16,f16,g16},
101                               {f16,g16,e16,g16},
102                               {g16,e16,g16,f16},
103                               {g16,f16,g16,e16} };
104   static double twt16[] = {a16,a16,a16,a16,
105                           b16,b16,b16,b16,
106                           b16,b16,b16,b16,
107                           b16,b16,b16,b16 };
108 
109   static double rstw17[][4] = {{a25,a25,a25,a25},
110                               {p17,q17,q17,q17},
111                               {q17,p17,q17,q17},
112                               {q17,q17,p17,q17},
113                               {q17,q17,q17,p17},
114 
115                               {e17,f17,g17,g17},
116                               {f17,e17,g17,g17},
117                               {e17,g17,g17,f17},
118                               {f17,g17,g17,e17},
119                               {g17,g17,e17,f17},
120                               {g17,g17,f17,e17},
121                               {g17,e17,f17,g17},
122                               {g17,f17,e17,g17},
123                               {e17,g17,f17,g17},
124                               {f17,g17,e17,g17},
125                               {g17,e17,g17,f17},
126                               {g17,f17,g17,e17} };
127 
128   static double twt17[] = {a17,b17,b17,b17,b17,
129                           c17,c17,c17,c17,c17,c17,
130                           c17,c17,c17,c17,c17,c17};
131 
132   static double twt29[] = {a29,b29,b29,b29,b29,
133                           c29,c29,c29,c29,c29,c29,
134                           c29,c29,c29,c29,c29,c29,
135                           d29,d29,d29,d29,d29,d29,
136                           d29,d29,d29,d29,d29,d29};
137 
138   static double rstw29[][4] = {{a25,a25,a25,a25},
139 
140                               {p29,q29,q29,q29},
141                               {q29,p29,q29,q29},
142                               {q29,q29,p29,q29},
143                               {q29,q29,q29,p29},
144 
145                               {e29,f29,g29,g29},
146                               {f29,e29,g29,g29},
147                               {e29,g29,g29,f29},
148                               {f29,g29,g29,e29},
149                               {g29,g29,e29,f29},
150                               {g29,g29,f29,e29},
151                               {g29,e29,f29,g29},
152                               {g29,f29,e29,g29},
153                               {e29,g29,f29,g29},
154                               {f29,g29,e29,g29},
155                               {g29,e29,g29,f29},
156                               {g29,f29,g29,e29},
157 
158                               {h29,i29,j29,j29},
159                               {i29,h29,j29,j29},
160                               {h29,j29,j29,i29},
161                               {i29,j29,j29,h29},
162                               {j29,j29,h29,i29},
163                               {j29,j29,i29,h29},
164                               {j29,h29,i29,j29},
165                               {j29,i29,h29,j29},
166                               {h29,j29,i29,j29},
167                               {i29,j29,h29,j29},
168                               {j29,h29,j29,i29},
169                               {j29,i29,j29,h29} };
170 
171 #ifdef __cplusplus
172 extern "C" {
173 #endif
174 
tetIntPnt(int nint,DARR4 ** bcord,double ** wt)175 int tetIntPnt(int nint, DARR4 **bcord, double **wt)
176 {
177   int retval = 1;
178 
179   if( nint == 1 ) {*bcord = rstw1 ; *wt = twt1; }
180   else if( nint == 4 ){*bcord = rstw4 ; *wt = twt4; }
181   else if( nint == 5 ){*bcord = rstw5 ; *wt = twt5; }
182   else if( nint == 16 ){*bcord = rstw16 ; *wt = twt16; }
183   else if( nint == 17 ){*bcord = rstw17 ; *wt = twt17; }
184   else if( nint == 29 ){*bcord = rstw29 ; *wt = twt29; }
185   else
186   {
187     fprintf(stderr,"\n%d integration points unsupported in symtet.c; give {1,4,5,16,17,29}\n",nint);
188     retval = 0;
189   }
190   return retval ;
191 }
192 
193 #ifdef __cplusplus
194 }
195 #endif
196