# -*- coding: utf-8 -*- """lane64 の GPU 実行台（limb 版）。thread 一本が lane 一本になる。 仕様は docs/common-compiler/arch/20。式 (2.1) の解釈器と同じ意味を、CUDA の kernel の中で lane ごとに走らせる。cuvcpu.py は語を bit 平面に持ち、命令一つが平面の演算の列になるので、 語が広いと launch の数で足が止まる。ここでは語を 64 bit の limb の並びで持ち、lane ごとの thread が block を跨いで命令を辿る。分岐は thread ごとに好きに飛ぶ（mask は要らない）。 並び 語 主。語 w の limb i、lane l は M[(w·NLIMB + i)·lanes + l] に置く。 同じ語を読む lane の番地が連続するので、warp の memory access が一度で済む。 lane 主に置くと lane ごとに words·NLIMB 語ずつ離れ、一命令が 32 回の転送になる 語 NLIMB 本の limb。二の補数、little endian。既定 16 本 = 1024 bit landmark R[0..7]。commonc の lane64 版は R1 R2 R3 だけを使う（back_lane64.cam の覆い） 作業領域 lane ごとに words 語。読む語（watch）だけを host へ戻す 語の幅を越えた乗算は報告 43 で止まる。lane64 の語は ℤ なので（20 §0.3）、この道は NLIMB·64 bit を越えた時点で core 意味論と離れる。 run(prog, seeds, words, watch, need=None) → lane ごとの (報告, {語: 値}) """ import numpy as np INSN_FIELDS = 7 NREG_BYTES = 8 * 8 # shared memory へ置く register。NREG 本 × 8 byte R_OK, R_KO, R_WIDE, R_BUDGET = 0, 42, 43, 44 SRC = r""" #define NLIMB %(nlimb)d #define NREG 8 typedef unsigned long long u64; typedef long long i64; struct Big { u64 v[NLIMB]; }; __device__ __forceinline__ void big_zero(Big &r) { for (int i = 0; i < NLIMB; i++) r.v[i] = 0; } __device__ __forceinline__ void big_imm(Big &r, i64 x) { r.v[0] = (u64)x; u64 f = (x < 0) ? ~0ULL : 0ULL; for (int i = 1; i < NLIMB; i++) r.v[i] = f; } __device__ __forceinline__ bool big_neg_p(const Big &a) { return (a.v[NLIMB - 1] >> 63) != 0; } __device__ __forceinline__ void big_add(Big &r, const Big &a, const Big &b) { u64 c = 0; for (int i = 0; i < NLIMB; i++) { u64 s = a.v[i] + b.v[i]; u64 c1 = (s < a.v[i]) ? 1 : 0; u64 t = s + c; u64 c2 = (t < s) ? 1 : 0; r.v[i] = t; c = c1 | c2; } } __device__ __forceinline__ void big_not(Big &r, const Big &a) { for (int i = 0; i < NLIMB; i++) r.v[i] = ~a.v[i]; } __device__ __forceinline__ void big_negate(Big &r, const Big &a) { Big n; big_not(n, a); Big one; big_imm(one, 1); big_add(r, n, one); } __device__ __forceinline__ void big_sub(Big &r, const Big &a, const Big &b) { Big nb; big_negate(nb, b); big_add(r, a, nb); } __device__ __forceinline__ bool big_eq(const Big &a, const Big &b) { for (int i = 0; i < NLIMB; i++) if (a.v[i] != b.v[i]) return false; return true; } __device__ __forceinline__ bool big_lt(const Big &a, const Big &b) { bool na = big_neg_p(a), nb = big_neg_p(b); if (na != nb) return na; for (int i = NLIMB - 1; i >= 0; i--) { if (a.v[i] != b.v[i]) return a.v[i] < b.v[i]; } return false; } __device__ __forceinline__ void big_shl(Big &r, const Big &a, int s) { if (s <= 0) { r = a; return; } int w = s / 64, bts = s %% 64; Big t; big_zero(t); for (int i = NLIMB - 1; i >= w; i--) { u64 hi = a.v[i - w] << bts; u64 lo = (bts && i - w - 1 >= 0) ? (a.v[i - w - 1] >> (64 - bts)) : 0; t.v[i] = hi | lo; } r = t; } // |a|·|b| を 2·NLIMB limb で作り、上半分が零でなければ幅を越えたとみなす __device__ bool big_mul(Big &r, const Big &a0, const Big &b0) { Big a, b; bool na = big_neg_p(a0), nb = big_neg_p(b0); if (na) big_negate(a, a0); else a = a0; if (nb) big_negate(b, b0); else b = b0; // 上位の零の limb は掛けない。小さい数どうしなら 1 回の乗算で済む int la = 0, lb = 0; for (int i = 0; i < NLIMB; i++) { if (a.v[i]) la = i + 1; if (b.v[i]) lb = i + 1; } u64 acc[2 * NLIMB]; for (int i = 0; i < 2 * NLIMB; i++) acc[i] = 0; for (int i = 0; i < la; i++) { u64 carry = 0; for (int j = 0; j < lb; j++) { u64 lo = a.v[i] * b.v[j]; u64 hi = __umul64hi(a.v[i], b.v[j]); u64 s = acc[i + j] + lo; u64 c1 = (s < lo) ? 1 : 0; u64 t = s + carry; u64 c2 = (t < s) ? 1 : 0; acc[i + j] = t; carry = hi + c1 + c2; } int k = i + lb; while (carry && k < 2 * NLIMB) { u64 s = acc[k] + carry; carry = (s < carry) ? 1 : 0; acc[k] = s; k++; } } bool wide = false; for (int i = NLIMB; i < 2 * NLIMB; i++) if (acc[i]) wide = true; if (acc[NLIMB - 1] >> 63) wide = true; for (int i = 0; i < NLIMB; i++) r.v[i] = acc[i]; if (na != nb) { Big t = r; big_negate(r, t); } return wide; } // 語の番地 → 詰めた並びの位置。work が 0 なら詰めない __device__ __forceinline__ i64 slot(i64 a, int need, int work, int words) { if (work == 0) return (a >= 0 && a < words) ? a : -1; if (a >= 0 && a < need) return a; if (a >= work && a - work + need < words) return a - work + need; return -1; } extern "C" __global__ void lane64_run(const i64 *prog, int ninsn, u64 *M, int words, int lanes, i64 budget, int need, int work, int *report, i64 *steps, const u64 *ROM, i64 romlo, i64 romn) { int lane = blockIdx.x * blockDim.x + threadIdx.x; if (lane >= lanes) return; Big R[NREG]; int pc = 0; i64 n = 0; int rep = %(r_budget)d; while (n < budget) { n++; if (pc < 0 || pc >= ninsn) { rep = %(r_end)d; break; } const i64 *ins = prog + (size_t)pc * %(fields)d; int op = (int)ins[0], a = (int)ins[1], b = (int)ins[2], c = (int)ins[3]; int sb = (int)ins[4], sc = (int)ins[5]; i64 im = ins[6]; if (op == 0) { rep = %(r_end)d; break; } // HALT else if (op == 1) { big_imm(R[a], im); pc++; } // MAT else if (op == 2) { Big x, y; big_shl(x, R[b], sb); big_shl(y, R[c], sc); big_add(R[a], x, y); pc++; } // ADD else if (op == 3) { Big t; if (big_mul(t, R[b], R[c])) { rep = %(r_wide)d; break; } R[a] = t; pc++; } // MUL else if (op == 4) { Big t; big_negate(t, R[b]); R[a] = t; pc++; } // NEG else if (op == 5) { Big x, y; big_shl(x, R[b], sb); big_shl(y, R[c], sc); pc += big_eq(x, y) ? (int)im : 1; } // EQ else if (op == 6) { Big x, y; big_shl(x, R[b], sb); big_shl(y, R[c], sc); pc += big_lt(x, y) ? (int)im : 1; } // LT else if (op == 7) { pc += (int)im; } // JMP else if (op == 8) { // LD i64 raw = (i64)R[b].v[0]; if (raw >= romlo) { // 共有帯。lane で割らない i64 o = raw - romlo; if (o >= romn) { rep = %(r_end)d; break; } const u64 *w = ROM + (size_t)o * NLIMB; for (int i = 0; i < NLIMB; i++) R[a].v[i] = w[i]; pc++; } else { i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; break; } const u64 *w = M + ((size_t)addr * NLIMB) * (size_t)lanes + lane; for (int i = 0; i < NLIMB; i++) R[a].v[i] = w[(size_t)i * lanes]; pc++; } } else if (op == 9) { // ST i64 raw = (i64)R[b].v[0]; if (raw >= romlo) { rep = %(r_end)d; break; } // 共有帯へは書けない i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; break; } u64 *w = M + ((size_t)addr * NLIMB) * (size_t)lanes + lane; for (int i = 0; i < NLIMB; i++) w[(size_t)i * lanes] = R[c].v[i]; pc++; } else if (op == 10) { rep = %(r_ok)d; break; } // OK else if (op == 11) { rep = %(r_ko)d; break; } // KO else if (op == 12) { pc++; } // PUT（GPU では捨てる） else { rep = %(r_end)d; break; } } report[lane] = rep; steps[lane] = n; } // ===== 語が 64 bit で足りるときの道 ===== // // Big R[NREG] は添字が命令から来るので、compiler が register に置けず local memory へ落ちる。 // 一命令ごとに local への往復が入り、thread 一本の歩みが桁で遅くなる。 // 語を 64 bit に決め打てる場合は、register を shared memory の // // R[k] = sreg[k * blockDim.x + threadIdx.x] // // へ置く。stride が blockDim.x なので bank が衝突せず、添字が動いても速い。 // 積が 64 bit を越えたら幅の報告で止める（意味は Big の道と同じ）。 extern "C" __global__ void lane64_run1(const i64 *prog, int ninsn, u64 *M, int words, int lanes, i64 budget, int need, int work, int *report, i64 *steps, const u64 *ROM, i64 romlo, i64 romn) { extern __shared__ i64 sreg[]; int lane = blockIdx.x * blockDim.x + threadIdx.x; if (lane >= lanes) return; i64 *R = sreg + threadIdx.x; const int rs = blockDim.x; int pc = 0; i64 n = 0; int rep = %(r_budget)d; while (n < budget) { n++; if (pc < 0 || pc >= ninsn) { rep = %(r_end)d; break; } const i64 *ins = prog + (size_t)pc * %(fields)d; int op = (int)ins[0], a = (int)ins[1], b = (int)ins[2], c = (int)ins[3]; int sb = (int)ins[4], sc = (int)ins[5]; i64 im = ins[6]; if (op == 0) { rep = %(r_end)d; break; } else if (op == 1) { R[a * rs] = im; pc++; } else if (op == 2) { R[a * rs] = (R[b * rs] << sb) + (R[c * rs] << sc); pc++; } else if (op == 3) { i64 x = R[b * rs], y = R[c * rs]; u64 ux = (u64)(x < 0 ? -x : x), uy = (u64)(y < 0 ? -y : y); if (__umul64hi(ux, uy) || ((ux * uy) >> 63)) { rep = %(r_wide)d; break; } R[a * rs] = x * y; pc++; } else if (op == 4) { R[a * rs] = -R[b * rs]; pc++; } else if (op == 5) { pc += ((R[b * rs] << sb) == (R[c * rs] << sc)) ? (int)im : 1; } else if (op == 6) { pc += ((R[b * rs] << sb) < (R[c * rs] << sc)) ? (int)im : 1; } else if (op == 7) { pc += (int)im; } else if (op == 8) { i64 raw = R[b * rs]; if (raw >= romlo) { // 共有帯。lane で割らない i64 o = raw - romlo; if (o >= romn) { rep = %(r_end)d; break; } R[a * rs] = (i64)ROM[(size_t)o]; pc++; } else { i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; break; } R[a * rs] = (i64)M[(size_t)addr * (size_t)lanes + lane]; pc++; } } else if (op == 9) { i64 raw = R[b * rs]; if (raw >= romlo) { rep = %(r_end)d; break; } // 共有帯へは書けない i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; break; } M[(size_t)addr * (size_t)lanes + lane] = (u64)R[c * rs]; pc++; } else if (op == 10) { rep = %(r_ok)d; break; } else if (op == 11) { rep = %(r_ko)d; break; } else if (op == 12) { pc++; } else { rep = %(r_end)d; break; } } report[lane] = rep; steps[lane] = n; } // ===== 多数の lane が一つの仕事を分ける道 ===== // // 段と段の間では、lane i が書いた値を lane j が読む。そのために // // 1. 共有帯へ書ける（既存の道は報告 45 で止める） // 2. 書いた後に grid 全体で足並みを揃える // // の二つが要る。合図は「共有帯への書き込み」そのものにする。境界へ書く所が // チェンバーの接合点なので、pc を外から指さなくて済む。全 lane が段ごとに // ちょうど一度ずつ書く形の平文でだけ使える。 // // 足並み揃えは世代を反転させる素朴な形。block が全部同時に載っていないと // 止まるので、grid は device に載る数までに抑える（host 側で決める）。 __device__ __forceinline__ void gbar(int *cnt, int *gen, int nblocks) { __syncthreads(); if (threadIdx.x == 0) { int g = *((volatile int *)gen); __threadfence(); if (atomicAdd(cnt, 1) == nblocks - 1) { *cnt = 0; __threadfence(); atomicExch(gen, g + 1); } else { while (*((volatile int *)gen) == g) __threadfence(); } } __syncthreads(); } extern "C" __global__ void lane64_runb(const i64 *prog, int ninsn, u64 *M, int words, int lanes, i64 budget, int need, int work, int *report, i64 *steps, u64 *ROM, i64 romlo, i64 romn, int *cnt, int *gen, int nblocks) { extern __shared__ i64 sregb[]; int lane = blockIdx.x * blockDim.x + threadIdx.x; i64 *R = sregb + threadIdx.x; const int rs = blockDim.x; // 足並みを揃える以上、休んでいる thread を作れない。lanes は grid x block に // ちょうど合わせる（余りの lane は自分の席で同じ文を走らせる）。 int pc = 0; i64 n = 0; int rep = %(r_budget)d; while (n < budget) { n++; if (pc < 0 || pc >= ninsn) { rep = %(r_end)d; break; } const i64 *ins = prog + (size_t)pc * %(fields)d; int op = (int)ins[0], a = (int)ins[1], b = (int)ins[2], c = (int)ins[3]; int sb = (int)ins[4], sc = (int)ins[5]; i64 im = ins[6]; if (op == 0) { rep = %(r_end)d; break; } else if (op == 1) { R[a * rs] = im; pc++; } else if (op == 2) { R[a * rs] = (R[b * rs] << sb) + (R[c * rs] << sc); pc++; } else if (op == 3) { i64 x = R[b * rs], y = R[c * rs]; u64 ux = (u64)(x < 0 ? -x : x), uy = (u64)(y < 0 ? -y : y); if (__umul64hi(ux, uy) || ((ux * uy) >> 63)) { rep = %(r_wide)d; break; } R[a * rs] = x * y; pc++; } else if (op == 4) { R[a * rs] = -R[b * rs]; pc++; } else if (op == 5) { pc += ((R[b * rs] << sb) == (R[c * rs] << sc)) ? (int)im : 1; } else if (op == 6) { pc += ((R[b * rs] << sb) < (R[c * rs] << sc)) ? (int)im : 1; } else if (op == 7) { pc += (int)im; } else if (op == 8) { i64 raw = R[b * rs]; if (raw >= romlo) { i64 o = raw - romlo; if (o >= romn) { rep = %(r_end)d; break; } R[a * rs] = (i64)ROM[(size_t)o]; pc++; } else { i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; break; } R[a * rs] = (i64)M[(size_t)addr * (size_t)lanes + lane]; pc++; } } else if (op == 9) { i64 raw = R[b * rs]; if (raw >= romlo) { // 共有帯へ書く。ここが接合点 i64 o = raw - romlo; if (o >= romn) { rep = %(r_end)d; break; } ROM[(size_t)o] = (u64)R[c * rs]; pc++; __threadfence(); gbar(cnt, gen, nblocks); continue; } i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; break; } M[(size_t)addr * (size_t)lanes + lane] = (u64)R[c * rs]; pc++; } else if (op == 10) { rep = %(r_ok)d; break; } else if (op == 11) { rep = %(r_ko)d; break; } else if (op == 12) { pc++; } else { rep = %(r_end)d; break; } } report[lane] = rep; steps[lane] = n; } // ===== warp を pc の小さい順にそろえて歩く道 ===== // // thread ごとに好きに飛ぶと、warp の中で pc が割れた瞬間から命令が直列になり、 // 語の読みも散って memory に縛られる（電力 50 W で頭打ち）。 // 20 §6 命題 6.1: 合流は最も小さい pc から取れば足りる。warp の中で最小の pc を // 一つ選び、そこに居る lane だけが一命令進む。分岐で割れても、後ろの lane が // 先に追いつくので、loop の出口で自然に合流する。 extern "C" __global__ void lane64_run1_sync(const i64 *prog, int ninsn, u64 *M, int words, int lanes, i64 budget, int need, int work, int *report, i64 *steps, const u64 *ROM, i64 romlo, i64 romn, int stride) { // stride は語の帯の中の lane の数（M2 の列数）。lanes はいま走らせる lane の数で、 // 先頭から lanes 本だけが歩く。段を長さの順に並べておけば、短い段が終わった後の // 歩は先頭の生きている lane だけで済む extern __shared__ i64 sreg[]; int lane = blockIdx.x * blockDim.x + threadIdx.x; bool alive = lane < lanes; i64 *R = sreg + threadIdx.x; const int rs = blockDim.x; int pc = 0; i64 n = 0; int rep = %(r_budget)d; bool done = !alive; while (true) { unsigned live = __ballot_sync(0xffffffffu, !done); if (live == 0) break; int mypc = done ? 0x7fffffff : pc; int pcmin = __reduce_min_sync(0xffffffffu, mypc); if (done || pc != pcmin) continue; n++; if (n > budget) { done = true; continue; } if (pc < 0 || pc >= ninsn) { rep = %(r_end)d; done = true; continue; } const i64 *ins = prog + (size_t)pc * %(fields)d; int op = (int)ins[0], a = (int)ins[1], b = (int)ins[2], c = (int)ins[3]; int sb = (int)ins[4], sc = (int)ins[5]; i64 im = ins[6]; if (op == 0) { rep = %(r_end)d; done = true; } else if (op == 1) { R[a * rs] = im; pc++; } else if (op == 2) { R[a * rs] = (R[b * rs] << sb) + (R[c * rs] << sc); pc++; } else if (op == 3) { i64 x = R[b * rs], y = R[c * rs]; u64 ux = (u64)(x < 0 ? -x : x), uy = (u64)(y < 0 ? -y : y); if (__umul64hi(ux, uy) || ((ux * uy) >> 63)) { rep = %(r_wide)d; done = true; } else { R[a * rs] = x * y; pc++; } } else if (op == 4) { R[a * rs] = -R[b * rs]; pc++; } else if (op == 5) { pc += ((R[b * rs] << sb) == (R[c * rs] << sc)) ? (int)im : 1; } else if (op == 6) { pc += ((R[b * rs] << sb) < (R[c * rs] << sc)) ? (int)im : 1; } else if (op == 7) { pc += (int)im; } else if (op == 8) { i64 raw = R[b * rs]; if (raw >= romlo) { i64 o = raw - romlo; if (o >= romn) { rep = %(r_end)d; done = true; } else { R[a * rs] = (i64)ROM[(size_t)o]; pc++; } } else { i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; done = true; } else { R[a * rs] = (i64)M[(size_t)addr * (size_t)stride + lane]; pc++; } } } else if (op == 9) { i64 raw = R[b * rs]; if (raw >= romlo) { rep = %(r_end)d; done = true; } else { i64 addr = slot(raw, need, work, words); if (addr < 0) { rep = %(r_end)d; done = true; } else { M[(size_t)addr * (size_t)stride + lane] = (u64)R[c * rs]; pc++; } } } else if (op == 10) { rep = %(r_ok)d; done = true; } else if (op == 11) { rep = %(r_ko)d; done = true; } else if (op == 12) { pc++; } else { rep = %(r_end)d; done = true; } } if (alive) { report[lane] = rep; steps[lane] = n; } } extern "C" __global__ void gather(const u64 *M, int words, int lanes, const int *watch, int nwatch, u64 *out) { int lane = blockIdx.x * blockDim.x + threadIdx.x; if (lane >= lanes) return; u64 *o = out + (size_t)lane * nwatch * NLIMB; for (int k = 0; k < nwatch; k++) { const u64 *w = M + ((size_t)watch[k] * NLIMB) * (size_t)lanes + lane; for (int i = 0; i < NLIMB; i++) o[k * NLIMB + i] = w[(size_t)i * lanes]; } } // 先頭 lanes 本だけを集める。stride は帯の中の lane の数 extern "C" __global__ void gather_s(const u64 *M, int stride, int lanes, const int *watch, int nwatch, u64 *out) { int lane = blockIdx.x * blockDim.x + threadIdx.x; if (lane >= lanes) return; u64 *o = out + (size_t)lane * nwatch; for (int k = 0; k < nwatch; k++) o[k] = M[(size_t)watch[k] * stride + lane]; } """ R_END = 45 _MOD = {} # ===== block 内の即値を畳む ===== # # commonc が出す lane64 の命令列は、番地を one の倍化と加算で組み立てる。 # CPU 側の vcpu.translate は block ごとに値を記号で追い、即値どうしの算術を # 一つの即値へ畳む（casm.lane64.vn）。GPU の解釈器は命令を一つずつ歩くので、 # 同じ畳み込みを先に済ませておかないと、番地を綴る鎖がそのまま実行費用になる。 # # R1 = 1 ; R1 = R1+R1 ; ... ; R2 = R1+R3 → R1 = 65554 の一命令 # # block の入口で値が分かっている register は無い、として始める。分岐は畳まない。 MAT, ADD, MUL, NEG, EQ, LT, JMP, LD, ST = 1, 2, 3, 4, 5, 6, 7, 8, 9 _BR = (EQ, LT, JMP) _TERM = (0, 10, 11) def leaders(prog): """block の先頭になる命令の番号""" out = set([0]) for i, ins in enumerate(prog): op, im = ins[0], ins[6] if op in _BR: out.add(i + int(im)) out.add(i + 1) elif op in _TERM: out.add(i + 1) return sorted(a for a in out if 0 <= a < len(prog)) def _reads(ins): op = ins[0] if op in (ADD, EQ, LT, MUL): return (ins[2], ins[3]) if op in (NEG, LD): return (ins[2],) if op == ST: return (ins[2], ins[3]) if op == 12: return (ins[2],) return () def _writes(ins): op = ins[0] return ins[1] if op in (MAT, ADD, MUL, NEG, LD) else None def _walk(prog, lo, hi, live, dead, rewrite): """block をひと通り追う。live には「自分の Store より前に読んだ定数番地」を集め、 rewrite が真なら即値の算術と語への往復を畳み、dead の番地への Store を落とす。""" val, mem = {}, {} for i in range(lo, hi): ins = prog[i] if ins is None: continue op = ins[0] if op == MAT: val[ins[1]] = int(ins[6]) elif op in (ADD, MUL, NEG): b, c, sb, sc = ins[2], ins[3], ins[4], ins[5] vb = val.get(b) vc = val.get(c) if op != NEG else 0 if vb is not None and vc is not None: v = (vb << sb) + (vc << sc) if op == ADD else (vb * vc if op == MUL else -vb) if rewrite: prog[i] = [MAT, ins[1], 0, 0, 0, 0, v] val[ins[1]] = v else: val.pop(ins[1], None) elif op == LD: a = val.get(ins[2]) if a is not None and a in mem: if rewrite: prog[i] = [MAT, ins[1], 0, 0, 0, 0, mem[a]] val[ins[1]] = mem[a] else: if a is not None: live.add(a) val.pop(ins[1], None) elif op == ST: a, v = val.get(ins[2]), val.get(ins[3]) if a is None: mem.clear() elif v is None: mem.pop(a, None) else: mem[a] = v # 同じ block の Load は先に繋いでから落とす。落とすのは block を跨ぐ分だけ if rewrite and a is not None and a in dead: prog[i] = None else: w = _writes(ins) if w is not None: val.pop(w, None) def fold(prog, wires=None): """block の中で値を記号で追い、即値どうしの算術と語への往復を畳む。意味は変わらない。 lane64 の一箱は MAT・ST・LD の往復に開くので、register だけを追っても畳めない。 定数番地への Store と、その番地からの Load を繋いで初めて鎖が一つの即値になる。 """ prog = [list(x) for x in prog] heads = leaders(prog) + [len(prog)] # 一巡目。block を跨いで読まれる定数番地を集める live = set() for k in range(len(heads) - 1): _walk(prog, heads[k], heads[k + 1], live, (), False) # commonc が出すワイヤの番地は定数だけで、本文の番地は WORK 以上にある。 # だから [64, need) のうち誰も跨いで読まない番地への Store は落としてよい dead = set() if wires: lo, hi = wires dead = set(a for a in range(lo, hi) if a not in live) # 二巡目。畳んで、死んだ Store を落とす for k in range(len(heads) - 1): lo, hi = heads[k], heads[k + 1] _walk(prog, lo, hi, set(), dead, True) # block の中だけの死んだ書き込みを落とす for k in range(len(heads) - 1): lo, hi = heads[k], heads[k + 1] last, drop = {}, set() for i in range(lo, hi): ins = prog[i] if ins is None: continue for r in _reads(ins): last.pop(r, None) w = _writes(ins) if w is not None: if ins[0] == MAT and w in last: drop.add(last[w]) last[w] = i if ins[0] == MAT else None if last[w] is None: last.pop(w) for i in drop: prog[i] = None # 落とした分だけ番号が詰まるので、分岐の変位を貼り直す keep = [i for i, x in enumerate(prog) if x is not None] newpos, n = {}, 0 for i in range(len(prog)): newpos[i] = n if prog[i] is not None: n += 1 newpos[len(prog)] = n out = [] for i in keep: ins = list(prog[i]) if ins[0] in _BR: ins[6] = newpos[min(i + int(ins[6]), len(prog))] - newpos[i] out.append(tuple(ins)) return out def compact(prog, wires, work=0): """ワイヤの番地を詰める。返りは (prog, 新しい need)。 commonc の一箱はワイヤ一本を語 64 + w に取るので、image が大きいと lane ごとの ワイヤの帯が数万語になる。畳んだ後に生き残る番地はその一部なので、番地として 使われている即値だけを 64 から順に付け直す。 即値が番地か値かは block の中の使われ方で決める。MAT で置いた register が、 上書きされるまでに LD / ST の番地としてしか使われなければ番地。値として 使われた即値は触らない。触らない即値がワイヤの帯に残れば、need は その最大までしか下がらない（返りの need がそれを示す）。 """ lo, hi = wires # 二値の八箱の退避語は定数番地の最大の上に置かれ、ワイヤの帯と本文の間の穴に落ちる # ことがある（16 §5.3）。穴は lane に無いので、そこも一緒に詰める hi = max(hi, work) if work else hi prog = [list(x) for x in prog] heads = leaders(prog) + [len(prog)] addr_mats, keep_hi = [], lo for k in range(len(heads) - 1): a, b = heads[k], heads[k + 1] for i in range(a, b): ins = prog[i] if ins[0] != MAT or not (lo <= ins[6] < hi): continue reg, ok, closed = ins[1], True, False for j in range(i + 1, b): jn = prog[j] op = jn[0] if op in (LD, ST) and jn[2] == reg and not (op == ST and jn[3] == reg): pass elif reg in _reads(jn): ok = False break if _writes(jn) == reg: closed = True break # block の外へ出る register は次の block でも番地として使われる（commonc は # ワイヤの番地を定数でしか作らない）。値として読まれた即値だけを残す if ok: addr_mats.append(i) else: keep_hi = max(keep_hi, ins[6] + 1) used = sorted(set(prog[i][6] for i in addr_mats)) base = keep_hi remap = {a: base + k for k, a in enumerate(used)} for i in addr_mats: prog[i][6] = remap[prog[i][6]] return [tuple(x) for x in prog], base + len(used) def _module(cp, nlimb): if nlimb not in _MOD: src = SRC % {"nlimb": nlimb, "fields": INSN_FIELDS, "r_ok": R_OK, "r_ko": R_KO, "r_wide": R_WIDE, "r_budget": R_BUDGET, "r_end": R_END} _MOD[nlimb] = cp.RawModule(code=src, options=("-std=c++11",)) return _MOD[nlimb] def to_limbs(v, nlimb): """整数を limb の並びへ（二の補数）""" mask = (1 << (64 * nlimb)) - 1 x = v & mask return [(x >> (64 * i)) & 0xFFFFFFFFFFFFFFFF for i in range(nlimb)] def from_limbs(limbs, nlimb): x = 0 for i, l in enumerate(limbs): x |= int(l) << (64 * i) if x >> (64 * nlimb - 1): x -= 1 << (64 * nlimb) return x _NOROM_BUF = {} def _NOROM(cp): """共有帯を使わないときに渡す一語の置き場。kernel は romlo = 2^62 で触らない""" if "b" not in _NOROM_BUF: _NOROM_BUF["b"] = cp.zeros((1,), dtype=cp.uint64) return _NOROM_BUF["b"] def pack(addr, need, work): """語番号 → 詰めた並びの位置。kernel の slot() と同じ式（式は device 側 §slot）。""" if work == 0: return addr if 0 <= addr < need: return addr return addr - work + need def run(prog, seeds, words, watch, nlimb=16, budget=1 << 26, chunk=64, need=0, work=0, base=None): """image（decode 済み）を seeds の lane 分だけ走らせる。 seeds は lane ごとの {語: 値}、watch は host へ戻す語の並び。 返りは lane ごとの (報告, {語: 値})。報告は 0 = tu、42 = te、43 = 幅を越えた、44 = 歩数の上限。 base は全 lane へ同じ値を置く語（{語: 値}）。模型の表のように lane ごとに同じ物を 渡すとき、seeds へ毎回入れると host の側だけで lane × 語 の手間になる。ここでは 語ごとに一度だけ device の帯へ書くので、手間は語の数で済む。 need と work を渡すと、ワイヤ（0 .. need）と本文の作業領域（work 以上）の間の穴を詰める。 commonc の lane64 版は本文を WORK = 2^16 から置くので、詰めないと lane ごとに 65536 語ぶんの空きを確保することになる（20 §1.2 の lane の作業領域は交わらない）。 words = need + (本文が使う語の数) 穴を詰めると lane あたりの領域が減り、同じ memory でより多くの lane を並べられる。 """ import cupy as cp mod = _module(cp, nlimb) k_run = mod.get_function("lane64_run1" if nlimb == 1 else "lane64_run") k_gather = mod.get_function("gather") prog = fold(prog) table = np.array([list(p) for p in prog], dtype=np.int64).reshape(-1) d_prog = cp.asarray(table) d_watch = cp.asarray(np.array([pack(w, need, work) for w in watch], dtype=np.int32)) out = [] for start in range(0, len(seeds), chunk): part = seeds[start:start + chunk] lanes = len(part) M = cp.zeros((lanes * words * nlimb,), dtype=cp.uint64) if base: for addr, v in base.items(): a = pack(addr, need, work) v = int(v) if 0 <= v < (1 << 63): M[(a * nlimb) * lanes:(a * nlimb) * lanes + lanes] = np.uint64(v) else: for i, limb in enumerate(to_limbs(v, nlimb)): if limb: o = (a * nlimb + i) * lanes M[o:o + lanes] = np.uint64(limb) idx, vals = [], [] for l, seed in enumerate(part): for addr, v in seed.items(): a = pack(addr, need, work) v = int(v) if 0 <= v < (1 << 63): idx.append((a * nlimb) * lanes + l) vals.append(v) else: for i, limb in enumerate(to_limbs(v, nlimb)): if limb: idx.append((a * nlimb + i) * lanes + l) vals.append(limb) if idx: M[cp.asarray(np.array(idx, dtype=np.int64))] = cp.asarray(np.array(vals, dtype=np.uint64)) report = cp.zeros((lanes,), dtype=cp.int32) steps = cp.zeros((lanes,), dtype=cp.int64) threads = 128 grid = (lanes + threads - 1) // threads sh = threads * NREG_BYTES if nlimb == 1 else 0 k_run((grid,), (threads,), (d_prog, np.int32(len(prog)), M, np.int32(words), np.int32(lanes), np.int64(budget), np.int32(need), np.int32(work), report, steps, _NOROM(cp), np.int64(1 << 62), np.int64(0)), shared_mem=sh) got = cp.zeros((lanes * len(watch) * nlimb,), dtype=cp.uint64) if watch: k_gather((grid,), (threads,), (M, np.int32(words), np.int32(lanes), d_watch, np.int32(len(watch)), got)) rep = report.get() raw = got.get() for l in range(lanes): vals_l = {} for k, w in enumerate(watch): o = (l * len(watch) + k) * nlimb vals_l[w] = from_limbs(raw[o:o + nlimb], nlimb) out.append((int(rep[l]), vals_l)) del M, got return out def run_shared(prog, base, lanes, words, watch, need, work, nlimb=4, budget=1 << 26, chunk=2048, romlo=None, romn=0): """lane が共有する種 base に、lane ごとの少しの語 lanes[k] を重ねて走らせる。 語の並びは詰める: ワイヤ [0, need) をそのまま、作業領域 [work, ...) を need の後ろに置く。 words は詰めた後の語数（need + 作業領域の幅）。平文の側の語は全部 work 以上にあるので、 間の番地は使われない。返りは lane ごとの (報告, {語: 値})。 romlo を渡すと、番地 [romlo, romlo + romn) は**全 lane で一枚だけ**持つ。 読むだけの表（値段の表など）をここへ置くと、lane あたりの重さが桁で落ちる。 lane あたりの語数 = need + (romlo - work) 共有帯の語数 = romn この帯へ Store すると報告 45 で止まる。書く語を入れてはいけない。 """ import cupy as cp mod = _module(cp, nlimb) k_run = mod.get_function("lane64_run1" if nlimb == 1 else "lane64_run") k_gather = mod.get_function("gather") prog = fold(prog, wires=(64, need)) table = np.array([list(p) for p in prog], dtype=np.int64).reshape(-1) d_prog = cp.asarray(table) hi = romlo if romlo is not None else (1 << 62) def pos(a): return a if a < need else a - work + need d_watch = cp.asarray(np.array([pos(w) for w in watch], dtype=np.int32)) # 共有の種を一枚の並びにする。共有帯の分は別の並びへ分ける row = np.zeros((words * nlimb,), dtype=np.uint64) rom = np.zeros((max(1, romn) * nlimb,), dtype=np.uint64) for addr, v in base.items(): v = int(v) p0 = (addr - hi) * nlimb if addr >= hi else pos(addr) * nlimb dst = rom if addr >= hi else row if 0 <= v < (1 << 63): dst[p0] = v else: for i, limb in enumerate(to_limbs(v, nlimb)): dst[p0 + i] = limb d_row = cp.asarray(row) d_rom = cp.asarray(rom) out = [] for start in range(0, len(lanes), chunk): part = lanes[start:start + chunk] n = len(part) M = cp.repeat(d_row, n) idx, vals = [], [] for l, seed in enumerate(part): for addr, v in seed.items(): idx.append((pos(addr) * nlimb) * n + l) vals.append(int(v)) if idx: M[cp.asarray(np.array(idx, dtype=np.int64))] = cp.asarray(np.array(vals, dtype=np.uint64)) report = cp.zeros((n,), dtype=cp.int32) steps = cp.zeros((n,), dtype=cp.int64) threads = 128 grid = (n + threads - 1) // threads sh = threads * NREG_BYTES if nlimb == 1 else 0 k_run((grid,), (threads,), (d_prog, np.int32(len(prog)), M, np.int32(words), np.int32(n), np.int64(budget), np.int32(need), np.int32(work), report, steps, d_rom, np.int64(hi), np.int64(romn)), shared_mem=sh) got = cp.zeros((n * len(watch) * nlimb,), dtype=cp.uint64) k_gather((grid,), (threads,), (M, np.int32(words), np.int32(n), d_watch, np.int32(len(watch)), got)) rep = report.get() raw = got.get().reshape(n, len(watch), nlimb) for l in range(n): out.append((int(rep[l]), dict((w, from_limbs(raw[l, k], nlimb)) for k, w in enumerate(watch)))) del M, got return out class Lanes(object): """lane の作業領域と共有帯を device に置いたまま、image を何度も流す。 host は語を一つずつ触らない。種は (lane, 語, 値) の並びを一度に散らし、読みは (lane × 語) の並びを一度に集める。どちらも cupy の一手で済むので、lane の数にも 語の数にも Python の反復が付かない。 L = Lanes(prog, need, work, romlo, romn, nlanes, base={語: 値}) L.put(lanes, addrs, vals) 種を散らす（三つとも同じ長さの並び） L.rom_put(offset, values) 共有帯の一部を書く rep, W = L.run(watch) rep[lane] が報告、W[lane, k] が watch[k] の値（device） L.M2 語 × lane の二次元の眺め（device） L.reset() 作業領域を base の並びへ戻す 語の並びは run_shared と同じ: ワイヤ [0, need) をそのまま、作業領域 [work, romlo) を need の後ろに詰める。M は語主で M2[pos(語), lane]。 共有帯 rom は cupy の並び一枚。別の Lanes の M をそのまま rom に渡せる（rom=）ので、 lane の畳みのように「全 lane の語を読む」kernel は、その M を共有帯として読めば書ける。 """ def __init__(self, prog, need, work, romlo, romn, nlanes, base=None, nlimb=1, rom=None, sync=True): """sync が真なら warp を pc の小さい順にそろえて歩く（lane64_run1_sync）。""" import cupy as cp self.cp = cp if nlimb != 1: raise ValueError("Lanes は 64 bit の語だけを扱う（nlimb=1）") self.need, self.work, self.romlo = need, work, romlo self.words = need + (romlo - work) self.n = nlanes mod = _module(cp, 1) self.sync = sync self.k_run = mod.get_function("lane64_run1_sync" if sync else "lane64_run1") self.k_gather = mod.get_function("gather_s") prog = fold(prog, wires=(64, need)) prog, need = compact(prog, (64, need), work=work) self.need = need self.words = need + (romlo - work) self.nprog = len(prog) self.d_prog = cp.asarray(np.array([list(p) for p in prog], dtype=np.int64).reshape(-1)) if rom is not None: self.d_rom = rom self.romn = int(rom.size) else: self.romn = romn self.d_rom = cp.zeros((max(1, romn),), dtype=cp.uint64) row = np.zeros((self.words,), dtype=np.uint64) for addr, v in (base or {}).items(): row[self.pos(addr)] = int(v) self.d_row = cp.asarray(row) self.M = cp.zeros((self.words * nlanes,), dtype=cp.uint64) self.reset() def pos(self, a): return a if a < self.need else a - self.work + self.need @property def M2(self): return self.M.reshape(self.words, self.n) def reset(self): self.M2[:, :] = self.d_row.reshape(-1, 1) def put(self, lanes, addrs, vals): """(lane, 語, 値) の並びを一度に散らす。addrs は語番号（pos は中で当てる）。""" cp = self.cp lanes = cp.asarray(np.asarray(lanes, dtype=np.int64)) addrs = np.asarray(addrs, dtype=np.int64) addrs = np.where(addrs < self.need, addrs, addrs - self.work + self.need) idx = cp.asarray(addrs) * self.n + lanes self.M[idx] = cp.asarray(np.asarray(vals, dtype=np.uint64)) def put_block(self, addr, block): """語 addr から続く語の帯へ、(語数 × lane) の並びを一度に置く。block は device でも host でもよい。""" cp = self.cp b = cp.asarray(block, dtype=cp.uint64) self.M2[self.pos(addr):self.pos(addr) + b.shape[0], :] = b def get_block(self, addr, n): """語 addr から n 語の帯を (n × lane) で返す（device）。""" return self.M2[self.pos(addr):self.pos(addr) + n, :] def rom_put(self, offset, values): self.d_rom[offset:offset + len(values)] = self.cp.asarray(np.asarray(values, dtype=np.uint64)) def run(self, watch=(), budget=1 << 36, active=None): """先頭 active 本の lane を走らせる（省くと全部）。返りは (報告, watch の値) で、どちらも active 本ぶん。""" cp = self.cp n = self.n if active is None else min(active, self.n) report = cp.zeros((n,), dtype=cp.int32) steps = cp.zeros((n,), dtype=cp.int64) threads = 128 grid = (n + threads - 1) // threads if self.sync: args = (self.d_prog, np.int32(self.nprog), self.M, np.int32(self.words), np.int32(n), np.int64(budget), np.int32(self.need), np.int32(self.work), report, steps, self.d_rom, np.int64(self.romlo), np.int64(self.romn), np.int32(self.n)) else: if n != self.n: raise ValueError("先頭だけを走らせるには sync の道が要る") args = (self.d_prog, np.int32(self.nprog), self.M, np.int32(self.words), np.int32(n), np.int64(budget), np.int32(self.need), np.int32(self.work), report, steps, self.d_rom, np.int64(self.romlo), np.int64(self.romn)) self.k_run((grid,), (threads,), args, shared_mem=threads * NREG_BYTES) W = None if watch: d_watch = cp.asarray(np.array([self.pos(w) for w in watch], dtype=np.int32)) got = cp.zeros((n * len(watch),), dtype=cp.uint64) self.k_gather((grid,), (threads,), (self.M, np.int32(self.n), np.int32(n), d_watch, np.int32(len(watch)), got)) W = got.reshape(n, len(watch)).astype(cp.int64) return report.get(), W