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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1169 by gross, Thu Mar 29 06:44:30 2007 UTC revision 1170 by btully, Fri May 25 05:24:57 2007 UTC
# Line 36  void Finley_Quad_getNodesTri(int numQuad Line 36  void Finley_Quad_getNodesTri(int numQuad
36    double Q1,Q2;    double Q1,Q2;
37    #define DIM 2    #define DIM 2
38        
39    /*  the easy case: */    /*  the easy cases: */
40        
41    if (numQuadNodes==1) {    if (numQuadNodes==1) {
42      QUADNODES(0,0)=1./3.;      QUADNODES(0,0)=1./3.;
43      QUADNODES(1,0)=1./3.;      QUADNODES(1,0)=1./3.;
44      QUADWEIGHTS(0)= .5;      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    } else {    } else {
69            
70      /*  get scheme on [0.1]^2 */      /*  get scheme on [0.1]^2 */
       
71      Finley_Quad_getNodesRec(numQuadNodes,quadNodes,quadWeights);      Finley_Quad_getNodesRec(numQuadNodes,quadNodes,quadWeights);
72      if (! Finley_noError()) return;      if (! Finley_noError()) return;
73            
# Line 72  void Finley_Quad_getNodesTet(int numQuad Line 94  void Finley_Quad_getNodesTet(int numQuad
94    double Q1,Q2,Q3,JA11,JA12,JA13,JA21,JA22,JA23,JA31,JA32,JA33,DET;    double Q1,Q2,Q3,JA11,JA12,JA13,JA21,JA22,JA23,JA31,JA32,JA33,DET;
95    #define DIM 3    #define DIM 3
96        
97    /*  the easy case: */    /*  the easy cases: */
     
98    if (numQuadNodes==1) {    if (numQuadNodes==1) {
99      QUADNODES(0,0)= .25;      QUADNODES(0,0)=0.25;
100      QUADNODES(1,0)= .25;      QUADNODES(1,0)=0.25;
101      QUADNODES(2,0)= .25;      QUADNODES(2,0)=0.25;
102      QUADWEIGHTS(0)=1./6.;      QUADWEIGHTS(0)=1./6.;
103      } 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    } else {    } else {
144            
145      /*  get scheme on [0.1]^3 */      /*  get scheme on [0.1]^3 */
# Line 459  int Finley_Quad_getNumNodesTri(int order Line 520  int Finley_Quad_getNumNodesTri(int order
520    int numQuadNodesLine;    int numQuadNodesLine;
521    if (order<=1) {    if (order<=1) {
522        return 1;        return 1;
523      } else if (order==2){
524          return 3;
525      } else if (order==3){
526          return 4;
527    } else {    } else {
528        numQuadNodesLine=Finley_Quad_getNumNodesLine(order+1);        numQuadNodesLine=Finley_Quad_getNumNodesLine(order+1);
529        if (Finley_noError()) {        if (Finley_noError()) {
# Line 483  int Finley_Quad_getNumNodesTet(int order Line 548  int Finley_Quad_getNumNodesTet(int order
548    int numQuadNodesLine;    int numQuadNodesLine;
549    if (order<=1) {    if (order<=1) {
550        return 1;        return 1;
551      } else if (order==2){
552          return 4;
553      } else if (order==3){
554          return 5;
555    } else {    } else {
556       numQuadNodesLine=Finley_Quad_getNumNodesLine(order+2);       numQuadNodesLine=Finley_Quad_getNumNodesLine(order+2);
557       if (Finley_noError()) {       if (Finley_noError()) {
Legend:
Removed from v.1169  
changed lines
  Added in v.1170
  ViewVC Help
Powered by ViewVC 1.1.26