Back to index

salome-med  6.5.0
ParaMEDMEMTest_OverlapDEC.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include "ParaMEDMEMTest.hxx"
00021 #include <cppunit/TestAssert.h>
00022 
00023 #include "CommInterface.hxx"
00024 #include "ProcessorGroup.hxx"
00025 #include "MPIProcessorGroup.hxx"
00026 #include "Topology.hxx"
00027 #include "OverlapDEC.hxx"
00028 #include "ParaMESH.hxx"
00029 #include "ParaFIELD.hxx"
00030 #include "ComponentTopology.hxx"
00031 
00032 #include "MEDCouplingUMesh.hxx"
00033 
00034 #include <set>
00035 
00036 void ParaMEDMEMTest::testOverlapDEC1()
00037 {
00038   std::string srcM("P0");
00039   std::string targetM("P0");
00040   int size;
00041   int rank;
00042   MPI_Comm_size(MPI_COMM_WORLD,&size);
00043   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00044 
00045   if (size != 3) return ;
00046    
00047   int nproc = 3;
00048   std::set<int> procs;
00049   
00050   for (int i=0; i<nproc; i++)
00051     procs.insert(i);
00052   
00053   ParaMEDMEM::CommInterface interface;
00054 
00055   ParaMEDMEM::OverlapDEC dec(procs);
00056 
00057   ParaMEDMEM::MEDCouplingUMesh* meshS=0;
00058   ParaMEDMEM::MEDCouplingUMesh* meshT=0;
00059   ParaMEDMEM::ParaMESH* parameshS=0;
00060   ParaMEDMEM::ParaMESH* parameshT=0;
00061   ParaMEDMEM::ParaFIELD* parafieldS=0;
00062   ParaMEDMEM::ParaFIELD* parafieldT=0;
00063   
00064   MPI_Barrier(MPI_COMM_WORLD);
00065   if(rank==0)
00066     {
00067       const double coordsS[10]={0.,0.,0.5,0.,1.,0.,0.,0.5,0.5,0.5};
00068       const double coordsT[6]={0.,0.,1.,0.,1.,1.};
00069       meshS=ParaMEDMEM::MEDCouplingUMesh::New();
00070       meshS->setMeshDimension(2);
00071       ParaMEDMEM::DataArrayDouble *myCoords=ParaMEDMEM::DataArrayDouble::New();
00072       myCoords->alloc(5,2);
00073       std::copy(coordsS,coordsS+10,myCoords->getPointer());
00074       meshS->setCoords(myCoords);
00075       myCoords->decrRef();
00076       int connS[7]={0,3,4,1, 1,4,2};
00077       meshS->allocateCells(2);
00078       meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS);
00079       meshS->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connS+4);
00080       meshS->finishInsertingCells();
00081       ParaMEDMEM::ComponentTopology comptopo;
00082       parameshS=new ParaMEDMEM::ParaMESH(meshS,*dec.getGrp(),"source mesh");
00083       parafieldS=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshS,comptopo);
00084       parafieldS->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
00085       double *valsS=parafieldS->getField()->getArray()->getPointer();
00086       valsS[0]=7.; valsS[1]=8.;
00087       //
00088       meshT=ParaMEDMEM::MEDCouplingUMesh::New();
00089       meshT->setMeshDimension(2);
00090       myCoords=ParaMEDMEM::DataArrayDouble::New();
00091       myCoords->alloc(3,2);
00092       std::copy(coordsT,coordsT+6,myCoords->getPointer());
00093       meshT->setCoords(myCoords);
00094       myCoords->decrRef();
00095       int connT[3]={0,2,1};
00096       meshT->allocateCells(1);
00097       meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
00098       meshT->finishInsertingCells();
00099       parameshT=new ParaMEDMEM::ParaMESH(meshT,*dec.getGrp(),"target mesh");
00100       parafieldT=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshT,comptopo);
00101       parafieldT->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
00102       double *valsT=parafieldT->getField()->getArray()->getPointer();
00103       valsT[0]=7.;
00104     }
00105   //
00106   if(rank==1)
00107     {
00108       const double coordsS[10]={1.,0.,0.5,0.5,1.,0.5,0.5,1.,1.,1.};
00109       const double coordsT[6]={0.,0.,0.5,0.5,0.,1.};
00110       meshS=ParaMEDMEM::MEDCouplingUMesh::New();
00111       meshS->setMeshDimension(2);
00112       ParaMEDMEM::DataArrayDouble *myCoords=ParaMEDMEM::DataArrayDouble::New();
00113       myCoords->alloc(5,2);
00114       std::copy(coordsS,coordsS+10,myCoords->getPointer());
00115       meshS->setCoords(myCoords);
00116       myCoords->decrRef();
00117       int connS[7]={0,1,2, 1,3,4,2};
00118       meshS->allocateCells(2);
00119       meshS->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connS);
00120       meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS+3);
00121       meshS->finishInsertingCells();
00122       ParaMEDMEM::ComponentTopology comptopo;
00123       parameshS=new ParaMEDMEM::ParaMESH(meshS,*dec.getGrp(),"source mesh");
00124       parafieldS=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshS,comptopo);
00125       parafieldS->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
00126       double *valsS=parafieldS->getField()->getArray()->getPointer();
00127       valsS[0]=9.; valsS[1]=11.;
00128       //
00129       meshT=ParaMEDMEM::MEDCouplingUMesh::New();
00130       meshT->setMeshDimension(2);
00131       myCoords=ParaMEDMEM::DataArrayDouble::New();
00132       myCoords->alloc(3,2);
00133       std::copy(coordsT,coordsT+6,myCoords->getPointer());
00134       meshT->setCoords(myCoords);
00135       myCoords->decrRef();
00136       int connT[3]={0,2,1};
00137       meshT->allocateCells(1);
00138       meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
00139       meshT->finishInsertingCells();
00140       parameshT=new ParaMEDMEM::ParaMESH(meshT,*dec.getGrp(),"target mesh");
00141       parafieldT=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshT,comptopo);
00142       parafieldT->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
00143       double *valsT=parafieldT->getField()->getArray()->getPointer();
00144       valsT[0]=8.;
00145     }
00146   //
00147   if(rank==2)
00148     {
00149       const double coordsS[8]={0.,0.5, 0.5,0.5, 0.,1., 0.5,1.};
00150       const double coordsT[6]={0.5,0.5,0.,1.,1.,1.};
00151       meshS=ParaMEDMEM::MEDCouplingUMesh::New();
00152       meshS->setMeshDimension(2);
00153       ParaMEDMEM::DataArrayDouble *myCoords=ParaMEDMEM::DataArrayDouble::New();
00154       myCoords->alloc(4,2);
00155       std::copy(coordsS,coordsS+8,myCoords->getPointer());
00156       meshS->setCoords(myCoords);
00157       myCoords->decrRef();
00158       int connS[4]={0,2,3,1};
00159       meshS->allocateCells(1);
00160       meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS);
00161       meshS->finishInsertingCells();
00162       ParaMEDMEM::ComponentTopology comptopo;
00163       parameshS=new ParaMEDMEM::ParaMESH(meshS,*dec.getGrp(),"source mesh");
00164       parafieldS=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshS,comptopo);
00165       parafieldS->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
00166       double *valsS=parafieldS->getField()->getArray()->getPointer();
00167       valsS[0]=10.;
00168       //
00169       meshT=ParaMEDMEM::MEDCouplingUMesh::New();
00170       meshT->setMeshDimension(2);
00171       myCoords=ParaMEDMEM::DataArrayDouble::New();
00172       myCoords->alloc(3,2);
00173       std::copy(coordsT,coordsT+6,myCoords->getPointer());
00174       meshT->setCoords(myCoords);
00175       myCoords->decrRef();
00176       int connT[3]={0,1,2};
00177       meshT->allocateCells(1);
00178       meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
00179       meshT->finishInsertingCells();
00180       parameshT=new ParaMEDMEM::ParaMESH(meshT,*dec.getGrp(),"target mesh");
00181       parafieldT=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshT,comptopo);
00182       parafieldT->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
00183       double *valsT=parafieldT->getField()->getArray()->getPointer();
00184       valsT[0]=9.;
00185     }
00186   dec.attachSourceLocalField(parafieldS);
00187   dec.attachTargetLocalField(parafieldT);
00188   dec.synchronize();
00189   dec.sendRecvData(true);
00190   //
00191   if(rank==0)
00192     {
00193       CPPUNIT_ASSERT_DOUBLES_EQUAL(8.75,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
00194     }
00195   if(rank==1)
00196     {
00197       CPPUNIT_ASSERT_DOUBLES_EQUAL(8.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
00198     }
00199   if(rank==2)
00200     {
00201       CPPUNIT_ASSERT_DOUBLES_EQUAL(10.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
00202     }
00203   delete parafieldS;
00204   delete parafieldT;
00205   delete parameshS;
00206   delete parameshT;
00207   meshS->decrRef();
00208   meshT->decrRef();
00209 
00210   MPI_Barrier(MPI_COMM_WORLD);
00211 }
00212