In mathematics you don't understand things, you just get used to them.

标签 maths 下的文章

link。

Lagrange Interpolation。

朴素的 dp 即设 $f_i(j)$ 表示前 $i$ 个位置,最大值为 $j$,位置 $i$ 可选可不选的方案数,转移即 $\displaystyle f_i(j) = f_{i-1}(j)+ [l_i \leqslant j \leqslant r_i]\sum_{k < j} f_{i-1}(k)$。把 $j$ 当作自变量,若没有 $l$,$r$ 的限制,则可以发现这是个多项式。加上 $l$,$r$,则是一个分段多项式。比如 $f_1(x) = [l_1 \leqslant x \leqslant r_1]$,$\displaystyle \sum_x f_1(x)$ 是分段多项式 $g_1(x)$,$f_2(x) = f_1(x) + [l_2 \leqslant x \leqslant r_2] g_1(x-1)$,每次前缀和加一次次数。

考虑如何削去 $l$,$r$ 的限制,我们可以像 codeforces - 1295f 一样,把区间们化成左闭右开并且把每一个端点当成“标兵”放到数轴上,对相邻标兵构成的区间段(同样是左闭右开)进行考察。容易发现,在同一个段中的多项式并没有分段,我们直接代入 $n+1$ 个点值,就可以确定这个最多 $n$ 次的多项式。

注意,当我们把区间端点离散化后,我们 dp 的状态就变了,$f_i(j)$ 的意义成为了前 $i$ 个位置,前 $j$ 个区间段(注意是值域)拉通了来看的方案数,即多项式在末尾处的取值。

$$ \texttt{|*****|*****|*****|*****|}\underbrace{{\color{blue}{\texttt {*****}}}}_{\text{cureent }j}{\texttt{|*****|*****|}} $$

看到我们研究的这个 $j$,段 $j$ 本身就是一个连续的多项式函数,我们需要从段 $j-1$ 的末尾转移过来,这就是为什么我们要求的是段 $j$ 末尾的函数值,即后面的转移需要用。

感觉说得不是很清楚,建议自己再画画图……

const int md = 1e9 + 7;
il int add(int x, int y) {
  if ((x += y) >= md)
    x -= md;
  return x;
}
il int sub(int x, int y) {
  if ((x -= y) < 0)
    x += md;
  return x;
}
il int mul(int x, int y) { return static_cast<long long>(x) * y % md; }
template <class T, class... P> il T mul(T x, T y, P... args) {
  return mul(x, mul(y, args...));
}
il void adds(int &x, int y) {
  if ((x += y) >= md)
    x -= md;
}
il void muls(int &x, int y) { x = static_cast<long long>(x) * y % md; }
int pow(int x, int y) {
  int z(1);
  for (; y; y >>= 1, muls(x, x))
    if (y & 1)
      muls(z, x);
  return z;
}
il int sgn(int n) { return n % 2 ? md - 1 : 1; }
int n, l[1 << 9], r[1 << 9], dat[2 << 9], tot, rec[2 << 9][1 << 9], sum[2 << 9],
    ifct[1 << 9], dp[2 << 9];
int p[1 << 9], q[1 << 9];
int fit(int t, int y[]) {
  int res = 0;
  p[0] = q[n + 2] = 1;
  rep(i, 1, n + 2) p[i] = mul(p[i - 1], sub(t, i));
  drep(i, n + 2, 1) q[i] = mul(q[i + 1], sub(t, i));
  rep(i, 1, n + 2) adds(res, mul(ifct[i - 1], ifct[n + 1 - i], sgn(n - i + 1),
                                 p[i - 1], q[i + 1], y[i]));
  return res;
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  cin >> n;
  ifct[0] = 1;
  rep(i, 1, n + 2) ifct[i] = mul(ifct[i - 1], i);
  rep(i, 1, n + 2) ifct[i] = pow(ifct[i], md - 2);
  rep(i, 1, n + 1) cin >> l[i] >> r[i];
  rep(i, 1, n + 1) dat[++tot] = l[i], dat[++tot] = ++r[i];
  sort(dat + 1, dat + tot + 1);
  tot = unique(dat + 1, dat + tot + 1) - dat - 1;
  rep(i, 1, n + 1) l[i] = lower_bound(dat + 1, dat + tot + 1, l[i]) - dat;
  rep(i, 1, n + 1) r[i] = lower_bound(dat + 1, dat + tot + 1, r[i]) - dat;
  rep(i, tot + 1) sum[i] = 1;
  rep(i, 1, n + 1) {
    rep(j, l[i], r[i]) {
      rep(k, 1, n + 2) adds(rec[j][k], sum[j - 1]);
      rep(k, 2, n + 2) adds(rec[j][k], rec[j][k - 1]);
      dp[j] = fit(dat[j + 1] - dat[j], rec[j]);
    }
    rep(j, 1, tot + 1) sum[j] = add(sum[j - 1], dp[j]);
  }
  cout << sum[tot] - 1 << "\n";
}

link。

关于 Lagrange Interpolation。

插值即是对于给出的向量组 $V = (v_0, ..., v_{n-1})$,其中 $v_i = (x_i, y_i)$,构造一个函数 $f(x)$ 使得 $\forall i$,有 $f(x_i) = y_i$。

朴素做法即对方程组消元,$O(n^3)$。

Lagrange 插值则是选择了类 CRT 的构造方法,即构造 $n$ 个函数 $f_i(x)$,使得 $\forall i$,有 $f_i(x_i) = y_i$,$\forall j \neq i$,有 $f_i(x_j) = 0$,那么构造结果 $f(x) = \sum \limits_{i=0}^{n-1} f_i(x)$。对于一个 $i$,我们有 $f_i(x) = a_i \times \prod \limits_{i \neq j} x-x_j$,此时若 $a_i \neq = 0$,则我们已经满足了条件二;考虑条件一,令 $a_i = \frac{y_i}{\prod \limits_{i \neq j} x-x_j}$ 即可,即有 $f(x) = \sum \limits_{i=0}^{n-1}y_i\prod\limits_{i\neq j}\frac{x-x_j}{x_i-x_j}$,$O(n^2)$。

看回这个题,朴素的 dp 做法即是令 $dp[u][d]$ 表示考虑子树 $u$ 的方案数,满足 $u$ 的赋权为 $d$,转移即各个子树的 dp 和的乘积,容易看出这是个多项式,即 $dp[u][d]$ 是关于 $d$ 的 $s_u-1$ 次多项式,其中 $s_u$ 为 $u$ 的子树大小。考虑差分 $dp[u][d]-dp[u][d-1] = \prod dp[v][d]$,又因为叶子结点满足条件即证。

那么求出 $dp[0][0 \dots n]$,插值求出其曲线方程,带入 $d$ 即可。

const int mod = 1e9+7;
il int sub(int x,int y) {
    if ((x-=y)<0) x+=mod;
    return x;
}
il int mul(int x,int y) {
    return static_cast<long long>(x)*y%mod;
}
il void adds(int& x,int y) {
    if ((x+=y)>=mod) x-=mod;
}
il void muls(int& x,int y) {
    x=static_cast<long long>(x)*y%mod;
}
il int qpow(int x,int y) {
    int z(1);
    for (; y; y>>=1,muls(x,x))
        if (y&1) muls(z,x);
    return z;
}
il int inv(int m) {
    assert(m>0);
    return qpow(m,mod-2);
}
int n, dp[3<<10][3<<10], d, Y[3<<10];
bsi<int> adj[3<<10];
void dfs(int x) {
    rep(i,1,n+1) dp[x][i] = 1;
    for (auto y : adj[x]) {
        dfs(y);
        rep(i,1,n+1) muls(dp[x][i], dp[y][i]);
    }
    rep(i,1,n+1) adds(dp[x][i], dp[x][i-1]);
}
signed main() {
    ios::sync_with_stdio(0);
    cin.tie(0);
    int x;
    cin >> n >> d;
    rep(i,1,n) {
        cin >> x;
        x--;
        adj[x] += i;
    }
    dfs(0);
    rep(i,n+1) Y[i] = dp[0][i];
    int ans=0;
    rep(i,n+1) {
        int u=Y[i],v=1;
        rep(j,n+1) if (i!=j) muls(u,sub(d,j)),muls(v,sub(i,j));
        adds(ans,mul(u,inv(v)));
    }
    cout<<ans<<"\n";
}

题都是好题,会的没几道。

「codeforces - 367E」Sereja and Intervals:注意到 $l_x < l_y$ 且 $r_x < r_y$,以及 $n \leqslant \sqrt{10^5}$,然后分配左右端点 dp。

using mint = modint1000000007;
/*
n <= sqrt
sort all ranges
构造结果长相: l_i < l_j, r_i < r_j
dp[i][l][r] [1, i], l ok, r ok l >= r
if i != x:
  dp[i][l][r] = dp[i-1][l-1][r]
*/
simple<mint> sip;
int n, m, X;
mint dp[2][320][320];
signed main() {
    ios::sync_with_stdio(0);
    cin.tie(0);
    cin >> n >> m >> X;
    if (n > m) {
        cout << "0\n";
        return 0;
    }
    dp[0][0][0] = 1;
    for (int i=1,turn=1; i<=m; ++i,turn^=1) {
        for (int u=0; u<=n+1; ++u) {
            for (int v=0; v<=n+1; ++v) {
                dp[turn][u][v] = 0;
            }
        }
        for (int u=0; u<=n; ++u) {
            for (int v=0; v<=u; ++v) {
                if (u > v && i != X) dp[turn][u][v+1] += dp[turn^1][u][v];
                dp[turn][u+1][v] += dp[turn^1][u][v];
                dp[turn][u+1][v+1] += dp[turn^1][u][v];
                if (i != X) dp[turn][u][v] += dp[turn^1][u][v];
            }
        }
    }
    cout << (sip.gfac(n)*dp[m&1][n][n]).val() << "\n";
}

「codeforces - 140E」New Year Garland:我觉得挺难想到的,但是被某些 alpha 人喷水题了 233。考虑没有行间限制的情况,你写出了 $m \times (m-1)^{\sum l_i}$,然后发现没用。写出一行的 dp,设 $f[i][j]$ 表示摆 $i$ 个球,用恰好 $j$ 种颜色,可以写出 $g[i][j] = g[i-1][j-1]+g[i-1][j] \times (j-1)$。然后考虑行间限制,设 $f[i][j]$ 表示前 $i$ 行的方案数,其中第 $i$ 行用了恰好 $j$ 种颜色,同样可以写出 $\displaystyle f[i][j] = \binom{m}{j} \times j! \times g[l[i]][j] \times \left(\sum f[i-1][k]\right) - g[l[i]][j] \times f[i-1][j]$。

using mint = dynamic_modint<-1>;
int n, m, p, a[1000100];
mint P[1000100], rnm[5100], tuiq[5100], sum = mint::raw(1), ans, fct[5100], sl[5100][5100];
signed main() {
    // ios::sync_with_stdio(0);
    // cin.tie(0);
    cin >> n >> m >> p;
    mint::set_mod(p);
    P[0] = 1;
    for (int i=1; i<=m; ++i) P[i] = P[i-1]*(m-i+1);
    for (int i=1; i<=n; ++i) {
        cin >> a[i];
    }
    sl[0][0] = 1;
    for (int i=1; i<=5000; ++i) {
        for (int j=1,up=min(i, m); j<=up; ++j) {
            sl[i][j] = sl[i-1][j-1]+sl[i-1][j]*(j-1);
        }
    }
    fct[0] = 1;
    for (int i=1; i<=5000; ++i) fct[i] = fct[i-1]*i;
    for (int i=1,turn=1; i<=n; ++i,turn^=1) {
        mint tmp(0);
        for (int j=1; j<=a[i]; ++j) {
            rnm[j] = P[j]*sum*sl[a[i]][j]-sl[a[i]][j]*tuiq[j]*fct[j];
            tmp += rnm[j];
        }
        sum = tmp;
        for (int j=1; j<=a[i-1]; ++j) tuiq[j] = 0;
        for (int j=1; j<=a[i]; ++j) tuiq[j] = rnm[j];
    }
    cout << sum.val() << "\n";
}

「codeforces - 1523E」Crypto Lights:trick:设概率。但是还是没有理解什么时候需要设概率。设出 $f_i$ 表示放了恰好 i 台灯停止的概率,则期望答案为 $\sum f_i \times i$,这是个经典形式,可以写成 f_i 的后缀和,那么期望答案为 $\sum s_i$,其中 s_i 是 f_i 的后缀和。考虑 s_i 有什么意义,一个比较容易想到的结论是 s_i 是放了 i-1 台灯依然合法 的概率。比如 n = 7,s_5 = f_5 + f_6 + f_7。注意到概率可加等价于各事件独立,这里再考察 f_i 的组合意义:放了恰好 i 台灯的停止的概率,所以取 f_{i+1} 的条件就是没有取到 f_i($1-f_i$,即放了 i 台灯依然合法,这和 f_i 的定义完全矛盾)因此各事件是独立的。

这就引出了另一个 trick:如果有概率的和式,可以考虑事件是否独立从而构造组合意义。

那么可以写出 $\displaystyle s_i = \frac{\binom{n-(i-2)\times (k-1)}{i-1}}{\binom{n}{i-1}}$。

using mint = modint1000000007;
int n, k;
mint fct[100100], ifct[100100];
void init(int maxn) {
  fct[0] = 1;
  for (int i=1; i<=maxn; ++i) fct[i] = fct[i-1]*i;
  ifct[maxn] = fct[maxn].inv();
  for (int i=maxn-1; i>=0; --i) ifct[i] = ifct[i+1]*(i+1);
}
mint gc(int x, int y) {
  if (x < y) return 0;
  return fct[x]*ifct[x-y]*ifct[y];
}
void mcase() {
  cin >> n >> k;
  mint ans = 1;
  for (int i=2; i<=n; ++i) {
    if (n-(i-2)*(k-1) <= 0) break;
    ans += gc(n-(i-2)*(k-1), i-1)/gc(n, i-1);
  }
  cout << ans.val() << "\n";
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  init(1e5);
  int tt;
  for (cin>>tt; tt--;) mcase();
}

「codeforces - 1237F」Balanced Domino Placements:第一步就不会了,太弱小了。分 $n$ 的奇偶讨论可以玩出来最优决策就是随便走,后手必胜,证明还没研究。然后发现终局一定是在 ABABA... 或者 BABAB... 的中间填一些空格,注意连续的空格不能超过一个,否则就是非法终局。

那么枚举 $x$ 表示一共走了多少步。

using mint = modint1000000007;
int n;
mint fct[1000100], ifct[1000100];
void init(int const maxn) {
    fct[0] = 1;
    for (int i=1; i<=maxn; ++i) fct[i] = fct[i-1]*i;
    ifct[maxn] = fct[maxn].inv();
    for (int i=maxn-1; i>=0; --i) ifct[i] = ifct[i+1]*(i+1);
}
mint gc(int x, int y) {
    if (x < y) return 0;
    return fct[x]*ifct[x-y]*ifct[y];
}
signed main() {
    ios::sync_with_stdio(0);
    cin.tie(0);
    init(1e6);
    cin >> n;
    mint ans = 0;
    for (int i=0; i<n; ++i) {
        if (i%2) {
            ans += gc(i+1, n-i-1)*fct[i];
        }
    }
    ans *= n*2;
    cout << ans.val() << "\n";
}

「codeforces - 1295F」 Good Contest:写成计数,缩减值域,用给出的区间把值域划分成 $O(n)$ 个部分(通常在通过取出区间左右端点混着排序划分值域或者下标啥的的时候,需要将区间变成左闭右开,避免一个元素同时属于多个划分后的部分),把第 j 个部分记作 $I_j = [d_j, d_{j+1})$。如果你设 $f[i][j]$ 表示决策了 i 个元素,最后一个元素放在 $I_j$,就会发现 $f[i+1][j+t]$ 好转移,然而 $I_j$ 内部的转移很难办。这里涉及一个 trick:整体转移,即一次决策把当前决策层做完。对于这个题,我们就改 $f[i][j]$ 的描述为「决策了 i 个元素,下一个元素放到 $I_j$ 后面的部分去的方案数」,那么转移即 $f[i+k][j+t] \gets f[i+k][j+t] + f[i][j] \times \binom{|I_{j+t}|+k-1}{k}$,后面那个组合数即在 $I_{j+t}$ 这个区间中选出一个由 k 个元素构成的不严格递增子序列的方案数(插板)。

但是这个转移我写了过不了啊,也不知道为什么,总之不求甚解抄了下飘飘蛋的转移…… 周末记得回档。

/*
有些题你两年前做不来两年后还是做不来
肯定写成计数
用区间把值域划成 O(n) 段
要把区间搞成 [)
f[i][j] 决策了 i 个元素, last element 在段 j
f[i+1][j+t] += f[i][j] * len[j+t]
分开 dp
f[i][j] 是后继不在段 j 的方案数
g[i][j][k] 是后继在段 j, 第 i 个元素放在段 j 的位置 k 的方案数
g[i+1][j][k] += g[i][j][k] * (realR[j] - k)
不行
k 记录了就寄
trick: 把一个决策阶段一下考虑完
那么 f[i][j], 考虑了前缀 [1, i], 都在前 j 段, 并且后继不能属于 j 方案数
f[i+k][j+t] += f[i][j] * C(len[j+t]+k-1, k)
I_{j+t} \subseteq Q_{i+k}
满足 i+k 的范围包括了 j
预处理组合数
f[i][j] += f[k][j-t] * C(len[j]+i-k-1, i-k), k from i-1 to 0
*/
using mint = modint998244353;
int n, xl[233], xr[233], len[233], dc[233], tot;
mint dp[233][233], c[233], iv[233];
mint gc(int x, int y) {
    mint res = 1;
    for (int i=x; i>=x-y+1; --i) res *= i;
    for (int i=y; i>=1; --i) res *= iv[i];
    return res;
}
signed main() {
    // ios::sync_with_stdio(0);
    // cin.tie(0);
    cin >> n;
    iv[0] = iv[1] = 1;
    for (int i=2; i<233; ++i) iv[i] = -998244353/i*iv[998244353%i];
    for (int i=1; i<=n; ++i) {
        cin >> xl[i] >> xr[i], xr[i]++;
        len[i] = xr[i]-xl[i];
        dc[++tot] = xl[i], dc[++tot] = xr[i];
    }
    sort(dc+1, dc+tot+1);
    tot = unique(dc+1, dc+tot+1)-dc-1;
    for (int i=1; i<=n; ++i) {
        xl[i] = lower_bound(dc+1, dc+tot+1, xl[i])-dc;
        xr[i] = lower_bound(dc+1, dc+tot+1, xr[i])-dc;
    }
    xl[0] = 1, xr[0] = tot+1;
    for (int i=1; i<=tot; ++i) dp[0][i] = 1;
    for (int i=1; i<=n; ++i) {
        for (int j=xl[i]; j<xr[i]; ++j) {
            for (int k=i; k>=1; --k) {
                if (j < xl[k] || j >= xr[k]) break;
                dp[i][j] += dp[k-1][j+1]*gc(dc[j+1]-dc[j]+i-k, i-k+1);
            }
        }
        for (int j=tot-1; j>=1; --j) dp[i][j] += dp[i][j+1];
    }
    mint dt = 1;
    for (int i=1; i<=n; ++i) dt *= len[i];
    cout << (dp[n][1]*dt.inv()).val() << "\n";
}

「codeforces - 1188C」Array Beauty:不容易一眼观察到,可以对数组排序(总之是子序列,根据美丽值的计算方式可知可以排序,其他贡献形式不知道可不可以,记得后面思考)。考虑反过来算贡献系数(这也是一个 trick,尤其是这种长得一脸坐标轴的贡献),即对每一个 $x$,算出 $f(x)$ 表示「美丽值至少为 x 的子序列的数量」,则答案为 $\sum f(x)$。设 $f[i][j]$ 表示「长度为 i 的子序列的数量,满足其以 a_j 结尾」,朴素转移比较平凡,注意到转移中「$a_j - a_u \geqslant x$」的限制,又结合对原数组的排序,可以知道这个有单调性,于是前缀和优化即可。最后还要注意到 $x$ 最大取 $\frac{10^5}{k-1}$。 证明考虑贡献形式是 min,比如我要让 10^5 - 1 作为 min,又要有恰好 k 个元素,必然不可能。

/*
算贡献
排序
钦定至少 x 的beauty
f[len][j] subsequences with length len, ending with a[j]
f[len][j] += f[len-1][u], a[j] - a[u] >= x, u < j
1e5/(k-1)
*/
using mint = modint998244353;
int n, kk, a[1100], ptr[1100];
mint ans, dp[1100][1100], cum[1100];
mint getans(int X) {
  for (int i = 0; i <= kk; ++i) {
    for (int j = 0; j <= n; ++j) dp[i][j] = 0;
  }
  ptr[0] = 0;
  for (int i=1; i<=n; ++i) {
    ptr[i] = ptr[i-1];
    while (a[i]-a[ptr[i]+1] >= X) ptr[i]++;
  }
  for (int i=0; i<=n; ++i) dp[1][i] = 1;
  for (int i=1; i<=n; ++i) cum[i] = cum[i-1]+dp[1][i];
  for (int i=2; i<=kk; ++i) {
    for (int j=1; j<=n; ++j) dp[i][j] = cum[ptr[j]];
    for (int j=1; j<=n; ++j) cum[j] = cum[j-1]+dp[i][j];
  }
  mint res = 0;
  for (int i=1; i<=n; ++i) res += dp[kk][i];
  return res;
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  cin >> n >> kk;
  for (int i=1; i<=n; ++i) cin >> a[i];
  sort(a+1, a+n+1);
  for (int i=1; i<=100000/(kk-1); ++i) ans += getans(i);
  cout << ans.val() << "\n";
}

「codeforces- 1237F」Balanced Domino Placements:很有启发性的一道题目。首先需要把问题转化为二分图匹配。如果不考虑题目中 k 个既有骨牌的限制,并把骨牌看作一个 $1 \times 1$ 的物体,我在 $(i, j)$ 上放一个骨牌就意味着行 i 和列 j 匹配了(which is a common trick)。现在把骨牌还原成 $1 \times 2$ 或 $2 \times 1$ 的物体,则问题这样转化:「把行看作一个长度为 $n$ 的纸带,列看作长度为 $m$ 的纸带,那么原题即在行和列纸带中选出一些长度为 2 的黑条(黑条之间不交不包含),纸带中的黑条与另一条纸带中的空位匹配」。可以看出这个二分图上的一个匹配就对应原问题的一个方案。那么写出答案的式子:$\displaystyle \sum_{i = 0}^{\lfloor\frac{n}{2}\rfloor} \sum_{j = 0}^{\lfloor\frac{m}{2}\rfloor} \binom{n-2i}{j} \binom{m-2j}{ i} f_i g_j i! j!$,其中 f_i,g_i 表示在纸带上选出 $i$ 个长度为 2 的黑条的方案数,直接 dp 即可。

/*
考虑成行和列匹配
两行, 每次占两个, 可以不占, 匹配另一行的空格
答案: \sum_{i=1}^{n/2} \sum_{j=1}^{m/2} C(n-2*i, j) * C(m-2*j, i) * f[i] * g[j]
f[i][j] = f[i-1][j]+f[i-2][j-1]
*/
using mint = modint998244353;
bool row[3700], col[3700];
int n, m, r, c, kk;
mint f[3700][3700], g[3700][3700], fct[3700], ifct[3700];
void init(int maxn) {
  fct[0] = 1;
  for (int i = 1; i <= maxn; ++i) fct[i] = fct[i-1]*i;
  ifct[maxn] = fct[maxn].inv();
  for (int i = maxn-1; i >= 0; --i) ifct[i] = ifct[i+1]*(i+1);
}
mint perm(int x, int y) {
  if (x < y) return 0;
  return fct[x]*ifct[x-y];
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  cin >> n >> m >> kk;
  init(max(n, m));
  for (int i = 1, r1, r2, c1, c2; i <= kk; ++i) {
    cin >> r1 >> c1 >> r2 >> c2;
    row[r1] = row[r2] = col[c1] = col[c2] = 1;
  }
  r = count(row+1, row+n+1, 0), c = count(col+1, col+m+1, 0);
  f[0][0] = 1;
  row[0] = col[0] = 1;
  for (int i = 1; i <= n; ++i) {
    for (int j = 0; j <= i/2; ++j) {
      f[i][j] = f[i-1][j];
      if (i >= 2 && j > 0 && !row[i] && !row[i-1]) f[i][j] += f[i-2][j-1];
    }
  }
  g[0][0] = 1;
  for (int i = 1; i <= m; ++i) {
    for (int j = 0; j <= i/2; ++j) {
      g[i][j] = g[i-1][j];
      if (i >= 2 && j > 0 && !col[i] && !col[i-1]) g[i][j] += g[i-2][j-1];
    }
  }
  mint ans = 0;
  for (int i = 0; i <= r/2; ++i) {
    for (int j = 0; j <= c/2; ++j) ans += perm(r-2*i, j)*perm(c-2*j, i)*f[n][i]*g[m][j];
  }
  cout << ans.val() << "\n";
}

「codeforces - 1111D」Destroy the Colony:背包还能删除…… 长见识了。在没有 x 和 y 的限制的时候,容易看出这就是个计数背包,考虑本质不同的询问个数只有 $|\Sigma|^2$ 个,考虑可持久化询问,那么时间复杂度就是 $O(|\Sigma|^3 \times n)$。计数 01 背包本质上是一个高维前缀和,那么如果钦定一个维度不选,那么我们完全可以删除这个温度的影响。$O(|\Sigma|^2 n)$

兔兔 2h 写了个复杂度一样但是跑满 T 了的 dp,默哀。

using mint = modint1000000007;
int n, q, freq[60];
char s[100100];
mint dp[100100], fct[100100], ifct[100100], ans[60][60];
void init(int mxn) {
  fct[0] = 1;
  for (int i = 1; i <= mxn; ++i) fct[i] = fct[i-1]*i;
  ifct[mxn] = fct[mxn].inv();
  for (int i = mxn-1; i >= 0; --i) ifct[i] = ifct[i+1]*(i+1);
}
int inline trs(char c) {
  if ('a' <= c && c <= 'z') return c-'a'+1;
  return 26+c-'A'+1;
}
void precompute() {
  dp[0] = 1;
  for (int i = 1; i <= 52; ++i) {
    if (!freq[i]) continue;
    for (int j = n; j >= freq[i]; --j) dp[j] += dp[j-freq[i]];
  } // dp[s]: 使用所有物品凑出 s 的方案数.
  mint prm = fct[n/2].pow(2);
  for (int i = 1; i <= 52; ++i) {
    if (freq[i]) prm *= ifct[freq[i]];
  }
  for (int i = 1; i <= 52; ++i) {
    if (!freq[i]) continue;
    for (int j = freq[i]; j <= n; ++j) dp[j] -= dp[j-freq[i]];
    for (int j = i+1; j <= 52; ++j) {
      if (!freq[j]) continue;
      for (int k = freq[j]; k <= n; ++k) dp[k] -= dp[k-freq[j]];
      ans[i][j] = ans[j][i] = 2*dp[n/2];
      for (int k = n; k >= freq[j]; --k) dp[k] += dp[k-freq[j]];
    }
    for (int j = n; j >= freq[i]; --j) dp[j] += dp[j-freq[i]];
  }
  for (int i = 1; i <= 52; ++i) ans[i][i] = dp[n/2];
  for (int i = 1; i <= 52; ++i) {
    for (int j = 1; j <= 52; ++j) ans[i][j] *= prm;
  }
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  cin >> s >> q;
  n = strlen(s);
  init(n);
  for (int i = 0; i < n; ++i) freq[trs(s[i])]++;
  precompute();
  for (int x, y; q--;) {
    cin >> x >> y;
    cout << ans[trs(s[x-1])][trs(s[y-1])].val() << "\n";
  }
}

「codeforces - 814E」An unavoidable detour for home:unweighted shortest paths 问题,考虑 bfs 树的性质,即非树边在同层或相邻层。因为这题说最短路唯一,所以非树边只能在同层,又因为编号大的点到 1 的最短路大于等于编号小的点到 1 的最短路,因此每一层的编号是有序的。于是 dp 即可,当我们加入点 $u$ 的时候,我们会关心上一层和当前层能够连边的结点数,这就是状态了。

/*
yuanlai unweighted, rnm, tq
bfs trees
non-tree arcs take places in the same levels
add vertexes from 1 to n
dp[u][i][j][k][l]: adding vertex u, i vertexes with 1 out-degree in the previous level,
  j vertexes with 2 out-degree in the previsou level, k ... 1 .. current, l ... 2 ... current
fuck
*/
using mint = modint1000000007;
int n, d[233];
mint dp[2][51][51][51][51];
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  cin >> n;
  for (int i = 1; i <= n; ++i) cin >> d[i];
  int cur = 1;
  dp[cur][d[1] == 2][d[1] == 3][d[2] == 2][d[2] == 3] = 1;
  for (int u = 2; u < n; ++u) {
    memset(dp[cur^1], 0, sizeof dp[cur^1]);
    for (int i = 0; i <= n; ++i) {
      for (int j = 0; i+j <= n; ++j) {
        for (int k = 0; i+j+k <= n; ++k) {
          for (int l = 0; i+j+k+l <= n; ++l) {
            if (dp[cur][i][j][k][l].val() == 0) continue;
            mint upd = dp[cur][i][j][k][l];
            if (i == 0 && j == 0) dp[cur][k][l][0][0] += upd;
            else {
              for (int qwq = 0; qwq <= 1; ++qwq) {
                int i_ = i, j_ = j, mc;
                if (qwq == 0) { // use 1-plug vertexes in the previous level
                  mc = i, i_--;
                  if (i == 0) continue;
                } else if (qwq == 1) { // use 2-plug vertexes in the previous level
                  mc = j, i_++, j_--;
                  if (j == 0) continue;
                }
                // i_: amount of vertexes with ... in the previous level
                // category: d[u+1]
                if (d[u+1] == 2) {
                  dp[cur^1][i_][j_][k+1][l] += upd*mc;
                  if (k > 0) dp[cur^1][i_][j_][k-1][l] += upd*mc*k;
                  if (l > 0) dp[cur^1][i_][j_][k+1][l-1] += upd*mc*l;
                } else {
                  dp[cur^1][i_][j_][k][l+1] += upd*mc;
                  if (k > 0) dp[cur^1][i_][j_][k][l] += upd*mc*k;
                  if (l > 0) dp[cur^1][i_][j_][k+2][l-1] += upd*mc*l;
                  if (k > 1) dp[cur^1][i_][j_][k-2][l] += upd*mc*k*(k-1)/2;
                  if (l > 1) dp[cur^1][i_][j_][k+2][l-2] += upd*mc*l*(l-1)/2;
                  if (k > 0 && l > 0) dp[cur^1][i_][j_][k][l-1] += upd*mc*k*l;
                }
              }
            }
          }
        }
      }
    }
    cur ^= 1;
  }
  cout << dp[cur][0][0][0][0].val() << "\n";
}

「hdu - 7060」Separated Number

我就会你妈个这

$$ C+\sum_{k=2}^K \left(\sum_{i=1}^n s_i \left( 10^{n-i} \times \binom{i-1}{k-1} + \sum_{j=0}^{n-i} 10^j\times \binom{n-j-2}{k-2} \right) \right) $$

其中 $C$ 是 k=1 时的答案,推导过程就考虑一个位置 i 作为 10^k 时的贡献系数即可。注意到如果我们能求出组合数前缀和,这个问题就解决了。

推导 $\displaystyle f_n(m) = \sum_{i=0}^m \binom{n}{i}$。考虑 $f_n$ 和 $f_{n-1}$ 的关系,放到杨辉三角上,有 $f_n(m) = f_{n-1}(m)-\binom{n-1}{m}$。

/*
Ans for k=1+\sum_{k=2}^K \left(\sum_{i=1}^n s_i \left( 10^{n-i} \times \binom{i-1}{k-1} + \sum_{j=0}^{n-i} 10^j\times \binom{n-j-2}{k-2} \right) \right) \\
我就会你妈个这
*/
using mint = modint998244353;
mint operator ""_m(ull x) {
    return mint(x);
}
mint fct[1000100], ifct[1000100], pw[1000100], pre[1000100], rep[1000100];
int n, kk;
char s[1000100];
void init(int up) {
  pw[0] = 1;
  for (int i = 1; i <= up; ++i) pw[i] = pw[i-1]*10;
  fct[0] = 1;
  for (int i = 1; i <= up; ++i) fct[i] = fct[i-1]*i;
  ifct[up] = fct[up].inv();
  for (int i = up-1; i >= 0; --i) ifct[i] = ifct[i+1]*(i+1);
}
mint comb(int x, int y) {
  if (x < 0 || y < 0 || x < y) return 0;
  return fct[x]*ifct[x-y]*ifct[y];
}
void mcase() {
  cin >> kk >> s;
  n = strlen(s);
  pre[0] = kk-1 >= 0;
  rep[0] = kk-2 >= 0;
  for (int i = 1; i <= n; ++i) rep[i] = rep[i-1]*2-comb(i-1, kk-2);
  for (int i = 1; i <= n; ++i) pre[i] = pre[i-1]*2-comb(i-1, kk-1);
  mint ans, sum;
  for (int i = n; i >= 1; --i) {
    if (i < n) sum += rep[i-1]*pw[n-i-1];
    ans += (s[i-1]-'0')*(pre[i-1]*pw[n-i]+sum);
  }
  cout << ans.val() << "\n";
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  init(1e6);
  int tt;
  for (cin >> tt; tt--;) mcase();
}

「hdu - 6822」Paperfolding:有意思。显然 [L R] 和 [U D] 对答案的影响是独立的,在现实中模拟可知答案是 $\frac{\sum (2^i+1)(2^{n-i}+1)\binom{n}{i}}{2^n}$,等差数列求和即可。

/*
有意思
至少四片, 考虑怎样的切痕会断
影响独立
分开算
*/
using mint = modint998244353;
mint operator ""_m(ull x) {
    return mint(x);
}
void mcase() {
  ll n; cin >> n;
  cout << ((2_m).pow(n).inv()*2*(3_m).pow(n)+(2_m).pow(n)+1).val() << "\n";
}
signed main() {
  ios::sync_with_stdio(0);
  cin.tie(0);
  int tt;
  for (cin >> tt; tt--;) mcase();
}

平面旋转。应该是比较好理解的版本。

我们对一个平面(逆时针)旋转 $\beta$ 度,无非就是对每一个有意义的向量 $\boldsymbol a = (x, y)$ 进行旋转。不妨考察单位向量 $\boldsymbol e = (\cos \alpha, \sin \alpha)$,令 $\boldsymbol a = k \cdot \boldsymbol e$,则由三角恒等变换得 $\boldsymbol e' = (\cos(\alpha + \beta), \sin (\alpha + \beta)) = (\cos \alpha \cos \beta - \sin \alpha \sin \beta, \sin \alpha \cos \beta + \sin \beta \cos \alpha)$,则 $\boldsymbol a' = k \cdot e'$。

令 $\beta = \frac{\pi}{4}$ 再把 $k$ 除以 $\sqrt 2$ 就得到切比雪夫距离转曼哈顿距离的式子了。

高维前缀和

给一个 intrada 性质的问题:

求 $\displaystyle F[mask] = \sum_{i \subseteq mask} A[i]$

这个形式看起来会很像一个 and-convolution,虽然并不完全是但这很重要。有个经典的朴素做法是以 $O(3^n)$ 枚举子集,从这个做法可以看出,$A[x]$,其中 $x$ 有 $k$ 个 off-bits,会被 $2^k$ 次计算。

到了这个时候网上的脑瘫博客就会开始告诉你二维怎么办?容斥!三维不想推怎么办?分维前缀和!$k$ 维怎么办?看代码。

事实上分维前缀和的思想的确非常重要,如果有人看这篇博客推荐先理解。

我们首先对 sum over subsets 问题和高维前缀和之间的关联建立认知。一个 multidimensional partial sums 具体是在一个 $k$-d 空间里求出向量 $((1, \dots, 1), (a_1, \dots, a_k))$ 的所有元素的运算结果,也就是说我们用 $(k, (n_1, \dots n_k))$ 这样一个向量来刻画问题的背景;而 sum over subsets 的一个 alternative ver. 则是给定 $k$ 个物品及其价值,询问全集的某个子集的所有子集的价值和的价值和,以及我们只需要 $k$ 就可刻画问题的背景(所有维度的有意义长度均为 $2$,即选与不选),而传统的 sum over subsets 则是给出全集所有子集的价值,询问全集的某个子集的所有子集价值和(即将某个子集的价值从由子集中的元素的价值决定到和子集中的元素的价值没关系)。

然后开始想象,typical sum over subsets problem 是在一个 $k$-d 空间中做高维前缀和做的事情,所以我们可以类比得出其解法,并且我们得出结论:typical sum over subsets problem 是 multidimensional partial sums problem 的弱化,区别在于每一维度的模长是 $\boldsymbol 2$ 还是任意

应该比网上大部分博客都要清楚,或者更接近本质吧。不过我觉得这个洞见显然不是 SOS 的本质。

for (int i=0; i<k; ++i) {
    for (int j=0; j<(1<<k); ++j) {
        if ((j>>i)&1) f[j] += f[j^(1<<i)]; // 相当于从于第 @i 维的 0 转移到 1
    }
}

dirichlet partial sums

这个东西就是上述问题的一个演绎应用,问题背景是:

已知数论函数 $g(T)$,求 $\displaystyle f(T) = \sum_{d \mid T} g(d)$。

$O(n \ln n)$ 的做法应该都会。

会发现这个的形式和上面的很像,不过直接套会有一个 trap。如果你想的是把数分解为 $\displaystyle \prod^k_{i=1} p_i^{c_i}$ 的形式后考虑一个 $k$-d 空间是不行的,因为 SOS 的维度模长为 $2$,所以应该考虑一个 $\left( \sum c_i+1 \right)$-d 的空间。复杂度是 $O(n \log \log n)$,dyh 有个高妙的证法,但是我忘了,就这。

for (int i=1; i<=tot; ++i) {
    for (int j=/*lower bound, depends on the specific problem*/; j<=n/prime[i]; ++j) {
        f[j*prime[i]] += f[j];
    }
}

当然这个有四个形式,不过很搞笑罢了,就是因数倍数乘上正着求和逆着求,本质都一样。

e.g.

看两道题。

第一个就是 codeforces - 585e,这个题求 $f$ 的部分我是直接上的莫反,结果只能调和级数,翻了下题解发现式子里面的 $\mu$ 直接没了,感觉很神奇,有空再想想。

第二个是 codeforces - 449d,一样的套路。记 $f_i$ 为子集 AND 和为 $i$ 的方案数,根据莫反的套路设出 $g_i$ 为子集 AND 和是 $i$ 超集的方案数。如果我们能求出 $h_i$ 表示 $i$ 的超集数量,则 $g_i = 2^{h_i}-1$,然后通过逆 sos 跑差分就行了,注意这里是高维后缀和。

从这两道题我们可以洞见 sos 和莫反有着莫大的关联,之前陀螺仪在 oj 指出过,现在才理解一点点,果然还是要看太阳神的。

之前也看到有人写出莫反不关于莫比乌斯函数的形式,有空 cancan。

gugugu……

[TODO]:莫反不关于莫比乌斯函数的形式;FMT;gcd-convolution。