Back to index

salome-med  6.5.0
ParaMEDMEMTest_StructuredCoincidentDEC.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 "DEC.hxx"
00028 #include "StructuredCoincidentDEC.hxx"
00029 #include "ParaMESH.hxx"
00030 #include "ParaFIELD.hxx"
00031 #include "ComponentTopology.hxx"
00032 #include "ICoCoMEDField.hxx"
00033 #include "MEDLoader.hxx"
00034 
00035 #include <string>
00036 
00037 // use this define to enable lines, execution of which leads to Segmentation Fault
00038 #define ENABLE_FAULTS
00039 
00040 // use this define to enable CPPUNIT asserts and fails, showing bugs
00041 #define ENABLE_FORCED_FAILURES
00042 
00043 using namespace std;
00044 using namespace ParaMEDMEM;
00045 
00046 /*
00047  * Check methods defined in StructuredCoincidentDEC.hxx
00048  *
00049  StructuredCoincidentDEC();
00050  StructuredCoincidentDEC(ProcessorGroup& local_group, ProcessorGroup& distant_group);
00051  virtual ~StructuredCoincidentDEC();
00052  void synchronize();
00053  void recvData();
00054  void sendData();
00055 */
00056 
00057 void ParaMEDMEMTest::testStructuredCoincidentDEC() {
00058   string testname="ParaMEDMEM - testStructured CoincidentDEC";
00059   //  MPI_Init(&argc, &argv); 
00060   int size;
00061   int rank;
00062   MPI_Comm_size(MPI_COMM_WORLD, &size);
00063   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00064   if (size<4) {
00065     return;
00066   }
00067 
00068   ParaMEDMEM::CommInterface interface;
00069 
00070   ParaMEDMEM::MPIProcessorGroup self_group (interface,rank,rank);
00071   ParaMEDMEM::MPIProcessorGroup target_group(interface,3,size-1);
00072   ParaMEDMEM::MPIProcessorGroup source_group (interface,0,2);
00073 
00074   ParaMEDMEM::MEDCouplingUMesh* mesh;
00075   ParaMEDMEM::ParaMESH* paramesh;
00076   ParaMEDMEM::ParaFIELD* parafield;
00077 
00078   string filename_xml1 = getResourceFile("square1_split");
00079   string filename_2    = getResourceFile("square1.med");
00080   //string filename_seq_wr  = makeTmpFile("");
00081   //string filename_seq_med = makeTmpFile("myWrField_seq_pointe221.med");
00082 
00083   // To remove tmp files from disk
00084   ParaMEDMEMTest_TmpFilesRemover aRemover;
00085 
00086   //loading the geometry for the source group
00087 
00088   ParaMEDMEM::StructuredCoincidentDEC dec(source_group, target_group);
00089 
00090   MPI_Barrier(MPI_COMM_WORLD);
00091   if (source_group.containsMyRank()) {
00092     string master = filename_xml1;
00093 
00094     ostringstream strstream;
00095     strstream <<master<<rank+1<<".med";
00096     ostringstream meshname;
00097     meshname<< "Mesh_2_"<< rank+1;
00098 
00099     mesh=MEDLoader::ReadUMeshFromFile(strstream.str().c_str(),meshname.str().c_str(),0);
00100     
00101 
00102     paramesh=new ParaMESH (mesh,source_group,"source mesh");
00103 
00104     ParaMEDMEM::ComponentTopology comptopo(6);
00105     parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00106 
00107     int nb_local=mesh->getNumberOfCells();
00108     const int* global_numbering = paramesh->getGlobalNumberingCell();
00109     
00110     double *value=parafield->getField()->getArray()->getPointer();
00111     for(int ielem=0; ielem<nb_local;ielem++)
00112       for (int icomp=0; icomp<6; icomp++)
00113         value[ielem*6+icomp]=global_numbering[ielem]*6+icomp;
00114 
00115     //ICoCo::Field* icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00116 
00117     dec.attachLocalField(parafield);
00118     dec.synchronize();
00119     dec.sendData();
00120     //delete icocofield;
00121   }
00122 
00123   //loading the geometry for the target group
00124   if (target_group.containsMyRank()) {
00125 
00126     string meshname2("Mesh_2");
00127     mesh = MEDLoader::ReadUMeshFromFile(filename_2.c_str(),meshname2.c_str(),0);
00128     
00129     paramesh=new ParaMESH (mesh,self_group,"target mesh");
00130     ParaMEDMEM::ComponentTopology comptopo(6, &target_group);
00131 
00132     parafield = new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
00133 
00134     int nb_local=mesh->getNumberOfCells();
00135     double *value=parafield->getField()->getArray()->getPointer();
00136     for (int ielem=0; ielem<nb_local; ielem++)
00137       for (int icomp=0; icomp<comptopo.nbLocalComponents(); icomp++)
00138         value[ielem*comptopo.nbLocalComponents()+icomp]=0.0;
00139     //ICoCo::Field* icocofield=new ICoCo::MEDField((MEDCouplingUMesh *)paramesh->getCellMesh(),parafield->getField());
00140 
00141     dec.attachLocalField(parafield);
00142     dec.synchronize();
00143     dec.recvData();
00144 
00145     //checking validity of field
00146     const double* recv_value = parafield->getField()->getArray()->getPointer();
00147     for (int i=0; i< nb_local; i++) {
00148       int first = comptopo.firstLocalComponent();
00149       for (int icomp = 0; icomp < comptopo.nbLocalComponents(); icomp++)
00150         CPPUNIT_ASSERT_DOUBLES_EQUAL(recv_value[i*comptopo.nbLocalComponents()+icomp],(double)(i*6+icomp+first),1e-12);
00151     }
00152     //delete icocofield;
00153   }
00154   delete parafield;
00155   delete paramesh;
00156   mesh->decrRef();
00157 
00158   //  MPI_Barrier(MPI_COMM_WORLD);
00159 
00160 }