Skip to main content

FPU_V

pub const FPU_V: &str = "// fpu.v \u{2014} a synthesizable IEEE-754 double-precision (binary64 / f64) floating-point unit,\n// proven BIT-EXACT to Rust `f64` in iverilog (tb_fpu_*.v, generated FROM Rust f64 vectors).\n//\n// The four float benchmarks (nbody, mandelbrot, spectral_norm, pi_leibniz) need f64: a constant\n// like 9.54791938424326609 has 17 significant figures and only survives in binary64. Every op here\n// is round-to-nearest-ties-to-even with correct \u{b1}0 / \u{b1}Inf / quiet-NaN / subnormal handling, matching\n// Rust\'s `f64` operator semantics exactly (Rust uses the platform\'s IEEE-754 hardware, which is what\n// we reproduce here). NO `real` / `$bitstoreal` appears in this module \u{2014} those are testbench-only;\n// this DUT is pure synthesizable logic. div/sqrt are multi-cycle digit-recurrence; add/sub/mul/convert\n// are ALSO multi-cycle (an FSM registers each stage) \u{2014} a start/done handshake, correctly rounded. No\n// stage holds more than one carry chain, barrel shift, leading-one search or multiplier: carries are\n// parallel-prefix, leading-one searches halve their word level by level, stickies test the shifted-out\n// bits with no (1 << s) - 1 mask. Every register boundary is between IDENTICAL arithmetic expressions,\n// so the RESULT stays bit-for-bit what the combinational form produced.\n//\n// Interface \u{2014} a single fpu module selects the op by a 3-bit code:\n//   OP_ADD=0 OP_SUB=1 OP_MUL=2 OP_DIV=3 OP_SQRT=4 OP_LT=5 OP_EQ=6\n// Arithmetic ops (add/sub/mul/div/sqrt) run a bounded FSM: `start` pulses to issue, `done` rises when\n// `y` is valid. Comparisons (lt/eq) stay combinational (no shifters). For comparisons the numeric\n// result low bit is the boolean; `result[0]` and `cmp` both carry it.\n\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\n// fp_addsub \u{2014} multi-cycle IEEE-754 binary64 add/sub, round-to-nearest-ties-even. start/done handshake.\n//   sub=0 -> a+b ; sub=1 -> a-b  (a - b == a + (-b), so we just flip b\'s sign)\n// The datapath is the SAME as the classic single-cycle form, split across the FSM below so each wide\n// shifter, the 128-bit add, the leading-one search and the rounding increment is its own small cone instead\n// of one 128-bit-wide nest. Bit-exact to the combinational reference.\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\nmodule fp_addsub(\n  input         clk,\n  input         rst,\n  input         start,\n  input  [63:0] a,\n  input  [63:0] b_in,\n  input         sub,      // 0: add, 1: subtract\n  output reg [63:0] y,\n  output reg    done\n);\n  localparam [63:0] QNAN = 64\'h7FF8000000000000;\n\n  // One cone per state, none holding more than one carry chain, barrel shift or leading-one search:\n  // DECODE (specials, the magnitude compare), ALIGN (the smaller operand\'s shift and its sticky), ADD (the 128-bit\n  // add or subtract), LZC (the leading one), NORM (the normalizing shift and its sticky), ROUND (round to nearest\n  // even and pack a normal result), and for a subnormal result SUBN (canonicalize), SUBR (denormalize) and SUBP\n  // (round and pack).\n  reg  [3:0] state;\n  localparam IDLE = 4\'d0, DECODE = 4\'d1, ALIGN = 4\'d2, ADD = 4\'d3, LZC = 4\'d4, NORM = 4\'d5, ROUND = 4\'d6,\n             SUBN = 4\'d7, SUBR = 4\'d8, SUBP = 4\'d9;\n\n  // latched inputs (captured at issue; the shared operand ports fav/fbv stay stable while we stall)\n  reg  [63:0] a_l, b_l;\n\n  // flip the sign of the second operand for subtraction\n  wire [63:0] b = sub ? {~b_l[63], b_l[62:0]} : b_l;\n\n  wire        sa = a_l[63], sb = b[63];\n  wire [10:0] ea = a_l[62:52], eb = b[62:52];\n  wire [51:0] ma = a_l[51:0],  mb = b[51:0];\n\n  wire a_is_zero  = (ea == 11\'d0)    && (ma == 52\'d0);\n  wire b_is_zero  = (eb == 11\'d0)    && (mb == 52\'d0);\n  wire a_is_inf   = (ea == 11\'h7FF)  && (ma == 52\'d0);\n  wire b_is_inf   = (eb == 11\'h7FF)  && (mb == 52\'d0);\n  wire a_is_nan   = (ea == 11\'h7FF)  && (ma != 52\'d0);\n  wire b_is_nan   = (eb == 11\'h7FF)  && (mb != 52\'d0);\n  wire a_is_sub   = (ea == 11\'d0)    && (ma != 52\'d0);  // subnormal (nonzero, exp field 0)\n  wire b_is_sub   = (eb == 11\'d0)    && (mb != 52\'d0);\n\n  // effective (unbiased-friendly) exponent: subnormals use exp=1, hidden bit 0; normals hidden bit 1\n  wire [11:0] xa = a_is_sub ? 12\'d1 : {1\'b0, ea};\n  wire [11:0] xb = b_is_sub ? 12\'d1 : {1\'b0, eb};\n  wire        ha = (ea != 11\'d0);\n  wire        hb = (eb != 11\'d0);\n  wire [52:0] siga = {ha, ma};\n  wire [52:0] sigb = {hb, mb};\n\n  reg  [11:0] exp_l;        // exponent of the larger-magnitude operand\n  reg  [52:0] sig_big, sig_small;\n  reg         sgn_big, sgn_small;\n  reg  [11:0] shamt;\n  reg  [127:0] small_ext;   // the smaller operand aligned, its sticky folded into bit 0\n  reg  [127:0] sum_ext;     // the exact aligned sum\n  reg          res_sign;\n  reg          op_is_sub;\n  reg  [11:0] norm_exp;\n  reg  [11:0] lead;         // position of the leading 1\n  reg  [127:0] norm_sig_r;  // leading 1 at bit 55: [55:3] the significand, [2:0] guard/round/sticky\n  reg  [11:0] rexp_r;       // the rounded exponent, for a subnormal result\n\n  // A word\'s bits below `s`, as one flag: the bits shifted to the top and tested, so no (1 << s) - 1 mask and no\n  // borrow chain. Zero when s is zero.\n  function below(input [127:0] x, input [11:0] s);\n    begin\n      below = (s == 12\'d0) ? 1\'b0 : ((x << (12\'d128 - s)) != 128\'d0);\n    end\n  endfunction\n\n  // The highest set bit of a nonzero 128-bit word, by halving: seven levels, not a 128-long priority chain.\n  function [11:0] lead_one(input [127:0] x);\n    reg [63:0] h6; reg [31:0] h5; reg [15:0] h4; reg [7:0] h3; reg [3:0] h2; reg [1:0] h1; reg [6:0] z;\n    begin\n      z[6] = ~|x[127:64]; h6 = z[6] ? x[63:0]   : x[127:64];\n      z[5] = ~|h6[63:32]; h5 = z[5] ? h6[31:0]  : h6[63:32];\n      z[4] = ~|h5[31:16]; h4 = z[4] ? h5[15:0]  : h5[31:16];\n      z[3] = ~|h4[15:8];  h3 = z[3] ? h4[7:0]   : h4[15:8];\n      z[2] = ~|h3[7:4];   h2 = z[2] ? h3[3:0]   : h3[7:4];\n      z[1] = ~|h2[3:2];   h1 = z[1] ? h2[1:0]   : h2[3:2];\n      z[0] = ~h1[1];\n      lead_one = {5\'d0, ~z};\n    end\n  endfunction\n\n  // a \u{b1} b over 128 bits as a Kogge-Stone prefix add, the subtract\'s +1 entering as bit 0\'s carry.\n  function [127:0] addsub128(input [127:0] x, input [127:0] yy, input s);\n    reg [127:0] yb, g, p, p0; integer d;\n    begin\n      yb = s ? ~yy : yy;\n      p0 = x ^ yb;\n      g = (x & yb) | (p0 & {127\'d0, s});\n      p = p0;\n      for (d = 1; d < 128; d = d * 2) begin\n        g = g | (p & (g << d));\n        p = p & (p << d);\n      end\n      addsub128 = p0 ^ {g[126:0], s};\n    end\n  endfunction\n\n  // v + cin through a parallel-prefix AND: log2(64) levels. cin is the scan\'s bit 0, so bit i flips on cin and\n  // every bit below it, and cin reaches the bits through the scan\'s tree rather than one net to all 64.\n  function [63:0] inc64(input [63:0] v, input cin);\n    reg [63:0] p; integer d;\n    begin\n      p = {v[62:0], cin};\n      for (d = 1; d < 64; d = d * 2) p = p & ((p << d) | ((64\'d1 << d) - 64\'d1));\n      inc64 = v ^ p;\n    end\n  endfunction\n\n  integer i;\n  reg [127:0] big_ext, small_w, norm_sig, win;\n  reg [63:0] rounded, rr;\n  reg [52:0] final_sig;\n  reg [11:0] rexp, final_exp, shr, rshift;\n  reg guard, roundb, sticky, lsb, round_up, g2, r2, s2, l2, up2;\n\n  initial begin state = IDLE; done = 1\'b0; y = 64\'d0; end\n\n  always @(posedge clk) begin\n    if (rst) begin state <= IDLE; done <= 1\'b0; y <= 64\'d0; end\n    else case (state)\n      IDLE: if (start) begin\n        done <= 1\'b0; a_l <= a; b_l <= b_in; state <= DECODE;\n      end\n      DECODE: begin\n        if (a_is_nan || b_is_nan) begin\n          y <= QNAN; done <= 1\'b1; state <= IDLE;\n        end\n        else if (a_is_inf && b_is_inf) begin\n          if (sa == sb) y <= {sa, 11\'h7FF, 52\'d0}; else y <= QNAN;\n          done <= 1\'b1; state <= IDLE;\n        end\n        else if (a_is_inf) begin y <= {sa, 11\'h7FF, 52\'d0}; done <= 1\'b1; state <= IDLE; end\n        else if (b_is_inf) begin y <= {sb, 11\'h7FF, 52\'d0}; done <= 1\'b1; state <= IDLE; end\n        else if (a_is_zero && b_is_zero) begin\n          if (sa == sb) y <= {sa, 63\'d0}; else y <= 64\'d0;\n          done <= 1\'b1; state <= IDLE;\n        end\n        else if (a_is_zero) begin y <= b; done <= 1\'b1; state <= IDLE; end\n        else if (b_is_zero) begin y <= a_l; done <= 1\'b1; state <= IDLE; end\n        else begin\n          if (xa > xb || (xa == xb && siga >= sigb)) begin\n            exp_l     <= xa;\n            sig_big   <= siga;   sgn_big   <= sa;\n            sig_small <= sigb;   sgn_small <= sb;\n            shamt     <= xa - xb;\n          end else begin\n            exp_l     <= xb;\n            sig_big   <= sigb;   sgn_big   <= sb;\n            sig_small <= siga;   sgn_small <= sa;\n            shamt     <= xb - xa;\n          end\n          state <= ALIGN;\n        end\n      end\n      ALIGN: begin\n        // The small operand, its LSB at bit 64, shifted right by shamt in a 128-bit lane; with 64 bits of room a\n        // shift up to 64 keeps every bit, and what falls off folds into a sticky bit 0.\n        small_w = {11\'d0, sig_small, 64\'d0};\n        if (shamt >= 12\'d128) begin\n          small_ext <= {127\'d0, |sig_small};\n        end else begin\n          small_ext <= (small_w >> shamt) | {127\'d0, below(small_w, shamt)};\n        end\n        op_is_sub <= (sgn_big != sgn_small);\n        res_sign  <= sgn_big;\n        state <= ADD;\n      end\n      ADD: begin\n        big_ext = {11\'d0, sig_big, 64\'d0};\n        sum_ext <= addsub128(big_ext, small_ext, op_is_sub);\n        state <= LZC;\n      end\n      LZC: begin\n        if (sum_ext == 128\'d0) begin\n          y <= 64\'d0; done <= 1\'b1; state <= IDLE;      // exact cancellation -> +0 (RNE)\n        end else begin\n          lead  <= lead_one(sum_ext);\n          state <= NORM;\n        end\n      end\n      NORM: begin\n        norm_exp <= exp_l + lead - 12\'d116;              // biased exponent of the leading 1\n        if (lead >= 12\'d55)\n          norm_sig_r <= (sum_ext >> (lead - 12\'d55)) | {127\'d0, below(sum_ext, lead - 12\'d55)};\n        else\n          norm_sig_r <= sum_ext << (12\'d55 - lead);\n        state <= ROUND;\n      end\n      ROUND: begin\n        norm_sig = norm_sig_r;\n        guard  = norm_sig[2];\n        roundb = norm_sig[1];\n        sticky = norm_sig[0];\n        lsb    = norm_sig[3];\n        round_up = guard & (roundb | sticky | lsb);\n        rounded = inc64({11\'d0, norm_sig[55:3]}, round_up);\n        rexp = norm_exp;\n        if (rounded[53]) begin\n          final_sig = rounded[53:1];\n          rexp      = rexp + 12\'d1;\n        end else begin\n          final_sig = rounded[52:0];\n        end\n        if ($signed(rexp) >= 12\'sd2047) begin\n          y <= {res_sign, 11\'h7FF, 52\'d0};              // overflow -> inf\n          done <= 1\'b1; state <= IDLE;\n        end else if ($signed(rexp) <= 12\'sd0) begin\n          rexp_r <= rexp;\n          state <= SUBN;                                // a subnormal result rounds from the exact sum\n        end else begin\n          final_exp = rexp;\n          y <= {res_sign, final_exp[10:0], final_sig[51:0]};\n          done <= 1\'b1; state <= IDLE;\n        end\n      end\n      SUBN: begin\n        // The exact sum canonicalized so its leading 1 sits at bit 55, bits below the window folded into sticky.\n        // The canonical window goes to norm_sig_r, free once ROUND has read it.\n        if (lead >= 12\'d55) begin\n          shr = lead - 12\'d55;\n          norm_sig_r <= (sum_ext >> shr) | {127\'d0, below(sum_ext, shr)};\n        end else begin\n          norm_sig_r <= sum_ext << (12\'d55 - lead);\n        end\n        state <= SUBR;\n      end\n      SUBR: begin\n        // Shift the window right by 1 - exponent, folding the discarded bits into sticky.\n        // The denormalized window goes to sum_ext, free once SUBN has read it.\n        rshift = 12\'d1 - rexp_r;\n        if (rshift >= 12\'d128)\n          sum_ext <= {127\'d0, |norm_sig_r};\n        else if (rshift != 12\'d0)\n          sum_ext <= (norm_sig_r >> rshift) | {127\'d0, below(norm_sig_r, rshift)};\n        else\n          sum_ext <= norm_sig_r;\n        state <= SUBP;\n      end\n      SUBP: begin\n        win = sum_ext;\n        g2 = win[2]; r2 = win[1]; s2 = win[0]; l2 = win[3];\n        up2 = g2 & (r2 | s2 | l2);\n        rr = inc64({11\'d0, win[55:3]}, up2);\n        if (rr[52]) y <= {res_sign, 11\'d1, rr[51:0]};\n        else        y <= {res_sign, 11\'d0, rr[51:0]};\n        done <= 1\'b1; state <= IDLE;\n      end\n      default: state <= IDLE;\n    endcase\n  end\nendmodule\n\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\n// fp_mul \u{2014} multi-cycle IEEE-754 binary64 multiply, round-to-nearest-ties-even. start/done handshake.\n// Same datapath as the single-cycle form, split across the FSM (DECODE: specials + leading-zero counts;\n// PRENORM: the pre-normalizing shifts; MUL and MUL2: the 53\u{d7}53 product\'s compressor tree (fp_csa_mul), its first\n// four levels and then the rest; PRODUCT: its final add + top-bit normalize; ROUND: normal RNE round + pack, or the\n// subnormal grid shift; SUBP: the subnormal round + pack). Bit-exact to the combinational reference.\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\nmodule fp_mul(\n  input         clk,\n  input         rst,\n  input         start,\n  input  [63:0] a,\n  input  [63:0] b,\n  output reg [63:0] y,\n  output reg    done\n);\n  localparam [63:0] QNAN = 64\'h7FF8000000000000;\n  parameter integer MUL_TREE = 0;   // 1: the significand product as a compressor tree (an ASIC build)\n\n  reg  [3:0] state;\n  localparam IDLE = 4\'d0, DECODE = 4\'d1, MUL = 4\'d2, ROUND = 4\'d3, PRENORM = 4\'d4, PRODUCT = 4\'d5, SUBP = 4\'d6,\n             MUL2 = 4\'d7;\n\n  // latched inputs (captured at issue; the shared operand ports stay stable while we stall)\n  reg  [63:0] a_l, b_l;\n\n  wire        sa = a_l[63], sb = b_l[63];\n  wire [10:0] ea = a_l[62:52], eb = b_l[62:52];\n  wire [51:0] ma = a_l[51:0],  mb = b_l[51:0];\n  wire        rs = sa ^ sb;                 // result sign\n\n  wire a_is_zero = (ea == 11\'d0)   && (ma == 52\'d0);\n  wire b_is_zero = (eb == 11\'d0)   && (mb == 52\'d0);\n  wire a_is_inf  = (ea == 11\'h7FF) && (ma == 52\'d0);\n  wire b_is_inf  = (eb == 11\'h7FF) && (mb == 52\'d0);\n  wire a_is_nan  = (ea == 11\'h7FF) && (ma != 52\'d0);\n  wire b_is_nan  = (eb == 11\'h7FF) && (mb != 52\'d0);\n\n  // significands with implicit bit\n  wire        ha = (ea != 11\'d0);\n  wire        hb = (eb != 11\'d0);\n  wire [52:0] siga = {ha, ma};\n  wire [52:0] sigb = {hb, mb};\n\n  reg  [105:0] prod;               // 53 x 53 = up to 106 bits\n  reg  signed [13:0] exp_sum;      // combined unbiased exponent (wide, signed)\n  reg  [13:0] lza, lzb;            // leading-zero adjustments for subnormals\n  reg  [52:0] sa_sh, sb_sh;        // pre-normalized significands (leading 1 at bit 52), registered\n  reg  [105:0] norm;\n  reg  signed [13:0] e;\n  reg  [55:0] mant;                // 53-bit + GRS\n  reg  g, r, s, l, up;\n  reg  [63:0] rnd;\n  reg  [52:0] fsig;\n  reg  signed [13:0] fe;\n  reg  [13:0] rsh;\n  reg  [127:0] wide;\n  reg  g2,r2,s2,l2,up2;\n  reg  [63:0] rr;\n  reg  signed [13:0] bexp;\n  reg          rs_l;               // result sign, registered for the ROUND stage\n  reg  [13:0]  lza_r, lzb_r;       // leading-zero counts, for PRENORM\n  reg  [63:0]  wide_r;             // a subnormal product on its grid (at most 56 bits), for SUBP\n\n  // The 53\u{d7}53 significand product, two clocks deep: sa_sh and sb_sh, written in PRENORM and held through PRODUCT,\n  // enter it at MUL and its product is read in PRODUCT. MUL_TREE (an ASIC build) makes it a compressor tree, its\n  // first four levels registered at the end of MUL and the rest at the end of MUL2 before a Kogge-Stone add; otherwise\n  // (an FPGA build) it is the `*` a DSP block takes, registered twice.\n  wire [105:0] tree_p;\n  generate\n    if (MUL_TREE) begin : tree\n      fp_csa_mul #(.AW(53), .BW(53), .CUT(4), .FINAL(1)) product_tree(.clk(clk), .a(sa_sh), .b(sb_sh), .p(tree_p));\n    end else begin : dsp\n      // One register: the operands hold through MUL2, so the product is stable when PRODUCT reads it.\n      reg [105:0] p1;\n      always @(posedge clk) p1 <= sa_sh * sb_sh;\n      assign tree_p = p1;\n    end\n  endgenerate\n\n  // The highest set bit of a nonzero word, by halving: six levels, not a priority chain.\n  function [5:0] lead_one64(input [63:0] m);\n    reg [31:0] h5; reg [15:0] h4; reg [7:0] h3; reg [3:0] h2; reg [1:0] h1; reg [5:0] z;\n    begin\n      z[5] = ~|m[63:32];  h5 = z[5] ? m[31:0]  : m[63:32];\n      z[4] = ~|h5[31:16]; h4 = z[4] ? h5[15:0] : h5[31:16];\n      z[3] = ~|h4[15:8];  h3 = z[3] ? h4[7:0]  : h4[15:8];\n      z[2] = ~|h3[7:4];   h2 = z[2] ? h3[3:0]  : h3[7:4];\n      z[1] = ~|h2[3:2];   h1 = z[1] ? h2[1:0]  : h2[3:2];\n      z[0] = ~h1[1];\n      lead_one64 = ~z;\n    end\n  endfunction\n\n  // Leading zeros of a 53-bit significand, 53 when it is zero.\n  function [13:0] clz53(input [52:0] v);\n    begin\n      clz53 = (v == 53\'d0) ? 14\'d53 : 14\'d52 - {8\'d0, lead_one64({11\'d0, v})};\n    end\n  endfunction\n\n  // v + cin through a parallel-prefix AND: log2(64) levels. cin is the scan\'s bit 0, so bit i flips on cin and\n  // every bit below it, and cin reaches the bits through the scan\'s tree rather than one net to all 64.\n  function [63:0] inc64(input [63:0] v, input cin);\n    reg [63:0] p; integer d;\n    begin\n      p = {v[62:0], cin};\n      for (d = 1; d < 64; d = d * 2) p = p & ((p << d) | ((64\'d1 << d) - 64\'d1));\n      inc64 = v ^ p;\n    end\n  endfunction\n\n  // A word\'s bits below `s` as one flag: shifted to the top and tested, no (1 << s) - 1 mask, no borrow chain.\n  function below(input [127:0] x, input [13:0] s);\n    begin\n      below = (s == 14\'d0) ? 1\'b0 : ((x << (14\'d128 - s)) != 128\'d0);\n    end\n  endfunction\n\n\n  initial begin state = IDLE; done = 1\'b0; y = 64\'d0; end\n\n  always @(posedge clk) begin\n    if (rst) begin state <= IDLE; done <= 1\'b0; y <= 64\'d0; end\n    else case (state)\n      IDLE: if (start) begin\n        done <= 1\'b0; a_l <= a; b_l <= b; state <= DECODE;\n      end\n      DECODE: begin\n        rs_l <= rs;\n        if (a_is_nan || b_is_nan) begin y <= QNAN; done <= 1\'b1; state <= IDLE; end\n        else if ((a_is_inf && b_is_zero) || (b_is_inf && a_is_zero)) begin y <= QNAN; done <= 1\'b1; state <= IDLE; end\n        else if (a_is_inf || b_is_inf) begin y <= {rs, 11\'h7FF, 52\'d0}; done <= 1\'b1; state <= IDLE; end\n        else if (a_is_zero || b_is_zero) begin y <= {rs, 63\'d0}; done <= 1\'b1; state <= IDLE; end\n        else begin\n          // unbiased exponents; leading-zero counts (this stage isolates the two clz53 cones + pre-shift).\n          lza_r <= clz53(siga);\n          lzb_r <= clz53(sigb);\n          state <= PRENORM;\n        end\n      end\n      PRENORM: begin\n        lza = lza_r;\n        lzb = lzb_r;\n        exp_sum <= (ea != 0 ? $signed({3\'d0, ea}) - 14\'sd1023 : -14\'sd1022)\n                 + (eb != 0 ? $signed({3\'d0, eb}) - 14\'sd1023 : -14\'sd1022)\n                 - $signed(lza) - $signed(lzb);\n        sa_sh <= (siga << lza);   // leading 1 at bit 52\n        sb_sh <= (sigb << lzb);\n        state <= MUL;\n      end\n      MUL: begin\n        // \u{2500}\u{2500} the 53\u{d7}53 product\'s compressor tree, its first four levels (registered inside product_tree) \u{2500}\u{2500}\n        state <= MUL2;\n      end\n      MUL2: begin\n        // \u{2500}\u{2500} its remaining levels (registered inside product_tree) \u{2500}\u{2500}\n        state <= PRODUCT;\n      end\n      PRODUCT: begin\n        // \u{2500}\u{2500} the tree\'s last two words added + top-bit normalize \u{2500}\u{2500}\n        prod = tree_p;\n        if (prod[105]) begin\n          norm <= prod;\n          e    <= exp_sum + 14\'sd1; // leading at 105 means value in [2,4) relative -> exp+1\n        end else begin\n          norm <= prod << 1;        // shift so leading 1 sits at bit 105\n          e    <= exp_sum;\n        end\n        state <= ROUND;\n      end\n      ROUND: begin\n        // Now norm[105] == 1 and its (would-be normal) biased exponent is `e`. Decide underflow BEFORE\n        // rounding \u{2014} rounding the exact product exactly once avoids the double-rounding bug.\n        if (e < -14\'sd1022) begin\n          // \u{2500}\u{2500} subnormal / underflow: the EXACT product shifted onto the subnormal grid, SUBP rounds it once \u{2500}\u{2500}\n          rsh = (-14\'sd1022) - e;                    // extra right-shift into the subnormal grid (>= 1)\n          wide = (norm >> (14\'d105 - 14\'d55)) | {127\'d0, below(norm, 14\'d105 - 14\'d55)};  // norm[105]->bit55\n          if (rsh >= 14\'d128)       wide_r <= {63\'d0, |wide};\n          else if (rsh != 14\'d0)    wide_r <= (wide >> rsh) | {127\'d0, below(wide, rsh)};\n          else                      wide_r <= wide;\n          state <= SUBP;\n        end else begin\n          // \u{2500}\u{2500} normal path: round the 53-bit significand \u{2500}\u{2500}\n          mant = {norm[105:53], (norm[52]), (norm[51]), (|norm[50:0])};  // 53 + g + r + s = 56 bits\n          g = mant[2]; r = mant[1]; s = mant[0]; l = mant[3];\n          up = g & (r | s | l);\n          rnd = inc64({11\'d0, mant[55:3]}, up);\n          if (rnd[53]) begin fsig = rnd[53:1]; fe = e + 14\'sd1; end\n          else         begin fsig = rnd[52:0]; fe = e; end\n          if (fe > 14\'sd1023) begin\n            y <= {rs_l, 11\'h7FF, 52\'d0};                          // overflow -> inf\n          end else if (fe < -14\'sd1022) begin\n            y <= {rs_l, 11\'d1, fsig[51:0]};\n          end else begin\n            bexp = fe + 14\'sd1023;\n            y <= {rs_l, bexp[10:0], fsig[51:0]};\n          end\n          done <= 1\'b1; state <= IDLE;\n        end\n      end\n      SUBP: begin\n        g2 = wide_r[2]; r2 = wide_r[1]; s2 = wide_r[0]; l2 = wide_r[3];\n        up2 = g2 & (r2 | s2 | l2);\n        rr = inc64({11\'d0, wide_r[55:3]}, up2);\n        if (rr[52]) y <= {rs_l, 11\'d1, rr[51:0]};   // rounded up to the smallest normal\n        else        y <= {rs_l, 11\'d0, rr[51:0]};\n        done <= 1\'b1; state <= IDLE;\n      end\n    endcase\n  end\nendmodule\n\n// fp_csa_mul \u{2014} an unsigned AW\u{d7}BW multiplier as a tree of 3:2 compressors over its BW partial-product rows. The\n// compressor outputs register after CUT levels (CUT = 0 or beyond the tree\'s depth: no register inside it), the last\n// two words register again when FINAL is set, and those two words resolve through a Kogge-Stone add. `p` is the\n// product of the operands sampled as many clocks earlier as there are registers on its path.\nmodule fp_csa_mul #(parameter integer AW = 53, parameter integer BW = 53, parameter integer CUT = 4,\n                 parameter integer FINAL = 1) (\n  input                clk,\n  input  [AW-1:0]      a,\n  input  [BW-1:0]      b,\n  output [AW+BW-1:0]   p\n);\n  localparam integer W = AW + BW;\n\n  // The words entering tree level `level`: each level turns every three words into two and passes the rest.\n  function integer words(input integer level);\n    integer n, l;\n    begin\n      n = BW;\n      for (l = 0; l < level; l = l + 1) n = 2 * (n / 3) + (n % 3);\n      words = n;\n    end\n  endfunction\n\n  function integer depth(input integer unused);\n    integer n, l;\n    begin\n      n = BW; l = 0;\n      while (n > 2) begin n = 2 * (n / 3) + (n % 3); l = l + 1; end\n      depth = l;\n    end\n  endfunction\n\n  localparam integer L = depth(0);\n\n  // Level-major bus: level `l` holds its words at word slots l*BW .. l*BW + words(l) - 1.\n  wire [W*BW*(L+1)-1:0] lv;\n\n  genvar i, l, g;\n  generate\n    for (i = 0; i < BW; i = i + 1) begin : row\n      assign lv[i*W +: W] = b[i] ? ({{BW{1\'b0}}, a} << i) : {W{1\'b0}};\n    end\n    for (l = 0; l < L; l = l + 1) begin : level\n      localparam integer N = words(l);\n      localparam integer G = N / 3;\n      for (g = 0; g < G; g = g + 1) begin : compress\n        wire [W-1:0] x = lv[((l*BW) + 3*g) * W +: W];\n        wire [W-1:0] y = lv[((l*BW) + 3*g + 1) * W +: W];\n        wire [W-1:0] z = lv[((l*BW) + 3*g + 2) * W +: W];\n        wire [W-1:0] xy = x ^ y;\n        wire [W-1:0] s = xy ^ z;\n        wire [W-1:0] c = ((x & y) | (z & xy)) << 1;\n        if (l + 1 == CUT) begin : registered\n          reg [W-1:0] sr, cr;\n          always @(posedge clk) begin sr <= s; cr <= c; end\n          assign lv[(((l+1)*BW) + 2*g) * W +: W] = sr;\n          assign lv[(((l+1)*BW) + 2*g + 1) * W +: W] = cr;\n        end else begin : wired\n          assign lv[(((l+1)*BW) + 2*g) * W +: W] = s;\n          assign lv[(((l+1)*BW) + 2*g + 1) * W +: W] = c;\n        end\n      end\n      for (g = 0; g < N % 3; g = g + 1) begin : carry_over\n        if (l + 1 == CUT) begin : registered\n          reg [W-1:0] r;\n          always @(posedge clk) r <= lv[((l*BW) + 3*G + g) * W +: W];\n          assign lv[(((l+1)*BW) + 2*G + g) * W +: W] = r;\n        end else begin : wired\n          assign lv[(((l+1)*BW) + 2*G + g) * W +: W] = lv[((l*BW) + 3*G + g) * W +: W];\n        end\n      end\n    end\n  endgenerate\n\n  wire [W-1:0] s_word = lv[(L*BW) * W +: W];\n  wire [W-1:0] c_word;\n  generate\n    if (BW == 1) begin : one_row\n      assign c_word = {W{1\'b0}};\n    end else begin : two_rows\n      assign c_word = lv[(L*BW + 1) * W +: W];\n    end\n  endgenerate\n  reg  [W-1:0] s_r, c_r;\n  always @(posedge clk) begin s_r <= s_word; c_r <= c_word; end\n  wire [W-1:0] s_in = FINAL ? s_r : s_word;\n  wire [W-1:0] c_in = FINAL ? c_r : c_word;\n\n  // The Kogge-Stone add of the last two words: group generate and propagate doubled in span each level.\n  function [W-1:0] prefix_add(input [W-1:0] u, input [W-1:0] v);\n    reg [W-1:0] gg, pp, p0; integer d;\n    begin\n      gg = u & v; p0 = u ^ v; pp = p0;\n      for (d = 1; d < W; d = d * 2) begin\n        gg = gg | (pp & (gg << d));\n        pp = pp & (pp << d);\n      end\n      prefix_add = p0 ^ (gg << 1);\n    end\n  endfunction\n\n  assign p = prefix_add(s_in, c_in);\nendmodule\n\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\n// fp_cmp \u{2014} combinational IEEE-754 binary64 less-than / equal (returns 1 bit each).\n//   lt: a < b (false if either is NaN, -0 == +0)\n//   eq: a == b (false if either is NaN, -0 == +0)\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\nmodule fp_cmp(\n  input  [63:0] a,\n  input  [63:0] b,\n  output        lt,\n  output        eq\n);\n  wire [10:0] ea = a[62:52], eb = b[62:52];\n  wire [51:0] ma = a[51:0],  mb = b[51:0];\n  wire a_is_nan = (ea == 11\'h7FF) && (ma != 52\'d0);\n  wire b_is_nan = (eb == 11\'h7FF) && (mb != 52\'d0);\n  wire any_nan  = a_is_nan || b_is_nan;\n\n  wire a_is_zero = (a[62:0] == 63\'d0);   // \u{b1}0 : exp and frac all zero\n  wire b_is_zero = (b[62:0] == 63\'d0);\n  wire both_zero = a_is_zero && b_is_zero;   // -0 == +0\n\n  wire sa = a[63], sb = b[63];\n\n  // Map to a monotonic ordering key (the classic \"flip\" trick): for a two\'s-complement view of the\n  // bit pattern, positives keep order; negatives are reversed. Convert to an unsigned key where the\n  // numeric order matches: if sign bit set -> ~bits ; else -> bits | 0x8000...0000.\n  wire [63:0] ka = a[63] ? ~a : (a | 64\'h8000000000000000);\n  wire [63:0] kb = b[63] ? ~b : (b | 64\'h8000000000000000);\n\n  assign eq = any_nan ? 1\'b0 : (both_zero ? 1\'b1 : (a == b));\n  assign lt = any_nan ? 1\'b0 : (both_zero ? 1\'b0 : (ka < kb));\nendmodule\n\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\n// fp_div \u{2014} multi-cycle IEEE-754 binary64 divide, digit-recurrence (restoring), correctly rounded.\n//   start/done handshake. Computes a/b with a 55-bit quotient (53 + guard/round) plus a sticky\n//   from the final nonzero remainder, then RNE-rounds; each recurrence step\'s subtract is a parallel-prefix\n//   subtract, and a subnormal quotient takes its grid shift and its rounding in stages of their own.\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\nmodule fp_div(\n  input         clk,\n  input         rst,\n  input         start,\n  input  [63:0] a,\n  input  [63:0] b,\n  output reg [63:0] y,\n  output reg    done\n);\n  localparam [63:0] QNAN = 64\'h7FF8000000000000;\n\n  reg  [3:0] state;\n  // The exponent datapath is PIPELINED (SETUP\u{2192}SETUP2\u{2192}SETUP3 and FIN\u{2192}FIN2) so no single cycle carries a\n  // long combinational exponent tail. On silicon that tail was the last-settling node of a deep cone (clz +\n  // barrel-shift + wide signed add in SETUP; 128-bit shifters + round + biased-exponent add in FIN), and it\n  // dropped a carry \u{2014} the mantissa (settling earlier in the same cone) stayed exact while the exponent came\n  // out 2^32 small. Registering between the clz/shift/add sub-steps keeps every cone short, matching the\n  // already-silicon-correct fp_mul, at the cost of a few extra cycles (the unit already runs ~70).\n  localparam IDLE = 4\'d0, SETUP = 4\'d1, SETUP2 = 4\'d2, SETUP3 = 4\'d3, RUN = 4\'d4, FIN = 4\'d5, FIN2 = 4\'d6,\n             FINS = 4\'d7, FINP = 4\'d8;\n\n  // latched operand decode\n  reg  sa, sb, rs;\n  reg  [10:0] ea, eb;\n  reg  [51:0] ma, mb;\n  reg  special;                 // special-case result already computed\n  reg  [63:0] special_y;\n  reg  signed [13:0] exp_res;   // unbiased result exponent (of the leading quotient bit)\n\n  // digit-recurrence datapath: divide two normalized 53-bit significands.\n  reg  [55:0] num;              // remainder / working numerator, gets shifted left each step\n  reg  [52:0] den;              // divisor significand (leading 1 at bit 52)\n  reg  [55:0] quo;              // quotient bits accumulate here (56 = 53 + 3 GRS)\n  reg  [6:0]  cnt;              // step counter\n  reg  [72:0] diff;\n  reg  [105:0] shl_a, shl_b;    // normalized significands (leading 1 at bit 52) \u{2014} temporaries\n  reg  signed [13:0] bexp;      // biased exponent temporary\n\n  reg  signed [13:0] lza, lzb;\n\n  // The highest set bit of a nonzero word, by halving: six levels, not a priority chain.\n  function [5:0] lead_one64(input [63:0] m);\n    reg [31:0] h5; reg [15:0] h4; reg [7:0] h3; reg [3:0] h2; reg [1:0] h1; reg [5:0] z;\n    begin\n      z[5] = ~|m[63:32];  h5 = z[5] ? m[31:0]  : m[63:32];\n      z[4] = ~|h5[31:16]; h4 = z[4] ? h5[15:0] : h5[31:16];\n      z[3] = ~|h4[15:8];  h3 = z[3] ? h4[7:0]  : h4[15:8];\n      z[2] = ~|h3[7:4];   h2 = z[2] ? h3[3:0]  : h3[7:4];\n      z[1] = ~|h2[3:2];   h1 = z[1] ? h2[1:0]  : h2[3:2];\n      z[0] = ~h1[1];\n      lead_one64 = ~z;\n    end\n  endfunction\n\n  // Leading zeros of a 53-bit significand, 53 when it is zero.\n  function [13:0] clz53(input [52:0] v);\n    begin\n      clz53 = (v == 53\'d0) ? 14\'d53 : 14\'d52 - {8\'d0, lead_one64({11\'d0, v})};\n    end\n  endfunction\n\n  // v + cin through a parallel-prefix AND: log2(64) levels. cin is the scan\'s bit 0, so bit i flips on cin and\n  // every bit below it, and cin reaches the bits through the scan\'s tree rather than one net to all 64.\n  function [63:0] inc64(input [63:0] v, input cin);\n    reg [63:0] p; integer d;\n    begin\n      p = {v[62:0], cin};\n      for (d = 1; d < 64; d = d * 2) p = p & ((p << d) | ((64\'d1 << d) - 64\'d1));\n      inc64 = v ^ p;\n    end\n  endfunction\n\n  // A word\'s bits below `s` as one flag: shifted to the top and tested, no (1 << s) - 1 mask, no borrow chain.\n  function below(input [127:0] x, input [13:0] s);\n    begin\n      below = (s == 14\'d0) ? 1\'b0 : ((x << (14\'d128 - s)) != 128\'d0);\n    end\n  endfunction\n\n  // x - y over 72 bits as a Kogge-Stone prefix subtract; bit 72 of the result is the borrow (1 when x < y).\n  function [72:0] sub72(input [71:0] x, input [71:0] yy);\n    reg [71:0] yb, g, p, p0; integer d;\n    begin\n      yb = ~yy;\n      p0 = x ^ yb;\n      g = (x & yb) | (p0 & 72\'d1);\n      p = p0;\n      for (d = 1; d < 72; d = d * 2) begin\n        g = g | (p & (g << d));\n        p = p & (p << d);\n      end\n      sub72 = {~g[71], p0 ^ {g[70:0], 1\'b1}};\n    end\n  endfunction\n\n  // rounding / packing registers\n  reg  g,r,s,l,up;\n  reg  [63:0] rnd;\n  reg  [52:0] fsig;\n  reg  signed [13:0] fe;\n  reg  [13:0] rsh;\n  reg  [127:0] wide;\n  reg  g2,r2,s2,l2,up2;\n  reg  [63:0] rr;\n\n  // FIN pipeline: the final biased-exponent add (fe_r + 1023) lives alone in FIN2, off a register, so the\n  // exponent output is never the tail of the deep FIN cone. normal_path picks the FIN2 assembly; the\n  // subnormal / overflow / special results (constant exponent fields, no add) are pre-packed into y_reg.\n  reg  [63:0] y_reg;\n  reg          normal_path;\n  reg  signed [13:0] fe_r;\n\n  initial begin state = IDLE; done = 1\'b0; y = 64\'d0; end\n\n  wire a_is_zero = (ea == 11\'d0)   && (ma == 52\'d0);\n  wire b_is_zero = (eb == 11\'d0)   && (mb == 52\'d0);\n  wire a_is_inf  = (ea == 11\'h7FF) && (ma == 52\'d0);\n  wire b_is_inf  = (eb == 11\'h7FF) && (mb == 52\'d0);\n  wire a_is_nan  = (ea == 11\'h7FF) && (ma != 52\'d0);\n  wire b_is_nan  = (eb == 11\'h7FF) && (mb != 52\'d0);\n\n  reg  [127:0] qwin;   // canonical quotient window for the subnormal path\n  reg  [63:0]  qwin_r; // that window (at most 57 bits), held for the subnormal stages\n  reg  [63:0]  den_win;// that window on the subnormal grid, its dropped bits folded into sticky\n  reg  [13:0]  rsh_r;  // the subnormal grid\'s extra right shift\n  reg          to_subnormal;\n  reg  [11:0]  qlead;\n\n  always @(posedge clk) begin\n    if (rst) begin state <= IDLE; done <= 1\'b0; y <= 64\'d0; end\n    else case (state)\n      IDLE: if (start) begin\n        done <= 1\'b0;\n        sa <= a[63]; sb <= b[63]; rs <= a[63]^b[63];\n        ea <= a[62:52]; eb <= b[62:52]; ma <= a[51:0]; mb <= b[51:0];\n        state <= SETUP;\n      end\n      SETUP: begin\n        // \u{2500}\u{2500} cycle 1: special-case detection + count-leading-zeros only (register lza/lzb) \u{2500}\u{2500}\n        if (a_is_nan || b_is_nan)        begin special <= 1\'b1; special_y <= QNAN; state <= FIN; end\n        else if (a_is_inf && b_is_inf)   begin special <= 1\'b1; special_y <= QNAN; state <= FIN; end\n        else if (a_is_inf)               begin special <= 1\'b1; special_y <= {rs, 11\'h7FF, 52\'d0}; state <= FIN; end\n        else if (b_is_inf)               begin special <= 1\'b1; special_y <= {rs, 63\'d0}; state <= FIN; end  // finite/inf=0\n        else if (a_is_zero && b_is_zero) begin special <= 1\'b1; special_y <= QNAN; state <= FIN; end          // 0/0=NaN\n        else if (b_is_zero)              begin special <= 1\'b1; special_y <= {rs, 11\'h7FF, 52\'d0}; state <= FIN; end // x/0=inf\n        else if (a_is_zero)              begin special <= 1\'b1; special_y <= {rs, 63\'d0}; state <= FIN; end    // 0/x=0\n        else begin\n          special <= 1\'b0;\n          lza <= clz53({(ea!=0), ma});\n          lzb <= clz53({(eb!=0), mb});\n          state <= SETUP2;\n        end\n      end\n      SETUP2: begin\n        // \u{2500}\u{2500} cycle 2: barrel-shift the significands and form the unbiased exponent \u{2014} each its own short\n        // cone off the registered lza/lzb (no clz\u{2192}shift\u{2192}wide-add chained together in one cycle) \u{2500}\u{2500}\n        shl_a <= ({53\'d0, (ea!=0), ma}) << lza;\n        shl_b <= ({53\'d0, (eb!=0), mb}) << lzb;\n        exp_res <= (ea != 0 ? $signed({3\'d0, ea}) - 14\'sd1023 : -14\'sd1022)\n                 - (eb != 0 ? $signed({3\'d0, eb}) - 14\'sd1023 : -14\'sd1022)\n                 - $signed(lza) + $signed(lzb);\n        state <= SETUP3;\n      end\n      SETUP3: begin\n        // \u{2500}\u{2500} cycle 3: seed the digit-recurrence from the registered normalized significands \u{2500}\u{2500}\n        den <= shl_b[52:0];\n        num <= {3\'d0, shl_a[52:0]};   // dividend significand (leading 1 at bit 52), room above to shift\n        quo <= 56\'d0;\n        cnt <= 7\'d56;                 // produce 56 quotient bits (53 + 3 GRS), first bit at position 55\n        state <= RUN;\n      end\n      RUN: begin\n        if (cnt != 7\'d0) begin\n          // \u{2500}\u{2500} one restoring-division step: compare THEN shift, so the first emitted bit is the 2^0\n          // digit of the significand ratio (both operands in [1,2), so the ratio is in (0.5,2)). \u{2500}\u{2500}\n          diff = sub72({16\'d0, num}, {19\'d0, den});\n          if (!diff[72]) begin                   // num >= den -> quotient bit 1, subtract\n            num <= {diff[54:0], 1\'b0};            // (num - den) << 1\n            quo <= {quo[54:0], 1\'b1};\n          end else begin\n            num <= {num[54:0], 1\'b0};             // num << 1\n            quo <= {quo[54:0], 1\'b0};\n          end\n          cnt <= cnt - 7\'d1;\n        end else begin\n          state <= FIN;\n        end\n      end\n      FIN: begin\n        // \u{2500}\u{2500} stage 1: normalize + round the quotient; register the significand (fsig) and unbiased\n        // exponent (fe_r). Overflow and special results carry constant exponent fields and are pre-packed into\n        // y_reg; a subnormal result goes through FINS and FINP. The biased-exponent ADD is deferred to FIN2. \u{2500}\u{2500}\n        if (special) begin\n          y_reg <= special_y; normal_path <= 1\'b0;\n        end else begin\n          // quotient of two [1,2) significands is in (0.5,2). Its leading 1 is at bit 55 (if >=1) or\n          // bit 54. Normalize so the leading 1 sits at bit 55 (window [55:3]=53-bit sig, [2:0]=GRS),\n          // with a non-empty division remainder `num` OR\'d into sticky.\n          if (quo[55]) begin\n            g = quo[2]; r = quo[1]; s = (quo[0] | (num != 56\'d0)); l = quo[3];\n            up = g & (r | s | l);\n            rnd = inc64({11\'d0, quo[55:3]}, up);\n            fe  = exp_res;\n            qwin = {72\'d0, quo};                 // leading 1 at bit 55 already\n            if (num != 56\'d0) qwin[0] = 1\'b1;\n          end else begin\n            g = quo[1]; r = quo[0]; s = (num != 56\'d0); l = quo[2];\n            up = g & (r | s | l);\n            rnd = inc64({11\'d0, quo[54:2]}, up);\n            fe  = exp_res - 14\'sd1;\n            qwin = {72\'d0, quo} << 1;            // shift leading 1 from bit 54 up to bit 55\n            if (num != 56\'d0) qwin[0] = 1\'b1;\n          end\n          if (rnd[53]) begin fsig <= rnd[53:1]; fe = fe + 14\'sd1; end\n          else         begin fsig <= rnd[52:0]; end\n\n          // pre-round exponent (of the leading quotient bit) drives the underflow decision so the\n          // subnormal path rounds the exact quotient window exactly ONCE (no double rounding).\n          if (quo[55]) bexp = exp_res; else bexp = exp_res - 14\'sd1;\n\n          if (bexp < -14\'sd1022) begin\n            // A subnormal quotient: the un-rounded window (leading 1 at bit 55) goes to FINS and FINP, which shift\n            // it onto the subnormal grid and round it once.\n            rsh_r  <= (-14\'sd1022) - bexp;\n            qwin_r <= qwin;\n            normal_path <= 1\'b0;\n            to_subnormal <= 1\'b1;\n          end else if (fe > 14\'sd1023) begin\n            y_reg <= {rs, 11\'h7FF, 52\'d0};                // overflow -> inf\n            normal_path <= 1\'b0;\n            to_subnormal <= 1\'b0;\n          end else begin\n            fe_r <= fe;                                   // registered exponent for FIN2\'s biased add\n            normal_path <= 1\'b1;\n            to_subnormal <= 1\'b0;\n          end\n        end\n        state <= (!special && bexp < -14\'sd1022) ? FINS : FIN2;\n      end\n      FINS: begin\n        if (rsh_r >= 14\'d128)    den_win <= {63\'d0, |qwin_r};\n        else if (rsh_r != 14\'d0) den_win <= (qwin_r >> rsh_r) | {127\'d0, below({64\'d0, qwin_r}, rsh_r)};\n        else                     den_win <= qwin_r;\n        state <= FINP;\n      end\n      FINP: begin\n        g2 = den_win[2]; r2 = den_win[1]; s2 = den_win[0]; l2 = den_win[3];\n        up2 = g2 & (r2 | s2 | l2);\n        rr = inc64({11\'d0, den_win[55:3]}, up2);\n        if (rr[52]) y_reg <= {rs, 11\'d1, rr[51:0]};       // rounded up to the smallest normal\n        else        y_reg <= {rs, 11\'d0, rr[51:0]};\n        state <= FIN2;\n      end\n      FIN2: begin\n        // \u{2500}\u{2500} stage 2: the sole remaining exponent op \u{2014} a short (fe_r + 1023) biased-exponent add off a\n        // register, then pack. This is what keeps the exponent out of the deep FIN cone on silicon. \u{2500}\u{2500}\n        if (normal_path) begin\n          bexp = fe_r + 14\'sd1023;\n          y <= {rs, bexp[10:0], fsig[51:0]};\n        end else begin\n          y <= y_reg;\n        end\n        done <= 1\'b1; state <= IDLE;\n      end\n    endcase\n  end\nendmodule\n\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\n// fp_sqrt \u{2014} multi-cycle IEEE-754 binary64 square root, digit-recurrence (restoring, non-restoring\n//   style bit-at-a-time), correctly rounded (RNE). start/done handshake.\n//   sqrt(neg) -> NaN, sqrt(+inf)=+inf, sqrt(\u{b1}0)=\u{b1}0, sqrt(NaN)=NaN. Each step\'s trial subtract is a\n//   parallel-prefix subtract; the leading-zero count, the normalizing shift and the round each take a stage.\n// \u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\u{2500}\nmodule fp_sqrt(\n  input         clk,\n  input         rst,\n  input         start,\n  input  [63:0] a,\n  output reg [63:0] y,\n  output reg    done\n);\n  localparam [63:0] QNAN = 64\'h7FF8000000000000;\n\n  reg  [2:0] state;\n  localparam IDLE = 3\'d0, SETUP = 3\'d1, RUN = 3\'d2, FIN = 3\'d3, NORM = 3\'d4, PACK = 3\'d5;\n\n  reg  [10:0] ea;\n  reg  [51:0] ma;\n  reg  asign;                   // input sign (only matters for -0)\n  reg  special;\n  reg  [63:0] special_y;\n\n  // Digit-recurrence sqrt. The radicand M (= significand, or 2\u{d7}significand for an odd input exponent)\n  // is in [1,4). We place M in a wide fixed-point lane with its binary point just below the top two\n  // integer bits, then peel 2 radicand bits per root bit. `root` accumulates the root; the leading root\n  // bit lands at a fixed position so we can extract sig/guard/round/sticky directly.\n  reg  signed [13:0] exp_in;    // unbiased exponent of the input\'s leading significand bit\n  reg  signed [13:0] exp_res;   // result exponent (floor(exp_in/2))\n  reg  [131:0] rad;             // radicand lane (2 bits consumed per step)\n  reg  [66:0]  root;            // accumulated root bits (leading root bit ends at a fixed position)\n  reg  [67:0]  rem;             // running remainder\n  reg  [7:0]   cnt;\n\n  reg  [67:0] trial, test;\n  reg  signed [13:0] lz;\n  reg  [13:0]  lz_r;        // leading zeros of a subnormal input, for NORM\n  reg  [72:0]  dt;\n  reg  [63:0]  rnd_r;       // the rounded significand, for PACK\n  // The highest set bit of a nonzero word, by halving: six levels, not a priority chain.\n  function [5:0] lead_one64(input [63:0] m);\n    reg [31:0] h5; reg [15:0] h4; reg [7:0] h3; reg [3:0] h2; reg [1:0] h1; reg [5:0] z;\n    begin\n      z[5] = ~|m[63:32];  h5 = z[5] ? m[31:0]  : m[63:32];\n      z[4] = ~|h5[31:16]; h4 = z[4] ? h5[15:0] : h5[31:16];\n      z[3] = ~|h4[15:8];  h3 = z[3] ? h4[7:0]  : h4[15:8];\n      z[2] = ~|h3[7:4];   h2 = z[2] ? h3[3:0]  : h3[7:4];\n      z[1] = ~|h2[3:2];   h1 = z[1] ? h2[1:0]  : h2[3:2];\n      z[0] = ~h1[1];\n      lead_one64 = ~z;\n    end\n  endfunction\n\n  // Leading zeros of a 53-bit significand, 53 when it is zero.\n  function [13:0] clz53(input [52:0] v);\n    begin\n      clz53 = (v == 53\'d0) ? 14\'d53 : 14\'d52 - {8\'d0, lead_one64({11\'d0, v})};\n    end\n  endfunction\n\n  // v + cin through a parallel-prefix AND: log2(64) levels. cin is the scan\'s bit 0, so bit i flips on cin and\n  // every bit below it, and cin reaches the bits through the scan\'s tree rather than one net to all 64.\n  function [63:0] inc64(input [63:0] v, input cin);\n    reg [63:0] p; integer d;\n    begin\n      p = {v[62:0], cin};\n      for (d = 1; d < 64; d = d * 2) p = p & ((p << d) | ((64\'d1 << d) - 64\'d1));\n      inc64 = v ^ p;\n    end\n  endfunction\n\n  // A word\'s bits below `s` as one flag: shifted to the top and tested, no (1 << s) - 1 mask, no borrow chain.\n  function below(input [127:0] x, input [13:0] s);\n    begin\n      below = (s == 14\'d0) ? 1\'b0 : ((x << (14\'d128 - s)) != 128\'d0);\n    end\n  endfunction\n\n  // x - y over 72 bits as a Kogge-Stone prefix subtract; bit 72 of the result is the borrow (1 when x < y).\n  function [72:0] sub72(input [71:0] x, input [71:0] yy);\n    reg [71:0] yb, g, p, p0; integer d;\n    begin\n      yb = ~yy;\n      p0 = x ^ yb;\n      g = (x & yb) | (p0 & 72\'d1);\n      p = p0;\n      for (d = 1; d < 72; d = d * 2) begin\n        g = g | (p & (g << d));\n        p = p & (p << d);\n      end\n      sub72 = {~g[71], p0 ^ {g[70:0], 1\'b1}};\n    end\n  endfunction\n\n  reg  g,r,rr_s,l,up;\n  reg  [63:0] rnd;\n  reg  [52:0] fsig;\n  reg  signed [13:0] fe;\n  reg  signed [13:0] bexp;\n  reg  [52:0] msig;             // normalized significand (leading 1 at bit 52)\n\n  // Number of root bits we generate. We need 53 significand + guard + round + a sticky, so 56 useful\n  // bits; produce a couple extra for the sticky tail. NROOT steps consume 2*NROOT radicand bits.\n  localparam integer NROOT = 59;\n\n  initial begin state = IDLE; done = 1\'b0; y = 64\'d0; end\n\n  wire a_is_zero = (ea == 11\'d0)   && (ma == 52\'d0);\n  wire a_is_inf  = (ea == 11\'h7FF) && (ma == 52\'d0);\n  wire a_is_nan  = (ea == 11\'h7FF) && (ma != 52\'d0);\n\n  always @(posedge clk) begin\n    if (rst) begin state <= IDLE; done <= 1\'b0; y <= 64\'d0; end\n    else case (state)\n      IDLE: if (start) begin\n        done <= 1\'b0; asign <= a[63]; ea <= a[62:52]; ma <= a[51:0]; state <= SETUP;\n      end\n      SETUP: begin\n        if (a_is_nan)                    begin special <= 1\'b1; special_y <= QNAN; state <= FIN; end\n        else if (asign && !a_is_zero)    begin special <= 1\'b1; special_y <= QNAN; state <= FIN; end  // sqrt(neg)=NaN\n        else if (a_is_inf)               begin special <= 1\'b1; special_y <= 64\'h7FF0000000000000; state <= FIN; end\n        else if (a_is_zero)              begin special <= 1\'b1; special_y <= {asign, 63\'d0}; state <= FIN; end // sqrt(\u{b1}0)=\u{b1}0\n        else begin\n          special <= 1\'b0;\n          lz_r <= clz53({(ea!=0), ma});\n          state <= NORM;\n        end\n      end\n      NORM: begin\n        begin\n          lz = lz_r;\n          msig = ({(ea!=0), ma} << lz);      // leading 1 at bit 52, value m in [1,2)\n          exp_in = (ea != 0 ? $signed({3\'d0, ea}) - 14\'sd1023 : -14\'sd1022) - $signed(lz);\n          // Radicand M = m (even exp) or 2m (odd exp), in [1,4). Place M\'s binary point so the top two\n          // lane bits [131:130] hold the integer part. m has its leading 1 at bit 52; shift it up to\n          // bit 130 (for even) or 131 (for odd, = 2m). Then process pairs from the top down.\n          if (exp_in[0]) begin\n            rad <= {msig, 79\'d0};             // 2m : leading 1 at bit 131 (52+79)\n            exp_res <= (exp_in - 14\'sd1) >>> 1;\n          end else begin\n            rad <= {1\'b0, msig, 78\'d0};       // m  : leading 1 at bit 130 (52+78)\n            exp_res <= exp_in >>> 1;\n          end\n          root <= 67\'d0; rem <= 68\'d0; cnt <= NROOT[7:0];\n          state <= RUN;\n        end\n      end\n      RUN: begin\n        if (cnt != 8\'d0) begin\n          // \u{2500}\u{2500} one restoring sqrt step \u{2500}\u{2500}\n          test  = {rem[65:0], rad[131:130]};   // bring down the top 2 radicand bits\n          trial = {root[65:0], 2\'b01};          // 2*root + 1  at this scale\n          dt = sub72({4\'d0, test}, {4\'d0, trial});\n          if (!dt[72]) begin                    // test >= trial\n            rem  <= dt[67:0];\n            root <= {root[65:0], 1\'b1};\n          end else begin\n            rem  <= test;\n            root <= {root[65:0], 1\'b0};\n          end\n          rad <= {rad[129:0], 2\'b00};           // consume the two processed bits\n          cnt <= cnt - 8\'d1;\n        end else state <= FIN;\n      end\n      FIN: begin\n        if (special) begin\n          y <= special_y; done <= 1\'b1; state <= IDLE;\n        end else begin\n          // After NROOT steps the root is an NROOT-bit integer approximating sqrt(M) scaled up. The root\n          // value is in [1,2) (since M in [1,4)), so its leading 1 is at bit (NROOT-1). We want the\n          // 53-bit significand at the top: significand = root[NROOT-1 : NROOT-53], guard = root[NROOT-54],\n          // round = root[NROOT-55], sticky = OR of the rest and the nonzero remainder.\n          // With NROOT=59: leading at bit 58; sig=root[58:6]; guard=root[5]; round=root[4];\n          // sticky = |root[3:0] | (rem!=0).\n          g    = root[5];\n          r    = root[4];\n          rr_s = (|root[3:0]) | (rem != 68\'d0);\n          l    = root[6];\n          up   = g & (r | rr_s | l);\n          rnd_r <= inc64({11\'d0, root[58:6]}, up);\n          state <= PACK;\n        end\n      end\n      PACK: begin\n        rnd = rnd_r;\n        if (rnd[53]) begin fsig = rnd[53:1]; fe = exp_res + 14\'sd1; end\n        else         begin fsig = rnd[52:0]; fe = exp_res; end\n        // sqrt of any finite positive f64 lands in the normal range (the smallest positive subnormal\n        // has sqrt ~ 2.2e-162, far above the subnormal threshold), so assemble directly.\n        bexp = fe + 14\'sd1023;\n        y <= {1\'b0, bexp[10:0], fsig[51:0]};\n        done <= 1\'b1; state <= IDLE;\n      end\n    endcase\n  end\nendmodule\n\n// fp_conv \u{2014} the f64 <-> i64 converters, FuncUnit::FFloorF (code 10, f64 -> i64 toward -inf) and FuncUnit::IntToF\n// (code 13, i64 -> f64 round-to-nearest-ties-even), as a five-stage pipeline: no 64-bit carry chain, leading-one\n// search or 64-bit shift shares a stage with another. `start` samples `x` and the direction into the unit\'s own\n// registers, as every unit samples its operands, so no stage\'s logic reads a port another unit shares; `y` and `done`\n// follow on the fifth clock. Bit-identical to the single-cycle conversions it replaces (== Rust `as` casts == the VM).\nmodule fp_conv(input clk, input rst, input start, input to_float, input [63:0] x, output reg [63:0] y, output reg done);\n  // Stage 0: the operand and the direction, as presented.\n  reg        v0, f0;\n  reg [63:0] x0;\n  // Stage 1: the sign and magnitude (IntToF), or the decoded fields (FFloor).\n  reg        v1, f1, s1, z1, small1, negsmall1;\n  reg [63:0] m1;\n  reg [10:0] ue1;\n  // Stage 2: the leading one (IntToF), or the significand shifted to the binary point (FFloor).\n  reg        v2, f2, s2, z2, small2, negsmall2, frac2;\n  reg [63:0] m2;\n  reg [5:0]  lead2;\n  // Stage 3: the 53-bit significand and its round-up (IntToF), or |floor| before the sign (FFloor).\n  reg        v3, f3, s3, z3, small3, negsmall3, up3;\n  reg [52:0] sig3;\n  reg [10:0] be3;\n  reg [63:0] mag3;\n\n  // The highest set bit of a nonzero word, by halving: each level asks whether the upper half is empty and keeps\n  // the half that holds the leading one, six levels deep rather than a 64-long priority chain.\n  function [5:0] lead_one(input [63:0] m);\n    reg [31:0] h5; reg [15:0] h4; reg [7:0] h3; reg [3:0] h2; reg [1:0] h1; reg [5:0] z;\n    begin\n      z[5] = ~|m[63:32];  h5 = z[5] ? m[31:0]  : m[63:32];\n      z[4] = ~|h5[31:16]; h4 = z[4] ? h5[15:0] : h5[31:16];\n      z[3] = ~|h4[15:8];  h3 = z[3] ? h4[7:0]  : h4[15:8];\n      z[2] = ~|h3[7:4];   h2 = z[2] ? h3[3:0]  : h3[7:4];\n      z[1] = ~|h2[3:2];   h1 = z[1] ? h2[1:0]  : h2[3:2];\n      z[0] = ~h1[1];\n      lead_one = ~z;\n    end\n  endfunction\n\n  // v + cin through a parallel-prefix AND of cin and v\'s low bits: log2(64) levels, never a ripple. cin is the\n  // scan\'s bit 0, so it reaches the bits through the scan\'s tree rather than one net to all 64.\n  function [63:0] inc64(input [63:0] v, input cin);\n    reg [63:0] p; integer d;\n    begin\n      p = {v[62:0], cin};\n      for (d = 1; d < 64; d = d * 2) p = p & ((p << d) | ((64\'d1 << d) - 64\'d1));\n      inc64 = v ^ p;\n    end\n  endfunction\n\n  wire [63:0] fsig = {11\'d0, 1\'b1, m1[51:0]};\n  wire        fbig = ue1 >= 11\'d52;\n  wire [63:0] fmag = fbig ? (fsig << (ue1 - 11\'d52)) : (fsig >> (11\'d52 - ue1));\n  wire        ffrac = !fbig && ((fsig << (ue1 + 11\'d12)) != 64\'d0);\n\n  wire        wide = lead2 > 6\'d52;\n  wire [5:0]  sh = lead2 - 6\'d52;\n  wire [63:0] kept = m2 >> sh;\n  wire [63:0] dropped = m2 << (7\'d64 - {1\'b0, sh});\n  wire [63:0] aligned = m2 << (6\'d52 - lead2);\n\n  wire [63:0] rsig = inc64({11\'d0, sig3}, up3);\n\n  always @(posedge clk) begin\n    if (rst) begin\n      v0 <= 1\'b0; v1 <= 1\'b0; v2 <= 1\'b0; v3 <= 1\'b0; done <= 1\'b0;\n    end else begin\n      v0 <= start;\n      if (start) begin\n        f0 <= to_float; x0 <= x;\n      end\n\n      v1 <= v0;\n      if (v0) begin\n        f1 <= f0;\n        s1 <= x0[63];\n        z1 <= x0 == 64\'d0;\n        m1 <= (f0 && x0[63]) ? inc64(~x0, 1\'b1) : x0;\n        small1 <= x0[62:52] < 11\'d1023;\n        negsmall1 <= x0[63] && !(x0[62:52] == 11\'d0 && x0[51:0] == 52\'d0);\n        ue1 <= x0[62:52] - 11\'d1023;\n      end\n\n      v2 <= v1; f2 <= f1; s2 <= s1; z2 <= z1; small2 <= small1; negsmall2 <= negsmall1;\n      if (f1) begin\n        m2 <= m1; lead2 <= lead_one(m1);\n      end else begin\n        m2 <= fmag; frac2 <= ffrac;\n      end\n\n      v3 <= v2; f3 <= f2; s3 <= s2; z3 <= z2; small3 <= small2; negsmall3 <= negsmall2;\n      if (f2) begin\n        be3 <= 11\'d1023 + {5\'d0, lead2};\n        if (wide) begin\n          sig3 <= kept[52:0];\n          up3 <= dropped[63] && ((|dropped[62:0]) || kept[0]);\n        end else begin\n          sig3 <= aligned[52:0];\n          up3 <= 1\'b0;\n        end\n      end else begin\n        mag3 <= inc64(m2, s2 && frac2);\n      end\n\n      done <= v3;\n      if (f3) begin\n        if (z3)            y <= 64\'d0;\n        else if (rsig[53]) y <= {s3, be3 + 11\'d1, 52\'d0};\n        else               y <= {s3, be3, rsig[51:0]};\n      end else begin\n        if (small3)        y <= negsmall3 ? 64\'hFFFF_FFFF_FFFF_FFFF : 64\'d0;\n        else if (s3)       y <= inc64(~mag3, 1\'b1);\n        else               y <= mag3;\n      end\n    end\n  end\nendmodule\n";
Expand description

The proven bit-exact sequential IEEE-754 f64 FPU (fpu.v beside this module; benchmarks/silicon/fpu.v is its golden copy for the FPGA flows) — fp_addsub/fp_mul/fp_div/ fp_sqrt (start/done handshake) + combinational fp_cmp. Appended verbatim after the engine module so the over-the-wire SoC can instantiate the float units; the engine’s 64-bit register file / heap already carry f64 bit patterns, so the FPU is the only float hardware the engine needs. Bit-exact vs Rust f64 over 41,521 vectors (its own verify-fpu.sh), and the exponent path is timing-safe (see the fp_div pipeline lock).