1 // Copyright (C) 2007-2016 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "ParaMEDMEMTest.hxx"
21 #include "CommInterface.hxx"
22 #include "ProcessorGroup.hxx"
23 #include "MPIProcessorGroup.hxx"
24 #include "InterpKernelDEC.hxx"
25 #include "MEDCouplingUMesh.hxx"
26 #include "ParaMESH.hxx"
27 #include "ParaFIELD.hxx"
28 #include "ComponentTopology.hxx"
32 using namespace MEDCoupling;
34 void ParaMEDMEMTest::testFabienAPI1()
38 MPI_Comm_size(MPI_COMM_WORLD,&size);
39 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
43 int procs_source_c[1]={0};
44 std::set<int> procs_source(procs_source_c,procs_source_c+1);
45 int procs_target_c[1]={1};
46 std::set<int> procs_target(procs_target_c,procs_target_c+1);
48 MEDCoupling::MEDCouplingUMesh *mesh=0;
49 MEDCoupling::ParaMESH *paramesh=0;
50 MEDCoupling::ParaFIELD *parafield=0;
52 MEDCoupling::CommInterface interface;
54 MPI_Barrier(MPI_COMM_WORLD);
55 double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
58 MEDCoupling::InterpKernelDEC *dec=new MEDCoupling::InterpKernelDEC(procs_source,procs_target);
59 if(dec->isInSourceSide())
61 mesh=MEDCouplingUMesh::New();
62 mesh->setMeshDimension(2);
63 DataArrayDouble *myCoords=DataArrayDouble::New();
65 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
66 mesh->setCoords(myCoords);
68 int targetConn[4]={0,2,3,1};
69 mesh->allocateCells(1);
70 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
71 mesh->finishInsertingCells();
72 MEDCoupling::ComponentTopology comptopo;
73 paramesh=new ParaMESH(mesh,*dec->getSourceGrp(),"source mesh");
74 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
75 parafield->getField()->setNature(IntensiveMaximum);
76 double *vals=parafield->getField()->getArray()->getPointer();
79 if(dec->isInTargetSide())
81 mesh=MEDCouplingUMesh::New();
82 mesh->setMeshDimension(2);
83 DataArrayDouble *myCoords=DataArrayDouble::New();
85 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
86 mesh->setCoords(myCoords);
88 int targetConn[6]={0,2,1,2,3,1};
89 mesh->allocateCells(2);
90 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
91 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
92 mesh->finishInsertingCells();
93 MEDCoupling::ComponentTopology comptopo;
94 paramesh=new ParaMESH(mesh,*dec->getTargetGrp(),"target mesh");
95 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
96 parafield->getField()->setNature(IntensiveMaximum);
98 dec->attachLocalField(parafield);
101 if(dec->isInTargetSide())
103 const double *valsToTest=parafield->getField()->getArray()->getConstPointer();
104 CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[0],7.,1e-14);
105 CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[1],7.,1e-14);
113 MPI_Barrier(MPI_COMM_WORLD);
117 * Idem testFabienAPI1 except that procs are shuffled. Test of the good management of group translation in newly created communicator.
119 void ParaMEDMEMTest::testFabienAPI2()
123 MPI_Comm_size(MPI_COMM_WORLD,&size);
124 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
128 int procs_source_c[1]={2};//difference with testFabienAPI1
129 std::set<int> procs_source(procs_source_c,procs_source_c+1);
130 int procs_target_c[1]={1};
131 std::set<int> procs_target(procs_target_c,procs_target_c+1);
133 MEDCoupling::MEDCouplingUMesh *mesh=0;
134 MEDCoupling::ParaMESH *paramesh=0;
135 MEDCoupling::ParaFIELD *parafield=0;
137 MEDCoupling::CommInterface interface;
139 MPI_Barrier(MPI_COMM_WORLD);
140 double targetCoords[8]={ 0.,0., 1., 0., 0., 1., 1., 1. };
143 MEDCoupling::InterpKernelDEC *dec=new MEDCoupling::InterpKernelDEC(procs_source,procs_target);
144 if(dec->isInSourceSide())
146 mesh=MEDCouplingUMesh::New();
147 mesh->setMeshDimension(2);
148 DataArrayDouble *myCoords=DataArrayDouble::New();
149 myCoords->alloc(4,2);
150 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
151 mesh->setCoords(myCoords);
153 int targetConn[4]={0,2,3,1};
154 mesh->allocateCells(1);
155 mesh->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,targetConn);
156 mesh->finishInsertingCells();
157 MEDCoupling::ComponentTopology comptopo;
158 paramesh=new ParaMESH(mesh,*dec->getSourceGrp(),"source mesh");
159 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
160 parafield->getField()->setNature(IntensiveMaximum);
161 double *vals=parafield->getField()->getArray()->getPointer();
164 if(dec->isInTargetSide())
166 mesh=MEDCouplingUMesh::New();
167 mesh->setMeshDimension(2);
168 DataArrayDouble *myCoords=DataArrayDouble::New();
169 myCoords->alloc(4,2);
170 std::copy(targetCoords,targetCoords+8,myCoords->getPointer());
171 mesh->setCoords(myCoords);
173 int targetConn[6]={0,2,1,2,3,1};
174 mesh->allocateCells(2);
175 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn);
176 mesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,targetConn+3);
177 mesh->finishInsertingCells();
178 MEDCoupling::ComponentTopology comptopo;
179 paramesh=new ParaMESH(mesh,*dec->getTargetGrp(),"target mesh");
180 parafield=new ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo);
181 parafield->getField()->setNature(IntensiveMaximum);
183 dec->attachLocalField(parafield);
186 if(dec->isInTargetSide())
188 const double *valsToTest=parafield->getField()->getArray()->getConstPointer();
189 CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[0],7.,1e-14);
190 CPPUNIT_ASSERT_DOUBLES_EQUAL(valsToTest[1],7.,1e-14);
198 MPI_Barrier(MPI_COMM_WORLD);