MFEM  v4.6.0
Finite element discretization library
communication.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #include "../config/config.hpp"
13 
14 #ifdef MFEM_USE_MPI
15 
16 #include <mpi.h>
17 #ifdef __bgq__
18 #include <mpix.h>
19 #endif
20 
21 #include "array.hpp"
22 #include "table.hpp"
23 #include "sets.hpp"
24 #include "communication.hpp"
25 #include "text.hpp"
26 #include "sort_pairs.hpp"
27 #include "globals.hpp"
28 
29 #include <iostream>
30 #include <map>
31 
32 using namespace std;
33 
34 namespace mfem
35 {
36 
37 GroupTopology::GroupTopology(const GroupTopology &gt)
38  : MyComm(gt.MyComm),
39  group_lproc(gt.group_lproc)
40 {
41  gt.groupmaster_lproc.Copy(groupmaster_lproc);
42  gt.lproc_proc.Copy(lproc_proc);
43  gt.group_mgroup.Copy(group_mgroup);
44 }
45 
46 void GroupTopology::ProcToLProc()
47 {
48  int NRanks;
49  MPI_Comm_size(MyComm, &NRanks);
50 
51  map<int, int> proc_lproc;
52 
53  // The local processor ids are assigned following the group order and within
54  // a group following their ordering in the group. In other words, the ids are
55  // assigned based on their order in the J array of group_lproc.
56  int lproc_counter = 0;
57  for (int i = 0; i < group_lproc.Size_of_connections(); i++)
58  {
59  const pair<const int, int> p(group_lproc.GetJ()[i], lproc_counter);
60  if (proc_lproc.insert(p).second)
61  {
62  lproc_counter++;
63  }
64  }
65  // Note: group_lproc.GetJ()[0] == MyRank --> proc_lproc[MyRank] == 0
66 
67  lproc_proc.SetSize(lproc_counter);
68  for (map<int, int>::iterator it = proc_lproc.begin();
69  it != proc_lproc.end(); ++it)
70  {
71  lproc_proc[it->second] = it->first;
72  }
73 
74  for (int i = 0; i < group_lproc.Size_of_connections(); i++)
75  {
76  group_lproc.GetJ()[i] = proc_lproc[group_lproc.GetJ()[i]];
77  }
78 
79  for (int i = 0; i < NGroups(); i++)
80  {
81  groupmaster_lproc[i] = proc_lproc[groupmaster_lproc[i]];
82  }
83 }
84 
85 void GroupTopology::Create(ListOfIntegerSets &groups, int mpitag)
86 {
87  groups.AsTable(group_lproc); // group_lproc = group_proc
88 
89  Table group_mgroupandproc;
90  group_mgroupandproc.SetDims(NGroups(),
91  group_lproc.Size_of_connections() + NGroups());
92  for (int i = 0; i < NGroups(); i++)
93  {
94  int j = group_mgroupandproc.GetI()[i];
95  group_mgroupandproc.GetI()[i+1] = j + group_lproc.RowSize(i) + 1;
96  group_mgroupandproc.GetJ()[j] = i;
97  j++;
98  for (int k = group_lproc.GetI()[i];
99  j < group_mgroupandproc.GetI()[i+1]; j++, k++)
100  {
101  group_mgroupandproc.GetJ()[j] = group_lproc.GetJ()[k];
102  }
103  }
104 
105  // build groupmaster_lproc with lproc = proc
106  groupmaster_lproc.SetSize(NGroups());
107 
108  // simplest choice of the group owner
109  for (int i = 0; i < NGroups(); i++)
110  {
111  groupmaster_lproc[i] = groups.PickElementInSet(i);
112  }
113 
114  // load-balanced choice of the group owner, which however can lead to
115  // isolated dofs
116  // for (i = 0; i < NGroups(); i++)
117  // groupmaster_lproc[i] = groups.PickRandomElementInSet(i);
118 
119  ProcToLProc();
120 
121  // Build 'group_mgroup':
122 
123  // Use aggregated neighbor communication: at most one send to and/or one
124  // receive from each neighbor.
125 
126  group_mgroup.SetSize(NGroups());
127  MFEM_DEBUG_DO(group_mgroup = -1);
128  for (int g = 0; g < NGroups(); g++)
129  {
130  if (IAmMaster(g)) { group_mgroup[g] = g; }
131  }
132 
133  // The Table 'lproc_cgroup': for each lproc, list the groups that are owned
134  // by this rank or by that lproc.
135  Table lproc_cgroup;
136  {
137  Array<Connection> lproc_cgroup_list;
138  for (int g = 1; g < NGroups(); g++)
139  {
140  if (IAmMaster(g))
141  {
142  const int gs = GetGroupSize(g);
143  const int *lprocs = GetGroup(g);
144  for (int i = 0; i < gs; i++)
145  {
146  if (lprocs[i])
147  {
148  lproc_cgroup_list.Append(Connection(lprocs[i],g));
149  }
150  }
151  }
152  else
153  {
154  lproc_cgroup_list.Append(Connection(GetGroupMaster(g),g));
155  }
156  }
157  lproc_cgroup_list.Sort();
158  lproc_cgroup_list.Unique();
159  lproc_cgroup.MakeFromList(GetNumNeighbors(), lproc_cgroup_list);
160  }
161 
162  // Determine size of the send-receive buffer. For each neighbor the buffer
163  // contains: <send-part><receive-part> with each part consisting of a list of
164  // groups. Each group, g, has group_lproc.RowSize(g)+2 integers: the first
165  // entry is group_lproc.RowSize(g) - the number of processors in the group,
166  // followed by the group-id in the master processor, followed by the ranks of
167  // the processors in the group.
168  Table buffer;
169  buffer.MakeI(2*lproc_cgroup.Size()-2); // excluding the "local" lproc, 0
170  for (int nbr = 1; nbr < lproc_cgroup.Size(); nbr++)
171  {
172  const int send_row = 2*(nbr-1);
173  const int recv_row = send_row+1;
174  const int ng = lproc_cgroup.RowSize(nbr);
175  const int *g = lproc_cgroup.GetRow(nbr);
176  for (int j = 0; j < ng; j++)
177  {
178  const int gs = group_lproc.RowSize(g[j]);
179  if (IAmMaster(g[j]))
180  {
181  buffer.AddColumnsInRow(send_row, gs+2);
182  }
183  else
184  {
185  MFEM_ASSERT(GetGroupMaster(g[j]) == nbr, "internal error");
186  buffer.AddColumnsInRow(recv_row, gs+2);
187  }
188  }
189  }
190  buffer.MakeJ();
191  for (int nbr = 1; nbr < lproc_cgroup.Size(); nbr++)
192  {
193  const int send_row = 2*(nbr-1);
194  const int recv_row = send_row+1;
195  const int ng = lproc_cgroup.RowSize(nbr);
196  const int *g = lproc_cgroup.GetRow(nbr);
197  for (int j = 0; j < ng; j++)
198  {
199  const int gs = group_lproc.RowSize(g[j]);
200  if (IAmMaster(g[j]))
201  {
202  buffer.AddConnection(send_row, gs);
203  buffer.AddConnections(
204  send_row, group_mgroupandproc.GetRow(g[j]), gs+1);
205  }
206  else
207  {
208  buffer.AddColumnsInRow(recv_row, gs+2);
209  }
210  }
211  }
212  buffer.ShiftUpI();
213  Array<MPI_Request> send_requests(lproc_cgroup.Size()-1);
214  Array<MPI_Request> recv_requests(lproc_cgroup.Size()-1);
215  send_requests = MPI_REQUEST_NULL;
216  recv_requests = MPI_REQUEST_NULL;
217  for (int nbr = 1; nbr < lproc_cgroup.Size(); nbr++)
218  {
219  const int send_row = 2*(nbr-1);
220  const int recv_row = send_row+1;
221  const int send_size = buffer.RowSize(send_row);
222  const int recv_size = buffer.RowSize(recv_row);
223  if (send_size > 0)
224  {
225  MPI_Isend(buffer.GetRow(send_row), send_size, MPI_INT, lproc_proc[nbr],
226  mpitag, MyComm, &send_requests[nbr-1]);
227  }
228  if (recv_size > 0)
229  {
230  MPI_Irecv(buffer.GetRow(recv_row), recv_size, MPI_INT, lproc_proc[nbr],
231  mpitag, MyComm, &recv_requests[nbr-1]);
232  }
233  }
234 
235  if (recv_requests.Size() > 0)
236  {
237  int idx;
238  IntegerSet group;
239  while (MPI_Waitany(recv_requests.Size(), recv_requests.GetData(), &idx,
240  MPI_STATUS_IGNORE),
241  idx != MPI_UNDEFINED)
242  {
243  const int recv_size = buffer.RowSize(2*idx+1);
244  const int *recv_buf = buffer.GetRow(2*idx+1);
245  for (int s = 0; s < recv_size; s += recv_buf[s]+2)
246  {
247  group.Recreate(recv_buf[s], recv_buf+s+2);
248  const int g = groups.Lookup(group);
249  MFEM_ASSERT(group_mgroup[g] == -1, "communication error");
250  group_mgroup[g] = recv_buf[s+1];
251  }
252  }
253  }
254 
255  MPI_Waitall(send_requests.Size(), send_requests.GetData(),
256  MPI_STATUSES_IGNORE);
257 
258  // debug barrier: MPI_Barrier(MyComm);
259 }
260 
261 void GroupTopology::Save(ostream &os) const
262 {
263  os << "\ncommunication_groups\n";
264  os << "number_of_groups " << NGroups() << "\n\n";
265 
266  os << "# number of entities in each group, followed by group ids in group\n";
267  for (int group_id = 0; group_id < NGroups(); ++group_id)
268  {
269  int group_size = GetGroupSize(group_id);
270  const int * group_ptr = GetGroup(group_id);
271  os << group_size;
272  for ( int group_member_index = 0; group_member_index < group_size;
273  ++group_member_index)
274  {
275  os << " " << GetNeighborRank( group_ptr[group_member_index] );
276  }
277  os << "\n";
278  }
279 
280  // For future use, optional ownership strategy.
281  // os << "# ownership";
282 }
283 
284 void GroupTopology::Load(istream &in)
285 {
286  // Load in group topology and create list of integer sets. Use constructor
287  // that uses list of integer sets.
288  std::string ident;
289 
290  // Read in number of groups
291  int number_of_groups = -1;
292  in >> ident;
293  MFEM_VERIFY(ident == "number_of_groups",
294  "GroupTopology::Load - expected 'number_of_groups' entry.");
295  in >> number_of_groups;
296 
297  // Skip number of entries in each group comment.
298  skip_comment_lines(in, '#');
299 
300  ListOfIntegerSets integer_sets;
301  for (int group_id = 0; group_id < number_of_groups; ++group_id)
302  {
303  IntegerSet integer_set;
304  Array<int>& array = integer_set;
305  int group_size;
306  in >> group_size;
307  array.Reserve(group_size);
308  for ( int index = 0; index < group_size; ++index )
309  {
310  int value;
311  in >> value;
312  array.Append(value);
313  }
314  integer_sets.Insert(integer_set);
315  }
316 
317  Create(integer_sets, 823);
318 }
319 
321 {
322  copy.SetComm(MyComm);
323  group_lproc.Copy(copy.group_lproc);
324  groupmaster_lproc.Copy(copy.groupmaster_lproc);
325  lproc_proc.Copy(copy.lproc_proc);
326  group_mgroup.Copy(copy.group_mgroup);
327 }
328 
330 {
331  mfem::Swap(MyComm, other.MyComm);
332  mfem::Swap(group_lproc, other.group_lproc);
333  mfem::Swap(groupmaster_lproc, other.groupmaster_lproc);
334  mfem::Swap(lproc_proc, other.lproc_proc);
335  mfem::Swap(group_mgroup, other.group_mgroup);
336 }
337 
338 // Initialize the static mpi_type for the specializations of MPITypeMap:
339 const MPI_Datatype MPITypeMap<int>::mpi_type = MPI_INT;
340 const MPI_Datatype MPITypeMap<double>::mpi_type = MPI_DOUBLE;
341 
342 
344  : gtopo(gt), mode(m)
345 {
346  group_buf_size = 0;
347  requests = NULL;
348  // statuses = NULL;
349  comm_lock = 0;
350  num_requests = 0;
351  request_marker = NULL;
352  buf_offsets = NULL;
353 }
354 
355 void GroupCommunicator::Create(const Array<int> &ldof_group)
356 {
358  for (int i = 0; i < ldof_group.Size(); i++)
359  {
360  int group = ldof_group[i];
361  if (group != 0)
362  {
364  }
365  }
366  group_ldof.MakeJ();
367 
368  for (int i = 0; i < ldof_group.Size(); i++)
369  {
370  int group = ldof_group[i];
371  if (group != 0)
372  {
373  group_ldof.AddConnection(group, i);
374  }
375  }
377 
378  Finalize();
379 }
380 
382 {
383  int request_counter = 0;
384 
385  // size buf_offsets = max(number of groups, number of neighbors)
386  buf_offsets = new int[max(group_ldof.Size(), gtopo.GetNumNeighbors())];
387  buf_offsets[0] = 0;
388  for (int gr = 1; gr < group_ldof.Size(); gr++)
389  {
390  if (group_ldof.RowSize(gr) != 0)
391  {
392  int gr_requests;
393  if (!gtopo.IAmMaster(gr)) // we are not the master
394  {
395  gr_requests = 1;
396  }
397  else
398  {
399  gr_requests = gtopo.GetGroupSize(gr)-1;
400  }
401 
402  request_counter += gr_requests;
403  group_buf_size += gr_requests * group_ldof.RowSize(gr);
404  }
405  }
406 
407  requests = new MPI_Request[request_counter];
408  // statuses = new MPI_Status[request_counter];
409  request_marker = new int[request_counter];
410 
411  // Construct nbr_send_groups and nbr_recv_groups: (nbr 0 = me)
414  for (int gr = 1; gr < group_ldof.Size(); gr++)
415  {
416  const int nldofs = group_ldof.RowSize(gr);
417  if (nldofs == 0) { continue; }
418 
419  if (!gtopo.IAmMaster(gr)) // we are not the master
420  {
422  }
423  else // we are the master
424  {
425  const int grp_size = gtopo.GetGroupSize(gr);
426  const int *grp_nbr_list = gtopo.GetGroup(gr);
427  for (int i = 0; i < grp_size; i++)
428  {
429  if (grp_nbr_list[i] != 0)
430  {
431  nbr_send_groups.AddAColumnInRow(grp_nbr_list[i]);
432  }
433  }
434  }
435  }
438  for (int gr = 1; gr < group_ldof.Size(); gr++)
439  {
440  const int nldofs = group_ldof.RowSize(gr);
441  if (nldofs == 0) { continue; }
442 
443  if (!gtopo.IAmMaster(gr)) // we are not the master
444  {
446  }
447  else // we are the master
448  {
449  const int grp_size = gtopo.GetGroupSize(gr);
450  const int *grp_nbr_list = gtopo.GetGroup(gr);
451  for (int i = 0; i < grp_size; i++)
452  {
453  if (grp_nbr_list[i] != 0)
454  {
455  nbr_send_groups.AddConnection(grp_nbr_list[i], gr);
456  }
457  }
458  }
459  }
462  // The above construction creates the Tables with the column indices
463  // sorted, i.e. the group lists are sorted. To coordinate this order between
464  // processors, we will sort the group lists in the nbr_recv_groups Table
465  // according to their indices in the master. This does not require any
466  // communication because we have access to the group indices in the master
467  // by calling: master_group_id = gtopo.GetGroupMasterGroup(my_group_id).
468  Array<Pair<int,int> > group_ids;
469  for (int nbr = 1; nbr < nbr_recv_groups.Size(); nbr++)
470  {
471  const int num_recv_groups = nbr_recv_groups.RowSize(nbr);
472  if (num_recv_groups > 0)
473  {
474  int *grp_list = nbr_recv_groups.GetRow(nbr);
475  group_ids.SetSize(num_recv_groups);
476  for (int i = 0; i < num_recv_groups; i++)
477  {
478  group_ids[i].one = gtopo.GetGroupMasterGroup(grp_list[i]);
479  group_ids[i].two = grp_list[i]; // my_group_id
480  }
481  group_ids.Sort();
482  for (int i = 0; i < num_recv_groups; i++)
483  {
484  grp_list[i] = group_ids[i].two;
485  }
486  }
487  }
488 }
489 
491 {
492  if (group_ltdof.Size() == group_ldof.Size()) { return; }
493 
495  for (int gr = 1; gr < group_ldof.Size(); gr++)
496  {
497  if (gtopo.IAmMaster(gr))
498  {
500  }
501  }
502  group_ltdof.MakeJ();
503  for (int gr = 1; gr < group_ldof.Size(); gr++)
504  {
505  if (gtopo.IAmMaster(gr))
506  {
507  const int *ldofs = group_ldof.GetRow(gr);
508  const int nldofs = group_ldof.RowSize(gr);
509  for (int i = 0; i < nldofs; i++)
510  {
511  group_ltdof.AddConnection(gr, ldof_ltdof[ldofs[i]]);
512  }
513  }
514  }
516 }
517 
519 {
520  nbr_ltdof.MakeI(nbr_send_groups.Size());
521  for (int nbr = 1; nbr < nbr_send_groups.Size(); nbr++)
522  {
523  const int num_send_groups = nbr_send_groups.RowSize(nbr);
524  if (num_send_groups > 0)
525  {
526  const int *grp_list = nbr_send_groups.GetRow(nbr);
527  for (int i = 0; i < num_send_groups; i++)
528  {
529  const int group = grp_list[i];
530  const int nltdofs = group_ltdof.RowSize(group);
531  nbr_ltdof.AddColumnsInRow(nbr, nltdofs);
532  }
533  }
534  }
535  nbr_ltdof.MakeJ();
536  for (int nbr = 1; nbr < nbr_send_groups.Size(); nbr++)
537  {
538  const int num_send_groups = nbr_send_groups.RowSize(nbr);
539  if (num_send_groups > 0)
540  {
541  const int *grp_list = nbr_send_groups.GetRow(nbr);
542  for (int i = 0; i < num_send_groups; i++)
543  {
544  const int group = grp_list[i];
545  const int nltdofs = group_ltdof.RowSize(group);
546  const int *ltdofs = group_ltdof.GetRow(group);
547  nbr_ltdof.AddConnections(nbr, ltdofs, nltdofs);
548  }
549  }
550  }
551  nbr_ltdof.ShiftUpI();
552 }
553 
555 {
556  nbr_ldof.MakeI(nbr_recv_groups.Size());
557  for (int nbr = 1; nbr < nbr_recv_groups.Size(); nbr++)
558  {
559  const int num_recv_groups = nbr_recv_groups.RowSize(nbr);
560  if (num_recv_groups > 0)
561  {
562  const int *grp_list = nbr_recv_groups.GetRow(nbr);
563  for (int i = 0; i < num_recv_groups; i++)
564  {
565  const int group = grp_list[i];
566  const int nldofs = group_ldof.RowSize(group);
567  nbr_ldof.AddColumnsInRow(nbr, nldofs);
568  }
569  }
570  }
571  nbr_ldof.MakeJ();
572  for (int nbr = 1; nbr < nbr_recv_groups.Size(); nbr++)
573  {
574  const int num_recv_groups = nbr_recv_groups.RowSize(nbr);
575  if (num_recv_groups > 0)
576  {
577  const int *grp_list = nbr_recv_groups.GetRow(nbr);
578  for (int i = 0; i < num_recv_groups; i++)
579  {
580  const int group = grp_list[i];
581  const int nldofs = group_ldof.RowSize(group);
582  const int *ldofs = group_ldof.GetRow(group);
583  nbr_ldof.AddConnections(nbr, ldofs, nldofs);
584  }
585  }
586  }
587  nbr_ldof.ShiftUpI();
588 }
589 
590 template <class T>
591 T *GroupCommunicator::CopyGroupToBuffer(const T *ldata, T *buf, int group,
592  int layout) const
593 {
594  switch (layout)
595  {
596  case 1:
597  {
598  return std::copy(ldata + group_ldof.GetI()[group],
599  ldata + group_ldof.GetI()[group+1],
600  buf);
601  }
602  case 2:
603  {
604  const int nltdofs = group_ltdof.RowSize(group);
605  const int *ltdofs = group_ltdof.GetRow(group);
606  for (int j = 0; j < nltdofs; j++)
607  {
608  buf[j] = ldata[ltdofs[j]];
609  }
610  return buf + nltdofs;
611  }
612  default:
613  {
614  const int nldofs = group_ldof.RowSize(group);
615  const int *ldofs = group_ldof.GetRow(group);
616  for (int j = 0; j < nldofs; j++)
617  {
618  buf[j] = ldata[ldofs[j]];
619  }
620  return buf + nldofs;
621  }
622  }
623 }
624 
625 template <class T>
626 const T *GroupCommunicator::CopyGroupFromBuffer(const T *buf, T *ldata,
627  int group, int layout) const
628 {
629  const int nldofs = group_ldof.RowSize(group);
630  switch (layout)
631  {
632  case 1:
633  {
634  std::copy(buf, buf + nldofs, ldata + group_ldof.GetI()[group]);
635  break;
636  }
637  case 2:
638  {
639  const int *ltdofs = group_ltdof.GetRow(group);
640  for (int j = 0; j < nldofs; j++)
641  {
642  ldata[ltdofs[j]] = buf[j];
643  }
644  break;
645  }
646  default:
647  {
648  const int *ldofs = group_ldof.GetRow(group);
649  for (int j = 0; j < nldofs; j++)
650  {
651  ldata[ldofs[j]] = buf[j];
652  }
653  break;
654  }
655  }
656  return buf + nldofs;
657 }
658 
659 template <class T>
660 const T *GroupCommunicator::ReduceGroupFromBuffer(const T *buf, T *ldata,
661  int group, int layout,
662  void (*Op)(OpData<T>)) const
663 {
664  OpData<T> opd;
665  opd.ldata = ldata;
666  opd.nldofs = group_ldof.RowSize(group);
667  opd.nb = 1;
668  opd.buf = const_cast<T*>(buf);
669 
670  switch (layout)
671  {
672  case 1:
673  {
674  MFEM_ABORT("layout 1 is not supported");
675  T *dest = ldata + group_ldof.GetI()[group];
676  for (int j = 0; j < opd.nldofs; j++)
677  {
678  dest[j] += buf[j];
679  }
680  break;
681  }
682  case 2:
683  {
684  opd.ldofs = const_cast<int*>(group_ltdof.GetRow(group));
685  Op(opd);
686  break;
687  }
688  default:
689  {
690  opd.ldofs = const_cast<int*>(group_ldof.GetRow(group));
691  Op(opd);
692  break;
693  }
694  }
695  return buf + opd.nldofs;
696 }
697 
698 template <class T>
699 void GroupCommunicator::BcastBegin(T *ldata, int layout) const
700 {
701  MFEM_VERIFY(comm_lock == 0, "object is already in use");
702 
703  if (group_buf_size == 0) { return; }
704 
705  int request_counter = 0;
706  switch (mode)
707  {
708  case byGroup: // ***** Communication by groups *****
709  {
710  T *buf;
711  if (layout != 1)
712  {
713  group_buf.SetSize(group_buf_size*sizeof(T));
714  buf = (T *)group_buf.GetData();
715  MFEM_VERIFY(layout != 2 || group_ltdof.Size() == group_ldof.Size(),
716  "'group_ltdof' is not set, use SetLTDofTable()");
717  }
718  else
719  {
720  buf = ldata;
721  }
722 
723  for (int gr = 1; gr < group_ldof.Size(); gr++)
724  {
725  const int nldofs = group_ldof.RowSize(gr);
726 
727  // ignore groups without dofs
728  if (nldofs == 0) { continue; }
729 
730  if (!gtopo.IAmMaster(gr)) // we are not the master
731  {
732  MPI_Irecv(buf,
733  nldofs,
736  40822 + gtopo.GetGroupMasterGroup(gr),
737  gtopo.GetComm(),
738  &requests[request_counter]);
739  request_marker[request_counter] = gr;
740  request_counter++;
741  }
742  else // we are the master
743  {
744  if (layout != 1)
745  {
746  CopyGroupToBuffer(ldata, buf, gr, layout);
747  }
748  const int gs = gtopo.GetGroupSize(gr);
749  const int *nbs = gtopo.GetGroup(gr);
750  for (int i = 0; i < gs; i++)
751  {
752  if (nbs[i] != 0)
753  {
754  MPI_Isend(buf,
755  nldofs,
757  gtopo.GetNeighborRank(nbs[i]),
758  40822 + gtopo.GetGroupMasterGroup(gr),
759  gtopo.GetComm(),
760  &requests[request_counter]);
761  request_marker[request_counter] = -1; // mark as send req.
762  request_counter++;
763  }
764  }
765  }
766  buf += nldofs;
767  }
768  break;
769  }
770 
771  case byNeighbor: // ***** Communication by neighbors *****
772  {
773  group_buf.SetSize(group_buf_size*sizeof(T));
774  T *buf = (T *)group_buf.GetData();
775  for (int nbr = 1; nbr < nbr_send_groups.Size(); nbr++)
776  {
777  const int num_send_groups = nbr_send_groups.RowSize(nbr);
778  if (num_send_groups > 0)
779  {
780  // Possible optimization:
781  // if (num_send_groups == 1) and (layout == 1) then we do not
782  // need to copy the data in order to send it.
783  T *buf_start = buf;
784  const int *grp_list = nbr_send_groups.GetRow(nbr);
785  for (int i = 0; i < num_send_groups; i++)
786  {
787  buf = CopyGroupToBuffer(ldata, buf, grp_list[i], layout);
788  }
789  MPI_Isend(buf_start,
790  buf - buf_start,
792  gtopo.GetNeighborRank(nbr),
793  40822,
794  gtopo.GetComm(),
795  &requests[request_counter]);
796  request_marker[request_counter] = -1; // mark as send request
797  request_counter++;
798  }
799 
800  const int num_recv_groups = nbr_recv_groups.RowSize(nbr);
801  if (num_recv_groups > 0)
802  {
803  // Possible optimization (requires interface change):
804  // if (num_recv_groups == 1) and the (output layout == 1) then
805  // we can receive directly in the output buffer; however, at
806  // this point we do not have that information.
807  const int *grp_list = nbr_recv_groups.GetRow(nbr);
808  int recv_size = 0;
809  for (int i = 0; i < num_recv_groups; i++)
810  {
811  recv_size += group_ldof.RowSize(grp_list[i]);
812  }
813  MPI_Irecv(buf,
814  recv_size,
816  gtopo.GetNeighborRank(nbr),
817  40822,
818  gtopo.GetComm(),
819  &requests[request_counter]);
820  request_marker[request_counter] = nbr;
821  request_counter++;
822  buf_offsets[nbr] = buf - (T*)group_buf.GetData();
823  buf += recv_size;
824  }
825  }
826  MFEM_ASSERT(buf - (T*)group_buf.GetData() == group_buf_size, "");
827  break;
828  }
829  }
830 
831  comm_lock = 1; // 1 - locked for Bcast
832  num_requests = request_counter;
833 }
834 
835 template <class T>
836 void GroupCommunicator::BcastEnd(T *ldata, int layout) const
837 {
838  if (comm_lock == 0) { return; }
839  // The above also handles the case (group_buf_size == 0).
840  MFEM_VERIFY(comm_lock == 1, "object is NOT locked for Bcast");
841 
842  switch (mode)
843  {
844  case byGroup: // ***** Communication by groups *****
845  {
846  if (layout == 1)
847  {
848  MPI_Waitall(num_requests, requests, MPI_STATUSES_IGNORE);
849  }
850  else if (layout == 0)
851  {
852  // copy the received data from the buffer to ldata, as it arrives
853  int idx;
854  while (MPI_Waitany(num_requests, requests, &idx, MPI_STATUS_IGNORE),
855  idx != MPI_UNDEFINED)
856  {
857  int gr = request_marker[idx];
858  if (gr == -1) { continue; } // skip send requests
859 
860  // groups without dofs are skipped, so here nldofs > 0.
861  T *buf = (T *)group_buf.GetData() + group_ldof.GetI()[gr];
862  CopyGroupFromBuffer(buf, ldata, gr, layout);
863  }
864  }
865  break;
866  }
867 
868  case byNeighbor: // ***** Communication by neighbors *****
869  {
870  // copy the received data from the buffer to ldata, as it arrives
871  int idx;
872  while (MPI_Waitany(num_requests, requests, &idx, MPI_STATUS_IGNORE),
873  idx != MPI_UNDEFINED)
874  {
875  int nbr = request_marker[idx];
876  if (nbr == -1) { continue; } // skip send requests
877 
878  const int num_recv_groups = nbr_recv_groups.RowSize(nbr);
879  if (num_recv_groups > 0)
880  {
881  const int *grp_list = nbr_recv_groups.GetRow(nbr);
882  const T *buf = (T*)group_buf.GetData() + buf_offsets[nbr];
883  for (int i = 0; i < num_recv_groups; i++)
884  {
885  buf = CopyGroupFromBuffer(buf, ldata, grp_list[i], layout);
886  }
887  }
888  }
889  break;
890  }
891  }
892 
893  comm_lock = 0; // 0 - no lock
894  num_requests = 0;
895 }
896 
897 template <class T>
898 void GroupCommunicator::ReduceBegin(const T *ldata) const
899 {
900  MFEM_VERIFY(comm_lock == 0, "object is already in use");
901 
902  if (group_buf_size == 0) { return; }
903 
904  int request_counter = 0;
905  group_buf.SetSize(group_buf_size*sizeof(T));
906  T *buf = (T *)group_buf.GetData();
907  switch (mode)
908  {
909  case byGroup: // ***** Communication by groups *****
910  {
911  for (int gr = 1; gr < group_ldof.Size(); gr++)
912  {
913  const int nldofs = group_ldof.RowSize(gr);
914  // ignore groups without dofs
915  if (nldofs == 0) { continue; }
916 
917  if (!gtopo.IAmMaster(gr)) // we are not the master
918  {
919  const int layout = 0;
920  CopyGroupToBuffer(ldata, buf, gr, layout);
921  MPI_Isend(buf,
922  nldofs,
925  43822 + gtopo.GetGroupMasterGroup(gr),
926  gtopo.GetComm(),
927  &requests[request_counter]);
928  request_marker[request_counter] = -1; // mark as send request
929  request_counter++;
930  buf += nldofs;
931  }
932  else // we are the master
933  {
934  const int gs = gtopo.GetGroupSize(gr);
935  const int *nbs = gtopo.GetGroup(gr);
936  buf_offsets[gr] = buf - (T *)group_buf.GetData();
937  for (int i = 0; i < gs; i++)
938  {
939  if (nbs[i] != 0)
940  {
941  MPI_Irecv(buf,
942  nldofs,
944  gtopo.GetNeighborRank(nbs[i]),
945  43822 + gtopo.GetGroupMasterGroup(gr),
946  gtopo.GetComm(),
947  &requests[request_counter]);
948  request_marker[request_counter] = gr;
949  request_counter++;
950  buf += nldofs;
951  }
952  }
953  }
954  }
955  break;
956  }
957 
958  case byNeighbor: // ***** Communication by neighbors *****
959  {
960  for (int nbr = 1; nbr < nbr_send_groups.Size(); nbr++)
961  {
962  // In Reduce operation: send_groups <--> recv_groups
963  const int num_send_groups = nbr_recv_groups.RowSize(nbr);
964  if (num_send_groups > 0)
965  {
966  T *buf_start = buf;
967  const int *grp_list = nbr_recv_groups.GetRow(nbr);
968  for (int i = 0; i < num_send_groups; i++)
969  {
970  const int layout = 0; // ldata is an array on all ldofs
971  buf = CopyGroupToBuffer(ldata, buf, grp_list[i], layout);
972  }
973  MPI_Isend(buf_start,
974  buf - buf_start,
976  gtopo.GetNeighborRank(nbr),
977  43822,
978  gtopo.GetComm(),
979  &requests[request_counter]);
980  request_marker[request_counter] = -1; // mark as send request
981  request_counter++;
982  }
983 
984  // In Reduce operation: send_groups <--> recv_groups
985  const int num_recv_groups = nbr_send_groups.RowSize(nbr);
986  if (num_recv_groups > 0)
987  {
988  const int *grp_list = nbr_send_groups.GetRow(nbr);
989  int recv_size = 0;
990  for (int i = 0; i < num_recv_groups; i++)
991  {
992  recv_size += group_ldof.RowSize(grp_list[i]);
993  }
994  MPI_Irecv(buf,
995  recv_size,
997  gtopo.GetNeighborRank(nbr),
998  43822,
999  gtopo.GetComm(),
1000  &requests[request_counter]);
1001  request_marker[request_counter] = nbr;
1002  request_counter++;
1003  buf_offsets[nbr] = buf - (T*)group_buf.GetData();
1004  buf += recv_size;
1005  }
1006  }
1007  MFEM_ASSERT(buf - (T*)group_buf.GetData() == group_buf_size, "");
1008  break;
1009  }
1010  }
1011 
1012  comm_lock = 2;
1013  num_requests = request_counter;
1014 }
1015 
1016 template <class T>
1017 void GroupCommunicator::ReduceEnd(T *ldata, int layout,
1018  void (*Op)(OpData<T>)) const
1019 {
1020  if (comm_lock == 0) { return; }
1021  // The above also handles the case (group_buf_size == 0).
1022  MFEM_VERIFY(comm_lock == 2, "object is NOT locked for Reduce");
1023 
1024  switch (mode)
1025  {
1026  case byGroup: // ***** Communication by groups *****
1027  {
1028  OpData<T> opd;
1029  opd.ldata = ldata;
1030  Array<int> group_num_req(group_ldof.Size());
1031  for (int gr = 1; gr < group_ldof.Size(); gr++)
1032  {
1033  group_num_req[gr] =
1034  gtopo.IAmMaster(gr) ? gtopo.GetGroupSize(gr)-1 : 0;
1035  }
1036  int idx;
1037  while (MPI_Waitany(num_requests, requests, &idx, MPI_STATUS_IGNORE),
1038  idx != MPI_UNDEFINED)
1039  {
1040  int gr = request_marker[idx];
1041  if (gr == -1) { continue; } // skip send requests
1042 
1043  // Delay the processing of a group until all receive requests, for
1044  // that group, are done:
1045  if ((--group_num_req[gr]) != 0) { continue; }
1046 
1047  opd.nldofs = group_ldof.RowSize(gr);
1048  // groups without dofs are skipped, so here nldofs > 0.
1049 
1050  opd.buf = (T *)group_buf.GetData() + buf_offsets[gr];
1051  opd.ldofs = (layout == 0) ?
1052  group_ldof.GetRow(gr) : group_ltdof.GetRow(gr);
1053  opd.nb = gtopo.GetGroupSize(gr)-1;
1054  Op(opd);
1055  }
1056  break;
1057  }
1058 
1059  case byNeighbor: // ***** Communication by neighbors *****
1060  {
1061  MPI_Waitall(num_requests, requests, MPI_STATUSES_IGNORE);
1062 
1063  for (int nbr = 1; nbr < nbr_send_groups.Size(); nbr++)
1064  {
1065  // In Reduce operation: send_groups <--> recv_groups
1066  const int num_recv_groups = nbr_send_groups.RowSize(nbr);
1067  if (num_recv_groups > 0)
1068  {
1069  const int *grp_list = nbr_send_groups.GetRow(nbr);
1070  const T *buf = (T*)group_buf.GetData() + buf_offsets[nbr];
1071  for (int i = 0; i < num_recv_groups; i++)
1072  {
1073  buf = ReduceGroupFromBuffer(buf, ldata, grp_list[i],
1074  layout, Op);
1075  }
1076  }
1077  }
1078  break;
1079  }
1080  }
1081 
1082  comm_lock = 0; // 0 - no lock
1083  num_requests = 0;
1084 }
1085 
1086 template <class T>
1088 {
1089  if (opd.nb == 1)
1090  {
1091  for (int i = 0; i < opd.nldofs; i++)
1092  {
1093  opd.ldata[opd.ldofs[i]] += opd.buf[i];
1094  }
1095  }
1096  else
1097  {
1098  for (int i = 0; i < opd.nldofs; i++)
1099  {
1100  T data = opd.ldata[opd.ldofs[i]];
1101  for (int j = 0; j < opd.nb; j++)
1102  {
1103  data += opd.buf[j*opd.nldofs+i];
1104  }
1105  opd.ldata[opd.ldofs[i]] = data;
1106  }
1107  }
1108 }
1109 
1110 template <class T>
1112 {
1113  for (int i = 0; i < opd.nldofs; i++)
1114  {
1115  T data = opd.ldata[opd.ldofs[i]];
1116  for (int j = 0; j < opd.nb; j++)
1117  {
1118  T b = opd.buf[j*opd.nldofs+i];
1119  if (data > b)
1120  {
1121  data = b;
1122  }
1123  }
1124  opd.ldata[opd.ldofs[i]] = data;
1125  }
1126 }
1127 
1128 template <class T>
1130 {
1131  for (int i = 0; i < opd.nldofs; i++)
1132  {
1133  T data = opd.ldata[opd.ldofs[i]];
1134  for (int j = 0; j < opd.nb; j++)
1135  {
1136  T b = opd.buf[j*opd.nldofs+i];
1137  if (data < b)
1138  {
1139  data = b;
1140  }
1141  }
1142  opd.ldata[opd.ldofs[i]] = data;
1143  }
1144 }
1145 
1146 template <class T>
1148 {
1149  for (int i = 0; i < opd.nldofs; i++)
1150  {
1151  T data = opd.ldata[opd.ldofs[i]];
1152  for (int j = 0; j < opd.nb; j++)
1153  {
1154  data |= opd.buf[j*opd.nldofs+i];
1155  }
1156  opd.ldata[opd.ldofs[i]] = data;
1157  }
1158 }
1159 
1160 void GroupCommunicator::PrintInfo(std::ostream &os) const
1161 {
1162  char c = '\0';
1163  const int tag = 46800;
1164  const int myid = gtopo.MyRank();
1165 
1166  int num_sends = 0, num_recvs = 0;
1167  size_t mem_sends = 0, mem_recvs = 0;
1168  int num_master_groups = 0, num_empty_groups = 0;
1169  int num_active_neighbors = 0; // for mode == byNeighbor
1170  switch (mode)
1171  {
1172  case byGroup:
1173  for (int gr = 1; gr < group_ldof.Size(); gr++)
1174  {
1175  const int nldofs = group_ldof.RowSize(gr);
1176  if (nldofs == 0)
1177  {
1178  num_empty_groups++;
1179  continue;
1180  }
1181  if (gtopo.IAmMaster(gr))
1182  {
1183  num_sends += (gtopo.GetGroupSize(gr)-1);
1184  mem_sends += sizeof(double)*nldofs*(gtopo.GetGroupSize(gr)-1);
1185  num_master_groups++;
1186  }
1187  else
1188  {
1189  num_recvs++;
1190  mem_recvs += sizeof(double)*nldofs;
1191  }
1192  }
1193  break;
1194 
1195  case byNeighbor:
1196  for (int gr = 1; gr < group_ldof.Size(); gr++)
1197  {
1198  const int nldofs = group_ldof.RowSize(gr);
1199  if (nldofs == 0)
1200  {
1201  num_empty_groups++;
1202  continue;
1203  }
1204  if (gtopo.IAmMaster(gr))
1205  {
1206  num_master_groups++;
1207  }
1208  }
1209  for (int nbr = 1; nbr < nbr_send_groups.Size(); nbr++)
1210  {
1211  const int num_send_groups = nbr_send_groups.RowSize(nbr);
1212  if (num_send_groups > 0)
1213  {
1214  const int *grp_list = nbr_send_groups.GetRow(nbr);
1215  for (int i = 0; i < num_send_groups; i++)
1216  {
1217  mem_sends += sizeof(double)*group_ldof.RowSize(grp_list[i]);
1218  }
1219  num_sends++;
1220  }
1221 
1222  const int num_recv_groups = nbr_recv_groups.RowSize(nbr);
1223  if (num_recv_groups > 0)
1224  {
1225  const int *grp_list = nbr_recv_groups.GetRow(nbr);
1226  for (int i = 0; i < num_recv_groups; i++)
1227  {
1228  mem_recvs += sizeof(double)*group_ldof.RowSize(grp_list[i]);
1229  }
1230  num_recvs++;
1231  }
1232  if (num_send_groups > 0 || num_recv_groups > 0)
1233  {
1234  num_active_neighbors++;
1235  }
1236  }
1237  break;
1238  }
1239  if (myid != 0)
1240  {
1241  MPI_Recv(&c, 1, MPI_CHAR, myid-1, tag, gtopo.GetComm(),
1242  MPI_STATUS_IGNORE);
1243  }
1244  else
1245  {
1246  os << "\nGroupCommunicator:\n";
1247  }
1248  os << "Rank " << myid << ":\n"
1249  " mode = " <<
1250  (mode == byGroup ? "byGroup" : "byNeighbor") << "\n"
1251  " number of sends = " << num_sends <<
1252  " (" << mem_sends << " bytes)\n"
1253  " number of recvs = " << num_recvs <<
1254  " (" << mem_recvs << " bytes)\n";
1255  os <<
1256  " num groups = " << group_ldof.Size() << " = " <<
1257  num_master_groups << " + " <<
1258  group_ldof.Size()-num_master_groups-num_empty_groups << " + " <<
1259  num_empty_groups << " (master + slave + empty)\n";
1260  if (mode == byNeighbor)
1261  {
1262  os <<
1263  " num neighbors = " << nbr_send_groups.Size() << " = " <<
1264  num_active_neighbors << " + " <<
1265  nbr_send_groups.Size()-num_active_neighbors <<
1266  " (active + inactive)\n";
1267  }
1268  if (myid != gtopo.NRanks()-1)
1269  {
1270  os << std::flush;
1271  MPI_Send(&c, 1, MPI_CHAR, myid+1, tag, gtopo.GetComm());
1272  }
1273  else
1274  {
1275  os << std::endl;
1276  }
1277  MPI_Barrier(gtopo.GetComm());
1278 }
1279 
1281 {
1282  delete [] buf_offsets;
1283  delete [] request_marker;
1284  // delete [] statuses;
1285  delete [] requests;
1286 }
1287 
1288 // @cond DOXYGEN_SKIP
1289 
1290 // instantiate GroupCommunicator::Bcast and Reduce for int and double
1291 template void GroupCommunicator::BcastBegin<int>(int *, int) const;
1292 template void GroupCommunicator::BcastEnd<int>(int *, int) const;
1293 template void GroupCommunicator::ReduceBegin<int>(const int *) const;
1294 template void GroupCommunicator::ReduceEnd<int>(
1295  int *, int, void (*)(OpData<int>)) const;
1296 
1297 template void GroupCommunicator::BcastBegin<double>(double *, int) const;
1298 template void GroupCommunicator::BcastEnd<double>(double *, int) const;
1299 template void GroupCommunicator::ReduceBegin<double>(const double *) const;
1300 template void GroupCommunicator::ReduceEnd<double>(
1301  double *, int, void (*)(OpData<double>)) const;
1302 
1303 // @endcond
1304 
1305 // instantiate reduce operators for int and double
1306 template void GroupCommunicator::Sum<int>(OpData<int>);
1307 template void GroupCommunicator::Min<int>(OpData<int>);
1308 template void GroupCommunicator::Max<int>(OpData<int>);
1309 template void GroupCommunicator::BitOR<int>(OpData<int>);
1310 
1311 template void GroupCommunicator::Sum<double>(OpData<double>);
1312 template void GroupCommunicator::Min<double>(OpData<double>);
1313 template void GroupCommunicator::Max<double>(OpData<double>);
1314 
1315 
1316 #ifdef __bgq__
1317 static void DebugRankCoords(int** coords, int dim, int size)
1318 {
1319  for (int i = 0; i < size; i++)
1320  {
1321  mfem::out << "Rank " << i << " coords: ";
1322  for (int j = 0; j < dim; j++)
1323  {
1324  mfem::out << coords[i][j] << " ";
1325  }
1326  mfem::out << endl;
1327  }
1328 }
1329 
1330 struct CompareCoords
1331 {
1332  CompareCoords(int coord) : coord(coord) {}
1333  int coord;
1334 
1335  bool operator()(int* const &a, int* const &b) const
1336  { return a[coord] < b[coord]; }
1337 };
1338 
1339 void KdTreeSort(int** coords, int d, int dim, int size)
1340 {
1341  if (size > 1)
1342  {
1343  bool all_same = true;
1344  for (int i = 1; i < size && all_same; i++)
1345  {
1346  for (int j = 0; j < dim; j++)
1347  {
1348  if (coords[i][j] != coords[0][j]) { all_same = false; break; }
1349  }
1350  }
1351  if (all_same) { return; }
1352 
1353  // sort by coordinate 'd'
1354  std::sort(coords, coords + size, CompareCoords(d));
1355  int next = (d + 1) % dim;
1356 
1357  if (coords[0][d] < coords[size-1][d])
1358  {
1359  KdTreeSort(coords, next, dim, size/2);
1360  KdTreeSort(coords + size/2, next, dim, size - size/2);
1361  }
1362  else
1363  {
1364  // skip constant dimension
1365  KdTreeSort(coords, next, dim, size);
1366  }
1367  }
1368 }
1369 
1370 MPI_Comm ReorderRanksZCurve(MPI_Comm comm)
1371 {
1372  MPI_Status status;
1373 
1374  int rank, size;
1375  MPI_Comm_rank(comm, &rank);
1376  MPI_Comm_size(comm, &size);
1377 
1378  int dim;
1379  MPIX_Torus_ndims(&dim);
1380 
1381  int* mycoords = new int[dim + 1];
1382  MPIX_Rank2torus(rank, mycoords);
1383 
1384  MPI_Send(mycoords, dim, MPI_INT, 0, 111, comm);
1385  delete [] mycoords;
1386 
1387  if (rank == 0)
1388  {
1389  int** coords = new int*[size];
1390  for (int i = 0; i < size; i++)
1391  {
1392  coords[i] = new int[dim + 1];
1393  coords[i][dim] = i;
1394  MPI_Recv(coords[i], dim, MPI_INT, i, 111, comm, &status);
1395  }
1396 
1397  KdTreeSort(coords, 0, dim, size);
1398 
1399  // DebugRankCoords(coords, dim, size);
1400 
1401  for (int i = 0; i < size; i++)
1402  {
1403  MPI_Send(&coords[i][dim], 1, MPI_INT, i, 112, comm);
1404  delete [] coords[i];
1405  }
1406  delete [] coords;
1407  }
1408 
1409  int new_rank;
1410  MPI_Recv(&new_rank, 1, MPI_INT, 0, 112, comm, &status);
1411 
1412  MPI_Comm new_comm;
1413  MPI_Comm_split(comm, 0, new_rank, &new_comm);
1414  return new_comm;
1415 }
1416 
1417 #else // __bgq__
1418 
1419 MPI_Comm ReorderRanksZCurve(MPI_Comm comm)
1420 {
1421  // pass
1422  return comm;
1423 }
1424 #endif // __bgq__
1425 
1426 } // namespace mfem
1427 
1428 #endif
int Lookup(IntegerSet &s)
Definition: sets.cpp:104
void Create(ListOfIntegerSets &groups, int mpitag)
Set up the group topology given the list of sets of shared entities.
int NRanks() const
Return the number of MPI ranks within this object&#39;s communicator.
int GetGroupMasterRank(int g) const
Return the rank of the group master for group &#39;g&#39;.
void Recreate(const int n, const int *p)
Create an integer set from C-array &#39;p&#39; of &#39;n&#39; integers. Overwrites any existing set data...
Definition: sets.cpp:68
int * GetJ()
Definition: table.hpp:114
void Unique()
Removes duplicities from a sorted array. This requires operator== to be defined for T...
Definition: array.hpp:259
void AsTable(Table &t)
Write the list of sets into table &#39;t&#39;.
Definition: sets.cpp:116
Helper struct to convert a C++ type to an MPI type.
GroupCommunicator(const GroupTopology &gt, Mode m=byNeighbor)
Construct a GroupCommunicator object.
MPI_Comm ReorderRanksZCurve(MPI_Comm comm)
void AddColumnsInRow(int r, int ncol)
Definition: table.hpp:78
void MakeI(int nrows)
Next 7 methods are used together with the default constructor.
Definition: table.cpp:81
T * CopyGroupToBuffer(const T *ldata, T *buf, int group, int layout) const
Copy the entries corresponding to the group group from the local array ldata to the buffer buf...
void Create(const Array< int > &ldof_group)
Initialize the communicator from a local-dof to group map. Finalize() is called internally.
void ReduceEnd(T *ldata, int layout, void(*Op)(OpData< T >)) const
Finalize reduction operation started with ReduceBegin().
void SetDims(int rows, int nnz)
Definition: table.cpp:140
T * GetData()
Returns the data.
Definition: array.hpp:115
int GetGroupMasterGroup(int g) const
Return the group number in the master for group &#39;g&#39;.
int NGroups() const
Return the number of groups.
int GetNeighborRank(int i) const
Return the MPI rank of neighbor &#39;i&#39;.
int MyRank() const
Return the MPI rank within this object&#39;s communicator.
int GetGroupMaster(int g) const
Return the neighbor index of the group master for a given group. Neighbor 0 is the local processor...
STL namespace.
void GetNeighborLDofTable(Table &nbr_ldof) const
Dofs to be received from communication neighbors.
void Swap(GroupTopology &other)
Swap the internal data with another GroupTopology object.
int RowSize(int i) const
Definition: table.hpp:108
void Save(std::ostream &out) const
Save the data in a stream.
void skip_comment_lines(std::istream &is, const char comment_char)
Check if the stream starts with comment_char. If so skip it.
Definition: text.hpp:31
void AddConnections(int r, const int *c, int nc)
Definition: table.cpp:104
int PickElementInSet(int i)
Return the value of the first element of the ith set.
Definition: sets.hpp:73
void MakeFromList(int nrows, const Array< Connection > &list)
Definition: table.cpp:279
static void Sum(OpData< T >)
Reduce operation Sum, instantiated for int and double.
int Append(const T &el)
Append element &#39;el&#39; to array, resize if necessary.
Definition: array.hpp:759
double b
Definition: lissajous.cpp:42
void Finalize()
Allocate internal buffers after the GroupLDofTable is defined.
const GroupTopology & gtopo
void ReduceBegin(const T *ldata) const
Begin reduction operation within each group where the master is the root.
void AddConnection(int r, int c)
Definition: table.hpp:80
void BcastBegin(T *ldata, int layout) const
Begin a broadcast within each group where the master is the root.
int Insert(IntegerSet &s)
Check to see if set &#39;s&#39; is in the list. If not append it to the end of the list. Returns the index of...
Definition: sets.cpp:91
void Reserve(int capacity)
Ensures that the allocated size is at least the given size.
Definition: array.hpp:159
const T * ReduceGroupFromBuffer(const T *buf, T *ldata, int group, int layout, void(*Op)(OpData< T >)) const
Perform the reduction operation Op to the entries of group group using the values from the buffer buf...
const int * GetGroup(int g) const
Return a pointer to a list of neighbors for a given group. Neighbor 0 is the local processor...
const T * CopyGroupFromBuffer(const T *buf, T *ldata, int group, int layout) const
Copy the entries corresponding to the group group from the buffer buf to the local array ldata...
static void Min(OpData< T >)
Reduce operation Min, instantiated for int and double.
void Sort()
Sorts the array in ascending order. This requires operator< to be defined for T.
Definition: array.hpp:251
void GetRow(int i, Array< int > &row) const
Return row i in array row (the Table must be finalized)
Definition: table.cpp:187
void Swap(Array< T > &, Array< T > &)
Definition: array.hpp:638
void GetNeighborLTDofTable(Table &nbr_ltdof) const
Dofs to be sent to communication neighbors.
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
static void Max(OpData< T >)
Reduce operation Max, instantiated for int and double.
Communications are performed one group at a time.
void AddAColumnInRow(int r)
Definition: table.hpp:77
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition: array.hpp:687
Helper struct for defining a connectivity table, see Table::MakeFromList.
Definition: table.hpp:27
void SetLTDofTable(const Array< int > &ldof_ltdof)
Initialize the internal group_ltdof Table.
void ShiftUpI()
Definition: table.cpp:115
Data structure on which we define reduce operations. The data is associated with (and the operation i...
double a
Definition: lissajous.cpp:41
~GroupCommunicator()
Destroy a GroupCommunicator object, deallocating internal data structures and buffers.
MPI_Comm GetComm() const
Return the MPI communicator.
int GetGroupSize(int g) const
Get the number of processors in a group.
int Size() const
Returns the number of TYPE I elements.
Definition: table.hpp:92
A set of integers.
Definition: sets.hpp:23
Mode
Communication mode.
void MakeJ()
Definition: table.cpp:91
int dim
Definition: ex24.cpp:53
void KdTreeSort(int **coords, int d, int dim, int size)
int index(int i, int j, int nx, int ny)
Definition: life.cpp:235
void Copy(Array &copy) const
Create a copy of the internal array to the provided copy.
Definition: array.hpp:864
int Size_of_connections() const
Definition: table.hpp:98
int Size() const
Return the logical size of the array.
Definition: array.hpp:141
void SetComm(MPI_Comm comm)
Set the MPI communicator to &#39;comm&#39;.
void Copy(Table &copy) const
Definition: table.cpp:389
int * GetI()
Definition: table.hpp:113
RefCoord s[3]
void BcastEnd(T *ldata, int layout) const
Finalize a broadcast started with BcastBegin().
List of integer sets.
Definition: sets.hpp:62
void Copy(GroupTopology &copy) const
Copy the internal data to the external &#39;copy&#39;.
void PrintInfo(std::ostream &out=mfem::out) const
Print information about the GroupCommunicator from all MPI ranks.
void Load(std::istream &in)
Load the data from a stream.
bool IAmMaster(int g) const
Return true if I am master for group &#39;g&#39;.
int GetNumNeighbors() const
Return the number of neighbors including the local processor.
static void BitOR(OpData< T >)
Reduce operation bitwise OR, instantiated for int only.