fast_photon 的博客

博客

技能离散对数

2026-05-04 10:57:19 By fast_photon

技能离散对数 $^{[1]}$

若无特殊说明,模数 $m$ 是质数,$g$ 是一个原根,$\log$ 默认以 $g$ 为底。$p$ 是质数。

传统😯的离散对数🧮,就是 $^{[1]}$

先用 BSGS $^{[2]}$ 预处理出 $\sqrt m + 1$ 以内质数的离散对数,时间复杂度 $\Theta(\dfrac {m^{0.75}} {\log m})$,然后根据 $\log ab = \log a + \log b$ 计算出 $\sqrt m + 1$ 以内所有数的离散对数。查询时,对于 $x<\sqrt m$,可以直接查表。对于 $x>\sqrt m$,令 $y=\lfloor \dfrac m x \rfloor, z_1 = m-xy, z_2 = x(y+1)-m$,由于 $z_1+z_2=x$,较小值必不超过 $\dfrac x 2$。以 $z_1$ 为例,对 $m$ 取模,有 $z_1\equiv -xy \pmod m$,两侧取对数得到 $\log z_1 \equiv \log x + \log y + \log (-1) \pmod {m-1}$。移项可得 $\log x \equiv \log z_1 - \log y - \log (-1) \pmod {m-1}$,于是递归计算 $\log z_1$。在 $\log_2 m$ 轮以内,$x$ 将减少至 $\sqrt m$ 以下,查询复杂度为 $q\log m$。

好无聊😩,好无趣😫。$^{[1]}$

由于 BSGS 固有的 $\sqrt m$ 复杂度,该做法难以在算法竞赛常见时间内处理 $10^{18}$ 左右的模数。

而技能离散对数🤓☝,就是在传统的离散对数🧐加入 $^{[1]}$

Index Calculus 算法 $^{[3]}$,其核心思想是通过随机 $1\le k\le m-1$,并分解 $g^k\bmod m$,得到关于 $\log_g p_i$ 的线性方程。然后通过解线性方程直接得到离散对数,绕开 BSGS 过程。

具体而言,$g^k\bmod m$ 可以认为是 $m$ 以内的随机数。若希望建立 $B$ 以内质数的线性方程组,则随机 $k$ 可以得到线性方程当且仅当 $g^k\bmod m$ 是 $B\ -$ 光滑数(所有质因子不超过 $B$)。令 $u=\log_B m$,则 $\rho \approx u^{-u}$ $^{[4]}$,通过选取合适的 $B$ 可以平衡到一个不错的复杂度。在 $10^{18}$ 左右的模数下,可以选取 $B=1000$。

好好玩🥰🥰🥰,要爆了💥💥💥$^{[1]}$

将上述两个做法结合起来,再加入一点神秘小巧思。在随机出 $x=g^k\bmod m$ 后,不立即将它分解。选取合适的 $B_2=10^7$,预处理出 $B_2$ 以内哪些数是 $B\ -$ 光滑数。执行类似基于值域预处理的快速离散对数中的过程。若 $y,y+1$ 中有至少一个是 $B_2$ 以内的 $B\ -$ 光滑数,则可将 $x$ 减小。若二者都是,贪心地选择较小的 $z$。若都不是,结束这一过程。限定这一过程最多进行 $32$ 轮。

然后通过试除检验 $x$ 的 $B\ -$ 光滑性。若是,则得到一组关于 $B$ 以内质数离散对数的线性方程。随机很多轮,直到线性方程组有唯一解。根据光滑数密度近似公式可以看出,将 $x$ 减小这一操作可以提高其是 $B\ -$ 光滑数的概率,进而减少随机轮数,降低光滑性检验代价。

维护一个对数表,初始包含 $B_2$ 以内的所有 $B\ -$ 光滑数。

接下来,再随机 $7\times 10^5$ 轮。在每一轮中,允许使用对数表里的任何数作为 $y/y+1$ 减小 $x$。无法继续减小时,试除以剔除 $x$ 在 $B$ 以内的质因子。若剩余的 $x$ 是一个 $B_2$ 以内的质数,则求出了其离散对数。若其不在表中,试图将其和其倍数插入对数表。这一部分的均摊复杂度为 $B_2\ln B_2$。

查询一个数 $t$ 时,随机 $k$ 并计算 $t\gets tg^k\bmod m$。使用对数表内的数减小 $t$。若得到的 $t$ 剔除小于 $B$ 的质因子后在对数表内,就可以得到答案。

技能离散对数,飞沙走石🌪️$^{[1]}$

这个做法实测跑得比较优秀,可以在 $1.5$ 秒内预处理一个模数并查询 $2^{17}$ 个离散对数。

技能离散对数,力拔山兮💪🏻$^{[1]}$

https://qoj.ac/submission/2300547

技能离散对数,飞沙走石🌪️$^{[1]}$

这个做法实测跑得比较优秀,可以在 $1.5$ 秒内预处理一个模数并查询 $2^{17}$ 个离散对数。

技能离散对数,时光倒流⏰$^{[1]}$

对不起,我暂时没有找到运行技能离散对数算法可能导致时光倒流的依据。要不要我为你整理离散对数在密码学中的常规应用或相关时间悖论的科普信息?

参考文献

$[1]$ 张兴朝, 李嘉诚, 张呈, 王广. 技能五子棋. 2025.

$[2]$ Daniel Shanks. Class number, a theory of factorization, and genera. 1971.

$[3]$ Leonard M. Adleman. A Subexponential Algorithm for the Discrete Logarithm Problem with Applications to Cryptography. 1979.

$[4]$ N.G. de Bruijn. The asymptotic behaviour of a function occurring in the theory of primes. 1951.

一种不太快速求一般图唯一完美匹配的确定性做法

2026-01-24 18:35:53 By fast_photon

对于无向图 $G=(V, E)$,记其邻接矩阵为 $A$,满足 $A_{i,j}=|\{(i,j)\}\cap E|$。将 $A$ 视作一个 $\mathbb F_2$ 下的矩阵。若 $\det A=0$,则 $G$ 无唯一完美匹配(可能有多组或没有)。否则求其逆 $A^{-1}$,令 $E'=\{(u,v)|A_{u,v}A^{-1}_{u,v}=1\}$。若 $E'$ 中所有边不构成一组完美匹配,则 $G$ 无唯一完美匹配。否则可以使用带花树在 $O(n^2)$ 时间内检查该组完美匹配是否唯一。

时间复杂度 $O(n^\omega)$。

另外,扔掉带花树部分,改为在 $E'$ 是一组完美匹配时直接报告唯一,可以通过 Judge Error 的所有数据。(当然,存在一组针对此做法的 hack)

co-first author: Max_s_xaM

一种不太快速计算多个少项式集合幂级数异或卷积的方法

2025-12-18 16:23:45 By fast_photon

首先求 ${\prod}_\oplus(a_ix^{\emptyset}+b_ix^{S_i})$,其中二项式个数为 $n$,$S_i\subseteq \{1,2,\cdots,m\}$。

考虑 FWT,要将 $\operatorname{FWT}(a_ix^{\emptyset}+b_ix^{S_i})$ 对位乘得到 $g(x)$ 之后 IFWT。

$\begin{aligned} & [x^T]g(x) \\ = & \prod [x^T]\operatorname{FWT}(a_i x^{\emptyset}+b_ix^{S_i}) \\ = & \prod (a_i+(-1)^{|S_i\cap T|}b_i) \\ = & \exp(\sum \ln (a_i+(-1)^{|S_i\cap T|}b_i)) \\ = & \exp(\sum [x^T] \operatorname{FWT}(\dfrac {\ln(a_i+b_i)+\ln(a_i-b_i)}{2} x^{\emptyset} + \dfrac {\ln(a_i+b_i)-\ln(a_i-b_i)}{2} x^{S_i})) \\ = & \exp([x^T]\operatorname{FWT}(\sum(\dfrac {\ln(a_i+b_i)+\ln(a_i-b_i)}{2} x^{\emptyset} + \dfrac {\ln(a_i+b_i)-\ln(a_i-b_i)}{2} x^{S_i})))\end{aligned}$

只要对数能求就能直接做。

如果有较多项呢?对于 $\sum a_{i,j}x^{S_{i,j}}$,对于每个 $i$,找一组 $F_2^m$ 中的基底,它们张成一个秩为 $2^{k_i}$ 的子空间。把 $S_{i,j}$ 变换到 $S'_{i,j}$,在子空间里做 FWT,取对数,做 IFWT,再把点值加到原空间的对应位置(具体来说,对于一个子空间的向量 $v$,取出所有 $v_j=1$ 的 $S_{i,j}$,将其 $\operatorname{xor}$ 起来),做 FWT,取指数,做 IFWT。

复杂度 $O(m2^m+\sum_i k_i2^{k_i}+L(\sum_i 2^{k_i}))$,其中 $L(x)$ 是计算 $x$ 个数的对数的复杂度,$k_i$ 是所有 $S_{i,j}$ 张成子空间的大小。

不能取对数的时候,考虑额外记录 $0$ 的指数。对于 $x\neq 0$,记 $\log x=(\ln x,0)$,记 $\ln 0=(0,1)$ 即可。做完原空间上的 FWT 后,若某个向量是 $(x,0)$,将其替换为数值 $\exp(x)$,否则替换为数值 $0$,再 IFWT 即可。

以下是示例代码。由于没有高消,其累加复杂度是 $O(2^{m_i})$,$m_i$ 是第 $i$ 个少项式的项数。可以在全局常量里调 $m$ 和模数,也可以改成输入模数。使用了较为朴素的方法求原根,在示例代码中这可能成为复杂度瓶颈。

#include<iostream>
#include<cstdio>
#include<fstream>
#include<queue>
#include<algorithm>
#include<cstring>
#include<random>
#include<map>
#include<cassert>
#include<stack>
#define fopen(x, y) freopen(x".in", "r", stdin); freopen(y".out", "w", stdout);
#define int long long
#define pii pair<int, int>
#define fi first
#define se second
#ifdef int
#define inf 0x3f3f3f3f3f3f3f3fll
#else
#define inf 0x3f3f3f3f
#endif

#define maxB 65536
#define maxm 1048576 

using namespace std;

//int m, M, mod, g, b, B, Mod;
const int m = 20, M = 1ll << m, mod = 998244353, g = 3, b = 16, B = 65536, Mod = 998244352ll * M;

inline int qpow(int x, int y) { int xum = 1; while(y) { if(y & 1) (xum *= x) %= mod; (x *= x) %= mod; y >>= 1; } return xum; }

int n, a, s[maxm], x[105], w[105], hfi[maxm], hse[maxm], ffi[maxm], fse[maxm];
int l3, l_1;

//void init_0();
void init_1();
int gxp(int x);
int gln(int x);

void fwt(int *f, int m, const int mod) {
    const int M = 1ll << m;
    for(int j = 0; j < m; j++) {
        int L = 1ll << (j + 1), bL = 1ll << j;
        for(int t = 0; t < M; t += L) {
            int *x = &f[t], *y = x + bL;
            for(int i = 0; i < bL; i++, x++, y++) {
                (*x += *y) %= mod;
                *y = (*x - (*y << 1)) % mod;
            }
        }
    }
}

void work() {
//    mod = 998244353; init_0();

    init_1();

    cin >> n;

    for(int i = 1; i <= n; i++) {
        int k;
        cin >> k;
        const int K = 1ll << k;
        for(int j = 0; j < K; j++) hfi[j] = hse[j] = 0;
        for(int j = 0; j < k; j++) {
            cin >> x[j] >> w[j];
            hfi[1ll << j] = w[j];
        }
        s[0] = 0;
        for(int j = 1; j < K; j++) {
            s[j] = s[j ^ (j & -j)] ^ x[__builtin_ctzll(j & -j)];
        }
        fwt(hfi, k, mod);
        for(int j = 0; j < K; j++) {
            (hfi[j] += mod) %= mod;
            if(hfi[j] == 0) hse[j] = 1ll << (m - k);
            else hfi[j] = (gln(hfi[j]) << (m - k)) % Mod;
        }
        fwt(hfi, k, Mod);
        fwt(hse, k, Mod);
        for(int j = 0; j < K; j++) {
            (ffi[s[j]] += hfi[j]) %= Mod;
            (fse[s[j]] += hse[j]) %= Mod;
        }
    }

    fwt(ffi, m, Mod);
    fwt(fse, m, Mod);

    int sum = 0, iv = qpow(M, mod - 2);
    for(int i = 0; i < M; i++) {
        if(fse[i]) {
            ffi[i] = 0;
            continue;
        }
        (ffi[i] >>= m) %= (mod - 1);
        if(ffi[i] < 0) ffi[i] += mod - 1;
        ffi[i] = gxp(ffi[i]);
    }
    fwt(ffi, m, mod);
    for(int i = 0; i < M; i++) {
        (ffi[i] *= iv) %= mod;
        cout << (ffi[i] + mod) % mod << ' ';
    }
    cout << '\n';
}

signed main() {
    ios::sync_with_stdio(false); cin.tie(0);
    int _ = 1;
    while(_--) work();
}

/*
void init_0() {
    Mod = (mod - 1) * maxn;
    vector<int> t;
    int s = mod - 1;
    B = 1, b = 0;
    while(B * B < mod) B <<= 1, b++;
    for(int i = 2; i * i <= s; i++) {
        if(s % i == 0) {
            while(s % i == 0) s /= i;
            t.push_back(i);
        }
    }
    for(int i = 2; i < mod; i++) {
        int flg = 1;
        for(int d : t) {
            if(qpow(i, (mod - 1) / d) == 1) {
                flg = 0;
                break;
            }
        }
        if(flg) {
            g = i;
            break;
        }
    }
}
//*/

int gb, gk[maxB], gbk[maxB];
pii gt[maxB + 1];
int find(int x) {
    pii *p = lower_bound(gt, gt + B, make_pair(x, 0ll));
    if(p->fi == x) return p->se;
    return -1;
}
int gln(int x) {
    int cnt = 0, p = find(x);
    while(p == -1) {
        (x *= gb) %= mod;
        cnt++;
        p = find(x);
    }
    return ((p - cnt * B) % (mod - 1) + (mod - 1)) % (mod - 1);
}
void init_1() {
    gb = 1;
    gk[0] = 1;
    for(int i = 1; i < B; i++) gk[i] = gk[i - 1] * g % mod;
    gb = gk[B - 1] * g % mod;
    gbk[0] = 1;
    for(int i = 1; i < B; i++) gbk[i] = gbk[i - 1] * gb % mod;

    for(int i = 0; i < B; i++) gt[i] = make_pair(gk[i], i);
    sort(gt, gt + B); 
    gt[B] = make_pair(mod, -1);
}
int gxp(int x) {
    return gk[x & (B - 1)] * gbk[x >> b] % mod;
}
共 3 篇博客