1//---------------------------------------------------------------------
2//---------------------------------------------------------------------
3#include <stdio.h>
4#include <stdlib.h>
5#include <math.h>
6#include "header.h"
7#include "mpinpb.h"
8
9#define mod(p,q) ((p)%(q))
10#define max(x,y)      ((x)>(y)? (x) : (y))
11#define min(x,y)      ((x)<(y)? (x) : (y))
12
13void make_set() {
14
15//---------------------------------------------------------------------
16//---------------------------------------------------------------------
17
18//---------------------------------------------------------------------
19//     This function allocates space for a set of cells and fills the set
20//     such that communication between cells on different nodes is only
21//     nearest neighbor
22//---------------------------------------------------------------------
23
24
25      int p, i, j, c, dir, size, excess, ierr,ierrcode;
26
27//---------------------------------------------------------------------
28//     compute square root; add small number to allow for roundoff
29//     (note: this is computed in setup_mpi.f also, but prefer to do
30//     it twice because of some include file problems).
31//---------------------------------------------------------------------
32      ncells = (int)(sqrt((double)(no_nodes) + 0.00001e0));
33
34//---------------------------------------------------------------------
35//     this makes coding easier
36//---------------------------------------------------------------------
37      p = ncells;
38
39//---------------------------------------------------------------------
40//     determine the location of the cell at the bottom of the 3D
41//     array of cells
42//---------------------------------------------------------------------
43      cell_coord(1,1) = mod(node,p) ;
44      cell_coord(2,1) = node/p ;
45      cell_coord(3,1) = 0;
46
47//---------------------------------------------------------------------
48//     set the cell_coords for cells in the rest of the z-layers;
49//     this comes down to a simple linear numbering in the z-direct-
50//     ion, and to the doubly-cyclic numbering in the other dirs
51//---------------------------------------------------------------------
52      for (c = 2; c <= p; c++) {
53         cell_coord(1,c) = mod(cell_coord(1,c-1)+1,p) ;
54         cell_coord(2,c) = mod(cell_coord(2,c-1)-1+p,p) ;
55         cell_coord(3,c) = c-1;
56      }
57
58//---------------------------------------------------------------------
59//     offset all the coordinates by 1 to adjust for Fortran arrays
60//---------------------------------------------------------------------
61      for (dir = 1; dir <= 3; dir++) {
62         for (c = 1; c <= p; c++) {
63            cell_coord(dir,c) = cell_coord(dir,c) + 1;
64         }
65      }
66
67//---------------------------------------------------------------------
68//     slice(dir,n) contains the sequence number of the cell that is in
69//     coordinate plane n in the dir direction
70//---------------------------------------------------------------------
71      for (dir = 1; dir <= 3; dir++) {
72         for (c = 1; c <= p; c++) {
73            slice(dir,cell_coord(dir,c)) = c;
74         }
75      }
76
77
78//---------------------------------------------------------------------
79//     fill the predecessor and successor entries, using the indices
80//     of the bottom cells (they are the same at each level of k
81//     anyway) acting as if full periodicity pertains; note that p is
82//     added to those arguments to the mod functions that might
83//     otherwise return wrong values when using the modulo function
84//---------------------------------------------------------------------
85      i = cell_coord(1,1)-1;
86      j = cell_coord(2,1)-1;
87
88      predecessor(1) = mod(i-1+p,p) + p*j;
89      predecessor(2) = i + p*mod(j-1+p,p);
90      predecessor(3) = mod(i+1,p) + p*mod(j-1+p,p);
91      successor(1)   = mod(i+1,p) + p*j;
92      successor(2)   = i + p*mod(j+1,p);
93      successor(3)   = mod(i-1+p,p) + p*mod(j+1,p);
94
95      int id = RCCE_ue();
96
97//---------------------------------------------------------------------
98//     now compute the sizes of the cells
99//---------------------------------------------------------------------
100      for (dir = 1; dir <= 3; dir++) {
101//---------------------------------------------------------------------
102//     set cell_coord range for each direction
103//---------------------------------------------------------------------
104         size   = grid_points(dir)/p;
105         excess = mod(grid_points(dir),p);
106         for (c = 1; c <= ncells; c++) {
107            if (cell_coord(dir,c) <= excess) {
108               cell_size(dir,c) = size+1;
109               cell_low(dir,c) = (cell_coord(dir,c)-1)*(size+1);
110               cell_high(dir,c) = cell_low(dir,c)+size;
111            } else {
112               cell_size(dir,c) = size;
113               cell_low(dir,c)  = excess*(size+1)+
114                    (cell_coord(dir,c)-excess-1)*size;
115               cell_high(dir,c) = cell_low(dir,c)+size-1;
116            }
117            if (cell_size(dir, c) <= 2) {
118               printf(" Error: Cell size too small. Min size is 3\n");
119               ierrcode = 1;
120               exit(1);
121            }
122         }
123      }
124
125      return;
126}
127
128//---------------------------------------------------------------------
129//---------------------------------------------------------------------
130
131
132void make_color() {
133
134//---------------------------------------------------------------------
135//---------------------------------------------------------------------
136
137//---------------------------------------------------------------------
138//     This function determines cycles in the communication graphs in
139//     the six coordinate directions, and colors the ranks so they know
140//     how to construct deadlock-free blocking communication schedules
141//---------------------------------------------------------------------
142
143      int p, i, j, dir, node_loc, comm_color, node_min, length, start_found;
144
145//---------------------------------------------------------------------
146//     compute square root; add small number to allow for roundoff
147//     (note: this is computed in setup_mpi.f also, but prefer to do
148//     it twice because of some include file problems).
149//---------------------------------------------------------------------
150      ncells = (int)(sqrt((double)(no_nodes) + 0.00001e0));
151
152//---------------------------------------------------------------------
153//     this makes coding easier
154//---------------------------------------------------------------------
155      p = ncells;
156
157      for (dir = 0; dir<6; dir++) {
158
159        node_loc = node_min = node; length = 1; start_found = 0;
160        while (!start_found) {
161          i = mod(node_loc,p) ;
162          j = node_loc/p ;
163
164          switch (dir) {
165            case (WESTDIR):   node_loc = mod(i-1+p,p) + p*j;          break;
166            case (EASTDIR):   node_loc = mod(i+1,p) + p*j;            break;
167            case (SOUTHDIR):  node_loc = i + p*mod(j-1+p,p);          break;
168            case (NORTHDIR):  node_loc = i + p*mod(j+1,p);            break;
169            case (BOTTOMDIR): node_loc = mod(i+1,p) + p*mod(j-1+p,p); break;
170            case (TOPDIR):    node_loc = mod(i-1+p,p) + p*mod(j+1,p); break;
171          }
172
173          // the next block ensures that the node with the lowest rank
174          // in this cycle is colored WHITE (=0), and that nodes an even
175          // number of jumps removed from that lowest-ranked member
176          // are also white. The others are RED (1).
177          if (node_loc <= node_min) {
178            node_min = node_loc;
179            comm_color = 0;
180          } else comm_color = !comm_color;
181          if (node_loc == node) start_found = 1;
182          else length++;
183        }
184        send_color[dir] = comm_color;
185        recv_color[dir] = !send_color[dir];
186        // if the number of nodes in this cycle is odd, we need to treat the
187        // last node before the "start" of the cycle differently
188        if (length%2) {
189          if (node == node_min) recv_color[dir] = 2;
190          i = mod(node,p) ;
191          j = node/p ;
192          switch (dir) {
193            case (WESTDIR):   node_loc = mod(i-1+p,p) + p*j;          break;
194            case (EASTDIR):   node_loc = mod(i+1,p) + p*j;            break;
195            case (SOUTHDIR):  node_loc = i + p*mod(j-1+p,p);          break;
196            case (NORTHDIR):  node_loc = i + p*mod(j+1,p);            break;
197            case (BOTTOMDIR): node_loc = mod(i+1,p) + p*mod(j-1+p,p); break;
198            case (TOPDIR):    node_loc = mod(i-1+p,p) + p*mod(j+1,p); break;
199          }
200          if (node_loc == node_min) send_color[dir] = 2;
201        }
202     }
203     return;
204}
205
206//---------------------------------------------------------------------
207//---------------------------------------------------------------------
208
209
210