1 // Copyright (C) 2007-2015 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 "BlockTopology.hxx"
21 #include "MEDCouplingMemArray.hxx"
22 #include "MEDCouplingCMesh.hxx"
23 #include "CommInterface.hxx"
24 #include "ProcessorGroup.hxx"
25 #include "MPIProcessorGroup.hxx"
26 #include "ComponentTopology.hxx"
27 #include "InterpKernelUtilities.hxx"
41 BlockTopology::BlockTopology() :
42 _dimension(0), _nb_procs_per_dim(0),
43 _local_array_indices(0), _cycle_type(0),
44 _proc_group(NULL),_nb_elems(0),
45 _owns_processor_group(false)
49 * Constructor of a block topology from a grid.
50 * This preliminary version simply splits along the first axis
51 * instead of making the best choice with respect to the
52 * values of the different axes.
54 BlockTopology::BlockTopology(const ProcessorGroup& group, MEDCouplingCMesh *grid):
55 _dimension(grid->getSpaceDimension()), _proc_group(&group), _owns_processor_group(false)
57 vector <int> axis_length(_dimension);
59 for (int idim=0; idim <_dimension; idim++)
61 DataArrayDouble *arr=grid->getCoordsAt(idim);
62 axis_length[idim]=arr->getNbOfElems();
63 _nb_elems*=axis_length[idim];
65 //default splitting along 1st dimension
66 _local_array_indices.resize(_dimension);
67 _nb_procs_per_dim.resize(_dimension);
69 _local_array_indices[0].resize(_proc_group->size()+1);
70 _local_array_indices[0][0]=0;
71 _nb_procs_per_dim[0]=_proc_group->size();
73 for (int i=1; i<=_proc_group->size(); i++)
75 _local_array_indices[0][i]=_local_array_indices[0][i-1]+
76 axis_length[0]/_proc_group->size();
77 if (i<= axis_length[0]%_proc_group->size())
78 _local_array_indices[0][i]+=1;
80 for (int i=1; i<_dimension; i++)
82 _local_array_indices[i].resize(2);
83 _local_array_indices[i][0]=0;
84 _local_array_indices[i][1]=axis_length[i];
85 _nb_procs_per_dim[i]=1;
87 _cycle_type.resize(_dimension);
88 for (int i=0; i<_dimension; i++)
89 _cycle_type[i]=ParaMEDMEM::Block;
93 * Creation of a block topology by composing
94 * a geometrical topology and a component topology.
95 * This constructor is intended for creating fields
96 * for which the parallel distribution is made on the
97 * components of the field rather than on the geometrical
98 * partitioning of the underlying mesh.
101 BlockTopology::BlockTopology(const BlockTopology& geom_topo, const ComponentTopology& comp_topo):_owns_processor_group(false)
103 // so far, the block topology can only be created if the proc group
104 // is either on geom_topo or on comp_topo
105 if (geom_topo.getProcGroup()->size()>1 && comp_topo.nbBlocks()>1)
106 throw INTERP_KERNEL::Exception(LOCALIZED("BlockTopology cannot yet be constructed with both complex geo and components topology"));
108 if (comp_topo.nbComponents()==1)
115 _dimension = geom_topo.getDimension()+1;
116 if (comp_topo.nbBlocks()>1)
117 _proc_group=comp_topo.getProcGroup();
119 _proc_group=geom_topo.getProcGroup();
120 _local_array_indices=geom_topo._local_array_indices;
121 vector<int> comp_indices = *(comp_topo.getBlockIndices());
122 _local_array_indices.push_back(comp_indices);
123 _nb_procs_per_dim=geom_topo._nb_procs_per_dim;
124 _nb_procs_per_dim.push_back(comp_topo.nbBlocks());
125 _cycle_type=geom_topo._cycle_type;
126 _cycle_type.push_back(Block);
127 _nb_elems=geom_topo.getNbElements()*comp_topo.nbComponents();
131 /*! Constructor for creating a one-dimensional
132 * topology from a processor group and a local
133 * number of elements on each processor
135 * The function must be called only by the processors belonging
136 * to group \a group. Calling it from a processor not belonging
137 * to \a group will cause an MPI error, while calling from a subset
138 * of \a group will result in a deadlock.
140 BlockTopology::BlockTopology(const ProcessorGroup& group, int nb_elem):_dimension(1),_proc_group(&group),_owns_processor_group(false)
142 int* nbelems_per_proc = new int[group.size()];
143 const MPIProcessorGroup* mpi_group=dynamic_cast<const MPIProcessorGroup*>(_proc_group);
144 const MPI_Comm* comm=mpi_group->getComm();
146 mpi_group->getCommInterface().allGather(&nbtemp, 1, MPI_INT,
147 nbelems_per_proc, 1, MPI_INT,
151 //splitting along only dimension
152 _local_array_indices.resize(1);
153 _nb_procs_per_dim.resize(1);
155 _local_array_indices[0].resize(_proc_group->size()+1);
156 _local_array_indices[0][0]=0;
157 _nb_procs_per_dim[0]=_proc_group->size();
159 for (int i=1; i<=_proc_group->size(); i++)
161 _local_array_indices[0][i]=_local_array_indices[0][i-1]+
162 nbelems_per_proc[i-1];
163 _nb_elems+=nbelems_per_proc[i-1];
165 _cycle_type.resize(1);
166 _cycle_type[0]=ParaMEDMEM::Block;
167 delete[] nbelems_per_proc;
170 BlockTopology::~BlockTopology()
172 if (_owns_processor_group)
176 //!converts a pair <subdomainid,local> to a global number
177 std::pair<int,int> BlockTopology::globalToLocal(const int global) const
182 int size_procs=_proc_group->size();
184 vector<int>axis_position(_dimension);
185 vector<int>axis_offset(_dimension);
186 for (int idim=0; idim<_dimension; idim++)
188 int axis_size=_local_array_indices[idim].size()-1;
189 int axis_nb_elem=_local_array_indices[idim][axis_size];
190 increment=increment/axis_nb_elem;
191 int proc_increment = size_procs/(axis_size);
192 int axis_pos=position/increment;
193 position=position%increment;
195 while (_local_array_indices[idim][iaxis]<=axis_pos)
197 subdomain_id+=proc_increment;
200 axis_position[idim]=axis_pos-_local_array_indices[idim][iaxis-1];
201 axis_offset[idim]=iaxis;
204 int local_increment=1;
205 for (int idim=_dimension-1; idim>=0; idim--)
207 local+=axis_position[idim]*local_increment;
208 local_increment*=_local_array_indices[idim][axis_offset[idim]]-_local_array_indices[idim][axis_offset[idim]-1];
210 return make_pair(subdomain_id,local);
213 //!converts local number to a global number
214 int BlockTopology::localToGlobal(const pair<int,int> local) const
217 int subdomain_id=local.first;
219 int loc=local.second;
220 int increment=_nb_elems;
221 int proc_increment=_proc_group->size();
222 int local_increment=getNbLocalElements();
223 for (int idim=0; idim < _dimension; idim++)
225 int axis_size=_local_array_indices[idim].size()-1;
226 int axis_nb_elem=_local_array_indices[idim][axis_size];
227 increment=axis_nb_elem==0?0:increment/axis_nb_elem;
228 proc_increment = proc_increment/(axis_size);
229 int proc_axis=subdomain_id/proc_increment;
230 subdomain_id=subdomain_id%proc_increment;
231 int local_axis_nb_elem=_local_array_indices[idim][proc_axis+1]-_local_array_indices[idim][proc_axis];
232 local_increment = (local_axis_nb_elem==0)?0:(local_increment/local_axis_nb_elem);
233 int iaxis=((local_increment==0)?0:(loc/local_increment))+_local_array_indices[idim][proc_axis];
234 global+=increment*iaxis;
235 loc = (local_increment==0)?0:(loc%local_increment);
240 //Retrieves the local number of elements
241 int BlockTopology::getNbLocalElements()const
243 int position=_proc_group->myRank();
246 for (int i=_dimension-1; i>=0; i--)
248 increment *=_nb_procs_per_dim[i];
249 int idim=position%increment;
250 position=position/increment;
251 int imin=_local_array_indices[i][idim];
252 int imax=_local_array_indices[i][idim+1];
253 nb_elem*=(imax-imin);
258 /*! Retrieves the min and max indices of the domain stored locally
259 * for each dimension. The output vector has the topology dimension
260 * as a size and each pair <int,int> contains min and max. Indices
261 * range from min to max-1.
263 std::vector<std::pair<int,int> > BlockTopology::getLocalArrayMinMax() const
265 vector<pair<int,int> > local_indices (_dimension);
266 int myrank=_proc_group->myRank();
268 for (int i=_dimension-1; i>=0; i--)
270 increment *=_nb_procs_per_dim[i];
271 int idim=myrank%increment;
272 local_indices[i].first=_local_array_indices[i][idim];
273 local_indices[i].second=_local_array_indices[i][idim+1];
274 cout << local_indices[i].first << " "<< local_indices[i].second<<endl;
276 return local_indices;
279 /*! Serializes the data contained in the Block Topology
280 * for communication purposes*/
281 void BlockTopology::serialize(int* & serializer, int& size) const
285 buffer.push_back(_dimension);
286 buffer.push_back(_nb_elems);
287 for (int i=0; i<_dimension; i++)
289 buffer.push_back(_nb_procs_per_dim[i]);
290 buffer.push_back(_cycle_type[i]);
291 buffer.push_back(_local_array_indices[i].size());
292 for (int j=0; j<(int)_local_array_indices[i].size(); j++)
293 buffer.push_back(_local_array_indices[i][j]);
296 //serializing the comm group
297 int size_comm=_proc_group->size();
298 buffer.push_back(size_comm);
299 MPIProcessorGroup world_group(_proc_group->getCommInterface());
300 for (int i=0; i<size_comm;i++)
302 int world_rank=world_group.translateRank(_proc_group, i);
303 buffer.push_back(world_rank);
306 serializer=new int[buffer.size()];
308 copy(buffer.begin(), buffer.end(), serializer);
313 * Unserializes the data contained in the Block Topology
314 * after communication. Uses the same structure as the one used for serialize()
317 void BlockTopology::unserialize(const int* serializer,const CommInterface& comm_interface)
319 const int* ptr_serializer=serializer;
320 cout << "unserialize..."<<endl;
321 _dimension=*(ptr_serializer++);
322 cout << "dimension "<<_dimension<<endl;
323 _nb_elems=*(ptr_serializer++);
324 cout << "nbelems "<<_nb_elems<<endl;
325 _nb_procs_per_dim.resize(_dimension);
326 _cycle_type.resize(_dimension);
327 _local_array_indices.resize(_dimension);
328 for (int i=0; i<_dimension; i++)
330 _nb_procs_per_dim[i]=*(ptr_serializer++);
331 _cycle_type[i]=(CYCLE_TYPE)*(ptr_serializer++);
332 _local_array_indices[i].resize(*(ptr_serializer++));
333 for (int j=0; j<(int)_local_array_indices[i].size(); j++)
334 _local_array_indices[i][j]=*(ptr_serializer++);
337 int size_comm=*(ptr_serializer++);
338 for (int i=0; i<size_comm; i++)
339 procs.insert(*(ptr_serializer++));
340 cout << "unserialize..."<<procs.size()<<endl;
341 _proc_group=new MPIProcessorGroup(comm_interface,procs);
342 _owns_processor_group=true;
343 //TODO manage memory ownership of _proc_group