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