Back to index

salome-med  6.5.0
ParaMEDMEMTest_FabienAPI.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 "CommInterface.hxx"
00022 #include "ProcessorGroup.hxx"
00023 #include "MPIProcessorGroup.hxx"
00024 #include "InterpKernelDEC.hxx"
00025 #include "MEDCouplingUMesh.hxx"
00026 #include "ParaMESH.hxx"
00027 #include "ParaFIELD.hxx"
00028 #include "ComponentTopology.hxx"
00029 
00030 #include <set>
00031 
00032 using namespace ParaMEDMEM;
00033 
00034 void ParaMEDMEMTest::testFabienAPI1()
00035 {
00036   int size;
00037   int rank;
00038   MPI_Comm_size(MPI_COMM_WORLD,&size);
00039   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00040   //
00041   if(size!=3)
00042     return ;
00043   int procs_source_c[1]={0};
00044   std::set<int> procs_source(procs_source_c,procs_source_c+1);
00045   int procs_target_c[1]={1};
00046   std::set<int> procs_target(procs_target_c,procs_target_c+1);
00047   //
00048   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
00049   ParaMEDMEM::ParaMESH *paramesh=0;
00050   ParaMEDMEM::ParaFIELD *parafield=0;
00051   //
00052   ParaMEDMEM::CommInterface interface;
00053   //
00054   MPI_Barrier(MPI_COMM_WORLD);
00055   double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
00056   CommInterface comm;
00057   //
00058   ParaMEDMEM::InterpKernelDEC *dec=new ParaMEDMEM::InterpKernelDEC(procs_source,procs_target);
00059   if(dec->isInSourceSide())
00060     {    
00061       mesh=MEDCouplingUMesh::New();
00062       mesh->setMeshDimension(2);
00063       DataArrayDouble *myCoords=DataArrayDouble::New();
00064       myCoords->alloc(4,2);
00065       std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
00066       mesh->setCoords(myCoords);
00067       myCoords->decrRef();
00068       int targetConn[4]={0,2,3,1};
00069       mesh->allocateCells(1);
00070       mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
00071       mesh->finishInsertingCells();
00072       ParaMEDMEM::ComponentTopology comptopo;
00073       paramesh=new ParaMESH(mesh,*dec->getSourceGrp(),"source mesh");
00074       parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00075       parafield->getField()->setNature(ConservativeVolumic);
00076       double *vals=parafield->getField()->getArray()->getPointer();
00077       vals[0]=7.;
00078     }
00079   if(dec->isInTargetSide())
00080     {
00081       mesh=MEDCouplingUMesh::New();
00082       mesh->setMeshDimension(2);
00083       DataArrayDouble *myCoords=DataArrayDouble::New();
00084       myCoords->alloc(4,2);
00085       std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
00086       mesh->setCoords(myCoords);
00087       myCoords->decrRef();
00088       int targetConn[6]={0,2,1,2,3,1};
00089       mesh->allocateCells(2);
00090       mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
00091       mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
00092       mesh->finishInsertingCells();
00093       ParaMEDMEM::ComponentTopology comptopo;
00094       paramesh=new ParaMESH(mesh,*dec->getTargetGrp(),"target mesh");
00095       parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00096       parafield->getField()->setNature(ConservativeVolumic);
00097     }
00098   dec->attachLocalField(parafield);
00099   dec->synchronize();
00100   dec->sendRecvData();
00101   if(dec->isInTargetSide())
00102     {
00103       const double *valsToTest=parafield->getField()->getArray()->getConstPointer();
00104       CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[0],7.,1e-14);
00105       CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[1],7.,1e-14);
00106     }
00107   //
00108   delete parafield;
00109   delete paramesh;
00110   if(mesh)
00111     mesh->decrRef();
00112   delete dec;
00113   MPI_Barrier(MPI_COMM_WORLD);
00114 }
00115 
00119 void ParaMEDMEMTest::testFabienAPI2()
00120 {
00121   int size;
00122   int rank;
00123   MPI_Comm_size(MPI_COMM_WORLD,&size);
00124   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00125   //
00126   if(size!=3)
00127     return ;
00128   int procs_source_c[1]={2};//difference with testFabienAPI1
00129   std::set<int> procs_source(procs_source_c,procs_source_c+1);
00130   int procs_target_c[1]={1};
00131   std::set<int> procs_target(procs_target_c,procs_target_c+1);
00132   //
00133   ParaMEDMEM::MEDCouplingUMesh *mesh=0;
00134   ParaMEDMEM::ParaMESH *paramesh=0;
00135   ParaMEDMEM::ParaFIELD *parafield=0;
00136   //
00137   ParaMEDMEM::CommInterface interface;
00138   //
00139   MPI_Barrier(MPI_COMM_WORLD);
00140   double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
00141   CommInterface comm;
00142   //
00143   ParaMEDMEM::InterpKernelDEC *dec=new ParaMEDMEM::InterpKernelDEC(procs_source,procs_target);
00144   if(dec->isInSourceSide())
00145     {    
00146       mesh=MEDCouplingUMesh::New();
00147       mesh->setMeshDimension(2);
00148       DataArrayDouble *myCoords=DataArrayDouble::New();
00149       myCoords->alloc(4,2);
00150       std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
00151       mesh->setCoords(myCoords);
00152       myCoords->decrRef();
00153       int targetConn[4]={0,2,3,1};
00154       mesh->allocateCells(1);
00155       mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
00156       mesh->finishInsertingCells();
00157       ParaMEDMEM::ComponentTopology comptopo;
00158       paramesh=new ParaMESH(mesh,*dec->getSourceGrp(),"source mesh");
00159       parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00160       parafield->getField()->setNature(ConservativeVolumic);
00161       double *vals=parafield->getField()->getArray()->getPointer();
00162       vals[0]=7.;
00163     }
00164   if(dec->isInTargetSide())
00165     {
00166       mesh=MEDCouplingUMesh::New();
00167       mesh->setMeshDimension(2);
00168       DataArrayDouble *myCoords=DataArrayDouble::New();
00169       myCoords->alloc(4,2);
00170       std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
00171       mesh->setCoords(myCoords);
00172       myCoords->decrRef();
00173       int targetConn[6]={0,2,1,2,3,1};
00174       mesh->allocateCells(2);
00175       mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
00176       mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
00177       mesh->finishInsertingCells();
00178       ParaMEDMEM::ComponentTopology comptopo;
00179       paramesh=new ParaMESH(mesh,*dec->getTargetGrp(),"target mesh");
00180       parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00181       parafield->getField()->setNature(ConservativeVolumic);
00182     }
00183   dec->attachLocalField(parafield);
00184   dec->synchronize();
00185   dec->sendRecvData();
00186   if(dec->isInTargetSide())
00187     {
00188       const double *valsToTest=parafield->getField()->getArray()->getConstPointer();
00189       CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[0],7.,1e-14);
00190       CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[1],7.,1e-14);
00191     }
00192   //
00193   delete parafield;
00194   delete paramesh;
00195   if(mesh)
00196     mesh->decrRef();
00197   delete dec;
00198   MPI_Barrier(MPI_COMM_WORLD);
00199 }