xref: /phasta/phSolver/common/symtri.c (revision 595995161822a203c8467e0e4a253d7bd7d6df32)
1*59599516SKenneth E. Jansen #include <FCMangle.h>
2*59599516SKenneth E. Jansen #define symtri FortranCInterface_GLOBAL_(symtri,SYMTRI)
3*59599516SKenneth E. Jansen 
4*59599516SKenneth E. Jansen typedef double DARR3[3];
5*59599516SKenneth E. Jansen 
6*59599516SKenneth E. Jansen int triIntPnt(int, DARR3**,double**);
7*59599516SKenneth E. Jansen 
8*59599516SKenneth E. Jansen symtri(int *n1, double pt[][4], double wt[], int *err)
9*59599516SKenneth E. Jansen {
10*59599516SKenneth E. Jansen   double *lwt;
11*59599516SKenneth E. Jansen   DARR3 *lpt;
12*59599516SKenneth E. Jansen   int i,j;
13*59599516SKenneth E. Jansen   *err = triIntPnt(*n1, &lpt, &lwt);
14*59599516SKenneth E. Jansen   for(i=0; i < *n1; i++){
15*59599516SKenneth E. Jansen     wt[i] = lwt[i];
16*59599516SKenneth E. Jansen     for(j=0; j < 3; j++)
17*59599516SKenneth E. Jansen       pt[i][j] = lpt[i][j];
18*59599516SKenneth E. Jansen   }
19*59599516SKenneth E. Jansen }
20*59599516SKenneth E. Jansen 
21*59599516SKenneth E. Jansen #include <stdio.h>
22*59599516SKenneth E. Jansen 
23*59599516SKenneth E. Jansen #define A 0.3333333333333330000000000000000000
24*59599516SKenneth E. Jansen #define B 0.666666666666667
25*59599516SKenneth E. Jansen #define C 0.166666666666667
26*59599516SKenneth E. Jansen #define D 0.562500000000000
27*59599516SKenneth E. Jansen #define E 0.520833333333333
28*59599516SKenneth E. Jansen 
29*59599516SKenneth E. Jansen #define F 0.109951743655322
30*59599516SKenneth E. Jansen #define G 0.223381589678011
31*59599516SKenneth E. Jansen #define H 0.816847572980459
32*59599516SKenneth E. Jansen #define I 0.091576213509771
33*59599516SKenneth E. Jansen #define J 0.108103018168070
34*59599516SKenneth E. Jansen #define K 0.445948490915965
35*59599516SKenneth E. Jansen 
36*59599516SKenneth E. Jansen #define L 0.225000000000000
37*59599516SKenneth E. Jansen #define M 0.125939180544827
38*59599516SKenneth E. Jansen #define N 0.132394152788506
39*59599516SKenneth E. Jansen 
40*59599516SKenneth E. Jansen #define O 0.797426985353087
41*59599516SKenneth E. Jansen #define P 0.101286507323456
42*59599516SKenneth E. Jansen #define Q 0.470142064105115
43*59599516SKenneth E. Jansen #define R 0.059715871789770
44*59599516SKenneth E. Jansen 
45*59599516SKenneth E. Jansen #define S 0.050844906370207
46*59599516SKenneth E. Jansen #define T 0.116786275726379
47*59599516SKenneth E. Jansen #define U 0.082851075618374
48*59599516SKenneth E. Jansen 
49*59599516SKenneth E. Jansen #define V 0.873821971016996
50*59599516SKenneth E. Jansen #define W 0.063089014491502
51*59599516SKenneth E. Jansen #define X 0.501426509658179
52*59599516SKenneth E. Jansen #define Y 0.249286745170910
53*59599516SKenneth E. Jansen #define Z 0.636502499121399
54*59599516SKenneth E. Jansen #define AA 0.310352451033785
55*59599516SKenneth E. Jansen #define BB 0.053145049844816
56*59599516SKenneth E. Jansen 
57*59599516SKenneth E. Jansen /* typedef double DARR3[3] ; */
58*59599516SKenneth E. Jansen 
59*59599516SKenneth E. Jansen static double rst1[][3] = {{A,A,A}};
60*59599516SKenneth E. Jansen static double wt1[] = {1.00000000000000000};
61*59599516SKenneth E. Jansen 
62*59599516SKenneth E. Jansen static double rst3[][3] = {{B,C,C},{C,B,C},{C,C,B}};
63*59599516SKenneth E. Jansen static double wt3[] = {A,A,A};
64*59599516SKenneth E. Jansen 
65*59599516SKenneth E. Jansen static double rst4[][3] = {{A,A,A},
66*59599516SKenneth E. Jansen                              {0.6,0.2,0.2},{0.2,0.6,0.2},{0.2,0.2,0.6}};
67*59599516SKenneth E. Jansen static double wt4[] = {-D,E,E,E};
68*59599516SKenneth E. Jansen 
69*59599516SKenneth E. Jansen static double rst6[][3] = {{H,I,I},{I,H,I},{I,I,H},
70*59599516SKenneth E. Jansen                              {J,K,K},{K,J,K},{K,K,J}};
71*59599516SKenneth E. Jansen static double wt6[] = {F,F,F,G,G,G};
72*59599516SKenneth E. Jansen 
73*59599516SKenneth E. Jansen static double rst7[][3] = {{A,A,A},
74*59599516SKenneth E. Jansen                              {O,P,P},{P,O,P},{P,P,O},
75*59599516SKenneth E. Jansen                              {Q,Q,R},{Q,R,Q},{R,Q,Q}};
76*59599516SKenneth E. Jansen static double wt7[] = {L,M,M,M,N,N,N};
77*59599516SKenneth E. Jansen 
78*59599516SKenneth E. Jansen static double rst12[][3] = {{V,W,1.0-V-W},{W,V,1.0-V-W},{W,W,1.0-W-W},
79*59599516SKenneth E. Jansen 			    {X,Y,1.0-X-Y},{Y,X,1.0-X-Y},{Y,Y,1.0-Y-Y},
80*59599516SKenneth E. Jansen 			    {Z,AA,1.0-Z-AA},{AA,Z,1.0-Z-AA},{Z,BB,1.0-Z-BB},
81*59599516SKenneth E. Jansen 			    {AA,BB,1.0-AA-BB},{BB,AA,1.0-BB-AA},{BB,Z,1.0-BB-Z}};
82*59599516SKenneth E. Jansen 
83*59599516SKenneth E. Jansen static double wt12[] = {S,S,S,T,T,T,U,U,U,U,U,U};
84*59599516SKenneth E. Jansen 
85*59599516SKenneth E. Jansen int triIntPnt(int nint, DARR3 **bcord, double **wt)
86*59599516SKenneth E. Jansen {
87*59599516SKenneth E. Jansen   int retval = 1 ;
88*59599516SKenneth E. Jansen   if( nint == 3 ) {*bcord = rst3 ; *wt = wt3; }
89*59599516SKenneth E. Jansen   else if( nint == 1 ){*bcord = rst1 ; *wt = wt1; }
90*59599516SKenneth E. Jansen   else if( nint == 4 ){*bcord = rst4 ; *wt = wt4; }
91*59599516SKenneth E. Jansen   else if( nint == 6 ){*bcord = rst6 ; *wt = wt6; }
92*59599516SKenneth E. Jansen   else if( nint == 7 ){*bcord = rst7 ; *wt = wt7; }
93*59599516SKenneth E. Jansen   else if( nint == 12){*bcord = rst12; *wt = wt12;}
94*59599516SKenneth E. Jansen   else
95*59599516SKenneth E. Jansen   {
96*59599516SKenneth E. Jansen     fprintf(stderr,"\n%d integration points unsupported in symtri.c; give {1,3,4,6,7,12}\n",nint);
97*59599516SKenneth E. Jansen     retval = 0;
98*59599516SKenneth E. Jansen   }
99*59599516SKenneth E. Jansen   return retval ;
100*59599516SKenneth E. Jansen }
101*59599516SKenneth E. Jansen 
102*59599516SKenneth E. Jansen 
103*59599516SKenneth E. Jansen 
104