]> SALOME platform Git repositories - tools/medcoupling.git/commitdiff
Salome HOME
[EDF26877] : Improvement of precision of integral computation of fields on Gauss...
authorAnthony Geay <anthony.geay@edf.fr>
Mon, 23 Jan 2023 16:29:31 +0000 (17:29 +0100)
committerAnthony Geay <anthony.geay@edf.fr>
Mon, 30 Jan 2023 09:05:09 +0000 (10:05 +0100)
src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx
src/INTERP_KERNEL/InterpKernelDenseMatrix.hxx
src/INTERP_KERNEL/InterpKernelDenseMatrix.txx
src/MEDCoupling/MEDCouplingFieldDiscretization.cxx
src/MEDCoupling/MEDCouplingGaussLocalization.cxx
src/MEDCoupling/MEDCouplingGaussLocalization.hxx
src/MEDCoupling_Swig/MEDCouplingBasicsTest4.py
src/MEDCoupling_Swig/MEDCouplingBasicsTest7.py
src/MEDCoupling_Swig/MEDCouplingCommon.i

index a74a5948fe2991b4ba2249fb919255055c21364c..05435efe39dfecdc8d835ce90c47784cb504d109 100644 (file)
@@ -916,6 +916,19 @@ void GaussInfo::tria3aInit()
   funValue[1] = -0.5*(gc[0] + gc[1]);
   funValue[2] = 0.5*(1.0 + gc[0]);
   SHAPE_FUN_MACRO_END;
+  
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = 0.0 ;
+  devFunValue[1] = 0.5 ;
+
+  devFunValue[2] = -0.5;
+  devFunValue[3] = -0.5;
+
+  devFunValue[4] = 0.5;
+  devFunValue[5] = 0.0;
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 /*!
@@ -944,6 +957,19 @@ void GaussInfo::tria3bInit()
   funValue[1] = gc[0];
   funValue[2] = gc[1];
   SHAPE_FUN_MACRO_END;
+
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = -1.0 ;
+  devFunValue[1] = -1.0 ;
+
+  devFunValue[2] = 1.0 ;
+  devFunValue[3] = 0.0 ;
+
+  devFunValue[4] = 0.0 ;
+  devFunValue[5] = 1.0 ;
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 /*!
@@ -987,6 +1013,28 @@ void GaussInfo::tria6aInit()
   funValue[4] = -1.0*(1.0 + gc[0])*(gc[0] + gc[1]);
   funValue[5] = (1.0 + gc[1])*(1.0 + gc[1]);
   SHAPE_FUN_MACRO_END;
+  
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = 0.0;
+  devFunValue[1] = 0.5*( 2*gc[1] + 1.0 );
+
+  devFunValue[2] = 0.5*( 2*gc[0] + 2.0*gc[1] + 1.0);
+  devFunValue[3] = 0.5*( 2*gc[1] + 2.0*gc[0] + 1.0);
+
+  devFunValue[4] = gc[0] + 0.5;
+  devFunValue[5] = 0.0;
+
+  devFunValue[6] = -1.0*(1.0 + gc[1]);
+  devFunValue[7] = -1.0*(2*gc[1]+gc[0]+1.0);
+
+  devFunValue[8] = -1.0*(2*gc[0]+gc[1]+1.0);
+  devFunValue[9] = -1.0*(1.0 + gc[0]);
+
+  devFunValue[10] = 0.0;
+  devFunValue[11] = (2*gc[1]+2.0);
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 /*!
@@ -1030,6 +1078,28 @@ void GaussInfo::tria6bInit()
   funValue[4] = 4.0*gc[0]*gc[1];
   funValue[5] = 4.0*gc[1]*(1.0 - gc[0] - gc[1]);
   SHAPE_FUN_MACRO_END;
+
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = 4*gc[0] + 4*gc[1] - 3.0 ;
+  devFunValue[1] = 4*gc[1] + 4*gc[0] - 3.0 ;
+
+  devFunValue[2] = 4*gc[0] - 1.0 ;
+  devFunValue[3] = 0.0 ;
+
+  devFunValue[4] = 0.0 ;
+  devFunValue[5] = 4*gc[1] - 1.0 ;
+
+  devFunValue[6] = -8.0*gc[0] - 4.0 * gc[1] + 4.0 ;
+  devFunValue[7] = -4.0*gc[0] ; 
+
+  devFunValue[8] = 4.0*gc[1];
+  devFunValue[9] = 4.0*gc[0];
+
+  devFunValue[10] = -4.0*gc[1] ;
+  devFunValue[11] = -8.0*gc[1] - 4.0*gc[0] + 4.0 ; 
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 void GaussInfo::tria7aInit()
@@ -1105,8 +1175,24 @@ void GaussInfo::quad4aInit()
   funValue[0] = 0.25*(1.0 + gc[1])*(1.0 - gc[0]);
   funValue[1] = 0.25*(1.0 - gc[1])*(1.0 - gc[0]);
   funValue[2] = 0.25*(1.0 - gc[1])*(1.0 + gc[0]);
-  funValue[3] = 0.25*(1.0 + gc[0])*(1.0 + gc[1]);
+  funValue[3] = 0.25*(1.0 + gc[1])*(1.0 + gc[0]);
   SHAPE_FUN_MACRO_END;
+  
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = -0.25*(1.0 + gc[1]);
+  devFunValue[1] = 0.25*(1.0 - gc[0]);
+
+  devFunValue[2] = -0.25*(1.0 - gc[1]);
+  devFunValue[3] = -0.25*(1.0 - gc[0]);
+
+  devFunValue[4] = 0.25*(1.0 - gc[1]);
+  devFunValue[5] = -0.25*(1.0 + gc[0]);
+
+  devFunValue[6] = 0.25*(1.0 + gc[1]);
+  devFunValue[7] = 0.25*(1.0 + gc[0]);
+  
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 /*!
@@ -1140,6 +1226,22 @@ void GaussInfo::quad4bInit()
   funValue[2] = 0.25*(1.0 + gc[0])*(1.0 + gc[1]);
   funValue[3] = 0.25*(1.0 - gc[0])*(1.0 + gc[1]);
   SHAPE_FUN_MACRO_END;
+  
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = -0.25*(1.0 - gc[1]);
+  devFunValue[1] = -0.25*(1.0 - gc[0]);
+
+  devFunValue[2] = 0.25*(1.0 - gc[1]);
+  devFunValue[3] = -0.25*(1.0 + gc[0]);
+
+  devFunValue[4] = 0.25*(1.0 + gc[1]);
+  devFunValue[5] = 0.25*(1.0 + gc[0]);
+
+  devFunValue[6] = -0.25*(1.0 + gc[1]);
+  devFunValue[7] = 0.25*(1.0 - gc[0]);
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 void GaussInfo::quad4cInit() 
@@ -1170,6 +1272,22 @@ void GaussInfo::quad4cInit()
    funValue[2] = 0.25*(1.0 + gc[0])*(1.0 + gc[1]);
    funValue[3] = 0.25*(1.0 + gc[0])*(1.0 - gc[1]);
    SHAPE_FUN_MACRO_END;
+   
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = -0.25*(1.0 - gc[1]);
+  devFunValue[1] = -0.25*(1.0 - gc[0]);
+
+  devFunValue[2] = -0.25*(1.0 + gc[1]);
+  devFunValue[3] = 0.25*(1.0 - gc[0]);
+
+  devFunValue[4] = 0.25*(1.0 + gc[0]);
+  devFunValue[5] = 0.25*(1.0 + gc[1]);
+
+  devFunValue[6] = 0.25*(1.0 - gc[1]);
+  devFunValue[7] = -0.25*(1.0 + gc[0]);
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 /*!
@@ -1255,6 +1373,34 @@ void GaussInfo::quad8aInit()
   funValue[6] = 0.5*(1.0 + gc[0])*(1.0 - gc[1])*(1.0 + gc[1]);
   funValue[7] = 0.5*(1.0 + gc[1])*(1.0 - gc[0])*(1.0 + gc[0]);
   SHAPE_FUN_MACRO_END;
+
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = 0.25*(1.0 + gc[1])*(2*gc[0]-gc[1]);
+  devFunValue[1] = 0.25*(1.0 - gc[0])*(2*gc[1]-gc[0]);
+
+  devFunValue[2] = 0.25*(1.0 - gc[1])*(2*gc[0]+gc[1]);
+  devFunValue[3] = 0.25*(1.0 - gc[0])*(2*gc[1]+gc[0]);
+
+  devFunValue[4] = 0.25*(1.0 - gc[1])*(2*gc[0]-gc[1]);
+  devFunValue[5] = 0.25*(1.0 + gc[0])*(2*gc[1]-gc[0]);
+
+  devFunValue[6] = 0.25*(1.0 + gc[1])*(2*gc[0]+gc[1]);
+  devFunValue[7] = 0.25*(1.0 + gc[0])*(2*gc[1]+gc[0]);
+
+  devFunValue[8] = -0.5*(1.0 - gc[1])*(1.0 + gc[1]);
+  devFunValue[9] = 0.5*(1.0 - gc[0])*(-2*gc[1]);
+
+  devFunValue[10] = 0.5*(1.0 - gc[1])*(-2*gc[0]);
+  devFunValue[11] = -0.5*(1.0 - gc[0])*(1.0 + gc[0]);
+
+  devFunValue[12] = 0.5*(1.0 - gc[1])*(1.0 + gc[1]);
+  devFunValue[13] = 0.5*(1.0 + gc[0])*(-2*gc[1]);
+
+  devFunValue[14] = 0.5*(1.0 + gc[1])*(-2*gc[0]);
+  devFunValue[15] = 0.5*(1.0 - gc[0])*(1.0 + gc[0]);
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 /*!
@@ -1308,6 +1454,34 @@ void GaussInfo::quad8bInit()
   funValue[6] = 0.5*(1.0 - gc[0]*gc[0])*(1.0 + gc[1]);
   funValue[7] = 0.5*(1.0 - gc[1]*gc[1])*(1.0 - gc[0]);
   SHAPE_FUN_MACRO_END;
+  
+  DEV_SHAPE_FUN_MACRO_BEGIN;
+
+  devFunValue[0] = 0.25*(1.0 - gc[1])*(2*gc[0] + gc[1] );
+  devFunValue[1] = 0.25*(1.0 - gc[0])*(2*gc[1] + gc[0] );
+
+  devFunValue[2] = 0.25*(1.0 - gc[1])*(2*gc[0] - gc[1] );
+  devFunValue[3] = 0.25*(1.0 + gc[0])*(2*gc[1] - gc[0] );
+
+  devFunValue[4] = 0.25*(1.0 + gc[1])*(2*gc[0] + gc[1]);
+  devFunValue[5] = 0.25*(1.0 + gc[0])*(2*gc[1] + gc[0]);
+
+  devFunValue[6] = 0.25*( 1.0 + gc[1])*(2*gc[0] - gc[1] );
+  devFunValue[7] = 0.25*(1.0 - gc[0])*(2*gc[1] - gc[0]);
+
+  devFunValue[8] = -(1.0 - gc[1])*gc[0];
+  devFunValue[9] = -0.5*(1.0 - gc[0]*gc[0]);
+
+  devFunValue[10] = 0.5*(1.0 - gc[1]*gc[1]);
+  devFunValue[11] = -(1.0 + gc[0])*gc[1];
+
+  devFunValue[12] = -(1.0 + gc[1])*gc[0];
+  devFunValue[13] = 0.5*(1.0 - gc[0]*gc[0]);
+
+  devFunValue[14] = -0.5*(1.0 - gc[1]*gc[1]);
+  devFunValue[15] = -(1.0 - gc[0])*gc[1];
+
+  DEV_SHAPE_FUN_MACRO_END;
 }
 
 void GaussInfo::quad9aInit()
index 38e4a61281f777a887ae4e7711ed0d5d5a9808a2..63080077f3c0d39c472f6f4e8d3255987177f3bc 100644 (file)
@@ -47,6 +47,7 @@ namespace INTERP_KERNEL
     void resize(mcIdType newn, mcIdType newm);
     void assign(mcIdType newn, mcIdType newm, const T &a);
     ~DenseMatrixT();
+    T determinant() const;
   };
 
   using DenseMatrix = DenseMatrixT<double>;
index 9705f4fdbacb7910a17fd04cf3371d36caaf9613..2a262c31b14648b9a9ba5aecb3f736040bb2f609 100644 (file)
@@ -20,6 +20,7 @@
 #pragma once
 
 #include "InterpKernelDenseMatrix.hxx"
+#include "InterpKernelException.hxx"
 
 namespace INTERP_KERNEL
 {
@@ -138,4 +139,22 @@ namespace INTERP_KERNEL
       delete[] (v);
     }
   }
+
+  template<class T>
+  T Determinant22(const T *m)
+  { return m[0]*m[3]-m[1]*m[2]; }
+
+  template<class T>
+  T Determinant33(const T *m)
+  { return m[0]*(m[4]*m[8]-m[7]*m[5])-m[1]*(m[3]*m[8]-m[6]*m[5])+m[2]*(m[3]*m[7]-m[6]*m[4]);}
+  
+  template <class T>
+  T DenseMatrixT<T>::determinant() const
+  {
+    if(nn==2 && mm==2)
+      return Determinant22(v[0]);
+    if(nn==3 && mm==3)
+      return Determinant33(v[0]);
+    THROW_IK_EXCEPTION("DenseMatrixT::determinant : only 2x2 and 3x3 implemented !");
+  }
 }
index d7cf86b66e9e3ed794092e7c3af5237daba5f18c..ebf44b23da8a165f829a53e1c80ec9a23096c851 100755 (executable)
@@ -30,6 +30,7 @@
 #include "InterpKernelAutoPtr.hxx"
 #include "InterpKernelGaussCoords.hxx"
 #include "InterpKernelMatrixTools.hxx"
+#include "InterpKernelDenseMatrix.hxx"
 
 #include <set>
 #include <list>
@@ -1711,8 +1712,11 @@ MEDCouplingFieldDouble *MEDCouplingFieldDiscretizationGauss::getMeasureField(con
 {
   if(!mesh)
     throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getMeasureField : mesh instance specified is NULL !");
-  MCAuto<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isAbs);
-  const double *volPtr=vol->getArray()->begin();
+  MCAuto<MEDCouplingUMesh> umesh(mesh->buildUnstructured());
+  const double *coordsOfMesh( umesh->getCoords()->begin() );
+  auto spaceDim(mesh->getSpaceDimension());
+  if( spaceDim != mesh->getMeshDimension())
+    THROW_IK_EXCEPTION("MEDCouplingFieldDiscretizationGauss::getMeasureField : only implemented when space dimension and mesh dimension are equal !");
   MCAuto<MEDCouplingFieldDouble> ret=MEDCouplingFieldDouble::New(ON_GAUSS_PT);
   ret->setMesh(mesh);
   ret->setDiscretization(const_cast<MEDCouplingFieldDiscretizationGauss *>(this));
@@ -1721,7 +1725,7 @@ MEDCouplingFieldDouble *MEDCouplingFieldDiscretizationGauss::getMeasureField(con
   _discr_per_cell->checkAllocated();
   if(_discr_per_cell->getNumberOfComponents()!=1)
     throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getMeasureField : no discr per cell array defined but with nb of components different from 1 !");
-  if(_discr_per_cell->getNumberOfTuples()!=vol->getNumberOfTuples())
+  if(_discr_per_cell->getNumberOfTuples()!=mesh->getNumberOfCells())
     throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getMeasureField : no discr per cell array defined but mismatch between nb of cells of mesh and size of spatial disr array !");
   MCAuto<DataArrayIdType> offset=getOffsetArr(mesh);
   MCAuto<DataArrayDouble> arr=DataArrayDouble::New(); arr->alloc(getNumberOfTuples(mesh),1);
@@ -1741,11 +1745,28 @@ MEDCouplingFieldDouble *MEDCouplingFieldDiscretizationGauss::getMeasureField(con
           const MEDCouplingGaussLocalization& loc=_loc[locId];
           mcIdType nbOfGaussPt=loc.getNumberOfGaussPt();
           INTERP_KERNEL::AutoPtr<double> weights=new double[nbOfGaussPt];
-          double sum=std::accumulate(loc.getWeights().begin(),loc.getWeights().end(),0.);
-          std::transform(loc.getWeights().begin(),loc.getWeights().end(),(double *)weights,std::bind(std::multiplies<double>(),std::placeholders::_1,1./sum));
           for(const mcIdType *cellId=curIds->begin();cellId!=curIds->end();cellId++)
-            for(mcIdType j=0;j<nbOfGaussPt;j++)
-              arrPtr[offsetPtr[*cellId]+j]=weights[j]*volPtr[*cellId];
+          {
+            std::vector<mcIdType> conn;
+            umesh->getNodeIdsOfCell(*cellId,conn);
+            std::vector<double> ptsInCell; ptsInCell.reserve(conn.size()*loc.getDimension());
+            std::for_each( conn.cbegin(), conn.cend(), [spaceDim,coordsOfMesh,&ptsInCell](mcIdType c) { ptsInCell.insert(ptsInCell.end(),coordsOfMesh+c*spaceDim,coordsOfMesh+(c+1)*spaceDim); } );
+            std::size_t nbPtsInCell(ptsInCell.size()/spaceDim);
+            INTERP_KERNEL::DenseMatrix jacobian(spaceDim,spaceDim);
+            MCAuto<DataArrayDouble> shapeFunc = loc.getDerivativeOfShapeFunctionValues();
+            for(mcIdType iGPt = 0 ; iGPt < nbOfGaussPt ; ++iGPt)
+            {
+              for(auto i = 0 ; i < spaceDim ; ++i)
+                for(auto j = 0 ; j < spaceDim ; ++j)
+                {
+                  double res = 0.0;
+                  for( std::size_t k = 0 ; k < nbPtsInCell ; ++k )
+                    res += ptsInCell[spaceDim*k+i] * shapeFunc->getIJ(iGPt,spaceDim*k+j);
+                  jacobian[ i ][ j ] = res;
+                }
+              arrPtr[offsetPtr[*cellId]+iGPt]=std::abs( jacobian.determinant() )*loc.getWeight(FromIdType<int>(iGPt));
+            }
+          }
         }
       else
         {
index e59f547417cd00c1e1a8f506c794c8da6adbfad7..5f1920d96dcefdedc82dffcafddb06bec53eca27 100644 (file)
@@ -143,7 +143,7 @@ double MEDCouplingGaussLocalization::getGaussCoord(int gaussPtIdInCell, int comp
   return _gauss_coord[gaussPtIdInCell*dim+comp];
 }
 
-double MEDCouplingGaussLocalization::getWeight(int gaussPtIdInCell, double newVal) const
+double MEDCouplingGaussLocalization::getWeight(int gaussPtIdInCell) const
 {
   checkCoherencyOfRequest(gaussPtIdInCell,0);
   return _weight[gaussPtIdInCell];
index 49cb19da3ba1fd83fa7aafbe1b10b3ed8d04d6a6..ccffabaf2989d4851b809af95f0962e16da61b76 100644 (file)
@@ -62,7 +62,7 @@ namespace MEDCoupling
     MEDCOUPLING_EXPORT const std::vector<double>& getGaussCoords() const { return _gauss_coord; }
     MEDCOUPLING_EXPORT double getGaussCoord(int gaussPtIdInCell, int comp) const;
     MEDCOUPLING_EXPORT const std::vector<double>& getWeights() const { return _weight; }
-    MEDCOUPLING_EXPORT double getWeight(int gaussPtIdInCell, double newVal) const;
+    MEDCOUPLING_EXPORT double getWeight(int gaussPtIdInCell) const;
     MEDCOUPLING_EXPORT void setRefCoord(int ptIdInCell, int comp, double newVal);
     MEDCOUPLING_EXPORT void setGaussCoord(int gaussPtIdInCell, int comp, double newVal);
     MEDCOUPLING_EXPORT void setWeight(int gaussPtIdInCell, double newVal);
index b4cbaace36ec95d957112a56f39b5e7dc5f870b0..b2850aca8ad9b0ef5c97cfb50047a241fdda832e 100644 (file)
@@ -2768,7 +2768,7 @@ class MEDCouplingBasicsTest4(unittest.TestCase):
         ft=MEDCouplingDataForTest.buildFieldOnGauss_1()
         mea=ft.buildMeasureField(False)
         mea.checkConsistencyLight()
-        self.assertTrue(mea.getArray().isEqual(DataArrayDouble([-0.08504076274779823,-0.06378057206084897,-0.08504076274779869,-0.10630095343474463,-0.12756114412169625,-0.10630095343474734,-0.0637805720608491,-0.0850407627477968,-0.1063009534347449,-0.0850407627477994,-0.10630095343474809,-0.1275611441216954,-0.037205333702161475,-0.037205333702161475,-0.037205333702161475,-0.037205333702161475,-0.047835429045636084,-0.047835429045636084,-0.047835429045636084,-0.047835429045636084,-0.05846552438911087,-0.05846552438911087,-0.05846552438911087,-0.05846552438911087,-0.037205333702161725,-0.037205333702161725,-0.037205333702161725,-0.037205333702161725,-0.047835429045635834,-0.047835429045635834,-0.047835429045635834,-0.047835429045635834,-0.05846552438911058,-0.05846552438911058,-0.05846552438911058,-0.05846552438911058,-0.03879154890291829,-0.03879154890291829,-0.03879154890291829,-0.04120270848015563,-0.04120270848015563,-0.04120270848015563,-0.03393028948486933,-0.03393028948486933,-0.03393028948486933,-0.03151955746491709,-0.03151955746491709,-0.03151955746491709,-0.02424752187358276,-0.02424752187358276,-0.02424752187358276,-0.026657914642918758,-0.026657914642918758,-0.026657914642918758,-0.04120270848015456,-0.04120270848015456,-0.04120270848015456,-0.03879154890291757,-0.03879154890291757,-0.03879154890291757,-0.031519557464916595,-0.031519557464916595,-0.031519557464916595,-0.03393028948487046,-0.03393028948487046,-0.03393028948487046,-0.0266579146429191,-0.0266579146429191,-0.0266579146429191,-0.024247521873582645,-0.024247521873582645,-0.024247521873582645,-0.01851718920904466,-0.01851718920904466,-0.01851718920904466,-0.01851718920904466,-0.029627502734471456,-0.029627502734471456,-0.029627502734471456,-0.029627502734471456,-0.04740400437515433,-0.015150427534672922,-0.015150427534672922,-0.015150427534672922,-0.015150427534672922,-0.024240684055476674,-0.024240684055476674,-0.024240684055476674,-0.024240684055476674,-0.038785094488762675,-0.011783665860301345,-0.011783665860301345,-0.011783665860301345,-0.011783665860301345,-0.018853865376482152,-0.018853865376482152,-0.018853865376482152,-0.018853865376482152,-0.030166184602371443,-0.018517189209044892,-0.018517189209044892,-0.018517189209044892,-0.018517189209044892,-0.029627502734471827,-0.029627502734471827,-0.029627502734471827,-0.029627502734471827,-0.04740400437515492,-0.015150427534672776,-0.015150427534672776,-0.015150427534672776,-0.015150427534672776,-0.02424068405547644,-0.02424068405547644,-0.02424068405547644,-0.02424068405547644,-0.03878509448876231,-0.011783665860301277,-0.011783665860301277,-0.011783665860301277,-0.011783665860301277,-0.01885386537648204,-0.01885386537648204,-0.01885386537648204,-0.01885386537648204,-0.030166184602371266]),1e-14))
+        self.assertTrue(mea.getArray().isEqual(DataArrayDouble([0.08504076274779823, 0.06378057206084897, 0.08504076274779869, 0.10630095343474463, 0.12756114412169625, 0.10630095343474734, 0.0637805720608491, 0.0850407627477968, 0.1063009534347449, 0.0850407627477994, 0.10630095343474809, 0.1275611441216954, 0.034136689498128064, 0.04027397790619449, 0.04027397790619487, 0.034136689498128446, 0.04476678484160338, 0.05090407324967005, 0.050904073249668745, 0.04476678484160208, 0.05539688018507733, 0.061534168593143805, 0.06153416859314442, 0.055396880185077955, 0.03413668949812847, 0.04027397790619517, 0.04027397790619497, 0.03413668949812827, 0.04476678484160235, 0.05090407324966857, 0.05090407324966929, 0.044766784841603076, 0.05539688018507761, 0.06153416859314382, 0.06153416859314359, 0.05539688018507738, 0.04109914601299856, 0.03748636227195687, 0.03771775769378482, 0.03872147624188224, 0.042334259982919814, 0.04261193448911228, 0.03507783977793162, 0.031465056036892154, 0.0353092351997554, 0.03041505840443063, 0.030229942066968697, 0.03384272580800827, 0.02658630560301859, 0.022973521861980706, 0.02311235911507622, 0.02420863583190326, 0.027821419572943488, 0.02800653591040381, 0.04233425998291915, 0.038721476241882075, 0.04261193448911155, 0.0377177576937831, 0.03748636227195874, 0.041099146012996766, 0.03384272580800756, 0.030229942066968624, 0.030415058404430012, 0.03146505603689226, 0.035077839777932385, 0.035309235199757406, 0.0278214195729437, 0.024208635831903073, 0.02800653591040439, 0.023112359115075753, 0.02297352186197992, 0.026586305603019535, 0.019921063519116766, 0.01730003285741078, 0.017300032857410663, 0.019921063519115587, 0.029776877101222114, 0.027364571036248495, 0.0297768771012211, 0.03151042326217079, 0.04709999543873555, 0.016537327484887662, 0.013916296823180364, 0.0139162968231802, 0.01653732748488678, 0.024362899446454346, 0.022012298827300398, 0.02436289944645354, 0.026158151053224943, 0.038536359904420205, 0.013153591450657661, 0.01053256078895069, 0.010532560788950497, 0.013153591450657937, 0.018948921791686592, 0.016660026618353183, 0.01894892179168666, 0.020805878844278196, 0.029972724370104942, 0.01992106351911645, 0.017300032857410823, 0.017300032857409463, 0.01992106351911769, 0.029776877101221757, 0.027364571036247978, 0.029776877101221573, 0.03151042326217261, 0.047099995438736386, 0.016537327484886743, 0.013916296823180036, 0.013916296823180596, 0.0165373274848873, 0.024362899446453434, 0.022012298827300106, 0.024362899446454343, 0.026158151053224426, 0.038536359904419636, 0.013153591450657875, 0.010532560788950653, 0.010532560788951054, 0.013153591450656942, 0.018948921791686835, 0.016660026618353606, 0.01894892179168645, 0.020805878844277627, 0.029972724370105074]),1e-8))
         f=MEDCouplingFieldDouble(ft)
         arr=DataArrayDouble(126,2)
         arr[:, 0] = list(range(126))
@@ -2776,10 +2776,10 @@ class MEDCouplingBasicsTest4(unittest.TestCase):
         arr[:,1]+=1000
         f.setArray(arr)
         f.checkConsistencyLight()
-        self.assertTrue(DataArrayDouble(f.integral(False)).isEqual(DataArrayDouble([-211.66121638700983,-4863.9563007698835]),1e-11))
-        self.assertTrue(DataArrayDouble(f.getWeightedAverageValue()).isEqual(DataArrayDouble([45.4960858131136,1045.496085813114]),1e-11))
-        self.assertTrue(DataArrayDouble(f.normL1()).isEqual(DataArrayDouble([45.49608581311362,1045.496085813114]),1e-11))
-        self.assertTrue(DataArrayDouble(f.normL2()).isEqual(DataArrayDouble([58.16846378340898,1046.1241521947334]),1e-11))
+        self.assertTrue(DataArrayDouble(f.integral(False)).isEqual(DataArrayDouble([211.67679879443182, 4863.855680512835]),1e-11))
+        self.assertTrue(DataArrayDouble(f.getWeightedAverageValue()).isEqual(DataArrayDouble([45.50057170549804, 1045.5005717054983]),1e-11))
+        self.assertTrue(DataArrayDouble(f.normL1()).isEqual(DataArrayDouble([45.50057170549804, 1045.5005717054983]),1e-11))
+        self.assertTrue(DataArrayDouble(f.normL2()).isEqual(DataArrayDouble([58.175073473810194, 1046.1288078361474]),1e-11))
         pass
 
     def testSwig2FieldDiscretizationComputeMeshRestrictionFromTupleIds1(self):
index ecdc04954c580c80e9ada812e4b0053a7346845a..d8815278c38319232d41d32e93bd11229c6e0feb 100644 (file)
@@ -1119,6 +1119,41 @@ class MEDCouplingBasicsTest7(unittest.TestCase):
             delta_Z.abs()
             self.assertTrue(delta_Z.findIdsNotInRange(-1e-5,+1e-5).empty())
 
+        # 2D cells
+        vec = [0.64,0.2]
+
+        for gt in [NORM_QUAD8,NORM_QUAD4,NORM_TRI3,NORM_TRI6]:
+            ref_coord = [list(elt) for elt in MEDCouplingGaussLocalization.GetDefaultReferenceCoordinatesOf(gt).getValuesAsTuple()]
+
+            der_computed = GetDerivative(ref_coord,vec)
+            der_computed.rearrange(2)
+
+            der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps,vec[1]])-GetShapeFunc(ref_coord,vec) ) / eps
+            delta_X = der_computed[:,0]-der_deduced
+            delta_X.abs()
+            self.assertTrue(delta_X.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+            der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+            delta_Y = der_computed[:,1]-der_deduced
+            delta_Y.abs()
+            self.assertTrue(delta_Y.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+        # B version of TRI6, QUAD4 and QUAD8
+        for gt,ref_coord in [(NORM_TRI6,[[0., 0.], [1., 0.], [0., 1.], [0.5, 0.], [0.5, 0.5], [0., 0.5]]),
+            (NORM_QUAD4,[[-1., -1.], [1., -1.], [1., 1.], [-1., 1.]]),(NORM_QUAD8,[[-1., -1.], [1., -1.], [1., 1.], [-1., 1.], [0., -1.], [1., 0.], [0., 1.], [-1., 0.]])]:
+            der_computed = GetDerivative(ref_coord,vec)
+            der_computed.rearrange(2)
+
+            der_deduced = ( GetShapeFunc(ref_coord,[vec[0]+eps,vec[1]])-GetShapeFunc(ref_coord,vec) ) / eps
+            delta_X = der_computed[:,0]-der_deduced
+            delta_X.abs()
+            self.assertTrue(delta_X.findIdsNotInRange(-1e-5,+1e-5).empty())
+
+            der_deduced = ( GetShapeFunc(ref_coord,[vec[0],vec[1]+eps])-GetShapeFunc(ref_coord,vec) ) / eps
+            delta_Y = der_computed[:,1]-der_deduced
+            delta_Y.abs()
+            self.assertTrue(delta_Y.findIdsNotInRange(-1e-5,+1e-5).empty())
+
     def testComputeTriangleHeight0(self):
         arr = DataArrayDouble([0,1])
         m = MEDCouplingCMesh() ; m.setCoords(arr,arr)
index 1793504a52499ce337150f8dadaac754ad6f67d5..d3404b6a4a3764cd8a3c2af32be3ca44c7db0b8e 100644 (file)
@@ -1289,7 +1289,7 @@ namespace MEDCoupling
     const std::vector<double>& getGaussCoords() const;
     double getGaussCoord(int gaussPtIdInCell, int comp) const;
     const std::vector<double>& getWeights() const;
-    double getWeight(int gaussPtIdInCell, double newVal) const;
+    double getWeight(int gaussPtIdInCell) const;
     void setRefCoord(int ptIdInCell, int comp, double newVal);
     void setGaussCoord(int gaussPtIdInCell, int comp, double newVal);
     void setWeight(int gaussPtIdInCell, double newVal);