Back to index

salome-med  6.5.0
ParaMEDMEMTest_NonCoincidentDEC.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 #ifdef MED_ENABLE_FVM
00021 
00022 #include "ParaMEDMEMTest.hxx"
00023 #include <cppunit/TestAssert.h>
00024 
00025 #include "MEDMEM_Exception.hxx"
00026 #include "CommInterface.hxx"
00027 #include "ProcessorGroup.hxx"
00028 #include "MPIProcessorGroup.hxx"
00029 #include "Topology.hxx"
00030 #include "DEC.hxx"
00031 #include "NonCoincidentDEC.hxx"
00032 #include "ParaMESH.hxx"
00033 #include "ParaFIELD.hxx"
00034 #include "UnstructuredParaSUPPORT.hxx"
00035 #include "ICoCoMEDField.hxx"
00036  
00037 #include <string>
00038 
00039 // use this define to enable lines, execution of which leads to Segmentation Fault
00040 #define ENABLE_FAULTS
00041 
00042 // use this define to enable CPPUNIT asserts and fails, showing bugs
00043 #define ENABLE_FORCED_FAILURES
00044 
00045 
00046 using namespace std;
00047 using namespace ParaMEDMEM;
00048 using namespace MEDMEM;
00049  
00050 /*
00051  * Check methods defined in InterpKernelDEC.hxx
00052  *
00053  InterpKernelDEC();
00054  InterpKernelDEC(ProcessorGroup& local_group, ProcessorGroup& distant_group);
00055  virtual ~InterpKernelDEC();
00056  void synchronize();
00057  void recvData();
00058  void sendData();
00059 */
00060 
00061 void ParaMEDMEMTest::testNonCoincidentDEC_2D()
00062 {
00063 
00064   int size;
00065   MPI_Comm_size(MPI_COMM_WORLD,&size);
00066 
00067   //the test is meant to run on five processors
00068   if (size !=5) return ;
00069   
00070   testNonCoincidentDEC( "/share/salome/resources/med/square1_split",
00071                         "Mesh_2",
00072                         "/share/salome/resources/med/square2_split",
00073                         "Mesh_3",
00074                         3,
00075                         1e-6);
00076 } 
00077 
00078 void ParaMEDMEMTest::testNonCoincidentDEC_3D()
00079 {
00080   int size;
00081   MPI_Comm_size(MPI_COMM_WORLD,&size);
00082   
00083   //the test is meant to run on five processors
00084   if (size !=4) return ;
00085   
00086   testNonCoincidentDEC( "/share/salome/resources/med/blade_12000_split2",
00087                         "Mesh_1",
00088                         "/share/salome/resources/med/blade_3000_split2",
00089                         "Mesh_1",
00090                         2,
00091                         1e4);
00092 } 
00093 
00094 void ParaMEDMEMTest::testNonCoincidentDEC(const string& filename1,
00095                                           const string& meshname1,
00096                                           const string& filename2,
00097                                           const string& meshname2,
00098                                           int nproc_source,
00099                                           double epsilon)
00100 {
00101   int size;
00102   int rank;
00103   MPI_Comm_size(MPI_COMM_WORLD,&size);
00104   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00105    
00106   set<int> self_procs;
00107   set<int> procs_source;
00108   set<int> procs_target;
00109   
00110   for (int i=0; i<nproc_source; i++)
00111     procs_source.insert(i);
00112   for (int i=nproc_source; i<size; i++)
00113     procs_target.insert(i);
00114   self_procs.insert(rank);
00115   
00116   ParaMEDMEM::CommInterface interface;
00117     
00118   ParaMEDMEM::ProcessorGroup* self_group = new ParaMEDMEM::MPIProcessorGroup(interface,self_procs);
00119   ParaMEDMEM::ProcessorGroup* target_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_target);
00120   ParaMEDMEM::ProcessorGroup* source_group = new ParaMEDMEM::MPIProcessorGroup(interface,procs_source);
00121   
00122   ParaMEDMEM::ParaMESH* source_mesh=0;
00123   ParaMEDMEM::ParaMESH* target_mesh=0;
00124   ParaMEDMEM::ParaSUPPORT* parasupport=0;
00125   //loading the geometry for the source group
00126 
00127   ParaMEDMEM::NonCoincidentDEC dec (*source_group,*target_group);
00128 
00129   MEDMEM::MESH* mesh;
00130   MEDMEM::SUPPORT* support;
00131   MEDMEM::FIELD<double>* field;
00132   ParaMEDMEM::ParaMESH* paramesh;
00133   ParaMEDMEM::ParaFIELD* parafield;
00134   
00135   string filename_xml1              = getResourceFile(filename1);
00136   string filename_xml2              = getResourceFile(filename2); 
00137   //string filename_seq_wr            = makeTmpFile("");
00138   //string filename_seq_med           = makeTmpFile("myWrField_seq_pointe221.med");
00139   
00140   // To remove tmp files from disk
00141   ParaMEDMEMTest_TmpFilesRemover aRemover;
00142   //aRemover.Register(filename_seq_wr);
00143   //aRemover.Register(filename_seq_med);
00144   MPI_Barrier(MPI_COMM_WORLD);
00145   ICoCo::Field* icocofield;
00146   if (source_group->containsMyRank())
00147     {
00148       string master = filename_xml1;
00149       
00150       ostringstream strstream;
00151       strstream <<master<<rank+1<<".med";
00152       ostringstream meshname ;
00153       meshname<< meshname1<<"_"<< rank+1;
00154       
00155       CPPUNIT_ASSERT_NO_THROW(mesh = new MESH(MED_DRIVER,strstream.str(),meshname.str()));
00156       support=new MEDMEM::SUPPORT(mesh,"all elements",MED_EN::MED_CELL);
00157     
00158       paramesh=new ParaMESH (*mesh,*source_group,"source mesh");
00159     
00160       parasupport=new UnstructuredParaSUPPORT( support,*source_group);
00161       ParaMEDMEM::ComponentTopology comptopo;
00162       parafield = new ParaFIELD(parasupport, comptopo);
00163 
00164       
00165       int nb_local=support->getNumberOfElements(MED_EN::MED_ALL_ELEMENTS);
00166       double * value= new double[nb_local];
00167       for(int ielem=0; ielem<nb_local;ielem++)
00168         value[ielem]=1.0;
00169       parafield->getField()->setValue(value);
00170 
00171       icocofield=new ICoCo::MEDField(paramesh,parafield);
00172      
00173       dec.attachLocalField(icocofield);
00174       delete [] value;
00175     }
00176   
00177   //loading the geometry for the target group
00178   if (target_group->containsMyRank())
00179     {
00180       string master= filename_xml2;
00181       ostringstream strstream;
00182       strstream << master<<(rank-nproc_source+1)<<".med";
00183       ostringstream meshname ;
00184       meshname<< meshname2<<"_"<<rank-nproc_source+1;
00185       
00186       CPPUNIT_ASSERT_NO_THROW(mesh = new MESH(MED_DRIVER,strstream.str(),meshname.str()));
00187       support=new MEDMEM::SUPPORT(mesh,"all elements",MED_EN::MED_CELL);
00188       
00189       paramesh=new ParaMESH (*mesh,*target_group,"target mesh");
00190       parasupport=new UnstructuredParaSUPPORT(support,*target_group);
00191       ParaMEDMEM::ComponentTopology comptopo;
00192       parafield = new ParaFIELD(parasupport, comptopo);
00193 
00194       
00195       int nb_local=support->getNumberOfElements(MED_EN::MED_ALL_ELEMENTS);
00196       double * value= new double[nb_local];
00197       for(int ielem=0; ielem<nb_local;ielem++)
00198         value[ielem]=0.0;
00199       parafield->getField()->setValue(value);
00200       icocofield=new ICoCo::MEDField(paramesh,parafield);
00201       
00202       dec.attachLocalField(icocofield);
00203       delete [] value;
00204     }
00205     
00206   
00207   //attaching a DEC to the source group 
00208   double field_before_int;
00209   double field_after_int;
00210   
00211   if (source_group->containsMyRank())
00212     { 
00213       field_before_int = parafield->getVolumeIntegral(1);
00214       MPI_Bcast(&field_before_int, 1,MPI_DOUBLE, 0,MPI_COMM_WORLD);
00215       dec.synchronize();
00216       cout<<"DEC usage"<<endl;
00217       dec.setOption("ForcedRenormalization",false);
00218 
00219       dec.sendData();
00220       //      paramesh->write(MED_DRIVER,"./sourcesquarenc");
00221       //parafield->write(MED_DRIVER,"./sourcesquarenc","boundary");  
00222    
00223       
00224     }
00225   
00226   //attaching a DEC to the target group
00227   if (target_group->containsMyRank())
00228     {
00229       MPI_Bcast(&field_before_int, 1,MPI_DOUBLE, 0,MPI_COMM_WORLD);
00230      
00231       dec.synchronize();
00232       dec.setOption("ForcedRenormalization",false);
00233       dec.recvData();
00234       //paramesh->write(MED_DRIVER, "./targetsquarenc");
00235       //parafield->write(MED_DRIVER, "./targetsquarenc", "boundary");
00236       field_after_int = parafield->getVolumeIntegral(1);
00237       
00238     }
00239   MPI_Bcast(&field_before_int,1,MPI_DOUBLE,0,MPI_COMM_WORLD);
00240   MPI_Bcast(&field_after_int, 1,MPI_DOUBLE, size-1,MPI_COMM_WORLD);
00241      
00242   CPPUNIT_ASSERT_DOUBLES_EQUAL(field_before_int, field_after_int, epsilon);
00243     
00244   delete source_group;
00245   delete target_group;
00246   delete self_group;
00247   delete icocofield;
00248   delete paramesh;
00249   delete parafield;
00250   delete support;
00251   delete parasupport;
00252   delete mesh;
00253   MPI_Barrier(MPI_COMM_WORLD);
00254   
00255 }
00256 #endif