Changeset 17865

Timestamp:

04/28/14 15:49:38 (11 years ago)

Author:

Mathieu Morlighem

Message:

BUG: fixed tetra gaussian points

Location:

issm/trunk-jpl/src/c

Files:

• classes/gauss/GaussTetra.cpp (modified) (1 diff)
• shared/Numerics/GaussPoints.cpp (modified) (1 diff)

issm/trunk-jpl/src/c/classes/gauss/GaussTetra.cpp

@@ -33 +33 @@
 /*FUNCTION GaussTetra::GaussTetra(int order) {{{*/
 GaussTetra::GaussTetra(int order){
+
         /*Get gauss points*/
         GaussLegendreTetra(&numgauss,&coords1,&coords2,&coords3,&coords4,&weights,order);
-        for(int i=0;i<numgauss;i++) this->weights[i]=this->weights[i]/6.; //FIXME: double check that
+
+        /*Rescale weights if necessary*/
+        IssmDouble sumweights = 0.;
+        for(int i=0;i<numgauss;i++) sumweights += this->weights[i];
+        if(sumweights==1.){
+                for(int i=0;i<numgauss;i++) this->weights[i] = this->weights[i]/6.; /*rescale volume to 1/6*/
+        }
+
+        /*Check final weights in debugging mode*/
+        #ifdef _ISSM_DEBUG_
+        sumweights = 0.; for(int i=0;i<numgauss;i++) sumweights += this->weights[i];
+        _assert_(sumweights>1./6.-1e-10);
+        _assert_(sumweights<1./6.+1e-10);
+        #endif
 
         /*Initialize static fields as undefinite*/

issm/trunk-jpl/src/c/shared/Numerics/GaussPoints.cpp

@@ -1425 +1425 @@
                         (*pl3 )[i]=l3p [iord-1][i];
                         (*pl4 )[i]=l4p [iord-1][i];
-                        (*pwgt)[i]=wgtp[iord-1][i]*2./3.*SQRT2;
+                        (*pwgt)[i]=wgtp[iord-1][i];
                 }
         }