Salome HOME
A forgotten C++ test
[tools/medcoupling.git] / src / ParaMEDMEMTest / ParaMEDMEMTest_OverlapDEC.cxx
index 63bfc9f9be256db379f92d31d348609848ad66df..f6b14cf496b72f11bb4d509fb395a1532e9bd77a 100644 (file)
@@ -46,6 +46,7 @@ using namespace ParaMEDMEM;
 
 typedef  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> MUMesh;
 typedef  MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> MFDouble;
+typedef  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> DADouble;
 
 //void ParaMEDMEMTest::testOverlapDEC_LMEC_seq()
 //{
@@ -180,17 +181,15 @@ typedef  MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> MFDouble;
 //
 //  MPI_Barrier(MPI_COMM_WORLD);
 //}
-
-void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
-                                  MEDCouplingUMesh *& meshS, MEDCouplingUMesh *& meshT,
-                                  ParaMESH*& parameshS, ParaMESH*& parameshT,
-                                  ParaFIELD*& parafieldS, ParaFIELD*& parafieldT)
+//
+void prepareData1(int rank, NatureOfField nature,
+                  MEDCouplingFieldDouble *& fieldS, MEDCouplingFieldDouble *& fieldT)
 {
   if(rank==0)
     {
       const double coordsS[10]={0.,0.,0.5,0.,1.,0.,0.,0.5,0.5,0.5};
       const double coordsT[6]={0.,0.,1.,0.,1.,1.};
-      meshS=MEDCouplingUMesh::New();
+      MUMesh meshS=MEDCouplingUMesh::New();
       meshS->setMeshDimension(2);
       DataArrayDouble *myCoords=DataArrayDouble::New();
       myCoords->alloc(5,2);
@@ -202,14 +201,14 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
       meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS);
       meshS->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connS+4);
       meshS->finishInsertingCells();
-      ComponentTopology comptopo;
-      parameshS=new ParaMESH(meshS, *grp,"source mesh");
-      parafieldS=new ParaFIELD(ON_CELLS,NO_TIME,parameshS,comptopo);
-      parafieldS->getField()->setNature(nature);
-      double *valsS=parafieldS->getField()->getArray()->getPointer();
+      fieldS = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+      DADouble arr = DataArrayDouble::New(); arr->alloc(meshS->getNumberOfCells(), 1);
+      fieldS->setMesh(meshS); fieldS->setArray(arr);
+      fieldS->setNature(nature);
+      double *valsS=fieldS->getArray()->getPointer();
       valsS[0]=7.; valsS[1]=8.;
       //
-      meshT=MEDCouplingUMesh::New();
+      MUMesh meshT=MEDCouplingUMesh::New();
       meshT->setMeshDimension(2);
       myCoords=DataArrayDouble::New();
       myCoords->alloc(3,2);
@@ -220,10 +219,11 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
       meshT->allocateCells(1);
       meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
       meshT->finishInsertingCells();
-      parameshT=new ParaMESH(meshT,*grp,"target mesh");
-      parafieldT=new ParaFIELD(ON_CELLS,NO_TIME,parameshT,comptopo);
-      parafieldT->getField()->setNature(nature);
-      double *valsT=parafieldT->getField()->getArray()->getPointer();
+      fieldT = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+      DADouble arr2 = DataArrayDouble::New(); arr2->alloc(meshT->getNumberOfCells(), 1);
+      fieldT->setMesh(meshT);  fieldT->setArray(arr2);
+      fieldT->setNature(nature);
+      double *valsT=fieldT->getArray()->getPointer();
       valsT[0]=7.;
     }
   //
@@ -231,7 +231,7 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
     {
       const double coordsS[10]={1.,0.,0.5,0.5,1.,0.5,0.5,1.,1.,1.};
       const double coordsT[6]={0.,0.,0.5,0.5,0.,1.};
-      meshS=MEDCouplingUMesh::New();
+      MUMesh meshS=MEDCouplingUMesh::New();
       meshS->setMeshDimension(2);
       DataArrayDouble *myCoords=DataArrayDouble::New();
       myCoords->alloc(5,2);
@@ -243,15 +243,14 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
       meshS->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connS);
       meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS+3);
       meshS->finishInsertingCells();
-      ComponentTopology comptopo;
-      parameshS=new ParaMESH(meshS,*grp,"source mesh");
-      parafieldS=new ParaFIELD(ON_CELLS,NO_TIME,parameshS,comptopo);
-      parafieldS->getField()->setNature(nature);
-      double *valsS=parafieldS->getField()->getArray()->getPointer();
-      valsS[0]=9.;
-      valsS[1]=11.;
+      fieldS = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+      DADouble arr = DataArrayDouble::New(); arr->alloc(meshS->getNumberOfCells(), 1);
+      fieldS->setMesh(meshS); fieldS->setArray(arr);
+      fieldS->setNature(nature);
+      double *valsS=fieldS->getArray()->getPointer();
+      valsS[0]=9.; valsS[1]=11.;
       //
-      meshT=MEDCouplingUMesh::New();
+      MUMesh meshT=MEDCouplingUMesh::New();
       meshT->setMeshDimension(2);
       myCoords=DataArrayDouble::New();
       myCoords->alloc(3,2);
@@ -262,10 +261,11 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
       meshT->allocateCells(1);
       meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
       meshT->finishInsertingCells();
-      parameshT=new ParaMESH(meshT,*grp,"target mesh");
-      parafieldT=new ParaFIELD(ON_CELLS,NO_TIME,parameshT,comptopo);
-      parafieldT->getField()->setNature(nature);
-      double *valsT=parafieldT->getField()->getArray()->getPointer();
+      fieldT = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+      DADouble arr2 = DataArrayDouble::New(); arr2->alloc(meshT->getNumberOfCells(), 1);
+      fieldT->setMesh(meshT);  fieldT->setArray(arr2);
+      fieldT->setNature(nature);
+      double *valsT=fieldT->getArray()->getPointer();
       valsT[0]=8.;
     }
   //
@@ -273,7 +273,7 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
     {
       const double coordsS[8]={0.,0.5, 0.5,0.5, 0.,1., 0.5,1.};
       const double coordsT[6]={0.5,0.5,0.,1.,1.,1.};
-      meshS=MEDCouplingUMesh::New();
+      MUMesh meshS=MEDCouplingUMesh::New();
       meshS->setMeshDimension(2);
       DataArrayDouble *myCoords=DataArrayDouble::New();
       myCoords->alloc(4,2);
@@ -284,14 +284,14 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
       meshS->allocateCells(1);
       meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS);
       meshS->finishInsertingCells();
-      ComponentTopology comptopo;
-      parameshS=new ParaMESH(meshS,*grp,"source mesh");
-      parafieldS=new ParaFIELD(ON_CELLS,NO_TIME,parameshS,comptopo);
-      parafieldS->getField()->setNature(nature);
-      double *valsS=parafieldS->getField()->getArray()->getPointer();
+      fieldS = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+      DADouble arr = DataArrayDouble::New(); arr->alloc(meshS->getNumberOfCells(), 1);
+      fieldS->setMesh(meshS); fieldS->setArray(arr);
+      fieldS->setNature(nature);
+      double *valsS=fieldS->getArray()->getPointer();
       valsS[0]=10.;
       //
-      meshT=MEDCouplingUMesh::New();
+      MUMesh meshT=MEDCouplingUMesh::New();
       meshT->setMeshDimension(2);
       myCoords=DataArrayDouble::New();
       myCoords->alloc(3,2);
@@ -302,152 +302,15 @@ void prepareData1(int rank, ProcessorGroup * grp, NatureOfField nature,
       meshT->allocateCells(1);
       meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
       meshT->finishInsertingCells();
-      parameshT=new ParaMESH(meshT,*grp,"target mesh");
-      parafieldT=new ParaFIELD(ON_CELLS,NO_TIME,parameshT,comptopo);
-      parafieldT->getField()->setNature(nature);
-      double *valsT=parafieldT->getField()->getArray()->getPointer();
+      fieldT = MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+      DADouble arr2 = DataArrayDouble::New(); arr2->alloc(meshT->getNumberOfCells(), 1);
+      fieldT->setMesh(meshT); fieldT->setArray(arr2);
+      fieldT->setNature(nature);
+      double *valsT=fieldT->getArray()->getPointer();
       valsT[0]=9.;
     }
 }
 
-/*! Test case from the official doc of the OverlapDEC.
- *  WARNING: bounding boxes are tweaked here to make the case more interesting (i.e. to avoid an all to all exchange
- *  between all procs).
- */
-void ParaMEDMEMTest::testOverlapDEC1()
-{
-  int size, rank;
-  MPI_Comm_size(MPI_COMM_WORLD,&size);
-  MPI_Comm_rank(MPI_COMM_WORLD,&rank);
-  //  char hostname[256];
-  //  printf("(%d) PID %d on localhost ready for attach\n", rank, getpid());
-  //  fflush(stdout);
-
-//    if (rank == 1)
-//      {
-//        int i=1, j=0;
-//        while (i!=0)
-//          j=2;
-//      }
-
-  if (size != 3) return ;
-  int nproc = 3;
-  std::set<int> procs;
-  for (int i=0; i<nproc; i++)
-    procs.insert(i);
-
-  CommInterface interface;
-  OverlapDEC dec(procs);
-  ProcessorGroup * grp = dec.getGroup();
-  MEDCouplingUMesh* meshS=0, *meshT=0;
-  ParaMESH* parameshS=0, *parameshT=0;
-  ParaFIELD* parafieldS=0, *parafieldT=0;
-
-  MPI_Barrier(MPI_COMM_WORLD);
-  prepareData1(rank, grp, ConservativeVolumic, meshS, meshT, parameshS, parameshT, parafieldS, parafieldT);
-
-  /*!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-   *  HACK ON BOUNDING BOX TO MAKE THIS CASE SIMPLE AND USABLE IN DEBUG
-   * Bounding boxes are slightly smaller than should be, thus localizing the work to be done
-   * and avoiding every proc talking to everyone else.
-   * Obviously this is NOT a good idea to do this in production code :-)
-   * !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-   */
-  dec.setBoundingBoxAdjustmentAbs(-1.0e-12);
-
-  dec.attachSourceLocalField(parafieldS);
-  dec.attachTargetLocalField(parafieldT);
-  dec.synchronize();
-  dec.sendRecvData(true);
-  //
-  MPI_Barrier(MPI_COMM_WORLD);
-  if(rank==0)
-    {
-      CPPUNIT_ASSERT_DOUBLES_EQUAL(8.75,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
-    }
-  if(rank==1)
-    {
-      CPPUNIT_ASSERT_DOUBLES_EQUAL(8.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
-    }
-  if(rank==2)
-    {
-      CPPUNIT_ASSERT_DOUBLES_EQUAL(10.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
-    }
-
-  delete parafieldS;
-  delete parafieldT;
-  delete parameshS;
-  delete parameshT;
-  meshS->decrRef();
-  meshT->decrRef();
-
-  MPI_Barrier(MPI_COMM_WORLD);
-}
-
-/*!
- * Same as testOverlapDEC1() but with regular bounding boxes. If you're looking for a nice debug case,
- * testOverlapDEC1() is identical in terms of geometry and field values, and more appropriate.
- */
-void ParaMEDMEMTest::testOverlapDEC2()
-{
-  int size, rank;
-  MPI_Comm_size(MPI_COMM_WORLD,&size);
-  MPI_Comm_rank(MPI_COMM_WORLD,&rank);
-
-  if (size != 3) return ;
-  int nproc = 3;
-  std::set<int> procs;
-  for (int i=0; i<nproc; i++)
-    procs.insert(i);
-
-//      if (rank == 0)
-//        {
-//          int i=1, j=0;
-//          while (i!=0)
-//            j=2;
-//        }
-
-
-  CommInterface interface;
-  OverlapDEC dec(procs);
-  ProcessorGroup * grp = dec.getGroup();
-  MEDCouplingUMesh* meshS=0, *meshT=0;
-  ParaMESH* parameshS=0, *parameshT=0;
-  ParaFIELD* parafieldS=0, *parafieldT=0;
-
-  MPI_Barrier(MPI_COMM_WORLD);
-  prepareData1(rank, grp, ConservativeVolumic,meshS, meshT, parameshS, parameshT, parafieldS, parafieldT);
-
-  /* Normal bounding boxes here: */
-  dec.setBoundingBoxAdjustmentAbs(+1.0e-12);
-
-  dec.attachSourceLocalField(parafieldS);
-  dec.attachTargetLocalField(parafieldT);
-  dec.synchronize();
-  dec.sendRecvData(true);
-  //
-  if(rank==0)
-    {
-      CPPUNIT_ASSERT_DOUBLES_EQUAL(8.75,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
-    }
-  if(rank==1)
-    {
-      CPPUNIT_ASSERT_DOUBLES_EQUAL(8.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
-    }
-  if(rank==2)
-    {
-      CPPUNIT_ASSERT_DOUBLES_EQUAL(10.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
-    }
-  delete parafieldS;
-  delete parafieldT;
-  delete parameshS;
-  delete parameshT;
-  meshS->decrRef();
-  meshT->decrRef();
-
-  MPI_Barrier(MPI_COMM_WORLD);
-}
-
 void prepareData2_buildOneSquare(MEDCouplingUMesh* & meshS_0, MEDCouplingUMesh* & meshT_0)
 {
   const double coords[10] = {0.0,0.0,  0.0,1.0,  1.0,1.0,  1.0,0.0, 0.5,0.5};
@@ -472,7 +335,6 @@ void prepareData2_buildOneSquare(MEDCouplingUMesh* & meshS_0, MEDCouplingUMesh*
   meshT_0->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT+3);
   meshT_0->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT+6);
   meshT_0->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT+9);
-
 }
 
 /**
@@ -589,6 +451,113 @@ void prepareData2(int rank, ProcessorGroup * grp, NatureOfField nature,
   meshT_0->decrRef();
 }
 
+/*! Test case from the official doc of the OverlapDEC.
+ *  WARNING: bounding boxes might be tweaked here to make the case more interesting (i.e. to avoid an all to all exchange
+ *  between all procs).
+ */
+void testOverlapDEC_generic(int workSharingAlgo, double bbAdj)
+{
+  int size, rank;
+  MPI_Comm_size(MPI_COMM_WORLD,&size);
+  MPI_Comm_rank(MPI_COMM_WORLD,&rank);
+  //  char hostname[256];
+  //  printf("(%d) PID %d on localhost ready for attach\n", rank, getpid());
+  //  fflush(stdout);
+
+//    if (rank == 0)
+//      {
+//        int i=1, j=0;
+//        while (i!=0)
+//          j=2;
+//      }
+
+  if (size != 3) return ;
+  int nproc = 3;
+  std::set<int> procs;
+  for (int i=0; i<nproc; i++)
+    procs.insert(i);
+
+  CommInterface interface;
+  OverlapDEC dec(procs);
+  MEDCouplingFieldDouble * mcfieldS=0, *mcfieldT=0;
+
+  prepareData1(rank, ConservativeVolumic, mcfieldS, mcfieldT);
+
+  // See comment in the caller:
+  dec.setBoundingBoxAdjustmentAbs(bbAdj);
+  dec.setWorkSharingAlgo(workSharingAlgo);  // just to ease debugging
+
+  dec.attachSourceLocalField(mcfieldS);
+  dec.attachTargetLocalField(mcfieldT);
+  dec.synchronize();
+//  dec.debugPrintWorkSharing(std::cout);
+  dec.sendRecvData(true);
+  //
+  MPI_Barrier(MPI_COMM_WORLD);
+  if(rank==0)
+    {
+      CPPUNIT_ASSERT_DOUBLES_EQUAL(8.75,mcfieldT->getArray()->getIJ(0,0),1e-12);
+    }
+  if(rank==1)
+    {
+      CPPUNIT_ASSERT_DOUBLES_EQUAL(8.5,mcfieldT->getArray()->getIJ(0,0),1e-12);
+    }
+  if(rank==2)
+    {
+      CPPUNIT_ASSERT_DOUBLES_EQUAL(10.5,mcfieldT->getArray()->getIJ(0,0),1e-12);
+    }
+
+  mcfieldS->decrRef();
+  mcfieldT->decrRef();
+
+  MPI_Barrier(MPI_COMM_WORLD);
+}
+
+void ParaMEDMEMTest::testOverlapDEC1()
+{
+  /*!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   *  HACK ON BOUNDING BOX TO MAKE THIS CASE SIMPLE AND USABLE IN DEBUG
+   * Bounding boxes are slightly smaller than should be, thus localizing the work to be done
+   * and avoiding every proc talking to everyone else.
+   * Obviously this is NOT a good idea to do this in production code :-)
+   * !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   */
+  testOverlapDEC_generic(0,-1.0e-12);
+}
+
+void ParaMEDMEMTest::testOverlapDEC1_bis()
+{
+   // Same BB hack as above
+  testOverlapDEC_generic(1,-1.0e-12);
+}
+
+void ParaMEDMEMTest::testOverlapDEC1_ter()
+{
+   // Same BB hack as above
+  testOverlapDEC_generic(2, -1.0e-12);
+}
+
+
+/*!
+ * Same as testOverlapDEC1() but with regular bounding boxes. If you're looking for a nice debug case,
+ * testOverlapDEC1() is identical in terms of geometry and field values, and more appropriate.
+ */
+void ParaMEDMEMTest::testOverlapDEC2()
+{
+  testOverlapDEC_generic(0,1.0e-12);
+}
+
+void ParaMEDMEMTest::testOverlapDEC2_bis()
+{
+  testOverlapDEC_generic(1,1.0e-12);
+}
+
+void ParaMEDMEMTest::testOverlapDEC2_ter()
+{
+  testOverlapDEC_generic(2,1.0e-12);
+}
+
+
 /*! Test focused on the mapping of cell IDs.
  * (i.e. when only part of the source/target mesh is transmitted)
  */
@@ -731,6 +700,5 @@ void ParaMEDMEMTest::testOverlapDEC4()
   meshT->decrRef();
 
   MPI_Barrier(MPI_COMM_WORLD);
-
 }