eam — wafer-scale molecular dynamics with the EAM potential

A classical molecular-dynamics benchmark that integrates Newton's equations for a binary Cu/Zr metal under the embedded-atom method (EAM) potential, with every simulation cell mapped onto its own PE. There is no host tile: all 512 core PEs run the same program (eam-par-eng.w), hold their own atoms, and talk only to neighbours over hand-configured fabric colors. Initial coordinates come from a LAMMPS dump (cuzr-500.dump, 500 atoms in a 20.9 Å periodic cube), the tabulated potential comes from Fdata.bin / phidata.bin / rhodata.bin, and the run is checked by comparing the potential energy at every sample against a stored reference (data1.dat).

The configuration in test.rc runs 10 MD steps on a 32 × 16 grid of core PEs (plus a one-tile "moat" ring, 34 × 18 = 612 tiles total), which represents a logical 8 × 8 × 8 grid of simulation cells.

---

Algorithm

The EAM potential

EAM models a metal as pair repulsion plus a many-body "embedding" term that accounts for the cost of placing an atom into the local electron density:

E  =  Σ_i F_{t(i)}( ρ_i )  +  Σ_{i<j} φ_{t(i)t(j)}( r_ij )
ρ_i =  Σ_{j≠i} ρ_{t(i)t(j)}( r_ij )

The force on atom i from atom j is

f_ij = −[ ( F'(ρ_i) + F'(ρ_j) ) · ρ'(r_ij) + φ'(r_ij) ] · r̂_ij,
r̂_ij = (x_i − x_j) / r_ij

which is why a single force evaluation needs two communication passes: ρ_i has to be complete (and F'(ρ_i) known) before any pair force that involves atom i can be formed.

φ, ρ and F are evaluated from cubic-spline coefficient tables, not analytic formulas. rhotab[3][TABSIZE][4] and phitab[3][TABSIZE][4] are indexed by type pair (aa, ab, bb) and Ftab[2][TABSIZE][4] by type (a, b); each entry is four Horner coefficients for one interval. A table lookup is i = round(s·x − ½), f = s·x − i, then ((c0·f + c1)·f + c2)·f + c3. The scale factors are rscale = 1/(7.6e-4·400) for r and rhoscale = 1/(3e-2·400) for ρ, set in make-inputbin.py; with rcut = 7.6 Å that is exactly TABSIZE = 25 intervals.

The time loop

main in eam-par-eng.w is a two-level loop over floating-point counters (iter, iinner). Each inner step is split into a first half and a second half, which is the two-pass structure above:

for iter in [0, niter):
  for iinner in [0, ninner]:
    first half :  exchange positions   ->  neighbour list  ->  ρ_i  ->  reduce ρ
    second half:  exchange F'(ρ)       ->  pair forces     ->  reduce f
                  velocity-Verlet kick/drift,  atom migration
  every nsample iterations: computeeng(); recordtime()

With dtv = ½dt at iinner == 0 and iinner == ninner, dt in between, and the drift x += v·dt skipped at iinner == ninner, the inner loop is velocity Verlet: ninner MD steps per outer iteration costing ninner+1 force evaluations. make-inputbin.py sets niter = 10, ninner = 1, nsample = 1, so the benchmark is 10 MD steps with an energy sample after each.

Spatial decomposition and the NZX/NZY fold

The physical box is cut into nx × ny × nz cells (nz = nzx·nzy), one cell per PE, sized so that RANGE cells cover the cutoff (8 × 8 × 8 cells of 2.61 Å, RANGE = 3, rcut = 7.6 Å). The wafer is two-dimensional, so the third dimension is folded into the first two: the PE grid is nxtot × nytot = (nx·nzx) × (ny·nzy), and a block of nzx × nzy adjacent PEs holds the nz cells of one logical (x, y) column.

Two index transforms (duplicated in C++ in 3dmd_route3sdr.cc and in Python in layoutindex.py) define the mapping:

• jtox(nx, j) lays the periodic line of nx logical cells onto nx physical positions as a serpentine: physical 0,1,2,3,… → logical 0, nx−1, 1, nx−2, …. Logical neighbours are then never more than two physical blocks apart and the wrap-around 0↔nx−1 is a single short hop — no long wrap wire is needed. hix = (pix/nzx) % 2 says which leg of the fold a PE is on.
• chain_jtox(nx, nzx, pix) composes that with the z-fold and returns the chain index lix = nzx·lix0 + izx, the PE's position along the 1-D communication chain that runs through a physical row (the Python version returns the tuple (lix0, izx, hix) instead). A step of +1 in logical x is a step of +nzx along the chain.

For nx = 4, nzx = 2 (physical row of 8 PEs):

physical PE  :  0   1  |  2   3  |  4   5  |  6   7
block / hix  :   0 / 0 |   1 / 1 |   2 / 0 |   3 / 1
logical x    :     0   |     3   |     1   |     2
chain idx lix:  0   1  |  7   6  |  2   3  |  5   4

The hix = 0 PEs carry chain indices 0, 1, 2, 3 left-to-right; the hix = 1 PEs carry 4, 5, 6, 7 right-to-left. The two legs are joined at the array edges by the moat (below), closing the periodic ring.

The z direction is handled differently: makeZloop/makeZchain build a Hamiltonian cycle through the nzx × nzy sub-block, and zidx<nzx>x<nzy>.dat records each PE's position on that cycle — its logical z index. For 4 × 2 the cycle is the obvious boustrophedon:

zidx4x2.dat:   0  7  6  5
               1  2  3  4

A +1 step in z is one hop along this cycle. make-inputbin.py and check-exchange.py read the same file so the host-side atom distribution agrees with the fabric routing.

Halo exchange: rods, not a cube

A naive halo would need (2·RANGE+1)³ − 1 neighbour cells. This code uses a neutral-territory / "platelet" decomposition that only needs three one-dimensional rods. Received atoms are packed into one pos array:

 index                      contents                         used as
 ------------------------------------------------------------------------
 0        .. nlocal         owned atoms                       i (integrate)
 ymptr    .. +nym           halo from −y  (RANGE cells)   ┐
 y0ptr    .. +nlocal        copy of the owned atoms       ├─  y-rod
 ypptr    .. +nyp           halo from +y  (RANGE cells)   ┘
 xmptr(=yend) .. +nxm       halo from −x                  ┐
 x0ptr    .. +nlocal        copy of the owned atoms       ├─  x-rod
 xpptr    .. +nxp           halo from +x                  ┘
 zpptr(=xend) .. +nzp       RANGE whole x-rods, from +z
 zend

detect_platelet() then enumerates exactly two pair loops:

 (A)   i ∈ y-rod      ×  j ∈ ( +x halo  ∪  +z block )
 (B)   i ∈ owned copy ×  j ∈ ( later owned atoms  ∪  +y halo )

Writing d for the cell-to-cell offset of a pair, (A) covers every d with dz > 0, plus every d with dz = 0, dx > 0; (B) covers dz = dx = 0, dy ≥ 0. That is precisely the lexicographic half-space, so every cell pair inside the cutoff is owned by exactly one PE and no pair is double-counted. The +z block is a copy of the whole x-rod of each of the next RANGE cells in +z, which is why the −x halo is received even though it never appears as a j index locally: it has to be forwarded along the z chain.

Forces (and ρ) computed on halo copies are sent back by reduce(), which runs the same chains in reverse with accumulating receives (<+-): first +z (folding the z block back into the local x-rod), then +x and −x (folding the x-rod into the x0 copy), then +y and −y. main finishes by adding ff1[x0ptr+j] and ff1[y0ptr+j] into ff1[j].

Fabric topology and routing

Ten colors are declared in md-eam.paint and bound to six sockets:

socket	color(s)	direction
s_xm	0 / 1	chain flowing toward logical −x
s_xp	2 / 3	chain flowing toward logical +x
s_ym	4 / 5	chain flowing toward logical −y
s_yp	6 / 7	chain flowing toward logical +y
s_zm	8	Hamiltonian z cycle, −z direction
s_zp	9	Hamiltonian z cycle, +z direction

Each x/y direction uses a pair of colors, one per leg of the serpentine fold. A PE taps the color matching its own hix/hiy parity and acts as a transparent wire for the partner color. The code() arguments bind the even color of each pair; at startup main rewrites its own input- and output-queue color fields (set_iq_color / set_oq_color) to the odd color when params.hix / params.hiy is set, so the same binary works on both legs.

Note what md-eam.paint does not do: all ten paint(... [C>C ; tip]) lines are commented out. No route is painted at all. The entire router and switch configuration is generated offline by 3dmd_route3sdr and patched into the image with miset image=… symbol=router / symbol=switch. bench.map declares these two symbols for every tile: router at 0x7c00+r and switch at 0x7c80+4*r+n, for the 24 colors r.

The moat ring exists to turn the fold around. For colors 0–7 its tiles get reflecting routers that use the hardware color swap. main() in 3dmd_route3sdr.cc picks them by position, with ix/iy running from −1 (the moat column or row at Paint x = 0 or y = 0) to nxtot/nytot:

if(ix < 0) moat = RtoRswap;
else if(ix >= nxtot) moat = LtoLswap;
else if(iy < 0) moat = UtoUswap;
else if(iy >= nytot) moat = DtoDswap;

RtoRswap takes a color in from the right (+x) neighbour and sends it back there with the horizontal swap, and LtoLswap does the same on the left. The y = 0 row reflects toward +Y (UtoUswap, the generator's up, Paint U) and the y = nytot+1 row toward −Y (DtoDswap, Paint D), with the vertical swap. A wavelet arriving on color c is sent straight back on color c ^ 1 — i.e. the turn-around and the leg-change happen in one tile, which is exactly what the serpentine needs.

Within a chain, data is moved with a marching switch. The generator gives each PE's switch four positions and an initial position derived from its chain index; 3dmd_route3sdr.cc documents the protocol in its own words:

//   if chainlen > 1
//      us->ds , c->ds , c->c , us->c
//      to advance, send (adv,adv,rst) , (adv) after data
//   else
//      c->ds , c->c , us->c
//      to advance, send (adv,rst) , (adv) after data

(us = upstream, ds = downstream, c = core.) chainlen is nzx/nzy for the x/y chains — one logical cell is nzx hops — and 1 for the z cycle. The two advance sequences are exactly the control wavelets in the Tungsten source:

//                 (adv*,adv*,rst*), (adv*)
wv xycmds[n=2] = { {0x0001, 0x2b40}, {0x0000, 0x0140} };
//                 (adv*,rst*)     , (adv*)
wv zcmds[n=2]  = { {0x0000, 0x0d40}, {0x0000, 0x0140} };
xp const ncmd = 2;
// * = no-CE bit set.

Every data burst in getcount(), exchange(), reduce() and the migration code is followed by the two-wavelet sequence for its chain.

Finally, the generator sets the teardown bit (bit 15) in the router configuration of colors 0–7 on every core tile; the z colors 8 and 9 never get it. While the bit is set the color's inputs are flow-controlled, so clear_teardown() in main clears it for colors 0–7 before the first exchange.

Phases of one timestep

FIRST HALF                                   SECOND HALF
  getcount()      exchange nlocal              pos1save <- pos1
                  along x±,y± (RANGE hops)     computeFp()        F'(ρ_i)
                  then the x-rod total          copy F' into x0, y0 slots
                  along z (RANGE hops)
  build ymptr..zend pointers,                  exchange()         1 word/atom
  firstblock/countblock, ×1 ×3 ×4 variants
  copy owned atoms into x0, y0
  exchange()      4 words/atom (x,y,z,id)      eamforce()         φ' and the
  remap_positions()  type from sign(id),                          EAM pair force
                  recentre + minimum image     restore pos1
  detect_platelet()  neighbour list            reduce()           3 words/atom
  rsqrt_rinv()    1/r by Newton on 0x5f37      add x0,y0 forces into local
  zerorho(); computerho()   ρ_ij into ff1,     velocity Verlet kick (+ drift)
                  dρ/dr kept in rhopij         atom migration in x, y, z
  reduce()        1 word/atom
  add x0,y0 contributions into local ρ

Atom migration (rebuildnext) runs after every second half. For each of the three directions it partitions the owned atoms into stays / goes − / goes + by comparing the minimum-image offset from the cell centre against cell[k][0..1], exchanges counts and then 8 words per migrating atom (4 position

• 3 velocity, padded) with the two adjacent cells, and appends the arrivals.

---

Source files

File	Role
bench.paint	what paint is invoked on: a copy of md-eam.paint
md-eam.paint	floorplan: colors, core rect, moat ring, code() placement
eam-par-eng.w	the whole MD program run by every core PE (~1260 lines)
swapmoat.w	program placed on the moat tiles; it only copies eight words into pcopy
Makefile	PAINTVARS, route/input generation, miset/miget wiring
test.rc	test configurations (see Building and running)
3dmd_route3sdr.cc	C++ generator for router/switch configuration + z-index map
3dmd_route3sdr	the generator, compiled from the .cc; the Makefile runs it and has no rule to rebuild it. The experimental SDK ships it compiled from the .cc
make-inputbin.py	builds params.bin, fparams.bin, posin.bin, velin.bin
layoutindex.py	Python copy of jtox/xtoj/chain_jtox/chain_xtoj; its chain_jtox returns a tuple
readdump.py	LAMMPS dump reader
writedump.py	LAMMPS dump writer
check-exchange.py	oracle for the halo counts/offsets the device computed
ecalc.py	correctness oracle: potential energy vs data1.dat; prints the sample timestamps
posbin2dump.py	writes per-PE atom data as a LAMMPS dump (manual; it does not take this app's outputs, see below)
cuzr-500.dump	input configuration: 500 Cu/Zr atoms, 20.9 Å periodic cube
data1.dat	reference energies in eV, 10 rows of PE KE TE; ecalc.py checks the PE column
Fdata.bin	Ftab[2][25][4] — embedding-function spline coefficients, Cu and Zr
phidata.bin	phitab[3][25][4] — pair-potential coefficients, aa/ab/bb
rhodata.bin	rhotab[3][25][4] — electron-density coefficients, aa/ab/bb
Fpdata.bin, phipdata.bin	[2][25][3] and [3][25][3] coefficient tables beside Fdata.bin and phidata.bin; nothing reads them

eam-par-eng.w

Memory layout first (~lines 1–125): overlapping address(&…) declarations give several views of the same storage — pos / pos0 / pos1, ff / ff1, countblock with xmcount…zpcount as windows into it, Ftab / Ftab1. pos0, vel1 and ff1 are extern; pos0, the local-atom prefix of pos, is used only by miset/miget. varstart / varend / velend / tabend are one-word markers for reading memory usage off bench.map.

Communication:

• getcount() — the all-gather of nlocal over the neighbourhood. For r = 0 … RANGE−1 it ships xmcount[r] out on s_xp and receives into xmcount[r+1] (what arrives came from the −x side), and symmetrically for xp/ym/yp. Then it forms the x-rod total nlocal + Σnxm + Σnxp and shifts that along z to fill zpcount[].
• exchange() — the position / F' halo exchange. Each round is a parallel { receive block r+1 ; { send block r ; control } }, driven entirely by firstblock[]/countblock[], so the same routine moves 4 words/atom (positions), 1 word/atom (F') or anything else simply by loading the right first4/count4, first1/count1 table.
• reduce() — the mirror image, walking r downward with accumulating receives ff1[…] <+- s_zp[]. Each direction is behind a debug flag (zpflag…ymflag, all 1 in a real run).

Physics:

• remap_positions() — recovers the atom type from sign(pos[i][3]) (the LAMMPS id is sign-encoded by make-inputbin.py) and re-centres every coordinate on the cell centre with a round-to-nearest trick (copysign(2²³, x), add, subtract). After this the minimum image is already applied, so no PBC test is needed in any distance loop.
• detect_platelet() — builds the neighbour list for the two pair loops above, filling ilist, jlist, klist (type-pair index × rhotablen), dxij, dist2, and setting nint. It also converts atomtype[y0ptr+…] from ±1 to 0/1 and builds klocal (type × Fptablen). On overflow it does ptr -= 1, i.e. it overwrites the last entry rather than failing, and nothing checks for it (see Parameters).
• rsqrt_rinv() — vectorised 1/√dist2 over the whole neighbour list: the 0x5f37 bit trick on the high half-word (id2 aliases rinv) followed by three Newton iterations.
• computerho() — table-interpolates ρ(r_ij), accumulates it into both ff1[i] and ff1[j] (ff1 doubles as the ρ accumulator in the first half), and leaves dρ/dr in rhopij and the table index/fraction in itmp/frac for later reuse.
• computeFp() — F'(ρ_i) for owned atoms, written into pos1[0:nlocal] (hence the pos1save copy in main). Negative ρ is clamped to zero.
• eamforce() — φ'(r_ij) from phitab, Fpsum = (F'_i + F'_j)·ρ'_ij + φ', scaled by 1/r, then ±dxij·Fpsum accumulated into ff[i] and ff[j].
• computeeng() — kinetic energy plus the full cubic evaluations of F(ρ_i) and φ(r_ij) into peval/keval. Because detect_platelet enumerates each pair once, summing peval over all PEs gives the total potential energy in eV. keval = ½Σ m·v² uses the masses in amu and v in Å/ps, so it is in amu·Å²/ps². ecalc.py does not check it: keflag = False turns its kinetic-energy code off.
• timerinit() / recordtime() — timerinit() sets the timestamp counter to count every cycle (control register 0x7b34); recordtime() reads the 48-bit counter at 0x7b30 into tstamp[128][3] and snapshots peval/keval into evec/kevec.

main also maintains a diagnostics array ptrs[22]: nlocal, the five halo counts, all the buffer pointers, nint, the two candidate-pair counters icount1/icount2, the number of force evaluations (ptrs[17]), and high-water marks for nlocal, zend and nint (ptrs[18..20]). Nothing checks those against NLOCMAX/NMAX/NINTMAX.

swapmoat.w

Its main copies the first eight words of router_config_regs (declared at address(0x7f20)) into pcopy and ends; it sends and receives nothing. The moat tiles only route: the Makefile writes their router and switch configuration into the image with the rest of the fabric (miset symbol=router, symbol=switch).

md-eam.paint

51 lines, the last 12 of them the commented-out paint block. It declares colors 0–9, reads the nine paint variables, defines corerect = rect(1 1 nxtot nytot) named core running eam-par-eng with the six socket→color bindings and the five scalar parameters, and a moat built as the union of four rectangles forming a one-tile ring around the core, running swapmoat. The name() calls are what make symbols appear as core.pos0 etc. in bench.map.

3dmd_route3sdr.cc

Run as 3dmd_route3sdr NX NY NZX NZY. It emits three files:

• zidx<nzx>x<nzy>.dat — the logical z index of each PE in the nzx × nzy sub-block (the Hamiltonian cycle order), consumed by make-inputbin.py and check-exchange.py.
• rconfig<nxtot>x<nytot>.bin — (nytot+2)·(nxtot+2)·24 u16 router-config words.
• swconfig<nxtot>x<nytot>.bin — the same shape × 4 switch words.

The interesting functions are layout2Dpbc() (per-PE, per-color x/y chain routers plus the fwd pass-through for the partner color; computes hix/hiy, the chain index, and the initial switch position from chainidx % (nz{x,y}+1)), makeZloop() (recursive construction of a Hamiltonian cycle on an na × nb grid, transposing when na is odd), makeZchain() (turns that cycle into a single-hop marching chain and records zidx), lookup() (fetches a z-chain router for a PE and applies xflip/yflip so the pattern is mirrored correctly on the folded legs), and routerreg::set_router/set_switch (the actual bit packing). main() fills the whole (nytot+2) × (nxtot+2) array, putting reflecting color-swap routers on the moat ring and setting the teardown bit on core tiles, both for colors 0–7.

The file also carries a fair amount of dead scaffolding — routerregx, the SW/SWK/SWANY/CMD macros, encode_command2, colselect, router::print2 and the unused zloop variable are not reachable from main().

---

Parameters

Set per row of test.rc; the Makefile forwards them with PAINTVARS = NLOCMAX NMAX NINTMAX TABSIZE NX NY NZX NZY RANGE, which flute.mk turns into set=NAME=value on the paint command line, where they arrive as bare globals in md-eam.paint. ARCH is not a paint variable; it goes into PAINTFLAGS and PIPELINEFLAGS.

Variable	t512	Meaning
NLOCMAX	20	capacity of the owned-atom arrays (vel[NLOCMAX][3], pos0, Fitmp, Ffrac, klocal). Must be at least the largest nlocal any cell reaches, including after migration.
NMAX	80	capacity of pos[NMAX][4] / ff[NMAX][3] / atomtype / pos1save, i.e. owned atoms plus all three halo rods plus the duplicated copies — the zend high-water mark. Atom migration also stages atoms in pos: in each direction it needs 3·nlocal + 2·(atoms sent toward + plus atoms received) ≤ NMAX, which nothing records or checks. Keep NMAX ≤ NINTMAX: the NINTMAX-sized scratch tmp1/tmp2 is indexed up to zend.
NINTMAX	200	capacity of the neighbour list (ilist, jlist, klist, dxij, dist2, rinv, coeff, …). detect_platelet clamps nint at NINTMAX − 1 by overwriting the last entry, and builds its candidate distances in dxij up to index nint + zend + 3, so it needs nint ≤ NINTMAX − 2 and nint + zend + 4 ≤ NINTMAX.
TABSIZE	25	number of spline intervals in each potential table; must match the shapes of Fdata.bin, phidata.bin, rhodata.bin and the rscale/rhoscale constants in make-inputbin.py.
NX, NY	8, 8	logical simulation cells in x and y.
NZX, NZY	4, 2	the z fold: nz = NZX·NZY logical z cells, carried by an NZX × NZY block of PEs. The PE grid is (NX·NZX) × (NY·NZY) = 32 × 16 here.
RANGE	3	communication range in cells. Must satisfy RANGE · cell_size ≥ rcut; here 3 × 2.61 Å ≥ 7.6 Å. Sets countsize = firstsize = RANGE+1 and the number of hops in every chain loop. check-exchange.py assumes 3.
ARCH	sdr	sdr = WSE-3. Only sdr is exercised; several hardware register addresses in eam-par-eng.w are hard-coded // SDR constants.

Nothing checks these bounds. ptrs[18..20] in ptrs.bin record the high-water marks of nlocal, zend and nint to compare against them.

Physical parameters are not in test.rc — dt = 1 fs, rcut = 7.6 Å, niter = 10, ninner = 1, nsample = 1, the masses (in amu) and the input dump file are all literals in make-inputbin.py. ecalc.py hard-codes 10 samples (nframes = 10) to match the 10 rows of data1.dat, so changing the step count also means editing nframes and regenerating data1.dat.

---

Building and running

The experimental SDK installs this example beside the tutorial, as projects/eam, with the route generator already compiled. Run every row of its test.rc with tungsten-tutorials -f projects/eam, or work by hand in a copy of the directory, set up as for the tutorials, where <sdk> is the SDK's directory:

cp -r <sdk>/share/tutorial/projects/eam . && cd eam
make TEST=t512                 # build, simulate and check
make TEST=pr_smoke             # same row (alternate tag)
make TEST='*'                  # every row of test.rc
make TEST=t512 GOAL=bench.img  # build the image only, no simulation
make TEST=t512 GOAL=debugx     # also trace, and unpack the debug viewer's files
make clean                     # remove all test-*/ directories

test.rc has 1 row. That row, tagged t512 and pr_smoke, sets the values in the t512 column of Parameters. MODE (default sim) only names the outputs (out-$(MODE).bin); sim is the only mode the Makefile provides.

make reports PASS or FAIL per row. Everything else lands in test-<tag>/: make-log (the real transcript), bench.map, bench.obj/bench.img, bench-sim.log, bench-sim.core, sim_stats.json, and the assembly listings eam-par-eng_artefacts/eam-par-eng1.yeti and swapmoat_artefacts/swapmoat2.yeti.

What the build does, in order:

1. paint compiles bench.paint → bench.obj + bench.map.
2. generate_bins runs the route generator as ../3dmd_route3sdr NX NY NZX NZY (router/switch configs + zidx*.dat) — the Makefile does not compile 3dmd_route3sdr.cc, so after changing it, rebuild 3dmd_route3sdr yourself in the project directory with g++ -O2 -o 3dmd_route3sdr 3dmd_route3sdr.cc (a copy of the project that has only the .cc needs this before its first run) — then runs ../make-inputbin.py (params.bin, fparams.bin, posin.bin, velin.bin, idmap.npy).
3. bench.img is bench.obj patched by a series of miset calls: the fabric configuration (symbol=router, symbol=switch), then core.params, core.fparams, core.pos0, core.vel1, then the three potential tables. The order= string is scraped out of bench.map; for the tables the x and y dimensions are stripped (tr -d xy) so one copy of the table is broadcast to every PE.
4. $(SIMFABRIC) runs the image on the Cerebras SDK simulator and writes bench-sim.core.
5. out-sim.bin extracts results with miget — posout.bin and velout.bin from pos0 and vel1, then one .bin file each for ptrs, firstblock, countblock, rinv, dist2, tstamp, ilist, jlist, klist, dxij, ff1, itmp, Fitmp, frac, Ffrac, coeff, rhopij, phipij, evec and kevec.

Correctness is checked two ways:

• check-exchange.py re-derives the expected halo counts and buffer offsets purely from the per-PE nlocal values — periodic np.rolls over the logical 3-D grid, un-folded from the physical 2-D layout with remap23() — and compares them element-for-element against the device's countblock, firstblock and ptrs, and fails on any mismatch. It hard-codes the range as 3 (commrange), so firstblock.bin and countblock.bin must hold 5·4 words per PE. This validates the addressing and communication bookkeeping; it checks neither the buffer sizes nor the physics.
• ecalc.py sums evec over all PEs for each of the 10 samples after the first and compares that potential energy against the first column of data1.dat. It fails if the relative error max|d| / (1 + max|ref|) is non-finite or exceeds 1e-5. With keflag = False it reads neither kevec nor the final velocities, so kinetic and total energy are not checked. It also prints PE 0's timestamps from tstamp, in cycles since the first sample.

A passing run ends make-log with errflag = 0 from check-exchange.py and Ok! from ecalc.py.

posbin2dump.py NLOCMAX fparams.bin PARAMS POS VEL FORCES out.dump writes per-PE atoms as one id-sorted LAMMPS dump of id, type, position, velocity and force. The Makefile does not run it, and it does not take this app's outputs as they are: it reads each PE's atom count from the nlocal field of 13-word paramstruct records (the layout of params.bin, which holds the counts from before the run), not from ptrs.bin; it needs FORCES to hold NLOCMAX atoms per PE, while ff1.bin holds NMAX; and it takes the box from fparams.bin[8:11], where this app stores invmass[0], invmass[1] and box[0].

References

• Parameters
• Sockets
• Struct declarations
• Hardware register access (address())
• Parallel blocks
• Built-in functions (control())
• Paint code attachment (code())
• Paint naming (name())

Run it

The experimental SDK installs this project beside the tutorial, as projects/eam. Run it like a tutorial: tungsten-tutorials -f projects/eam.