1 // Copyright (C) 2007-2013 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.
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"
39 //!converts a pair <subdomainid,local> to a global number
40 std::pair<int,int> BlockTopology::globalToLocal(const int global) const
45 int size_procs=_proc_group->size();
47 vector<int>axis_position(_dimension);
48 vector<int>axis_offset(_dimension);
49 for (int idim=0; idim<_dimension; idim++)
51 int axis_size=_local_array_indices[idim].size()-1;
52 int axis_nb_elem=_local_array_indices[idim][axis_size];
53 increment=increment/axis_nb_elem;
54 int proc_increment = size_procs/(axis_size);
55 int axis_pos=position/increment;
56 position=position%increment;
58 while (_local_array_indices[idim][iaxis]<=axis_pos)
60 subdomain_id+=proc_increment;
63 axis_position[idim]=axis_pos-_local_array_indices[idim][iaxis-1];
64 axis_offset[idim]=iaxis;
67 int local_increment=1;
68 for (int idim=_dimension-1; idim>=0; idim--)
70 local+=axis_position[idim]*local_increment;
71 local_increment*=_local_array_indices[idim][axis_offset[idim]]-_local_array_indices[idim][axis_offset[idim]-1];
73 return make_pair(subdomain_id,local);
76 //!converts local number to a global number
77 int BlockTopology::localToGlobal(const pair<int,int> local) const
80 int subdomain_id=local.first;
83 int increment=_nb_elems;
84 int proc_increment=_proc_group->size();
85 int local_increment=getNbLocalElements();
86 for (int idim=0; idim < _dimension; idim++)
88 int axis_size=_local_array_indices[idim].size()-1;
89 int axis_nb_elem=_local_array_indices[idim][axis_size];
90 increment=axis_nb_elem==0?0:increment/axis_nb_elem;
91 proc_increment = proc_increment/(axis_size);
92 int proc_axis=subdomain_id/proc_increment;
93 subdomain_id=subdomain_id%proc_increment;
94 int local_axis_nb_elem=_local_array_indices[idim][proc_axis+1]-_local_array_indices[idim][proc_axis];
95 local_increment = (local_axis_nb_elem==0)?0:(local_increment/local_axis_nb_elem);
96 int iaxis=((local_increment==0)?0:(loc/local_increment))+_local_array_indices[idim][proc_axis];
97 global+=increment*iaxis;
98 loc = (local_increment==0)?0:(loc%local_increment);
103 //Retrieves the local number of elements
104 int BlockTopology::getNbLocalElements()const
106 int position=_proc_group->myRank();
109 for (int i=_dimension-1; i>=0; i--)
111 increment *=_nb_procs_per_dim[i];
112 int idim=position%increment;
113 position=position/increment;
114 int imin=_local_array_indices[i][idim];
115 int imax=_local_array_indices[i][idim+1];
116 nb_elem*=(imax-imin);
122 * Constructor of a block topology from a grid.
123 * This preliminary version simply splits along the first axis
124 * instead of making the best choice with respect to the
125 * values of the different axes.
127 BlockTopology::BlockTopology(const ProcessorGroup& group, MEDCouplingCMesh *grid):
128 _dimension(grid->getSpaceDimension()), _proc_group(&group), _owns_processor_group(false)
130 vector <int> axis_length(_dimension);
132 for (int idim=0; idim <_dimension; idim++)
134 DataArrayDouble *arr=grid->getCoordsAt(idim);
135 axis_length[idim]=arr->getNbOfElems();
136 _nb_elems*=axis_length[idim];
138 //default splitting along 1st dimension
139 _local_array_indices.resize(_dimension);
140 _nb_procs_per_dim.resize(_dimension);
142 _local_array_indices[0].resize(_proc_group->size()+1);
143 _local_array_indices[0][0]=0;
144 _nb_procs_per_dim[0]=_proc_group->size();
146 for (int i=1; i<=_proc_group->size(); i++)
148 _local_array_indices[0][i]=_local_array_indices[0][i-1]+
149 axis_length[0]/_proc_group->size();
150 if (i<= axis_length[0]%_proc_group->size())
151 _local_array_indices[0][i]+=1;
153 for (int i=1; i<_dimension; i++)
155 _local_array_indices[i].resize(2);
156 _local_array_indices[i][0]=0;
157 _local_array_indices[i][1]=axis_length[i];
158 _nb_procs_per_dim[i]=1;
160 _cycle_type.resize(_dimension);
161 for (int i=0; i<_dimension; i++)
162 _cycle_type[i]=ParaMEDMEM::Block;
166 * Creation of a block topology by composing
167 * a geometrical topology and a component topology.
168 * This constructor is intended for creating fields
169 * for which the parallel distribution is made on the
170 * components of the field rather than on the geometrical
171 * partitioning of the underlying mesh.
174 BlockTopology::BlockTopology(const BlockTopology& geom_topo, const ComponentTopology& comp_topo):_owns_processor_group(false)
176 // so far, the block topology can only be created if the proc group
177 // is either on geom_topo or on comp_topo
178 if (geom_topo.getProcGroup()->size()>1 && comp_topo.nbBlocks()>1)
179 throw INTERP_KERNEL::Exception(LOCALIZED("BlockTopology cannot yet be constructed with both complex geo and components topology"));
181 if (comp_topo.nbComponents()==1)
188 _dimension = geom_topo.getDimension()+1;
189 if (comp_topo.nbBlocks()>1)
190 _proc_group=comp_topo.getProcGroup();
192 _proc_group=geom_topo.getProcGroup();
193 _local_array_indices=geom_topo._local_array_indices;
194 vector<int> comp_indices = *(comp_topo.getBlockIndices());
195 _local_array_indices.push_back(comp_indices);
196 _nb_procs_per_dim=geom_topo._nb_procs_per_dim;
197 _nb_procs_per_dim.push_back(comp_topo.nbBlocks());
198 _cycle_type=geom_topo._cycle_type;
199 _cycle_type.push_back(Block);
200 _nb_elems=geom_topo.getNbElements()*comp_topo.nbComponents();
204 /*! Constructor for creating a one-dimensional
205 * topology from a processor group and a local
206 * number of elements on each processor
208 * The function must be called only by the processors belonging
209 * to group \a group. Calling it from a processor not belonging
210 * to \a group will cause an MPI error, while calling from a subset
211 * of \a group will result in a deadlock.
213 BlockTopology::BlockTopology(const ProcessorGroup& group, int nb_elem):_dimension(1),_proc_group(&group),_owns_processor_group(false)
215 int* nbelems_per_proc = new int[group.size()];
216 const MPIProcessorGroup* mpi_group=dynamic_cast<const MPIProcessorGroup*>(_proc_group);
217 const MPI_Comm* comm=mpi_group->getComm();
219 mpi_group->getCommInterface().allGather(&nbtemp, 1, MPI_INT,
220 nbelems_per_proc, 1, MPI_INT,
224 //splitting along only dimension
225 _local_array_indices.resize(1);
226 _nb_procs_per_dim.resize(1);
228 _local_array_indices[0].resize(_proc_group->size()+1);
229 _local_array_indices[0][0]=0;
230 _nb_procs_per_dim[0]=_proc_group->size();
232 for (int i=1; i<=_proc_group->size(); i++)
234 _local_array_indices[0][i]=_local_array_indices[0][i-1]+
235 nbelems_per_proc[i-1];
236 _nb_elems+=nbelems_per_proc[i-1];
238 _cycle_type.resize(1);
239 _cycle_type[0]=ParaMEDMEM::Block;
240 delete[] nbelems_per_proc;
243 BlockTopology::~BlockTopology()
245 if (_owns_processor_group)
249 /*! Retrieves the min and max indices of the domain stored locally
250 * for each dimension. The output vector has the topology dimension
251 * as a size and each pair <int,int> contains min and max. Indices
252 * range from min to max-1.
254 std::vector<std::pair<int,int> > BlockTopology::getLocalArrayMinMax() const
256 vector<pair<int,int> > local_indices (_dimension);
257 int myrank=_proc_group->myRank();
259 for (int i=_dimension-1; i>=0; i--)
261 increment *=_nb_procs_per_dim[i];
262 int idim=myrank%increment;
263 local_indices[i].first=_local_array_indices[i][idim];
264 local_indices[i].second=_local_array_indices[i][idim+1];
265 cout << local_indices[i].first << " "<< local_indices[i].second<<endl;
267 return local_indices;
270 /*! Serializes the data contained in the Block Topology
271 * for communication purposes*/
272 void BlockTopology::serialize(int* & serializer, int& size) const
276 buffer.push_back(_dimension);
277 buffer.push_back(_nb_elems);
278 for (int i=0; i<_dimension; i++)
280 buffer.push_back(_nb_procs_per_dim[i]);
281 buffer.push_back(_cycle_type[i]);
282 buffer.push_back(_local_array_indices[i].size());
283 for (int j=0; j<(int)_local_array_indices[i].size(); j++)
284 buffer.push_back(_local_array_indices[i][j]);
287 //serializing the comm group
288 int size_comm=_proc_group->size();
289 buffer.push_back(size_comm);
290 MPIProcessorGroup world_group(_proc_group->getCommInterface());
291 for (int i=0; i<size_comm;i++)
293 int world_rank=world_group.translateRank(_proc_group, i);
294 buffer.push_back(world_rank);
297 serializer=new int[buffer.size()];
299 copy(buffer.begin(), buffer.end(), serializer);
304 * Unserializes the data contained in the Block Topology
305 * after communication. Uses the same structure as the one used for serialize()
308 void BlockTopology::unserialize(const int* serializer,const CommInterface& comm_interface)
310 const int* ptr_serializer=serializer;
311 cout << "unserialize..."<<endl;
312 _dimension=*(ptr_serializer++);
313 cout << "dimension "<<_dimension<<endl;
314 _nb_elems=*(ptr_serializer++);
315 cout << "nbelems "<<_nb_elems<<endl;
316 _nb_procs_per_dim.resize(_dimension);
317 _cycle_type.resize(_dimension);
318 _local_array_indices.resize(_dimension);
319 for (int i=0; i<_dimension; i++)
321 _nb_procs_per_dim[i]=*(ptr_serializer++);
322 _cycle_type[i]=(CYCLE_TYPE)*(ptr_serializer++);
323 _local_array_indices[i].resize(*(ptr_serializer++));
324 for (int j=0; j<(int)_local_array_indices[i].size(); j++)
325 _local_array_indices[i][j]=*(ptr_serializer++);
328 int size_comm=*(ptr_serializer++);
329 for (int i=0; i<size_comm; i++)
330 procs.insert(*(ptr_serializer++));
331 cout << "unserialize..."<<procs.size()<<endl;
332 _proc_group=new MPIProcessorGroup(comm_interface,procs);
333 _owns_processor_group=true;
334 //TODO manage memory ownership of _proc_group