From f3a72b58ed76db4ef2265fdcba4dd51c8a910629 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Pawe=C5=82=20Bylica?= Date: Mon, 10 Aug 2026 12:08:23 +0200 Subject: [PATCH] crypto: Use a sliding window in modexp The fixed window precomputes every power b^1..b^(2^w - 1). That all-powers table is the shape constant-time modexp needs, because sliding windows leak through data-dependent access. This code is not constant-time and EVM exponents are public calldata, so the table was twice the size it had to be. Store only the odd powers and align each window to a set bit, trimming trailing zeros so the value stays odd. Zero runs then cost one squaring each, and the smaller table affords one more width within the same scratch. Read each window with a single two-byte access through Exponent::window() rather than one call per exponent bit, and trim it to odd with countr_zero, which also removes the scan for the lowest set bit. A few percent fewer instructions over the benchmark matrix and over a corpus of mainnet inputs, up to about 17% on small exponents. One case regresses slightly: a modulus that is almost a power of two, where the odd part is a single word, so the table build is not amortised. --- lib/evmone_precompiles/modexp.cpp | 106 ++++++++++++--------- test/unittests/precompiles_expmod_test.cpp | 55 +++++++---- 2 files changed, 96 insertions(+), 65 deletions(-) diff --git a/lib/evmone_precompiles/modexp.cpp b/lib/evmone_precompiles/modexp.cpp index 923600983b..d54b0a2368 100644 --- a/lib/evmone_precompiles/modexp.cpp +++ b/lib/evmone_precompiles/modexp.cpp @@ -292,6 +292,19 @@ class Exponent const auto bit = (byte >> bit_index) & 1; return bit != 0; } + + /// Returns bits [lo, hi] as an integer, the bit at hi being the most significant. + /// The range must span at most 8 bits, so it covers at most two adjacent bytes. + [[nodiscard]] size_t window(size_t lo, size_t hi) const noexcept + { + assert(lo <= hi && hi - lo < 8); + const auto exp_size = (bit_width_ + 7) / 8; + const auto byte_index = exp_size - 1 - lo / 8; + auto bytes = size_t{data_[byte_index]}; + if (byte_index != 0) // Prepend the next more significant byte if there is one. + bytes |= size_t{data_[byte_index - 1]} << 8; + return (bytes >> (lo % 8)) & ((size_t{1} << (hi + 1 - lo)) - 1); + } }; /// Performs the Almost Montgomery Multiplication (AMM). @@ -369,27 +382,26 @@ template <> } /// Maximum window width used by the windowed method in modexp_odd. -constexpr unsigned MAX_WINDOW_WIDTH = 4; +constexpr unsigned MAX_WINDOW_WIDTH = 5; +static_assert(MAX_WINDOW_WIDTH <= 8, "Exponent::window() covers at most two adjacent bytes"); -/// Number of precomputed values for the max width windowed method. -constexpr size_t MAX_PRECOMPUTED = (size_t{1} << MAX_WINDOW_WIDTH) - 1; +/// Number of precomputed base odd powers for the max width windowed method. +constexpr size_t MAX_PRECOMPUTED = size_t{1} << (MAX_WINDOW_WIDTH - 1); -/// Selects the fixed-window width from the exponent bit length. -/// -/// TODO: Switch to a sliding window: the table then holds only odd powers. Measured ~3-7%. -/// TODO: Tune for the densest exponent instead of the average, because gas is charged on -/// exponent bit length and ignores Hamming weight. The width then collapses to -/// min(MAX_WINDOW_WIDTH, (bit_width(exp_bits) + 1) / 2). Measured +10.7% worst case. +/// Selects the sliding-window width from the exponent bit length. constexpr unsigned window_width(size_t exp_bits) noexcept { - // Break-even points for a random exponent, where the 2^w extra table multiplies stop - // being repaid: 2^w / ((1-2^-w)/w - (1-2^-(w+1))/(w+1)) = 16, 48, 140 (rounded to 144). - if (exp_bits <= 16) + // Break-even points for a random exponent, where the table's extra multiply stops being + // repaid: 2^(w-1) / (1/(w+1) - 1/(w+2)) = 6, 24, 80, 240. Each narrower width is kept + // one bit longer, which measures better on the sparse small exponents seen in practice. + if (exp_bits <= 7) return 1; - if (exp_bits <= 48) + if (exp_bits <= 25) return 2; - if (exp_bits <= 144) + if (exp_bits <= 81) return 3; + if (exp_bits <= 241) + return 4; return MAX_WINDOW_WIDTH; } @@ -408,7 +420,7 @@ void modexp_odd(std::span result, std::span base, Expo const auto exp_bits = exp.bit_width(); const auto w = window_width(exp_bits); - const auto table_size = (size_t{1} << w) - 1; + const auto table_size = size_t{1} << (w - 1); // Layout: u[n + base.size()] | table[MAX_PRECOMPUTED*n] // | rem_scratch[2*n + 2*base.size() + 2]. @@ -433,46 +445,50 @@ void modexp_odd(std::span result, std::span base, Expo auto r_tmp = std::span{u.first(n)}; const auto m = std::span{mod}; - // base_mont^j, for j in 1..table_size. - const auto precomputed = [table, n](size_t j) noexcept { - return std::span{table.subspan((j - 1) * n, n)}; + // base_mont^v, for odd v. + const auto precomputed = [table, n](size_t v) noexcept { + return std::span{table.subspan((v / 2) * n, n)}; }; - // precomputed(1) = base_mont is already set. - for (size_t j = 2; j <= table_size; ++j) - mul_amm(precomputed(j), precomputed(j - 1), precomputed(1), m, mod_inv); - - // Reads the `width` exponent bits starting at index `lo`. - // TODO: A window spans at most two adjacent bytes, so it could be read with one - // two-byte load, a shift and a mask. Est. 1-3%, and only for a 4-word modulus - // with a very long exponent; measure before doing it. - const auto window = [&](size_t lo, size_t width) noexcept { - size_t v = 0; - for (size_t b = 0; b < width; ++b) - v |= size_t{exp[lo + b]} << b; - return v; + // Fill the precomputed table (precomputed(1) is already set). + if (table_size > 1) + { + mul_amm(r_tmp, precomputed(1), precomputed(1), m, mod_inv); // r_tmp = base_mont^2. + for (size_t v = 3; v < 2 * table_size; v += 2) + mul_amm(precomputed(v), precomputed(v - 2), r_tmp, m, mod_inv); + } + + // The widest window of at most w bits ending at bit `hi`, which must be set. Trailing + // zero bits are trimmed off, so the value is odd and only odd table entries are + // needed. Returns the value and the index of its lowest bit. + const auto window = [exp, w](size_t hi) noexcept { + const auto lo = hi + 1 >= w ? hi + 1 - w : size_t{0}; + const auto v = exp.window(lo, hi); + const auto tz = static_cast(std::countr_zero(v)); // v != 0: exp[hi] is set. + return std::pair{v >> tz, lo + tz}; }; - // Windows tile from the bottom, so the ragged one is processed first, at the top. - // The top bit is always set, so that first window is nonzero. - // TODO: Tiling from the top instead would save w - top_width squarings when - // exp_bits % w != 0 (up to 4%), at the cost of a special-cased final iteration. - const size_t top_width = (exp_bits - 1) % w + 1; - std::ranges::copy(precomputed(window(exp_bits - top_width, top_width)), r_cur.begin()); + // The top bit is always set, so the first window ends there and is loaded directly. + auto [v_top, pos] = window(exp_bits - 1); + std::ranges::copy(precomputed(v_top), r_cur.begin()); - for (size_t pos = exp_bits - top_width; pos != 0;) + while (pos != 0) { - pos -= w; - for (unsigned s = 0; s != w; ++s) // square w times + --pos; + mul_amm(r_tmp, r_cur, r_cur, m, mod_inv); // Square for this bit. + std::swap(r_cur, r_tmp); + if (!exp[pos]) + continue; + + const auto [v, lo] = window(pos); + for (auto b = lo; b != pos; ++b) // One more square for each remaining window bit. { mul_amm(r_tmp, r_cur, r_cur, m, mod_inv); std::swap(r_cur, r_tmp); } - if (const size_t v = window(pos, w); v != 0) // multiply by base_mont^v - { - mul_amm(r_tmp, r_cur, precomputed(v), m, mod_inv); - std::swap(r_cur, r_tmp); - } + mul_amm(r_tmp, r_cur, precomputed(v), m, mod_inv); + std::swap(r_cur, r_tmp); + pos = lo; } // Convert from Montgomery form: multiply by 1. Reuses precomputed(1) storage. diff --git a/test/unittests/precompiles_expmod_test.cpp b/test/unittests/precompiles_expmod_test.cpp index 6726ba77d5..b03fe58c6c 100644 --- a/test/unittests/precompiles_expmod_test.cpp +++ b/test/unittests/precompiles_expmod_test.cpp @@ -291,55 +291,70 @@ TEST_P(expmod, inputs) // 2^129 mod (7 * 2^128): carry propagates and is absorbed in nonzero word. {"02", "0081", "0700000000000000000000000000000000", "0200000000000000000000000000000000"}, - // Fixed-window exponentiation in modexp_odd. One case per window width w, and - // per width of the leading partial window ((exp_bits - 1) % w + 1), which is what - // aligns the remaining windows. The exponents are picked so that the windows - // consumed cover 0 (multiply skipped), 1 and 2^w - 1 (first and last precomputed - // power). Modulus is the secp256k1 field prime: odd, 4 words, so these also cover - // the mul_amm<4> specialization. - // exp_bits=16, w=1: plain binary square-and-multiply, no table. + // Sliding-window exponentiation in modexp_odd. One case per window width w=1..5. + // Each exponent is built as: top bit (1) | zero run of w+1 bits | one run of w bits + // | trailing zeros, so its windows exercise both the first (b^1) and last + // (b^(2^w-1)) precomputed odd powers, the zero run in between being wide enough to + // keep them in separate windows, traversed by squarings alone. Modulus is the + // secp256k1 field prime: odd, 4 words, so these also cover the mul_amm<4> + // specialization. + // exp_bits=6, w=1: plain binary square-and-multiply, no table. + {"03", "24", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", + "00000000000000000000000000000000000000000000000002153e468b91c6d1"}, + // exp_bits=10, w=2: windows hit b^1 and b^3. + {"03", "0230", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", + "e123f780b153ebd75b17a6e7a7133dba60d90a7dbc0f770f08af0055f8e2c7ed"}, + // exp_bits=30, w=3: windows hit b^1 and b^7. + {"03", "21c00000", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", + "6cc581d10c7d071216edf63238959949056d7cddf5a90711a7c7cdec6b3e861f"}, + // exp_bits=100, w=4: windows hit b^1 and b^15. + {"03", "083c0000000000000000000000", + "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", + "7ff2c68783b688439f7c43de4cbfe265f8875ec726564a442c2cbd1244f6d99e"}, + // exp_bits=254 (mainnet-typical size), w=5: windows hit b^1 and b^31. + {"03", "207c000000000000000000000000000000000000000000000000000000000000", + "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", + "40ea9ce0f6a2c94a7bec98114179d8e1a21287312a25c1fdd7bf46e3d723984a"}, + // Same exponent as the w=5 case above, with a 5-word modulus: the cases above only + // ever run through the mul_amm<4> specialization, this covers the generic + // std::dynamic_extent instantiation at w > 1. + {"03", "207c000000000000000000000000000000000000000000000000000000000000", + "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", + "aa50260a96f69a722fc965bbfec20c21195eda68068b20e9899976f80ed8f6d4f6816bec10fc4ee6"}, + + // Random exponents straddling the width thresholds of the fixed-window + // implementation this replaced, which no longer coincide with the bands above: + // exp_bits 16..18 (w=2), 48..51 (w=3), 144..148 (w=4). Same modulus as above, + // except for the last case, which repeats exp_bits=148 with a 5-word modulus. {"03", "8005", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "79c4559d064ab3615f6da729a1f67265b88ee2eaba22838109bea30fb7bee31b"}, - // exp_bits=17, w=2, leading window 1 bit. {"03", "01001b", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "a890a61d8d745fae67a345fb031b048c0cf8952b43622263de0fdc4391a6c6a9"}, - // exp_bits=18, w=2, leading window 2 bits. {"03", "0200c9", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "600614416289329cf72ef906cdfc1dea20339051ec80ed3ff692eb14ed33be81"}, - // exp_bits=48, w=2: last exponent size before w becomes 3. {"03", "80013b71b865", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "fd66fdbe1f0c43e6640c121c366b9061c7f13964a572828c8e3968a50dba847f"}, - // exp_bits=49, w=3, leading window 1 bit. {"03", "0100d2c92fc182", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "651aace134976d8456fcc35686a57cf12670b2e596dabecd0ddae9984ced96c4"}, - // exp_bits=50, w=3, leading window 2 bits. {"03", "0200a6a7ef231d", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "f722a91e1faa3b57f0a19af8d4506b395a0a342e9ee2cbe65cd7a63155d38537"}, - // exp_bits=51, w=3, leading window 3 bits. {"03", "04013929f7999c", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "06f41e370c4ef45a2bc5e1ade1504fbe35e5a42a8f8c2b17ad16a6c657900d48"}, - // exp_bits=144, w=3: last exponent size before w becomes 4. {"03", "8004cb3ff13151bb9f84a488a5d62e79a680", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "97265df41405de7f9b35c1037c349ef367cffd34ed6a86cb933fe14f84bb12d1"}, - // exp_bits=145, w=4, leading window 1 bit. {"03", "010014b0a1922289f0b19f56c6c373b0e5cd4a", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "3587c0d41ce1eb59ec2fa686877d8166aa9740f2410f9271592e5f283e3bd738"}, - // exp_bits=146, w=4, leading window 2 bits. {"03", "02008d61508c16734bdbe4a9578f4c8185d260", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "f65d573e0ba5bdc7cc0e31072eb946ffe5138d0cd4bc936cc1a714d17cdaf954"}, - // exp_bits=147, w=4, leading window 3 bits. {"03", "040160dce60c2531e93ae750b53938d5b04faf", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "0648a7caabfd3d4b972c034830faf933179ed038e1e6a6c4c3ad26f330fe1397"}, - // exp_bits=148, w=4, leading window 4 bits. {"03", "0802ae8d294c48793907af3e71b536ed84fa84", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "40c2770e749bcbf7949855252da0258cc5ae80658427a4af8ba3489a81182ee9"}, - // Same, with a 5-word modulus: the windowed loop above only ever runs through the - // mul_amm<4> specialization, this covers the generic instantiation. {"03", "08f83d563ebc382e09e4b8245edebc817af708", "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f", "8016137e4c542dd66f4ab5f668fc0ac76d43353a675f3d4616a56f23757e463ca1093164385ef006"},