Back to index

salome-med  6.5.0
ParaMEDMEMTestMPI2_2.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 <cppunit/extensions/HelperMacros.h>
00021 
00022 #include "MPI2Connector.hxx"
00023 #include "ParaMESH.hxx"
00024 #include "ParaFIELD.hxx"
00025 #include "MEDCouplingUMesh.hxx"
00026 #include "MEDCouplingFieldDouble.hxx"
00027 #include "InterpKernelDEC.hxx"
00028 #include "MPIProcessorGroup.hxx"
00029 #include "CommInterface.hxx"
00030 
00031 #include <mpi.h>
00032 #include <iostream>
00033 #include <stdlib.h>
00034 
00035 class MPI2ParaMEDMEMTest : public CppUnit::TestFixture
00036 {
00037   CPPUNIT_TEST_SUITE( MPI2ParaMEDMEMTest );
00038   CPPUNIT_TEST( testBasicMPI2_1 );
00039   CPPUNIT_TEST_SUITE_END();
00040 public:
00041   void testBasicMPI2_1();
00042 };
00043 
00044 using namespace ParaMEDMEM;
00045 
00046 void MPI2ParaMEDMEMTest::testBasicMPI2_1()
00047 {
00048   int lsize, lrank, gsize, grank;
00049   MPI_Comm gcom;
00050   std::string service = "SERVICE";
00051   std::ostringstream meshfilename, meshname;
00052   ParaMEDMEM::ParaMESH *paramesh=0;
00053   ParaMEDMEM::MEDCouplingUMesh* mesh;
00054   ParaMEDMEM::ParaFIELD *parafield=0;
00055   ParaMEDMEM::CommInterface* interface;
00056   ParaMEDMEM::MPIProcessorGroup* source, *target;
00057   
00058   MPI_Comm_size( MPI_COMM_WORLD, &lsize );
00059   MPI_Comm_rank( MPI_COMM_WORLD, &lrank );
00060   if(lsize!=3)
00061     {
00062       CPPUNIT_ASSERT(false);
00063       return;
00064     }
00065 
00066   /* Connection to remote programm */
00067   MPI2Connector *mpio = new MPI2Connector;
00068   gcom = mpio->remoteMPI2Connect(service);
00069   
00070   MPI_Comm_size( gcom, &gsize );
00071   MPI_Comm_rank( gcom, &grank );
00072   if(gsize!=5)
00073     {
00074       CPPUNIT_ASSERT(false);
00075       return;
00076     }
00077 
00078   interface = new ParaMEDMEM::CommInterface;
00079   source = new ParaMEDMEM::MPIProcessorGroup(*interface,0,gsize-lsize-1,gcom);
00080   target = new ParaMEDMEM::MPIProcessorGroup(*interface,gsize-lsize,gsize-1,gcom);
00081 
00082   const double targetCoordsAll[3][16]={{0.7,1.45,0.7,1.65,0.9,1.65,0.9,1.45,  1.1,1.4,1.1,1.6,1.3,1.6,1.3,1.4},
00083                                        {0.7,-0.6,0.7,0.7,0.9,0.7,0.9,-0.6,  1.1,-0.7,1.1,0.6,1.3,0.6,1.3,-0.7},
00084                                        {0.7,-1.55,0.7,-1.35,0.9,-1.35,0.9,-1.55,  1.1,-1.65,1.1,-1.45,1.3,-1.45,1.3,-1.65}};
00085   int conn4All[8]={0,1,2,3,4,5,6,7};
00086   double targetResults[3][2]={{34.,34.},{38.333333333333336,42.666666666666664},{47.,47.}};
00087 
00088   std::ostringstream stream; stream << "targetmesh2D proc " << grank-(gsize-lsize);
00089   mesh=MEDCouplingUMesh::New(stream.str().c_str(),2);
00090   mesh->allocateCells(2);
00091   mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All);
00092   mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,conn4All+4);
00093   mesh->finishInsertingCells();
00094   DataArrayDouble *myCoords=DataArrayDouble::New();
00095   myCoords->alloc(8,2);
00096   const double *targetCoords=targetCoordsAll[grank-(gsize-lsize)];
00097   std::copy(targetCoords,targetCoords+16,myCoords->getPointer());
00098   mesh->setCoords(myCoords);
00099   myCoords->decrRef();
00100   paramesh=new ParaMESH (mesh,*target,"target mesh");
00101   ParaMEDMEM::ComponentTopology comptopo;
00102   parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00103 
00104   ParaMEDMEM::InterpKernelDEC dec(*source,*target);
00105   parafield->getField()->setNature(ConservativeVolumic);
00106 
00107   dec.setMethod("P0");
00108   dec.attachLocalField(parafield);
00109   dec.synchronize();
00110   dec.setForcedRenormalization(false);
00111   dec.recvData();
00112   const double *res=parafield->getField()->getArray()->getConstPointer();
00113   const double *expected=targetResults[grank-(gsize-lsize)];
00114   CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[0],res[0],1e-13);
00115   CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[1],res[1],1e-13);
00116   /* Deconnection of remote programm */
00117   mpio->remoteMPI2Disconnect(service);
00118   /* clean-up */
00119   delete mpio;
00120   delete parafield;
00121   mesh->decrRef();
00122   delete paramesh;
00123   delete source;
00124   delete target;
00125   delete interface;
00126 }
00127 
00128 CPPUNIT_TEST_SUITE_REGISTRATION( MPI2ParaMEDMEMTest );
00129 
00130 #include "MPIMainTest.hxx"