1 // Copyright (C) 2007-2014 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 <cppunit/TestAssert.h>
23 #include "CommInterface.hxx"
24 #include "ProcessorGroup.hxx"
25 #include "MPIProcessorGroup.hxx"
26 #include "Topology.hxx"
27 #include "OverlapDEC.hxx"
28 #include "ParaMESH.hxx"
29 #include "ParaFIELD.hxx"
30 #include "ComponentTopology.hxx"
32 #include "MEDCouplingUMesh.hxx"
36 void ParaMEDMEMTest::testOverlapDEC1()
38 std::string srcM("P0");
39 std::string targetM("P0");
42 MPI_Comm_size(MPI_COMM_WORLD,&size);
43 MPI_Comm_rank(MPI_COMM_WORLD,&rank);
45 if (size != 3) return ;
50 for (int i=0; i<nproc; i++)
53 ParaMEDMEM::CommInterface interface;
55 ParaMEDMEM::OverlapDEC dec(procs);
57 ParaMEDMEM::MEDCouplingUMesh* meshS=0;
58 ParaMEDMEM::MEDCouplingUMesh* meshT=0;
59 ParaMEDMEM::ParaMESH* parameshS=0;
60 ParaMEDMEM::ParaMESH* parameshT=0;
61 ParaMEDMEM::ParaFIELD* parafieldS=0;
62 ParaMEDMEM::ParaFIELD* parafieldT=0;
64 MPI_Barrier(MPI_COMM_WORLD);
67 const double coordsS[10]={0.,0.,0.5,0.,1.,0.,0.,0.5,0.5,0.5};
68 const double coordsT[6]={0.,0.,1.,0.,1.,1.};
69 meshS=ParaMEDMEM::MEDCouplingUMesh::New();
70 meshS->setMeshDimension(2);
71 ParaMEDMEM::DataArrayDouble *myCoords=ParaMEDMEM::DataArrayDouble::New();
73 std::copy(coordsS,coordsS+10,myCoords->getPointer());
74 meshS->setCoords(myCoords);
76 int connS[7]={0,3,4,1, 1,4,2};
77 meshS->allocateCells(2);
78 meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS);
79 meshS->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connS+4);
80 meshS->finishInsertingCells();
81 ParaMEDMEM::ComponentTopology comptopo;
82 parameshS=new ParaMEDMEM::ParaMESH(meshS,*dec.getGrp(),"source mesh");
83 parafieldS=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshS,comptopo);
84 parafieldS->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
85 double *valsS=parafieldS->getField()->getArray()->getPointer();
86 valsS[0]=7.; valsS[1]=8.;
88 meshT=ParaMEDMEM::MEDCouplingUMesh::New();
89 meshT->setMeshDimension(2);
90 myCoords=ParaMEDMEM::DataArrayDouble::New();
92 std::copy(coordsT,coordsT+6,myCoords->getPointer());
93 meshT->setCoords(myCoords);
96 meshT->allocateCells(1);
97 meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
98 meshT->finishInsertingCells();
99 parameshT=new ParaMEDMEM::ParaMESH(meshT,*dec.getGrp(),"target mesh");
100 parafieldT=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshT,comptopo);
101 parafieldT->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
102 double *valsT=parafieldT->getField()->getArray()->getPointer();
108 const double coordsS[10]={1.,0.,0.5,0.5,1.,0.5,0.5,1.,1.,1.};
109 const double coordsT[6]={0.,0.,0.5,0.5,0.,1.};
110 meshS=ParaMEDMEM::MEDCouplingUMesh::New();
111 meshS->setMeshDimension(2);
112 ParaMEDMEM::DataArrayDouble *myCoords=ParaMEDMEM::DataArrayDouble::New();
113 myCoords->alloc(5,2);
114 std::copy(coordsS,coordsS+10,myCoords->getPointer());
115 meshS->setCoords(myCoords);
117 int connS[7]={0,1,2, 1,3,4,2};
118 meshS->allocateCells(2);
119 meshS->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connS);
120 meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS+3);
121 meshS->finishInsertingCells();
122 ParaMEDMEM::ComponentTopology comptopo;
123 parameshS=new ParaMEDMEM::ParaMESH(meshS,*dec.getGrp(),"source mesh");
124 parafieldS=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshS,comptopo);
125 parafieldS->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
126 double *valsS=parafieldS->getField()->getArray()->getPointer();
127 valsS[0]=9.; valsS[1]=11.;
129 meshT=ParaMEDMEM::MEDCouplingUMesh::New();
130 meshT->setMeshDimension(2);
131 myCoords=ParaMEDMEM::DataArrayDouble::New();
132 myCoords->alloc(3,2);
133 std::copy(coordsT,coordsT+6,myCoords->getPointer());
134 meshT->setCoords(myCoords);
136 int connT[3]={0,2,1};
137 meshT->allocateCells(1);
138 meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
139 meshT->finishInsertingCells();
140 parameshT=new ParaMEDMEM::ParaMESH(meshT,*dec.getGrp(),"target mesh");
141 parafieldT=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshT,comptopo);
142 parafieldT->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
143 double *valsT=parafieldT->getField()->getArray()->getPointer();
149 const double coordsS[8]={0.,0.5, 0.5,0.5, 0.,1., 0.5,1.};
150 const double coordsT[6]={0.5,0.5,0.,1.,1.,1.};
151 meshS=ParaMEDMEM::MEDCouplingUMesh::New();
152 meshS->setMeshDimension(2);
153 ParaMEDMEM::DataArrayDouble *myCoords=ParaMEDMEM::DataArrayDouble::New();
154 myCoords->alloc(4,2);
155 std::copy(coordsS,coordsS+8,myCoords->getPointer());
156 meshS->setCoords(myCoords);
158 int connS[4]={0,2,3,1};
159 meshS->allocateCells(1);
160 meshS->insertNextCell(INTERP_KERNEL::NORM_QUAD4,4,connS);
161 meshS->finishInsertingCells();
162 ParaMEDMEM::ComponentTopology comptopo;
163 parameshS=new ParaMEDMEM::ParaMESH(meshS,*dec.getGrp(),"source mesh");
164 parafieldS=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshS,comptopo);
165 parafieldS->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
166 double *valsS=parafieldS->getField()->getArray()->getPointer();
169 meshT=ParaMEDMEM::MEDCouplingUMesh::New();
170 meshT->setMeshDimension(2);
171 myCoords=ParaMEDMEM::DataArrayDouble::New();
172 myCoords->alloc(3,2);
173 std::copy(coordsT,coordsT+6,myCoords->getPointer());
174 meshT->setCoords(myCoords);
176 int connT[3]={0,1,2};
177 meshT->allocateCells(1);
178 meshT->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,connT);
179 meshT->finishInsertingCells();
180 parameshT=new ParaMEDMEM::ParaMESH(meshT,*dec.getGrp(),"target mesh");
181 parafieldT=new ParaMEDMEM::ParaFIELD(ParaMEDMEM::ON_CELLS,ParaMEDMEM::NO_TIME,parameshT,comptopo);
182 parafieldT->getField()->setNature(ParaMEDMEM::ConservativeVolumic);//IntegralGlobConstraint
183 double *valsT=parafieldT->getField()->getArray()->getPointer();
186 dec.attachSourceLocalField(parafieldS);
187 dec.attachTargetLocalField(parafieldT);
189 dec.sendRecvData(true);
193 CPPUNIT_ASSERT_DOUBLES_EQUAL(8.75,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
197 CPPUNIT_ASSERT_DOUBLES_EQUAL(8.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
201 CPPUNIT_ASSERT_DOUBLES_EQUAL(10.5,parafieldT->getField()->getArray()->getIJ(0,0),1e-12);
210 MPI_Barrier(MPI_COMM_WORLD);