Lesson 11 — A grid of blocks: halo exchange and what it costs

Lessons 8 to 10 used a row of PEs. This lesson uses a rectangle of them. It splits a 2-D field into square blocks, one block per PE, and runs a Jacobi stencil on it. Every sweep, each PE sends the four edges of its block to its four neighbours and receives theirs. This is a halo exchange, the communication pattern of most grid codes. The lesson also writes the Paint layout as a function of the grid size, and checks what the kernel uses against what a PE has.

What you will be able to do

• Split a 2-D field into B × B blocks, one per PE, and exchange the block edges every sweep.
• Route a grid of any width and height with four colorpairs, and close it with a ring of moat PEs.
• Write a layout as a Paint function of W and H with define, and repeat a PE across the grid with hstackrep and vstackrep.
• Check a kernel against a PE's budget of colors, queues, microthreads and memory.
• Measure the exchange against the computation, and predict how each grows with B.

The program

The field has H·B rows and W·B columns of sp values. One Jacobi sweep replaces every point by the mean of its four neighbours:

u_new[i][j] = 0.25 * (u[i-1][j] + u[i+1][j] + u[i][j-1] + u[i][j+1])

Points outside the field count as 0. The test starts from random integers between 0 and 63 and runs T sweeps.

W × H PEs hold the field, one B × B block each. A ring of moat PEs surrounds them. A moat PE answers every wavelet a grid PE sends it with a 0: the zero outside the field. Here is the layout for W = H = 2:

        x=0    x=1    x=2    x=3
y=3   [moat] [moat] [moat] [moat]
y=2   [moat] [blk ] [blk ] [moat]       blk:  block.w, one B x B block
y=1   [moat] [blk ] [blk ] [moat]       moat: moat.w, answers with 0
y=0   [moat] [moat] [moat] [moat]

The four corner moats never receive anything. The whole (W+2) × (H+2) rectangle is the beach of lesson 8: every PE in it gets the same route rules.

After T sweeps, the test compares every PE's u, ghost cells included (see below), with a NumPy computation in float32. They must be equal bit for bit. After the last sweep, the ghost cells hold the neighbours' edges from before the last update, and ref.py builds them that way. Each row also prints the cycles that the last sweep took.

Row	B	W × H	T
tutorial	4	2 × 2	3
b8	8	3 × 2	2
b16	16	2 × 2	1
weak	8	4 × 4	2

Why a block per PE

With one point per PE, a PE does three additions and one multiplication per sweep, and exchanges four values. With a B × B block it computes B² points and exchanges 4B values. The computation grows with the area of the block and the communication with its perimeter, so a larger block spends a smaller share of each sweep communicating. Memory sets the largest block (see The budget).

Why colorpairs

A 2 × 2 grid can be wired with one color per direction:

          x=0                 x=1
y=1     pe01  ---hfwd--->   pe11
              <--hrev----
        |  ^                |  ^
   vrev |  | vfwd      vrev |  | vfwd
        v  |                v  |
y=0     pe00  ---hfwd--->   pe10
              <--hrev----

Each color then has one sender and one receiver in each row or column. pe00 gets hfwd [C>R] and pe10 gets hfwd [L>C]. The four PEs need four code() calls with different socket bindings and sixteen route rules.

This does not extend to a third column. The middle PE of a row of three would receive hfwd from the left and also send its own value on hfwd to the right. That is two inputs for one color, and a router takes one input per color (lesson 4). Lesson 8 solved this for a row with a colorpair. This lesson uses one colorpair per direction, so every PE of any W × H grid gets the same eight rules.

The code

The block: block.w

socket a;    // sends right, receives from the left
socket b;    // sends left,  receives from the right
socket c;    // sends up,    receives from below
socket d;    // sends down,  receives from above

sp u[n=B+2][m=B+2];    // written by miset
sp v[n=B][m=B];

u holds the block in u[1..B][1..B] and a ring of ghost cells around it. The ghost cells hold copies of the neighbours' edge values, which the update of the block's edge points reads. The first index, n, grows with Paint's y, so n = B is the top row of the block; the second, m, grows with x. For B = 4:

n=5    .  d  d  d  d  .        ghost row:    from the PE above, on d
n=4    a  *  *  *  *  b
n=3    a  *  *  *  *  b        ghost column: from the left on a (m=0),
n=2    a  *  *  *  *  b                      from the right on b (m=5)
n=1    a  *  *  *  *  b
n=0    .  c  c  c  c  .        ghost row:    from the PE below, on c
      m=0 1  2  3  4  5

The four corner ghosts are never written and stay 0; the stencil does not read them. Each PE sends its own edges the other way: column m = B right on a, column m = 1 left on b, row n = B up on c and row n = 1 down on d.

The exchange

        parallel {
            ∀e a[] <- u[e][B];
            ∀e b[] <- u[e][1];
            ∀e c[] <- u[B][e];
            ∀e d[] <- u[1][e];
            ∀e u[e][0] <- a[];
            ∀e u[e][B+1] <- b[];
            ∀e u[0][e] <- c[];
            ∀e u[B+1][e] <- d[];
        }

e runs over [1, B+1), the B rows or columns of the block. Each arm moves one edge: four arms send, four receive. u[e][B] is a column of the block, a vector whose elements lie one row of u, B + 2 values, apart in memory. Every sp value travels as one 32-bit wavelet.

The eight arms are active at the same time: each is a microthread, as in lesson 5. They still share one core: its pipeline takes an instruction from one thread at a time, and its one port to the router moves one wavelet per cycle each way (see Exchange against update). A PE sends to a neighbour while that neighbour sends to it, so each PE must be ready to receive while it sends. In sequence, the sends would finish only while the buffers between the PEs could hold a whole edge (lesson 6); see Pitfalls.

The update

        ∀p ∀q v[p][q] <- u[p][q+1] + u[p+2][q+1];
        ∀p ∀q v[p][q] <+- u[p+1][q];
        ∀p ∀q v[p][q] <+- u[p+1][q+2];
        ∀p ∀q u[p+1][q+1] <- 0.25:sp * v[p][q];

p and q run over [0, B). Point (p, q) of the block is u[p+1][q+1], so u[p][q+1] is the point below it, u[p+2][q+1] the point above, and u[p+1][q] and u[p+1][q+2] its left and right neighbours. <+- adds to its destination, like +=. The sum goes into v first, because each new value of u must be computed from old values only.

The update is four statements because each one is a single vector instruction over the whole block (lesson 3). Written as one statement, v[p][q] <- 0.25:sp * (u[p][q+1] + u[p+2][q+1] + u[p+1][q] + u[p+1][q+2]), Carbide computes one point at a time. In a test kernel on a 16 × 16 block it took 3,596 cycles, 7.8 times as long as the four statements.

Timestamps

u16 tsc[1] address(0x7b30);
u16 tsc_control[1] address(0x7b34);
u16 stamp[k=3];        // timestamps of the last sweep

These are the timestamp counter of lesson 3. Each sweep writes stamp[0] before the exchange, stamp[1] between the exchange and the update, and stamp[2] after the update, so the last sweep's values remain. cycles.py subtracts them modulo 65536 and prints the smallest and largest value over the PEs. In a row with T = 1, such as b16, stamp[0] comes right after the enable and reads 0, as t[0] does in lesson 9; the times are still right in the simulator.

The moat: moat.w

parallel {

    forever dispatch(b) {
        wavelet() {
            a[] <- 0:sp;
        }
    }

    forever dispatch(a) {
        wavelet() {
            b[] <- 0:sp;
        }
    }

    forever dispatch(d) {
        wavelet() {
            c[] <- 0:sp;
        }
    }

    forever dispatch(c) {
        wavelet() {
            d[] <- 0:sp;
        }
    }

}

A moat PE at the left of the grid receives the grid's left edge on b, from the right. For each wavelet, the first handler sends a 0 on a, which the routes carry right into ghost column m = 0 of the PE that sent it. The other three handlers do the same for the other sides. Each is the handler with no arguments from lesson 7, and forever re-arms it, so it answers any number of wavelets. Every moat PE has all four handlers, so one kernel serves the whole ring; each moat PE uses only the one that faces the grid.

The four dispatches are the arms of one parallel block. A forever dispatch never finishes, so in sequence only the first handler would ever be enabled. With parallel { replaced by a plain {, Carbide gives no message, the moats answer only on b, and the tutorial row hangs in its first exchange.

The layout: define, hstackrep and vstackrep

bench.paint builds the layout from three Paint functions of the size, and combines them:

union(beach(W+2 H+2) shore(W+2 H+2) blocks(W H))

define name(arguments) { ... } makes a Paint function. Its body runs when the function is called, and its last expression is the value it returns. The body can use the names the file defined before the function, such as the colors and the variables B and T from test.rc.

define beach(w h) {
    let r = rect(0 0 w h)
    r : paint(ap [CXR|R])
    r : paint(bp [CXL|L])
    r : paint(cp [CXU|U])
    r : paint(dp [CXD|D])
    r : paint(a  [L|LXC])
    r : paint(b  [R|RXC])
    r : paint(c  [D|DXC])
    r : paint(d  [U|UXC])
    r : colorpair(a ap)
    r : colorpair(b bp)
    r : colorpair(c cp)
    r : colorpair(d dp)
    r
}

This is lesson 8's colorpair, once per direction. ap [CXR|R] and a [L|LXC] carry what a PE sends on a to the right, where the next PE receives it on a. bp and b do the same to the left, cp and c up, and dp and d down. Eight colors carry the four directions. All eight rules are painted on one rect, because a | term acts at the edges of each rectangle a region is built from (lesson 8).

define shore(w h) {
    let s = union( rect(0   0   w 1)
                   rect(0   h-1 w 1)
                   rect(0   1   1 h-2)
                   rect(w-1 1   1 h-2) )
    s : code('moat [a:a b:b c:c d:d])
    s
}

The shore is the outer ring of the beach: the bottom row, the top row, and the left and right columns between them.

define blocks(w h) {
    let pe = tile()
    pe : code('block [a:a b:b c:c d:d B:B T:T])
    pe : hstackrep(w 0 'x)
    pe : vstackrep(h 0 'y)
    pe : place(1 1)
    pe : name('grid)
    pe
}

blocks makes one PE with the kernel and repeats it. pe : f(...) stands for pe = f(pe ...), so each line replaces pe with the region f returns.

• hstackrep(pe w 0 'x) puts w copies of pe side by side along x, with 0 PEs between them, and names the new dimension x.
• vstackrep(pe h 0 'y) repeats that row h times along y.
• place(1 1) moves the region by (1, 1), inside the shore.
• name('grid) gives the symbols their prefix: grid.u, grid.stamp.

Every copy has the same parameters, so all of them run one compiled block.w. With a one-PE unit, rect(1 1 W H 'x 'y) gives the same PEs. Replication matters when the repeated unit is several PEs with different code or routes: the heat-equation application repeats a 2 × 2 unit whose four PEs bind their sockets differently.

Loading the field: Makefile and ref.py

	miset image=$@ out=$@ map=bench.map symbol=grid.u order="y x n m _word" data=init-u.bin

The map lists grid.u with five dimensions (here for B=4 W=2 H=2):

grid.u                       [ y<2 x<2 n<6 m<6 _word<2 ] ( 0x0+12*n+2*m+_word 1+x 1+y )

y and x come from vstackrep and hstackrep, n and m from the array, and _word from the two 16-bit words of each sp value. The order string lists them slowest first, so init-u.bin holds the PEs row by row from the map's y = 0, and each PE's u row by row from n = 0. The map's x and y count the copies from 0: its y = 0 is the grid row at Paint y = 1, the first row above the shore. ref.py writes the file in that order, and the test reads grid.u back in the same order and compares the files with cmp.

What the hardware does

Each statement of the update is one instruction

The listing of the b8 row, test-b8/block_artefacts/block2.yeti, describes the operands of the update with two-dimensional DSDs, which the listing calls tensor:

xr u_xr_2 register(X0) = { .shape={0, 8, 8}, .strides={0, 10, 1}, };
tensor u_src0_3  = {.base=(u+2), .shape={0, 8, 8}, .length={1, 8, 8}, .save=0, .xr=0, .vp=0, .ua=0, };

The descriptor covers 8 rows of 8 values, rows 10 values apart: one row of u. It starts at word 2 of u, which is u[0][1], the first point the first statement reads. The four statements are four instructions:

	fadds [T2c] = [T08], [T17];  //@../block.w:39
	fadds [T2d] = [T09], [T18];  //@../block.w:40
	fadds [T2e] = [T0a], [T19];  //@../block.w:41
	fscales [T2f] = [T0b], R6;  //@../block.w:42

Each fadds adds 64 pairs of sp values. Both of its sources are in memory, so it does two per cycle, the memory limit of lesson 3. Carbide turns the multiplication by 0.25 into fscales, which adds R6 = −2 to each value's exponent; it reads one array and does four per cycle. A point costs three fadds at half a cycle and one fscales at a quarter: 1.75 cycles.

The exchange's eight arms are eight microthreads (.utid=0 to .utid=7 in the listing), and the column send ∀e a[] <- u[e][B] is one instruction over a vector with stride 10:

vector u_src1_1  = {.base=(u+36), .length=0, .stride=10, .save=0, .vp=0, .ua=0, .rchoasjt=0, };

Each microthread, when it finishes, activates or unblocks one of seven small tasks. A task runs only when it is both activated and unblocked. Each of the seven, when it runs, unblocks the next, and the last one runs the update. So the update starts only after all eight transfers are done.

The budget

Resource	block.w uses	A PE has	Note
Colors	8	24	four colorpairs
Input queues	4, each doubled by Carbide	8	queues 0, 2, 4 and 6
Output queues	4	8
General microthreads	8	8	one per arm of the exchange
Memory	2(B+2)² + 2B² words for u and v	24,576 words (48 KB), shared with the code

An input queue holds 4 or 8 wavelets: queues 0 and 1 hold 8, queues 2 to 7 hold 4. The hardware can let an even-numbered queue use the storage of the odd-numbered queue after it, when that queue is not in use, and so double in size. This is set for each pair separately: queue 0 can grow to 16 wavelets, and queues 2, 4 and 6 to 8.

The listing shows what Carbide chose for block.w: the four receive sockets on queues 0, 2, 4 and 6, all doubled.

iqueue IQ6 register(IQ6) = { .color=17, .filter=0, .doubling=1, };
iqueue IQ4 register(IQ4) = { .color=19, .filter=0, .doubling=1, };
iqueue IQ2 register(IQ2) = { .color=21, .filter=0, .doubling=1, };
iqueue IQ0 register(IQ0) = { .color=23, .filter=0, .doubling=1, };

The recv lines give each socket's queue: a has .qid=0 (16 wavelets), b .qid=2, c .qid=4 and d .qid=6 (8 each). Together they take the storage of all eight queues.

With a fifth receive socket, Carbide doubles none of them: a test kernel with five receives got queues 0 to 4, all undoubled. That is the compiler's choice; the hardware could still double some of the pairs. Lesson 12 has seven receive sockets and runs without doubling.

The parallel block uses all eight general microthreads; the core's four other microthreads run only umov copies. With a ninth transfer arm, Carbide runs all the arms in order on the main thread, and says so only in a warning in the make-log:

warning: a parallel block has 9 fabric arms and there are 8 microthreads, so its arms run in order

The exchange would then depend on buffering, as the first pitfall shows.

The largest block that fits is B = 76: u and v take 23,720 of the 24,576 words, and the code and the task table the rest. With B = 77 Paint stops with:

../block.w: error: [WMC TASK ALLOC] Task table starting at [0x5f18] with size [0x200] exceeds the main memory limit [0x6000].

Exchange against update, measured

Cycles of the last sweep:

B	Grid	Values sent per sweep	Exchange	Update	Update per point
4	2 × 2	16	94-97	41	2.6
8	3 × 2	32	107-113	128	2.0
16	2 × 2	64	131-136	463	1.8
76	1 × 1	304	377	10,178	1.8

A range is the smallest and largest value over the PEs. To run B = 76 yourself, add the row b76 B=76 W=1 H=1 T=1 CYCLES=20000 ARCH=sdr to test.rc: its update alone is longer than the Makefile's default cap of 5,000 cycles.

• The exchange costs about 70 cycles plus 4 cycles for each value of an edge. A PE sends four edges, and all four leave through the one port between its core and its router, which takes one wavelet per cycle (lesson 4). The four edges it receives come in through the same port, in the other direction. At B = 76, 304 of the 377 cycles are this per-value cost. The links between the routers are not the limit: each carries one edge. A test version that exchanged only the left and right edges took 77-80 cycles at B = 16 and 173-176 at B = 64, against 131-136 and 323-328 for all four: 2 cycles per value of an edge instead of 4.
• The update grows with B². Above B = 16 it settles at about 1.8 cycles per point: the 1.75 of its four instructions, plus their start.
• The two are close at B = 8. At B = 16 the update takes 3.4 times as long as the exchange, and at B = 76, 27 times as long.

The weak row runs the b8 block on 4 × 4 PEs instead of 3 × 2. Its exchange takes 102-113 cycles and its update 128, and its last memory write happens at cycle 521, as in b8. This is weak scaling: with the same block per PE, a grid of more PEs takes the same time per sweep. It holds because each PE talks only to its four neighbours.

Run it

In your copy of the tutorial (see Lesson 1):

cd tutorial/main-4-grids/11-grid-of-blocks
make TEST=tutorial
cat test-tutorial/make-log

make TEST='*' runs all four rows. A passing row prints one line, such as 11-grid-of-blocks : 0:03.44 PASS tutorial [ B=4 W=2 H=2 T=3 ARCH=sdr ]. Near its end, the make-log shows the line that cycles.py prints:

B=4, last sweep, over all PEs: exchange 94-97 cycles, update 41-41 cycles (2.56 per point)

The test does not check the cycle counts.

If u differs, the row prints FAIL, and the make-log shows where cmp found the first difference. With ∀e a[] <- u[e][1]; in place of ∀e a[] <- u[e][B];, each PE sends its left column to the right, and the tutorial row fails with out-ref.bin out-sim.bin differ: byte 38, line 1. Each PE's u takes (B+2)² × 4 = 144 bytes, so byte 38 is in the first PE of the file, the one at x = 1, y = 1.

After you edit a source file, run make clean before the next test.

Pitfalls

• An exchange in sequence works only while the buffers last. With the parallel block replaced by a plain { }, the eight transfers run in order on the main thread: four sends, then four receives. That passed at B = 4 and B = 16, and hung at B = 32: at the cycle cap no PE had finished its exchange, and cmp failed. At B = 4 the sequence was even faster, 44-49 cycles against 94-97, because it starts no microthreads and no tasks. It is still wrong: whether it finishes depends on how much the fabric buffers.
• Paint the | rules on one rectangle. A | term acts at the edges of each rectangle a region is built from, not at the outline of the region. Painted on a union, such as shore or the union of the grid and the shore, the eight rules would break where the rectangles meet, and Paint gives no message.
• Name the replicated dimensions. Without a name, hstackrep(pe w) calls its dimension x and vstackrep(pe h) calls its dimension y, but a count of 1 leaves the dimension out of the map. With H = 1, order="y x n m _word" then fails with order: 'y' is not a name or number. hstackrep(pe w 0 'x) and vstackrep(pe h 0 'y) keep the dimension at any count.
• The moat answers once per wavelet. With an hp field and hp replies, two values can share a wavelet (lesson 4), so the moat answers fewer values than the grid PE waits for. On a 2 × 2 grid with B = 4, no PE finished its first exchange.
• Keep the index convention in one place. n grows upward, with Paint's y. If ref.py put row 0 at the top, the test would fail although the kernel is right.
• One statement per vector operation. Combining the update into one statement made it 7.8 times slower (see The update).

Exercises

1. From the table, predict the exchange and the update for B = 32. Then add the row b32 B=32 W=2 H=2 T=1 ARCH=sdr to test.rc, run make TEST=b32, and compare.
2. Hold the top edge of the field at 1 instead of 0. Write moat_hot.w, a copy of moat.w whose dispatch(c) handler sends 1.0:sp on d. In bench.paint, give the top row of the shore its own code('moat_hot ...) and the rest of the ring moat. In ref.py, set the top ring row of g to 1 before the sweeps. Corner ghosts still stay 0.

Next

Lesson 12 — A complete iterative solver runs these sweeps until the whole grid has converged, which needs a sum over all PEs every sweep.

References

• Function definitions (define)
• Replication (hstackrep, vstackrep)
• Placement (place)
• Boundary operators
• Color pairs
• Color swapping
• The CE's queues
• Microthreads
• Parallel blocks
• Multi-level loops