Salome HOME
Debug on GENERAL_24 and GENERAL_48
[tools/medcoupling.git] / src / INTERP_KERNELTest / MeshTestToolkit.txx
index 84bb2df8e325afcae121562acc181b607300717a..0c551bbfa83ef952cb09d86bdc6bed2f931219f0 100644 (file)
@@ -1,41 +1,42 @@
-//  Copyright (C) 2007-2008  CEA/DEN, EDF R&D
+// Copyright (C) 2007-2013  CEA/DEN, EDF R&D
 //
-//  This library is free software; you can redistribute it and/or
-//  modify it under the terms of the GNU Lesser General Public
-//  License as published by the Free Software Foundation; either
-//  version 2.1 of the License.
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
 //
-//  This library is distributed in the hope that it will be useful,
-//  but WITHOUT ANY WARRANTY; without even the implied warranty of
-//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-//  Lesser General Public License for more details.
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
 //
-//  You should have received a copy of the GNU Lesser General Public
-//  License along with this library; if not, write to the Free Software
-//  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
 //
-//  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
 //
-#include "MEDNormalizedUnstructuredMesh.hxx"
-#include "MEDNormalizedUnstructuredMesh.txx"
+#include "TestInterpKernelUtils.hxx"
 
 #include "MeshTestToolkit.hxx"
-#include "MEDMEM_Mesh.hxx"
 
-#include "Interpolation3DSurf.txx"
+#include "MEDFileMesh.hxx"
+
+#include "MEDCouplingNormalizedUnstructuredMesh.hxx"
+#include "MEDCouplingNormalizedUnstructuredMesh.txx"
+#include "MEDCouplingFieldDouble.hxx"
+
+#include "Interpolation3DSurf.hxx"
 #include "Interpolation2D.txx"
 #include "Interpolation3D.txx"
 
-#include <iostream>
 #include <map>
-#include <vector>
 #include <cmath>
+#include <vector>
+#include <cstring>
+#include <iostream>
 #include <algorithm>
 
-//#include "VectorUtils.hxx"
-
-#include "MEDMEM_Field.hxx"
-#include "MEDMEM_Support.hxx"
 
 // levels : 
 // 1 - titles and volume results
 #include <cppunit/extensions/HelperMacros.h>
 
 //#define VOL_PREC 1.0e-6
-
-using namespace MEDMEM;
-using namespace std;
-using namespace MED_EN;
+using namespace ParaMEDMEM;
 using namespace INTERP_KERNEL;
 
 namespace INTERP_TEST
 {
-
   /**
    * Calculates the sum of a row of an intersection matrix
    *
@@ -73,7 +70,7 @@ namespace INTERP_TEST
       {
         if(iter->count(i) != 0.0)
           {
-            map<int, double>::const_iterator iter2 = iter->find(i);
+            std::map<int, double>::const_iterator iter2 = iter->find(i);
             vol += fabs(iter2->second);
           }
       }
@@ -93,7 +90,7 @@ namespace INTERP_TEST
   {
     double vol = 0.0;
     const std::map<int, double>& col = m[i];
-    for(map<int, double>::const_iterator iter = col.begin() ; iter != col.end() ; ++iter)
+    for(std::map<int, double>::const_iterator iter = col.begin() ; iter != col.end() ; ++iter)
       {
         vol += fabs(iter->second);
       }
@@ -107,23 +104,10 @@ namespace INTERP_TEST
    * @param tab    pointer to double[no. elements of mesh] array in which to store the volumes
    */
   template <int SPACEDIM, int MESHDIM>
-  void MeshTestToolkit<SPACEDIM,MESHDIM>::getVolumes(MEDMEM::MESH& mesh, double* tab) const
+  void MeshTestToolkit<SPACEDIM,MESHDIM>::getVolumes(ParaMEDMEM::MEDCouplingUMesh& mesh, double *tab) const
   {
-    SUPPORT *sup=new SUPPORT(&mesh,"dummy",MED_CELL);
-    FIELD<double>* f;
-    switch (MESHDIM)
-      {
-      case 2:
-        f=mesh.getArea(sup);
-        break;
-      case 3:
-        f=mesh.getVolume(sup);
-        break;
-      }
-    const double *tabS = f->getValue();
-    std::copy(tabS,tabS+mesh.getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS),tab);
-    delete sup;
-    delete f;
+    MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=mesh.getMeasureField(true);
+    std::copy(vol->getArray()->begin(),vol->getArray()->end(),tab);
   }
 
   /**
@@ -136,21 +120,20 @@ namespace INTERP_TEST
   template <int SPACEDIM, int MESHDIM>
   double MeshTestToolkit<SPACEDIM,MESHDIM>::sumVolume(const IntersectionMatrix& m) const
   {
-    
-    vector<double> volumes;
+    std::vector<double> volumes;
     for(IntersectionMatrix::const_iterator iter = m.begin() ; iter != m.end() ; ++iter)
       {
-        for(map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
+        for(std::map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
           {
             volumes.push_back(fabs(iter2->second));
           }
       }
-    
+
     // sum in ascending order to avoid rounding errors
-    
+
     sort(volumes.begin(), volumes.end());
     const double vol = accumulate(volumes.begin(), volumes.end(), 0.0);
-    
+
     return vol;
   }
 
@@ -167,17 +150,17 @@ namespace INTERP_TEST
    * @return true if the condition is verified, false if not.
    */
   template <int SPACEDIM, int MESHDIM>
-  bool MeshTestToolkit<SPACEDIM,MESHDIM>::testVolumes(const IntersectionMatrix& m,  MEDMEM::MESH& sMesh,  MEDMEM::MESH& tMesh) const
+  bool MeshTestToolkit<SPACEDIM,MESHDIM>::testVolumes(const IntersectionMatrix& m,  ParaMEDMEM::MEDCouplingUMesh& sMesh,  ParaMEDMEM::MEDCouplingUMesh& tMesh) const
   {
     bool ok = true;
 
     // source elements
-    double* sVol = new double[sMesh.getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS)];
+    double* sVol = new double[sMesh.getNumberOfCells()];
     getVolumes(sMesh, sVol);
 
-    for(int i = 0; i < sMesh.getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS); ++i)
+    for(int i = 0; i < sMesh.getNumberOfCells(); ++i)
       {
-        const double sum_row = sumRow(m, i+1);
+        const double sum_row = sumRow(m, i);
         if(!epsilonEqualRelative(sum_row, fabs(sVol[i]), _precision))
           {
             LOG(1, "Source volume inconsistent : vol of cell " << i << " = " << sVol[i] << " but the row sum is " << sum_row );
@@ -187,9 +170,9 @@ namespace INTERP_TEST
       }
 
     // target elements
-    double* tVol = new double[tMesh.getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS)];
+    double* tVol = new double[tMesh.getNumberOfCells()];
     getVolumes(tMesh, tVol);
-    for(int i = 0; i < tMesh.getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS); ++i)
+    for(int i = 0; i < tMesh.getNumberOfCells(); ++i)
       {
         const double sum_col = sumCol(m, i);
         if(!epsilonEqualRelative(sum_col,fabs(tVol[i]), _precision))
@@ -220,14 +203,14 @@ namespace INTERP_TEST
     int i = 0;
     for(IntersectionMatrix::const_iterator iter = m1.begin() ; iter != m1.end() ; ++iter)
       {
-        for(map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
+        for(std::map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
           {
             int j = iter2->first;
-            if(m2.at(j-1).count(i+1) == 0)
+            if(m2.at(j).count(i) == 0)
               {
                 if(!epsilonEqual(iter2->second, 0.0, _precision))
                   {
-                    LOG(2, "V1( " << i << ", " << j << ") exists, but V2( " << j - 1 << ", " << i + 1 << ") " << " does not " );
+                    LOG(2, "V1( " << i << ", " << j << ") exists, but V2( " << j << ", " << i << ") " << " does not " );
                     LOG(2, "(" << i << ", " << j << ") fails");
                     compatitable = false;
                   }
@@ -241,7 +224,7 @@ namespace INTERP_TEST
       }
     return compatitable;
   }
-      
+
   /**
    * Tests if two intersection matrices are each others' transposes.
    *
@@ -252,7 +235,6 @@ namespace INTERP_TEST
   template <int SPACEDIM, int MESHDIM>
   bool MeshTestToolkit<SPACEDIM,MESHDIM>::testTranspose(const IntersectionMatrix& m1, const IntersectionMatrix& m2) const
   {
-
     int i = 0;
     bool isSymmetric = true;
 
@@ -263,17 +245,17 @@ namespace INTERP_TEST
 
     for(IntersectionMatrix::const_iterator iter = m1.begin() ; iter != m1.end() ; ++iter)
       {
-        for(map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
+        for(std::map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
           {
             int j = iter2->first;
             const double v1 = fabs(iter2->second);
             //if(m2[j - 1].count(i+1) > 0)
             //  {
-            map<int, double> theMap =  m2.at(j-1);
-            const double v2 = fabs(theMap[i + 1]); 
+            std::map<int, double> theMap =  m2.at(j);
+            const double v2 = fabs(theMap[i]); 
             if(v1 != v2)
               {
-                LOG(2, "V1( " << i << ", " << j << ") = " << v1 << " which is different from V2( " << j - 1 << ", " << i + 1 << ") = " << v2 << " | diff = " << v1 - v2 );
+                LOG(2, "V1( " << i << ", " << j << ") = " << v1 << " which is different from V2( " << j << ", " << i << ") = " << v2 << " | diff = " << v1 - v2 );
                 if(!epsilonEqualRelative(v1, v2, _precision))
                   {
                     LOG(2, "(" << i << ", " << j << ") fails");
@@ -301,17 +283,17 @@ namespace INTERP_TEST
   bool MeshTestToolkit<SPACEDIM,MESHDIM>::testDiagonal(const IntersectionMatrix& m) const
   {
     LOG(1, "Checking if matrix is diagonal" );
-    int i = 1;
+    int i = 0;
     bool isDiagonal = true;
     for(IntersectionMatrix::const_iterator iter = m.begin() ; iter != m.end() ; ++iter)
       {
-        for(map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
+        for(std::map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
           {
             int j = iter2->first;
             const double vol = iter2->second;
             if(vol != 0.0 && (i != j))
               {
-                LOG(2, "V( " << i - 1 << ", " << j << ") = " << vol << " which is not zero" );
+                LOG(2, "V( " << i << ", " << j << ") = " << vol << " which is not zero" );
                 if(!epsilonEqual(vol, 0.0, _precision))
                   {
                     LOG(2, "(" << i << ", " << j << ") fails");
@@ -337,14 +319,12 @@ namespace INTERP_TEST
   void MeshTestToolkit<SPACEDIM,MESHDIM>::dumpIntersectionMatrix(const IntersectionMatrix& m) const
   {
     int i = 0;
-    std::cout << "Intersection matrix is " << endl;
+    std::cout << "Intersection matrix is " << std::endl;
     for(IntersectionMatrix::const_iterator iter = m.begin() ; iter != m.end() ; ++iter)
       {
-        for(map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
+        for(std::map<int, double>::const_iterator iter2 = iter->begin() ; iter2 != iter->end() ; ++iter2)
           {
-    
-            std::cout << "V(" << i << ", " << iter2->first << ") = " << iter2->second << endl;
-    
+            std::cout << "V(" << i << ", " << iter2->first << ") = " << iter2->second << std::endl;
           }
         ++i;
       }
@@ -364,19 +344,18 @@ namespace INTERP_TEST
   template <int SPACEDIM, int MESHDIM>
   void MeshTestToolkit<SPACEDIM,MESHDIM>::calcIntersectionMatrix(const char* mesh1path, const char* mesh1, const char* mesh2path, const char* mesh2, IntersectionMatrix& m) const
   {
-    const string dataBaseDir = getenv("MED_ROOT_DIR");
-    const string dataDir = dataBaseDir + string("/share/salome/resources/med/");
-
     LOG(1, std::endl << "=== -> intersecting src = " << mesh1path << ", target = " << mesh2path );
 
     LOG(5, "Loading " << mesh1 << " from " << mesh1path);
-    MESH sMesh(MED_DRIVER, dataDir+mesh1path, mesh1);
-  
+    MEDCouplingAutoRefCountObjectPtr<MEDFileUMesh> sMeshML=MEDFileUMesh::New(INTERP_TEST::getResourceFile(mesh1path).c_str(),mesh1);
+    MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> sMesh=sMeshML->getMeshAtLevel(0);
+
     LOG(5, "Loading " << mesh2 << " from " << mesh2path);
-    MESH tMesh(MED_DRIVER, dataDir+mesh2path, mesh2);
+    MEDCouplingAutoRefCountObjectPtr<MEDFileUMesh> tMeshML=MEDFileUMesh::New(INTERP_TEST::getResourceFile(mesh2path).c_str(),mesh2);
+    MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> tMesh=tMeshML->getMeshAtLevel(0);
 
-    MEDNormalizedUnstructuredMesh<SPACEDIM,MESHDIM> sMesh_wrapper(&sMesh);
-    MEDNormalizedUnstructuredMesh<SPACEDIM,MESHDIM> tMesh_wrapper(&tMesh);
+    MEDCouplingNormalizedUnstructuredMesh<SPACEDIM,MESHDIM> sMesh_wrapper(sMesh);
+    MEDCouplingNormalizedUnstructuredMesh<SPACEDIM,MESHDIM> tMesh_wrapper(tMesh);
 
     if (SPACEDIM==2 && MESHDIM==2)
       {
@@ -397,19 +376,18 @@ namespace INTERP_TEST
       }
     else
       {
-        throw MEDEXCEPTION("Wrong dimensions");
+        throw INTERP_KERNEL::Exception("Wrong dimensions");
       }
     // if reflexive, check volumes
     if(strcmp(mesh1path,mesh2path) == 0)
       {
-        const bool row_and_col_sums_ok = testVolumes(m, sMesh, tMesh);
+        const bool row_and_col_sums_ok = testVolumes(m, *sMesh, *tMesh);
         CPPUNIT_ASSERT_EQUAL_MESSAGE("Row or column sums incorrect", true, row_and_col_sums_ok);
         const bool is_diagonal =testDiagonal(m);
         CPPUNIT_ASSERT_EQUAL_MESSAGE("Self intersection matrix is not diagonal", true, is_diagonal);
       }
 
     LOG(1, "Intersection calculation done. " << std::endl );
-  
   }
 
   /**
@@ -462,14 +440,13 @@ namespace INTERP_TEST
       }
     else
       {
-      
         IntersectionMatrix matrix2;
         calcIntersectionMatrix(mesh2path, mesh2, mesh1path, mesh1, matrix2);    
 
 #if LOG_LEVEL >= 2
         dumpIntersectionMatrix(matrix2);
 #endif
-      
+
         const double vol2 = sumVolume(matrix2);
 
         LOG(1, "vol1 =  " << vol1 << ", vol2 = " << vol2 << ", correctVol = " << correctVol );
@@ -479,7 +456,6 @@ namespace INTERP_TEST
         CPPUNIT_ASSERT_DOUBLES_EQUAL(correctVol, vol2, prec * std::max(vol2, correctVol));
         CPPUNIT_ASSERT_DOUBLES_EQUAL(vol1, vol2, prec * std::max(vol1, vol2));
       }
-
   }
 
   /**
@@ -498,12 +474,10 @@ namespace INTERP_TEST
   template <int SPACEDIM, int MESHDIM>
   void MeshTestToolkit<SPACEDIM,MESHDIM>::intersectMeshes(const char* mesh1, const char* mesh2, const double correctVol, const double prec, bool doubleTest) const
   {
-    const string path1 = string(mesh1) + string(".med");
+    const std::string path1 = std::string(mesh1) + std::string(".med");
     std::cout << "here :" << path1 << std::endl;
-    const string path2 = string(mesh2) + string(".med");
+    const std::string path2 = std::string(mesh2) + std::string(".med");
 
     intersectMeshes(path1.c_str(), mesh1, path2.c_str(), mesh2, correctVol, prec, doubleTest);
   }
-
-
 }