XGC1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
my_mpi.hpp
Go to the documentation of this file.
1 #ifndef MY_MPI_HPP
2 #define MY_MPI_HPP
3 
4 #include <mpi.h>
5 #include <vector>
6 #include <cstdio>
7 
8 extern MPI_Comm SML_COMM_WORLD;
9 extern int SML_COMM_RANK;
10 extern int SML_COMM_SIZE;
11 
12 void create_comm_world(int color);
13 void destroy_comm_world();
14 
15 /* To illustrate the different communicators: Imagine SML_COMM_WORLD has 12 ranks:
16  * 0 1 2 3 4 5 6 7 8 9 10 11
17  */
18 
19 struct MyMPI{
20  // General communicator
21  MPI_Comm comm;
22  MPI_Group group;
23  int nranks;
24  int my_rank;
25  /* Reorders SML_COMM_WORLD e.g. (with 12 ranks, 4 planes):
26  * 0 3 6 9
27  * 1 4 7 10
28  * 2 5 8 11
29  * So e.g. processes with ranks 0 and 3 in SML_COMM_WORLD are reordered to be
30  * consecutive in comm
31  * Not sure why this transpose is done; maybe the order affects some MPI algorithms
32  * but not clear when/why this would be better.
33  * Also, in the Fortran code comm overwrites SML_COMM_WORLD in the setup (both
34  * called sml_comm), so something similar could be done here. - ALS
35  */
36 
37  // Intraplane communicator
38  MPI_Comm plane_comm;
41  /* Takes a column from comm, e.g. (with 12 ranks, 4 planes, this is rank 4):
42  * 3
43  * 4
44  * 5
45  */
46 
47  // Interplane communicator
48  MPI_Comm intpl_comm;
51  /* Takes a row from comm, e.g. (with 12 ranks, 4 planes, this is rank 4):
52  *
53  * 1 4 7 10
54  *
55  */
56 
57  MyMPI(){}
58 
59  // Spoofed MPI for dry run
60  MyMPI(int nranks_in)
61  : nranks(nranks_in),
62  my_rank(0) {}
63 
64  MyMPI(const MPI_Comm& comm_world, int nplanes){
65  // Get total number of processors
66  MPI_Comm_size( comm_world, &nranks );
67 
68  // Check that nranks is divisible by nplanes
69  if(nranks%nplanes != 0){
70  fflush(stdout);
71  if(SML_COMM_RANK==0){
72  printf("\nERROR: Invalid number of MPI ranks: must be divisible by sml_nphi_total\n");
73  printf("#MPI ranks=%d", nranks);
74  printf("#sml_nphi_total=%d", nplanes);
75  fflush(stdout);
76  }
77  MPI_Abort(SML_COMM_WORLD, 1);
78  }
79 
80  int ranks_per_plane = nranks/nplanes;
81 
82  if(false/*sml_plane_major*/){
83  // get the group underlying sml_comm
84  //mpi_comm_group(sml_comm, sml_comm_group, ierr)
85  }else{
86  // redefine sml_comm pe ordering from consecutive within planes
87  // (Plane_Major) to consecutive across planes (Interplane-Major)
88  std::vector<int> new_sml_comm_ranks(nranks);
89  int k = 0;
90  for(int j=0; j<ranks_per_plane; j++){
91  for(int i=0; i<nplanes; i++){
92  new_sml_comm_ranks[ranks_per_plane*i + j] = k;
93  k++;
94  }
95  }
96 
97  // get the group underlying sml_comm (== mpi_comm_world)
98  // Use sml_comm_world instead of MPI_COMM_WORLD (for MPMD/staging execution)
99  MPI_Group mpi_comm_world_group;
100  MPI_Comm_group(comm_world, &mpi_comm_world_group);
101 
102  // create new group permuting sml_comm pe ordering to Interplane-Major
103  MPI_Group_incl(mpi_comm_world_group, nranks, new_sml_comm_ranks.data(), &group);
104 
105  // Create the new communicator
106  MPI_Comm_create(comm_world, group, &comm);
107  MPI_Comm_rank( comm, &my_rank );
108  }
109 
110  // INTRA-PLANE MPI COMMUNICATOR
111  //checkpoint("\nPlane mpi communication");
112  int plane_0_pe=int(my_rank/ranks_per_plane)*ranks_per_plane;
113  std::vector<int> plane_ranks(ranks_per_plane);
114  for(int j=0; j<ranks_per_plane; j++){
115  plane_ranks[j]=plane_0_pe + j;
116  }
117 
118  // Create the new plane group
119  MPI_Group plane_group;
120  MPI_Group_incl(group, ranks_per_plane, plane_ranks.data(), &plane_group);
121 
122  // Create the new plane communicator
123  MPI_Comm_create(comm, plane_group, &plane_comm);
124 
125  //call mpi_comm_size(plane_comm, plane_nranks, ierr)
126  MPI_Comm_rank(plane_comm, &my_plane_rank);
127  MPI_Comm_size(plane_comm, &n_plane_ranks);
128 
129  // INTER-PLANE MPI COMMUNICATOR
130  //checkpoint("\nInter-plane mpi communication");
131  std::vector<int> intpl_ranks(nplanes);
132  for(int i=0; i<nplanes; i++){
133  intpl_ranks[i]=my_plane_rank + i*ranks_per_plane;
134  }
135 
136  // Create the new inter-plane group
137  MPI_Group intpl_group;
138  MPI_Group_incl(group, nplanes, intpl_ranks.data(), &intpl_group);
139 
140  // Create the new inter-plane communicator
141  MPI_Comm_create(comm, intpl_group, &intpl_comm);
142 
143  MPI_Comm_rank(intpl_comm, &my_intpl_rank);
144  MPI_Comm_size(intpl_comm, &n_intpl_ranks );
145 
146  /*
147  call check_point('adios mpi communication?')
148  do i=0, sml_intpl_nranks-1
149  adios_ranks(i)= i*ranks_per_plane
150  enddo
151  call mpi_group_incl(sml_comm_group, sml_intpl_nranks, adios_ranks, sml_adios_group, ierr)
152  call mpi_comm_create(sml_comm, sml_adios_group, sml_adios_comm, ierr)
153 
154 
155  //#ifndef NO_TASKMAP
156  // output task-to-node mapping
157  write(c_color,'(i8)') sml_comm_color
158  if (sml_mype .eq. 0) then
159  open(unit=14, file='TASKMAP_Color'//trim(adjustl(c_color))//'.txt', &
160  status='OLD', access='SEQUENTIAL', position='APPEND' )
161  endif
162  call taskmap_write(14, sml_comm, &
163  'SML_COMM COLOR #'//trim(adjustl(c_color)), sml_plane_mype, sml_intpl_mype, .false.)
164  if (sml_mype .eq. 0) close(14)
165  //#endif
166  //#ifdef XGC1
167  //jyc print mpi placement info
168  if (sml_verbose) then
169  call mpi_get_processor_name(nodename, len, ierr)
170  print *, sml_mype, 'Process (plane,node,color):', sml_plane_index, trim(nodename), ' ', sml_comm_color
171  endif
172  //#endif
173  */
174  }
175 };
176 
177 // Some useful debugging tools
178 template <typename F>
179 void execute_in_rank_order(const MPI_Comm& comm, F func){
180  int my_rank;
181  int n_ranks;
182  MPI_Comm_rank( comm, &my_rank );
183  MPI_Comm_size( comm, &n_ranks );
184 
185  // double to make sure *every* rank is flushing at *every* barrier
186  fflush(stdout);
187  MPI_Barrier(comm);
188  fflush(stdout);
189  MPI_Barrier(comm);
190  for(int i=0; i<n_ranks; i++){
191  if(my_rank==i) func(i);
192  fflush(stdout);
193  MPI_Barrier(comm);
194  fflush(stdout);
195  MPI_Barrier(comm);
196  }
197 };
198 
199 #endif
MPI_Comm SML_COMM_WORLD
Definition: my_mpi.cpp:4
MPI_Comm plane_comm
Definition: my_mpi.hpp:38
MyMPI(const MPI_Comm &comm_world, int nplanes)
Definition: my_mpi.hpp:64
int my_intpl_rank
Definition: my_mpi.hpp:49
int my_rank
Definition: my_mpi.hpp:24
Definition: my_mpi.hpp:19
void create_comm_world(int color)
Definition: my_mpi.cpp:13
int SML_COMM_SIZE
Definition: my_mpi.cpp:6
MPI_Comm comm
Definition: my_mpi.hpp:21
int my_plane_rank
Definition: my_mpi.hpp:39
int n_plane_ranks
Definition: my_mpi.hpp:40
int SML_COMM_RANK
Definition: my_mpi.cpp:5
void execute_in_rank_order(const MPI_Comm &comm, F func)
Definition: my_mpi.hpp:179
MyMPI()
Definition: my_mpi.hpp:57
MPI_Comm intpl_comm
Definition: my_mpi.hpp:48
void destroy_comm_world()
Definition: my_mpi.cpp:28
MPI_Group group
Definition: my_mpi.hpp:22
int nranks
Definition: my_mpi.hpp:23
int n_intpl_ranks
Definition: my_mpi.hpp:50
MyMPI(int nranks_in)
Definition: my_mpi.hpp:60