专栏文章

声速幂

___cplusplus
算法·理论参与者 5已保存评论 12

文章操作

快速查看文章及其快照的属性，并进行相关操作。

当前评论
8 条
当前快照
1 份
快照标识符
@mj4joaf2
此快照首次捕获于
2025/12/14 01:03
2 个月前
此快照最后确认于
2026/02/14 01:29
6 天前
查看原文
省流：固定底数；对于 k≤wk≤w，O(w2k/k)O(w2^k/k) 预处理，O(w/k)O(w/k) 询问。

前言

前置知识：快速幂 光速幂也可以在里面找到。
符号表示：在本文中，WW 为指数的数域，w=log⁡2Ww = \log_2 W。aa 为底数。
对于结果有模数的幂，注意到 ab≡a(b mod φ(p))+φ(p)(modp)a^b \equiv a^{(b \bmod \varphi(p)) + \varphi(p)} \pmod p，因此有 W=2φ(p)W = 2\varphi(p)。对于素模数有 W=p−1W = p - 1。
由于洛谷上没有讲光速幂的文章，这里简要讲一下。
光速幂
对于很小的 WW，我们可以 O(W)O(W) 地预处理出幂的值。
对于较大的 WW 则不行。
取 s≤Ws \le W，我们知道
ab=a⌊b/s⌋⋅ab mod sa^b = a^{\lfloor b/s\rfloor}\cdot a^{b \bmod s}
因此我们预处理 fi=ai  (0≤i<s)f_i = a^i\;(0 \le i < s) 和 gi=ais  (0≤i≤(W−1)/s)g_i = a^{is} \; (0 \le i \le (W-1)/s)，则 ab=gb/s⋅fb mod sa^b = g_{b/s}\cdot f_{b\bmod s}。
注意到取 s=O(W)s = O(\sqrt{W}) 时，复杂度是最佳的。

原理

同样解决固定底数的幂，对于 k≤wk\le w，可以在 O(w2k/k)O(w2^k/k) 的时间内预处理一个底数，在 O(w/k)O(w/k) 的时间内回答一个询问。
我们考虑处理一个 s=2ks = 2^k 进制的快速幂。很明显，这个 ss 进制的数一共有 w/kw/k 位。
比如 w=8,k=4w=8,k=4：
a100101102=a100100002⋅a000001102a^{10010110_2} = a^{1001\textcolor{lightgray}{0000}_2} \cdot a^{\textcolor{lightgray}{0000}0110_2}
ss 进制快速幂的每一位都有 ss 个值，这个值很大，不能立即算出，需要预处理。
对于这个例子，我们要预处理出 （〇）a000000002∼a000011112a^{\textcolor{lightgray}{0000}0000_2}\sim a^{\textcolor{lightgray}{0000}1111_2}， （一）a000000002∼a111100002a^{0000\textcolor{lightgray}{0000}_2}\sim a^{1111\textcolor{lightgray}{0000}_2}。
考虑缩小 kk，比如 w=8,k=2w = 8,k = 2：
a100101102=a000000102⋅a000001002⋅a000100002⋅a100000002. a^{10010110_2} = a^{\textcolor{lightgray}{00}\textcolor{lightgray}{0000}10_2} \cdot a^{\textcolor{lightgray}{0000}01\textcolor{lightgray}{00}_2} \cdot a^{\textcolor{lightgray}{00}01\textcolor{lightgray}{0000}_2} \cdot a^{10\textcolor{lightgray}{00}\textcolor{lightgray}{0000}_2} .
这时候，我们就需要预处理出
  • （〇）a000000002∼a000000112a^{\textcolor{lightgray}{00}\textcolor{lightgray}{0000}00_2} \sim a^{\textcolor{lightgray}{00}\textcolor{lightgray}{0000}11_2}，
  • （一）a000000002∼a000011002a^{\textcolor{lightgray}{0000}00\textcolor{lightgray}{00}_2} \sim a^{\textcolor{lightgray}{0000}11\textcolor{lightgray}{00}_2}，
  • （二）a000000002∼a001100002a^{\textcolor{lightgray}{00}00\textcolor{lightgray}{0000}_2} \sim a^{\textcolor{lightgray}{00}11\textcolor{lightgray}{0000}_2}，
  • （三）a000000002∼a110000002a^{00\textcolor{lightgray}{00}\textcolor{lightgray}{0000}_2} \sim a^{11\textcolor{lightgray}{00}\textcolor{lightgray}{0000}_2}。
考虑如何预处理。我们令 f(i,j)=aisj  (0≤i<s,0≤j<w/k)f(i,j) = a^{is^j}\;(0\le i<s,0\le j<w/k)。可以验证，这就是以上预处理的内容。换言之，f(i,j)=ai2jkf(i,j)=a^{i2^{jk}}。
  • 首先枚举计算 a2i,0≤i<wa^{2^i}, 0\le i<w。
  • 对于 i=0i=0，明显有 f(i,j)=1f(i,j)=1。
  • f(i+1,j)=a(i+1)sj=aisj+sj=aisj⋅asj=f(i,j)⋅asjf(i+1,j)=a^{(i+1)s^j}=a^{is^j+s^j}=a^{is^j}\cdot a^{s^j}=f(i,j)\cdot a^{s^j}。
形式化地，
f(i,j)=1(i=0)f(i,j)=f(i−1,j)⋅asj(i>0)\begin{aligned} f(i,j)&=1\quad&(i=0)\\ f(i,j)&=f(i-1,j)\cdot a^{s^j}\quad&(i>0) \end{aligned}
显然是 O(ws/k)O(ws/k) 的。
对于询问，只需要将这个 ss 进制的指数 bb 的每一位预处理出的值相乘。
形式化地，
ab=∏i=0w/k−1f((a/si) mod s,i)a^b = \prod_{i=0}^{w/k-1}f\left((a/s^i)\bmod s,i\right)
很明显是 O(w/k)O(w/k) 的。

例题

例题
令 W=264W = 2^{64}。
给定 a,b,na,b,n，求 ⨁i=1nabi mod W mod W\displaystyle \bigoplus_{i=1}^n a^{b^i\bmod W}\bmod W。
0≤a,b<W,1≤n≤1080 \le a,b < W, 1 \le n \le 10^8。
令 ϕ=φ(W)=263\phi = \varphi(W) = 2^{63}，则当 bi≥ϕb^i\ge \phi 时，abi≡a(bi mod ϕ)+ϕ(modW)a^{b^i}\equiv a^{(b^i \bmod \phi)+\phi}\pmod {W}，否则不能降幂。因为此时 bi<2ϕb^i<2\phi，故没有降幂。
此题中，w=64w = 64，我们取 k=16k = 16。
参考代码CPP
// C++14(GCC 9) O2 洛谷 IDE
#include <bits/stdc++.h>
using namespace std;
typedef size_t Z;
constexpr Z N = 5 + 1e9, M = 65535, S = 65536;

Z a, b, cur = 1, n, ax[64], ans;
Z f0[S], f1[S], f2[S], f3[S];

int main() {
  cin >> a >> b >> n;

  // 预处理
  // ax[i] = pow(a, pow(2, i));
  ax[0] = a; 
  for (Z i = 1; i < 64; ++i)
    ax[i] = ax[i - 1] * ax[i - 1];
  Z x0 = ax[0],  x1 = ax[16],
    x2 = ax[32], x3 = ax[48];
  // f##j[i] = pow(a, i * pow(65536, j));
  f0[0] = f1[0] = f2[0] = f3[0] = 1;
  for (Z i = 1; i < 65536; ++i) {
    f0[i] = f0[i - 1] * x0;
    f1[i] = f1[i - 1] * x1;
    f2[i] = f2[i - 1] * x2;
    f3[i] = f3[i - 1] * x3;
  }

  for (Z i = 1; i <= n; ++i) {
    cur *= b;
    ans ^= f0[cur & M] * f1[cur >> 16 & M]
           * f2[cur >> 32 & M] * f3[cur >> 48 & M]; // 询问
  }
  cout << ans;
  return 0;
}
运行时间约 500500 毫秒，无氧 900900 毫秒，空间 25602560 kb，正确性已经过对比验证。
一般快速幂应该过不了。

评论

共 12 条评论，欢迎与作者交流。

正在加载评论...