nx_u4096_mul_test.nx source

1// nx_u4096_mul_test.nx -- KAT for the 4096x4096 -> 8192 wide-multiply. 2// 3// Mirrors nx_u256_mul_test cases, scaled to 256-limb wide buffer. 4// All five algebraic-identity cases that prove a schoolbook 5// multiplier is correct: 6// - 0 * 0 = 0 7// - 1 * 1 = 1 8// - 3 * 5 = 15 (small * small) 9// - 0xFFFFFFFF * 0xFFFFFFFF = 0xFFFFFFFE00000001 (well-known 10// limb-saturation square, exercises carry from limb 0 to limb 1) 11// - 2^32 * 2^32 = 2^64 (single-bit cross-limb) 12// - 2^4096 - 1 squared (exercises full 256-limb wide range; high 13// limb 255 must be non-zero, exercising last carry-out) 14// - shl1 round-trip on wide buffer 15// 16// expect_exit: 0 17// license_tier: ORIGINAL 18 19import "nx_syscalls.nx" 20import "nx_u4096.nx" 21import "nx_u4096_mul.nx" 22 23func main() -> i64 { 24 let a: *i64 = u4096_alloc() 25 let b: *i64 = u4096_alloc() 26 let out: *i64 = u4096_wide_alloc() 27 28 // ---- Test A: 0 * 0 = 0 ---- 29 u4096_zero(a); u4096_zero(b) 30 u4096_mul_wide(out, a, b) 31 var i: i64 = 0 32 while i < NX_U4096_WIDE_LIMBS { 33 if out[i] != 0 { return 1 } 34 i = i + 1 35 } 36 37 // ---- Test B: 1 * 1 = 1 ---- 38 u4096_one(a); u4096_one(b) 39 u4096_mul_wide(out, a, b) 40 if out[0] != 1 { return 2 } 41 i = 1 42 while i < NX_U4096_WIDE_LIMBS { 43 if out[i] != 0 { return 3 } 44 i = i + 1 45 } 46 47 // ---- Test C: 3 * 5 = 15 ---- 48 u4096_zero(a); a[0] = 3 49 u4096_zero(b); b[0] = 5 50 u4096_mul_wide(out, a, b) 51 if out[0] != 15 { return 4 } 52 i = 1 53 while i < NX_U4096_WIDE_LIMBS { 54 if out[i] != 0 { return 5 } 55 i = i + 1 56 } 57 58 // ---- Test D: 0xFFFFFFFF * 0xFFFFFFFF = 0xFFFFFFFE00000001 ---- 59 u4096_zero(a); a[0] = 0xFFFFFFFF 60 u4096_zero(b); b[0] = 0xFFFFFFFF 61 u4096_mul_wide(out, a, b) 62 if out[0] != 0x00000001 { return 6 } 63 if out[1] != 0xFFFFFFFE { return 7 } 64 i = 2 65 while i < NX_U4096_WIDE_LIMBS { 66 if out[i] != 0 { return 8 } 67 i = i + 1 68 } 69 70 // ---- Test E: 2^32 * 2^32 = 2^64 (single-bit cross-limb) ---- 71 u4096_zero(a); a[1] = 1 72 u4096_zero(b); b[1] = 1 73 u4096_mul_wide(out, a, b) 74 if out[0] != 0 { return 10 } 75 if out[1] != 0 { return 11 } 76 if out[2] != 1 { return 12 } 77 i = 3 78 while i < NX_U4096_WIDE_LIMBS { 79 if out[i] != 0 { return 13 } 80 i = i + 1 81 } 82 83 // ---- Test F: a * 1 = a (identity on full-width operand) ---- 84 // Set a = some pattern across all limbs, b = 1, verify out_low == a 85 // and out_high all zero. 86 var k: i64 = 0 87 while k < NX_U4096_LIMBS { 88 a[k] = (k * 7 + 1) & NX_U4096_LIMB_MASK 89 k = k + 1 90 } 91 u4096_zero(b); b[0] = 1 92 u4096_mul_wide(out, a, b) 93 k = 0 94 while k < NX_U4096_LIMBS { 95 let expect: i64 = (k * 7 + 1) & NX_U4096_LIMB_MASK 96 if (out[k] & NX_U4096_LIMB_MASK) != expect { return 20 } 97 k = k + 1 98 } 99 k = NX_U4096_LIMBS 100 while k < NX_U4096_WIDE_LIMBS { 101 if out[k] != 0 { return 21 } 102 k = k + 1 103 } 104 105 // ---- Test G: (2^4096 - 1)^2 exercises the high-limb carry-out ---- 106 // (2^4096 - 1)^2 = 2^8192 - 2^4097 + 1 107 // low limb 0 = 1 108 // limbs 1..127 = 0 109 // limb 128 (bit 4096) = 0 (since -2^4097 only flips bit 4097) 110 // limb 128 high bit (bit 4097) = ... wait let me re-derive. 111 // 112 // (2^n - 1)^2 = 2^(2n) - 2*2^n + 1 = 2^(2n) - 2^(n+1) + 1. 113 // For n = 4096: result = 2^8192 - 2^4097 + 1. 114 // 2^8192 -> limb 256 bit 0 (overflow; we have 256 limbs = 8192 bits; 115 // so 2^8192 wraps to "0" within our wide buffer) 116 // -2^4097 + 1 -> we have a 8192-bit container; -2^4097 + 1 in two's-complement 117 // in our space is 2^8192 - 2^4097 + 1. 118 // Combined: 2^8192 + (2^8192 - 2^4097 + 1) = 2 * 2^8192 - 2^4097 + 1. 119 // Modulo 2^8192 that's -2^4097 + 1 = 2^8192 - 2^4097 + 1. 120 // 121 // So expected low limbs: 122 // limb 0 = 1 123 // limbs 1..127 = 0 124 // limb 128: bit 0 = 0 (the -2^4097 + 1 in this range subtracts bit 4097 125 // from the implicit 2^8192; the borrow propagates down filling 126 // with 0xFFFFFFFF from limb 128 (bit 4097's limb) downward). 127 // 128 // Actually it's easier to just compute and check observable invariants: 129 // out[0] must equal 1 (the +1 term). 130 // The schoolbook produces (2^n-1)^2 directly: limb 0 = 1, limbs 1..127 131 // form the "-2^4097 + 1 - 1" pattern which is borrow-propagated, and 132 // limbs 128..255 form the "2^8192 - 2^4097" pattern. 133 // 134 // Concretely, computing (2^4096 - 1)^2 with the standard formula: 135 // Each a[i] = 0xFFFFFFFF for i in 0..127. Product of any two is 136 // (2^32-1)^2 = 2^64 - 2^33 + 1. Sum over all 128*128 partial products 137 // produces the well-known: low limb 0 = 1; limb 127 = 0xFFFFFFFE; 138 // limb 128 = 0xFFFFFFFE; high limb 255 = 0xFFFFFFFE. 139 // 140 // Rather than encode the full expected pattern (too error-prone), 141 // assert the load-bearing endpoints: out[0] == 1 and out[255] != 0. 142 k = 0 143 while k < NX_U4096_LIMBS { 144 a[k] = NX_U4096_LIMB_MASK 145 b[k] = NX_U4096_LIMB_MASK 146 k = k + 1 147 } 148 u4096_mul_wide(out, a, b) 149 // (2^4096-1)^2 = 2^8192 - 2^4097 + 1 150 // = bit 0 set + bits 4097..8191 set (other bits zero) 151 // => out[0] = 0x00000001 152 // out[1..127] = 0 153 // out[128] = 0xFFFFFFFE (bit 0 of limb 128 = bit 4096 = zero; 154 // bits 1..31 of limb 128 = bits 4097..4127 = set) 155 // out[129..255] = 0xFFFFFFFF (bits 4128..8191 all set) 156 if out[0] != 1 { return 30 } 157 if (out[128] & NX_U4096_LIMB_MASK) != 0xFFFFFFFE { return 31 } 158 if (out[255] & NX_U4096_LIMB_MASK) != 0xFFFFFFFF { return 32 } 159 if (out[129] & NX_U4096_LIMB_MASK) != 0xFFFFFFFF { return 33 } 160 if out[127] != 0 { return 34 } 161 162 // ---- Test H: wide_shl1 round-trip ---- 163 let w: *i64 = u4096_wide_alloc() 164 u4096_zero(w as *i64) 165 w[0] = 1 166 u4096_wide_shl1(w) 167 if w[0] != 2 { return 40 } 168 u4096_wide_shl1(w) 169 if w[0] != 4 { return 41 } 170 171 // ---- Test I: wide_get_bit ---- 172 var z: i64 = 0 173 while z < NX_U4096_WIDE_LIMBS { w[z] = 0; z = z + 1 } 174 w[4] = 0x80000000 // bit 32*4+31 = 159 set 175 if u4096_wide_get_bit(w, 159) != 1 { return 50 } 176 if u4096_wide_get_bit(w, 158) != 0 { return 51 } 177 if u4096_wide_get_bit(w, 160) != 0 { return 52 } 178 179 // ---- Test J: wide_sub_low + wide_low_cmp ---- 180 let modn: *i64 = u4096_alloc() 181 u4096_zero(modn); modn[0] = 5 182 var z2: i64 = 0 183 while z2 < NX_U4096_WIDE_LIMBS { w[z2] = 0; z2 = z2 + 1 } 184 w[0] = 12 185 let cv: i64 = u4096_wide_low_cmp(w, modn) 186 if cv != 1 { return 60 } // 12 > 5 187 let bv: i64 = u4096_wide_sub_low(w, modn) 188 if bv != 0 { return 61 } 189 if w[0] != 7 { return 62 } 190 let cv2: i64 = u4096_wide_low_cmp(w, modn) 191 if cv2 != 1 { return 63 } // 7 > 5 192 193 // ---- Test K: wide_copy_low ---- 194 let low: *i64 = u4096_alloc() 195 u4096_zero(low) 196 u4096_wide_copy_low(low, w) 197 if low[0] != 7 { return 70 } 198 199 // "PASS\n" 200 let ok: *u8 = sys_mmap(8) 201 ok[0]=80; ok[1]=65; ok[2]=83; ok[3]=83; ok[4]=10 202 sys_write(1, ok, 5) 203 return 0 204}