Implementation of CurveLinear::getMeasureField

diff --git a/src/MEDCoupling/MEDCouplingCurveLinearMesh.cxx b/src/MEDCoupling/MEDCouplingCurveLinearMesh.cxx
index 99af63bc2c862b40ba84148406f8c501890d4deb..c5568f043ec52d7dabad3556307d3e1020353db6 100644 (file)
--- a/src/MEDCoupling/MEDCouplingCurveLinearMesh.cxx
+++ b/src/MEDCoupling/MEDCouplingCurveLinearMesh.cxx
@@ -23,6 +23,8 @@
 #include "MEDCouplingMemArray.hxx"
 #include "MEDCouplingFieldDouble.hxx"
 
+#include "VolSurfUser.txx"
+
 #include <functional>
 #include <algorithm>
 #include <sstream>
@@ -338,7 +340,104 @@ void MEDCouplingCurveLinearMesh::getBoundingBox(double *bbox) const
 
 MEDCouplingFieldDouble *MEDCouplingCurveLinearMesh::getMeasureField(bool isAbs) const
 {
-  throw INTERP_KERNEL::Exception("MEDCouplingCurveLinearMesh::getMeasureField : not implemented yet !");
+  checkCoherency();
+  int meshDim=getMeshDimension();
+  std::string name="MeasureOfMesh_"; name+=getName();
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> field=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
+  field->setName(name.c_str()); field->setMesh(const_cast<MEDCouplingCurveLinearMesh *>(this)); field->synchronizeTimeWithMesh();
+  switch(meshDim)
+    {
+    case 3:
+      { getMeasureFieldMeshDim3(isAbs,field); return field.retn(); }
+    case 2:
+      { getMeasureFieldMeshDim2(isAbs,field); return field.retn(); }
+    case 1:
+      { getMeasureFieldMeshDim1(isAbs,field); return field.retn(); }
+    default:
+      throw INTERP_KERNEL::Exception("MEDCouplingCurveLinearMesh::getMeasureField : space dimension must be in [1,2,3] !");
+    }
+}
+
+/*!
+ * \param [in,out] f field feeded with good values.
+ * \sa MEDCouplingCurveLinearMesh::getMeasureField
+ */
+void MEDCouplingCurveLinearMesh::getMeasureFieldMeshDim1(bool isAbs, MEDCouplingFieldDouble *field) const throw(INTERP_KERNEL::Exception)
+{
+  int nbnodes=getNumberOfNodes();
+  int spaceDim=getSpaceDimension();
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=DataArrayDouble::New(); field->setArray(arr);
+  if(nbnodes==0)
+    { arr->alloc(0,1); return; }
+  if(spaceDim==1)
+    {
+      arr->alloc(nbnodes-1,1);
+      std::transform(_coords->begin()+1,_coords->end(),_coords->begin(),arr->getPointer(),std::minus<double>());
+      if(isAbs)
+        arr->abs();
+    }
+  else
+    {
+      MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> tmp=DataArrayDouble::New(); tmp->alloc(nbnodes-1,spaceDim);
+      std::transform(_coords->begin()+spaceDim,_coords->end(),_coords->begin(),tmp->getPointer(),std::minus<double>());
+      MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> tmp2=tmp->magnitude(); field->setArray(tmp2);
+    }
+}
+
+/*!
+ * \param [in,out] f field feeded with good values.
+ * \sa MEDCouplingCurveLinearMesh::getMeasureField
+ */
+void MEDCouplingCurveLinearMesh::getMeasureFieldMeshDim2(bool isAbs, MEDCouplingFieldDouble *field) const throw(INTERP_KERNEL::Exception)
+{
+  int nbcells=getNumberOfCells();
+  int spaceDim=getSpaceDimension();
+  if(spaceDim!=2 && spaceDim!=3)
+    throw INTERP_KERNEL::Exception("MEDCouplingCurveLinearMesh::getMeasureFieldMeshDim2 : with meshDim 2 only space dimension 2 and 3 are possible !");
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=DataArrayDouble::New(); field->setArray(arr);
+  arr->alloc(nbcells,1);
+  double *pt=arr->getPointer();
+  const double *coords=_coords->begin();
+  int nX=_structure[0]-1;
+  int conn[4];
+  for(int i=0;i<nbcells;i++,pt++)
+    {
+      int cy=i/nX,cx=i-cy*nX;
+      conn[0]=cy*(nX+1)+cx; conn[1]=(cy+1)*(nX+1)+cx; conn[2]=(cy+1)*(nX+1)+1+cx; conn[3]=cy*(nX+1)+cx+1;
+      *pt=INTERP_KERNEL::computeVolSurfOfCell2<int,INTERP_KERNEL::ALL_C_MODE>(INTERP_KERNEL::NORM_QUAD4,conn,4,coords,spaceDim);
+    }
+  if(isAbs)
+    arr->abs();
+}
+
+/*!
+ * \param [in,out] f field feeded with good values.
+ * \sa MEDCouplingCurveLinearMesh::getMeasureField
+ */
+void MEDCouplingCurveLinearMesh::getMeasureFieldMeshDim3(bool isAbs, MEDCouplingFieldDouble *field) const throw(INTERP_KERNEL::Exception)
+{
+  int nbcells=getNumberOfCells();
+  int spaceDim=getSpaceDimension();
+  if(spaceDim!=3)
+    throw INTERP_KERNEL::Exception("MEDCouplingCurveLinearMesh::getMeasureFieldMeshDim3 : with meshDim 3 only space dimension 3 is possible !");
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=DataArrayDouble::New(); field->setArray(arr);
+  arr->alloc(nbcells,1);
+  double *pt=arr->getPointer();
+  const double *coords=_coords->begin();
+  int nX=_structure[0]-1,nY=(_structure[0]-1)*(_structure[1]-1);
+  int nY1=_structure[0]*_structure[1];
+  int conn[8];
+  for(int i=0;i<nbcells;i++,pt++)
+    {
+      int cz=i/nY;
+      int cy=(i-cz*nY)/nX;
+      int cx=(i-cz*nY)-nX*cy;
+      conn[0]=cz*nY1+cy*(nX+1)+cx; conn[1]=cz*nY1+(cy+1)*(nX+1)+cx; conn[2]=cz*nY1+(cy+1)*(nX+1)+1+cx; conn[3]=cz*nY1+cy*(nX+1)+cx+1;
+      conn[4]=(cz+1)*nY1+cy*(nX+1)+cx; conn[5]=(cz+1)*nY1+(cy+1)*(nX+1)+cx; conn[6]=(cz+1)*nY1+(cy+1)*(nX+1)+1+cx; conn[7]=(cz+1)*nY1+cy*(nX+1)+cx+1;
+      *pt=INTERP_KERNEL::computeVolSurfOfCell2<int,INTERP_KERNEL::ALL_C_MODE>(INTERP_KERNEL::NORM_HEXA8,conn,8,coords,3);
+    }
+  if(isAbs)
+    arr->abs();
 }
 
 /*!
@@ -421,8 +520,7 @@ void MEDCouplingCurveLinearMesh::scale(const double *point, double factor)
 
 MEDCouplingMesh *MEDCouplingCurveLinearMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
 {
-  //not implemented yet !
-  return 0;
+  throw INTERP_KERNEL::Exception("MEDCouplingCurveLinearMesh::mergeMyselfWith : not available for CurveLinear Mesh !");
 }
 
 DataArrayDouble *MEDCouplingCurveLinearMesh::getCoordinatesAndOwner() const

diff --git a/src/MEDCoupling/MEDCouplingCurveLinearMesh.hxx b/src/MEDCoupling/MEDCouplingCurveLinearMesh.hxx
index 8a883791bb70fe6426af97c901d106f34e103c4a..1370794944cb473390c4c2460737ec97ff01084f 100644 (file)
--- a/src/MEDCoupling/MEDCouplingCurveLinearMesh.hxx
+++ b/src/MEDCoupling/MEDCouplingCurveLinearMesh.hxx
@@ -85,6 +85,10 @@ namespace ParaMEDMEM
     void serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const;
     void unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
                          const std::vector<std::string>& littleStrings);
+  private:
+    void getMeasureFieldMeshDim1(bool isAbs, MEDCouplingFieldDouble *field) const throw(INTERP_KERNEL::Exception);
+    void getMeasureFieldMeshDim2(bool isAbs, MEDCouplingFieldDouble *field) const throw(INTERP_KERNEL::Exception);
+    void getMeasureFieldMeshDim3(bool isAbs, MEDCouplingFieldDouble *field) const throw(INTERP_KERNEL::Exception);
   private:
     MEDCouplingCurveLinearMesh();
     MEDCouplingCurveLinearMesh(const MEDCouplingCurveLinearMesh& other, bool deepCpy);

diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
index dac1bd32dc96fafe1f618cc9baecbd276fdb76bc..33b4ff2f5768ec6f3fe9a209473252040559bb66 100644 (file)
--- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
+++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
@@ -11044,6 +11044,50 @@ class MEDCouplingBasicsTest(unittest.TestCase):
         self.assertTrue(d2.isEqual(DataArrayInt([0,1,1,1,1,1,1,2,2,2,2,2,2,3])))
         pass
 
+    def testSwigCurveLinearMesh2(self):
+        c=MEDCouplingCMesh()
+        #2D
+        arr1=DataArrayDouble([0,1,3,7])
+        arr2=DataArrayDouble([0,1,1.5])
+        c.setCoords(arr1,arr2)
+        u=c.buildUnstructured()
+        coo=u.getCoords()
+        cl=MEDCouplingCurveLinearMesh()
+        cl.setCoords(coo)
+        cl.setNodeGridStructure([4,3])
+        cl.checkCoherency2()
+        li1=[1.,2.,4.,0.5,1.,2.]
+        self.assertTrue(cl.getMeasureField(False).getArray().isEqual(DataArrayDouble(li1),1e-14))
+        self.assertTrue(u.getMeasureField(False).getArray().isEqual(DataArrayDouble(li1),1e-14))
+        #3D
+        c.setCoords(arr1,arr2,arr2)
+        u=c.buildUnstructured()
+        coo=u.getCoords()
+        cl=MEDCouplingCurveLinearMesh()
+        cl.setCoords(coo)
+        cl.setNodeGridStructure([4,3,3])
+        cl.checkCoherency2()
+        li2=[1.,2.,4.,0.5, 1.,2.,0.5,1.,2.,0.25,0.5,1.]
+        self.assertTrue(cl.getMeasureField(False).getArray().isEqual(DataArrayDouble(li2),1e-14))
+        self.assertTrue(u.getMeasureField(False).getArray().isEqual(DataArrayDouble(li2),1e-14))
+        #1D spaceDim 1
+        coo=DataArrayDouble(5) ; coo.iota(0.)
+        coo=coo*coo
+        cl.setCoords(coo)
+        cl.setNodeGridStructure([5])
+        cl.checkCoherency2()
+        li3=[1.,3.,5.,7.]
+        self.assertTrue(cl.getMeasureField(False).getArray().isEqual(DataArrayDouble(li3),1e-14))
+        self.assertTrue(cl.buildUnstructured().getMeasureField(False).getArray().isEqual(DataArrayDouble(li3),1e-14))
+        #1D spaceDim 2
+        coo=DataArrayDouble.Meld(coo,coo)
+        cl.setCoords(coo)
+        cl.checkCoherency2()
+        li4=[sqrt(2.)*elt for elt in [1.,3.,5.,7.]]
+        self.assertTrue(cl.getMeasureField(False).getArray().isEqual(DataArrayDouble(li4),1e-14))
+        self.assertTrue(cl.buildUnstructured().getMeasureField(False).getArray().isEqual(DataArrayDouble(li4),1e-14))
+        pass
+
     def setUp(self):
         pass
     pass