s_r_f 的博客

s_r_f 的博客

这里是一个只会背板和fst的蒟蒻

[OI笔记]数位DP合集 & 对数位DP的一点理解

posted on 2020-04-28 10:28:25 | under OI笔记 |

$About$

几天前和 $pcf$一起口胡了一波典型的数位 $DP$题 $,$

感觉如果不写代码的话就要变成口胡选手了 $,$就来写写 $blog$ $+$ 刷一波题。

前置芝士: $DP$

本文仅针对对 $DP$有一定了解 $,$并有一定基础的选手 $.$

本文不会过多提及数位 $DP$本身 $,$主要还是讲应用 $.$

定义

数位 $DP,$就是一类和数位有关的 $DP,$一般是计数 $/$求和 类型的题 $.$

数位 $DP$有两种理解方式 $:$

第一种理解

数位 $DP$的状态只和 已经确定好的数的信息(比如说值/位数/前几个数等)没有确定好的数字的位数 有关 $.$

在状态比较多的题中 $,$推荐使用记忆化搜索来实现 $DP$的过程 $.$

数位 $DP$分为两个部分 $:$ $DP$ 和 对给你的数查询答案 $.$

其中 $,$ 对给你的数查询答案 这一部分的本质就是 $,$ $O($位数 × 进制 $)$次关于 $DP$ 状态的查询 $,$ 下文中除非特殊情况不会再赘述 $.$


第二种理解

在状态中计入 $0/1$表示是否达到上界 $,$就是另一种 $DP$方式 $,$但是对于每组数据需要重新 $DP$一遍 $.$


提示:

下文中复杂度分析中可能存在字母 $T,N,M,$ 在无特殊说明的情况下 $T$为询问组数 $,$ $N$为进制 $,$ $M$为最大位数 $.$

$Problems$

$Part.I$ 基础数位 $DP$


洛谷 P4317 花神的数论题

比较基础的数位 $DP.$

考虑记 $dp(n,m)$表示有 $n$个二进制位 $,$其中有 $m$个 $1$的二进制数个数 $.$

$($ 实际上这个 $dp(n,m)$ 就是 $C_{n}^{m},$但因为组合数学和本文内容没什么关系 $,$故省略 $.$ $)$

为什么这个题里面没有考虑有没有前导 $0$的情况呢 $?$因为前导 $0$对答案没有影响 $.$

然后就可以直接回答 $O(log)$次关于 $DP$状态的查询了 $.$

复杂度 $O(M^2+TNM),$其中 $T=1,M=60,N=2$

$DP$部分代码 $:$

LL dp[64][64]; bool vis[64][64];
inline LL DP(int n,int m){
    if (n < m || m < 0) return 0; if (n == 0) return 1;
    if (vis[n][m]) return dp[n][m]; vis[n][m] = 1;
    return dp[n][m] = DP(n-1,m) + DP(n-1,m-1);
}

SP1433 KPSUM - The Sum

这个题有两种思路 $:$

一种是枚举交错和然后求出每种交错和的数的个数 $,$并且还要考虑奇偶性 $,$比较繁琐 $.$

另一种是直接把这些数作为一个整体来 $dp$交错和 $.$

记 $dpsum(presum,prelen,n,1/0)$表示已经确定的前缀的交错和为 $presum$ $,$已经确定的前缀排除掉前导零的长度为 $prelen,$还有 $n$位没有确定 $,$目前的位数中是 $/$否有非零数 的交错和 $.$

为了 $dp$这个数值 $,$我们还需要知道 $dplen(prelen,n,1/0)$表示已经已经确定的前缀排除掉前导零的长度为 $prelen$ $,$ 还有 $n$位没有确定 $,$目前的位数中是 $/$ 否有非零数 的数字的长度总和 $.$

$dp$转移的时候考虑下位置的奇偶性即可 $.$

注意到 $len/prelen$本身并不重要 $,$只和奇偶性有关 $,$所以 $dp$过程中的 $len/prelen$可以对 $2$ 取模 $,$ 可以进一步压缩状态 $.$

复杂度 $O(N^2M+TNM),$其中 $T$的范围没给 $,$ $N=10,M=15$

代码见 题解链接


「BalticOI 2013」非回文数 Palindrome-Free Numbers

这题相比上一题的交错和而言 $,$思维难度有一点提高 $,$而代码难度却大幅降低。

考虑如果存在回文 $,$那么一定存在长度 $\leq3$的回文 $,$所以我们只要记录上一个数字和上上个数字就可以了 $.$

然后同样的 $,$很明显前导 $0$会影响答案 $,$所以记状态的时候再记一维 $0/1$即可 $.$

复杂度 $O(N^3M+TNM),$其中 $T=1,N=10,M=18$

LOJ评测链接

为啥我代码是LOJ上最优解啊(雾,不过70ms一共67个点,有没有人挑战一下卡到每个点只跑1ms啊

一个需要高精度 $(M=1000)$的双倍经验 $:$ P3413 SAC#1 - 萌数


LOJ#6274. 数字

考虑数位 $dp,$记 $dp(n,lx,rx,ly,ry)$ $($ 后4维都是 $0/1$ $)$表示目前在决策第 $n$位 $,$ 目前 $x/y$的值是否正好 $=$上 $/$下界

然后注意一下在转移 $x$ $and$ $y=0$ 的位数的时候 $,$要把得到的答案取 $max.$

代码见评测链接


$Part.II$ 对状态进行剪枝或合并

CF55D Beautiful numbers

双倍经验 $:$ SP8177 JZPEXT - Beautiful numbers EXTREME但是似乎出锅了不能 $remote$ $judge?$

首先不难设出状态 $:$ $dp(s,r,n)$ 表示目前已经确定下来的前缀 $mod$ $2520$ 为 $r,$ 目前的 $1..9$ 是否存在状压到 $s$ 里 $,$ 还有 $n$ 位没有确定 $.$

但如果直接设状态状态总数为 $2\times10^7,$再乘上每次转移是 $O(10)$ 显然会 $T$ 飞

由于我们只需要考虑集合里面数的 $lcm$大小 $,$所以有些状态是可以合并的 $,$就写个预处理把它合并一下就好了 $.$

$P.S:$ 我代码里面没有合并干净 $,$是 $73$种状态 $,$但是也能过 $.$

如果真正合并干净了那就是 $48$种状态 $.$

复杂度 $O(V1*V2*NM+TNM),$其中 $V1$为 $lcm_{i=1}^{10}i=2520,$ $V2$为合并过后的状态总数 $.$

SP8177题解链接

$DP$部分代码 $:$

LL pw[19];
int trans[1<<9],tid[1<<9]; bool used[1<<9];
int cnt,v[73+5],nxt[73+5][10];
inline void init(){
    int i,j,k,id,l = 1<<9,s;
    for (pw[0] = i = 1; i <= 18; ++i) pw[i] = pw[i-1] * 10;
    for (i = 0; i < l; ++i){
        trans[i] = i;
        for (j = 9; j >= 1; --j) if (trans[i] >> (j-1) & 1)
        for (k = j + j; k <= 9; k += j) if (trans[i] >> (k-1) & 1){
            trans[i] &= ((1<<9)-1)^(1<<(j-1)); break;
        }
        used[trans[i]] = 1;
    }
    for (i = 0; i < l; ++i) if (used[i]) ++cnt,v[cnt] = i,tid[i] = cnt;//cnt==73
    for (i = 0; i < l; ++i) if (id = tid[i]){
        for (v[id] = 1,j = 1; j <= 9; ++j) if (i>>(j-1)&1) v[id] = lcm(v[id],j);
        for (nxt[id][0] = id,j = 1; j <= 9; ++j) nxt[id][j] = tid[trans[i|(1<<j-1)]];
    }
}
inline bool check(int id,int r){ return r % v[id] == 0; }
LL dp[74][2520][20]; bool vis[74][2520][20];
inline LL DP(int s,int r,int n){
    if (!n) return check(s,r) ? 1 : 0;
    if (vis[s][r][n]) return dp[s][r][n]; vis[s][r][n] = 1;
    LL tot = 0;
    for (int i = 0; i <= 9; ++i) tot += DP(nxt[s][i],(r*10+i)%2520,n-1);
    return dp[s][r][n] = tot;
}

[AHOI2009]同类分布

首先枚举余数 $mod$ $,$然后对于每个余数 $DP.$

设 $dp(s,r,n)$表示目前已经确定的前缀的数字的总和为 $s,$目前已经确定的前缀对 $mod$取模之后的结果 $.$

然后就可以过了 $.$

但是如果直接 $dp$的话 $,$一个点在 $luogu$机器上开 $O2$要跑 $1.2s,$所以可以剪枝来优化常数 $.$

剪枝之后不开 $O2$能跑在 $1s$ 以内 $.$

因为 $mod$的值域为 $O(NM),$而对每个 $mod$ $dp$的状态数总数是 $O(N^3M^2)$的

所以复杂度 $O(N^5M^3+TNM),$其中 $N=10,M=18,T=1$

当然因为加了剪枝以及实际上的 $mod$只到 $M*9=162,$所以效率还是可以的 $.$

$DP$ 部分代码 $:$

const int S = 162+1;
LL dp[S][S][19]; bool vis[S][S][19]; int M;
inline void Clear(int mod){
    register int i,j,k;
    M = mod;
    for (i = 0; i <= mod; ++i) for (j = 0; j < mod; ++j) for (k = 0; k <= 18; ++k) vis[i][j][k] = dp[i][j][k] = 0;
}
inline LL DP(int s,int r,int n){
    if (s < 0 || s > n * 9) return 0; if (!n) return (s == 0 && r == 0) ? 1 : 0;
    if (vis[s][r][n]) return dp[s][r][n]; vis[s][r][n] = 1;
    LL tot = 0; for (int i = 0; i <= 9; ++i) tot += DP(s-i,(r*10+i)%M,n-1);
    return dp[s][r][n] = tot;
}
inline LL solve(LL n){
    int mod,i,j,x,s,r; LL tot = 0; bool flag = 0;
    for (mod = 1; mod <= 162; ++mod){
        Clear(mod); flag = 0;
        s = r = 0;
        for (i = 18; i >= 0; --i) if ((x = n/pw[i]%10) || flag){
            for (j = 0; j < x; ++j) tot += DP(mod-s-j,(r*10+j)%M,i);
            s += x,r = (r*10+x) % M; flag = 1;
        }
    }
    return tot;
}

[SDOI2013]淘金

首先写个爆搜 $,$发现对于 $1\leq i \leq10^{12}$的 $i$ $,$ $f(i)$ 只有 $8282$ 种 $.$ 所以考虑对于这 $8282$ 种数求出有多少数字满足 $f(i) = $这个数字 $.$

数位 $dp,$ 记 $dp(x,n,1/0)$ 表示目前还有 $n$位数字没有决定 $,$除去前导零之后剩下的数字乘积为 $x$ $($ 所有的 $x$应当事先搜出来放到数组里 $),$ 目前有 $/$无非零数的方案数 $.$

求出来之后 $,$令 $a_i$为第 $i$个数字被用到的次数 $,$

把 $a_i$从大到小排序 $,$用堆贪心就可以了 $.$

复杂度 $O((8282*NM+TNM) \times hash()).$

代码见我的题解


$Part.III$ 利用其它方式优化计算


$1.$算贡献

CF908G New Year and Original Order

这题也有两种思路 $,$不过交错和那道题两种思路不管优劣都能过 $,$但是在这个题目里较劣复杂度的做法就过不去了 $.$

一种思路是对于 $O(700*10)$个询问直接整体 $DP$求答案 $,$

记 $dp[i][j]$表示目前已经安排到了数字 $i,$已经放了 $j$个数位的方案数和答案 $($ 可以用一个 $struct$来存 $)$

$.$每次 $DP$的复杂度都是 $O(700^2*10)$的 $,$总复杂度 $O(700^3*10^2)$过不去 $,$并且不能在只 $dp$一次的情况下处理多组询问 $.$

这种思路的复杂度是 $O(TM^3N^2)$

另一种思路是 $,$我们考虑贡献 $.$

记 $sum(n) = \sum\limits_{i=0}^{n-1}10^i$

那么对于一个数字答案显然为 $\sum\limits_{i=0}^{9}$ $sum($ 满足 $k\leq i$ 的数字的个数 $).$

所以对这个 $($ 满足 $k\leq i$ 的数字的个数 $)$ $dp,$就可以降低复杂度 $.$

记 $dp(n,m,k)$表示目前还有 $n$位还没确定 $,$目前已经确定下来的前缀中有 $m$位数字 $\geq$ $k$ $,$直接 $dp$即可 $.$

处理一组询问的复杂度仍然是 $O(10*700)$即 $O(NM).$

这种做法的复杂度是 $O(M^2N+TNM),$ 可以通过本题 $.$

两种思路的代码见我的题解


$2.$结合位运算/反演等其它计数工具

P3791 普通数学题

考虑到这个题是 $xor,$所以容易想到是进制为 $2$的数位 $dp.$

不难发现 $\sum\limits_{i=1}^{n}d(i)=\sum\limits_{i=1}^n \lfloor \frac{n}{i} \rfloor,$可以 $O(\sqrt n)$计算 $.$

然后就 $O(M^2)$枚举两维的限制 $,$再计算即可 $.$

怎么计算呢 $?$

设两维分别有 $a/b$位没有确定 $.$

可以发现 $,$前面有一些位置的部分是被钦定的 $,$而后一部分 $($ 后 $max(a,b)$位 $)$是没有被确定的 $,$且可以任意取 $.$

然后答案就是两个 $\sum\limits_{i=1}^{n}d(i)$的差 $ \times 2^{min(a,b)}$ $.$

不难发现虽然询问有 $O(M^2)$个 $,$但是本质不同的只有 $O(M)$个 $,$ 因为关于 $\sum d(i)$的询问只和未确定数位较多的那个数字的已确定的前缀有关 $.$

所以复杂度为 $O(\sqrt{n} \times M)=O(\sqrt{N^M} \times M),$其中 $M=33.$

注意这个复杂度里的 $n$是值域 $,$其实际值为 $N^M$

所以这个题就可以拿来加强了 $,$比如 $d(i$ $xor$ $j$ $xor$ $x$ $)-> $ $φ(i$ $xor$ $j$ $xor$ $x$ $)$或者 $d(i$ $xor$ $j$ $xor$ $x$ $)$ $->$ $μ(i$ $xor$ $j$ $xor$ $x$ $)$ /cy

代码 $:$

  const int P = 998244353,MAX = 33;
map<LL,int>FF;
inline int F(LL n){
    if (n <= 0) return 0; if (FF.count(n)) return FF[n];
    LL sum = 0,l = 1,r;
    while (l <= n) r = n/(n/l),sum = (sum + (r-l+1) * (n/l)) % P,l = r+1;
    return FF[n] = sum;
}
inline LL get(int l,int r){ return (1ll<<r+1) - (1ll<<l); }
inline int calc(LL x,LL v1,LL v2,int a,int b){
    if (a < b) swap(a,b); static LL s;
    s = (x^v1^v2)&get(a,MAX);
    return (1ll<<b)%P * (F(s+(1ll<<a)-1)-F(s-1)+P) % P;
}
LL n,m,x,ans;
int main(){
    int i,j; LL v1,v2;
    cin >> n >> m >> x; ++n,++m;
    for (i = MAX,v1 = 0; i >= 0; --i) if (n>>i&1){
        for (j = MAX,v2 = 0; j >= 0; --j) if (m>>j&1) ans += calc(x,v1,v2,i,j),v2 += 1ll<<j;
        v1 += 1ll<<i;
    }
    cout << (ans%P+P)%P << '\n'; 
    return 0;
}

CF582D Number of Binominal Coefficients

数位 $DP.$

考虑 $lucas$定理的过程 $,$设 $a = k,b = n-k,$则在 $p$进制下 $a+b$有 $>=k$次进位 $,$且 $a+b<=A.$

设 $dp[n][m][1/0][1/0]$表示当前考虑到第 $n$位 $,$ 还需要进位 $m$次 $,$这一位的数值是 $/$否紧贴着上界 $,$这一位是 $/$否往上一位进位 $.$

复杂度 $O(3500^2),$其中 $3500$为 $A$在 $p$进制下最大位数 $.$

这个题目为什么不能用 $O(3500*p)$次对 $dp$状态的询问呢 $?$

因为 $p$的值域太大 $,$并且这个题没有多组询问 $.((($

$DP$部分代码 $:$

int dp[N][N][2][2];
bool vis[N][N][2][2];
inline LL s(int l,int r){ return l > r ? 0 : (l+r) * (r-l+1ll) / 2 % P;  }
inline LL calc(int n){
    if (n < 0) return 0; if (n < p) return 1ll * (n+1) * (n+2) / 2 % P;
    n = min(n,p*2-2);
    return (1ll * (n-p+2) * p % P + s(-p+n+2,p-1)) % P;
}
inline LL calc(int l,int r){ l = max(l,0),r = min(r,p*2-2); return (l > r) ? (0) : (calc(r) - calc(l-1) + P) % P; }
inline LL DP(int n,int m,int s1,int s2){
    m = max(m,0);
    if (!n) return (!s2 && !m) ? 1 : 0;
    if (vis[n][m][s1][s2]) return dp[n][m][s1][s2]; vis[n][m][s1][s2] = 1;
    int ans = 0;
    if (s1){
        if (!s2){
            upd(ans,calc(0,p-2) * DP(n-1,m-1,1,1) % P);
            upd(ans,calc(0,p-1) * DP(n-1,m,1,0) % P);
        }
        else{
            upd(ans,calc(p-1,p+p-2) * DP(n-1,m-1,1,1) % P);
            upd(ans,calc(p,p+p-1) * DP(n-1,m,1,0) % P);
        }
    }
    else{
        if (!s2){
            upd(ans,calc(0,A[n]-1) * DP(n-1,m,1,0) % P);
            upd(ans,calc(A[n],A[n]) * DP(n-1,m,0,0) % P);

            upd(ans,calc(0,A[n]-2) * DP(n-1,m-1,1,1) % P);
            upd(ans,calc(A[n]-1,A[n]-1) * DP(n-1,m-1,0,1) % P);
        }
        else{
            upd(ans,calc(p,p+A[n]-1) * DP(n-1,m,1,0) % P);
            upd(ans,calc(p+A[n],p+A[n]) * DP(n-1,m,0,0) % P);

            upd(ans,calc(p-1,p+A[n]-2) * DP(n-1,m-1,1,1) % P);
            upd(ans,calc(p+A[n]-1,p+A[n]-1) * DP(n-1,m-1,0,1) % P);
        }
    }
    dp[n][m][s1][s2] = ans;
    return ans;
}

CF585F Digits of Number Pi

可以利用 $AC$自动机 $,$也可以利用 $SAM$来表示状态 $.$

然后数位 $dp$就行了.

dp部分代码 $:$

pair<int,int> trans[V][D][10];
struct SuffixAutomation{
    int cnt,ed; int ch[V][10],fa[V],len[V];
    inline void init(){ cnt=ed=1; }
    inline void extend(int c){
        int p = ed,np = ++cnt; len[np] = len[p] + 1;
        while (p && !ch[p][c]) ch[p][c] = np,p = fa[p];
        if (!p) fa[np] = 1;
        else{
            int q = ch[p][c];
            if (len[q] == len[p] + 1) fa[np] = q;
            else{
                int nq = ++cnt;
                fa[nq] = fa[q],len[nq] = len[p]+1,memcpy(ch[nq],ch[q],40);
                fa[q] = fa[np] = nq;
                while (ch[p][c] == q) ch[p][c] = nq,p = fa[p];
            }
        }
        ed = np;
    }
    inline pair<int,int> getnxt(int now,int l,int c){
        while (fa[now] && !ch[now][c]){
            now = fa[now]; if (l > len[now]) l = len[now];
        }
        if (ch[now][c]) now = ch[now][c],++l;
        if (l >= (m>>1)) now = l = -1;
        return make_pair(now,l);
    }
    inline void get(){
        int i,j,k;
        for (i = 1; i <= cnt; ++i) for (j = 0; j <= (m>>1); ++j) for (k = 0; k <= 9; ++k)
            trans[i][j][k] = getnxt(i,j,k);
    }
}SAM;

int f[2][2][M]; bool visf[2][2][M];
inline int F(int tx,int ty,int i){
    if (i > m) return 1;
    if (visf[tx][ty][i]) return f[tx][ty][i]; visf[tx][ty][i] = 1;
    int r = 0,c;
    for (c = 0; c <= 9; ++c){
        if (!tx && c < x[i]) continue;
        if (!ty && c > y[i]) continue;
        upd(r,F(tx|(c>x[i]),ty|(c<y[i]),i+1));
    }
    return f[tx][ty][i] = r;
}
int g[2][2][M][V][D]; bool visg[2][2][M][V][D];

inline int G(int tx,int ty,int i,int now,int l){
    if (i > m) return 0;
    if (visg[tx][ty][i][now][l]) return g[tx][ty][i][now][l]; visg[tx][ty][i][now][l] = 1;
    int r = 0,c; pair<int,int>t;
    for (c = 0; c <= 9; ++c){
        if (!tx && c < x[i]) continue;
        if (!ty && c > y[i]) continue;
        t = trans[now][l][c];
        if (t.first != -1) upd(r,G(tx|(c>x[i]),ty|(c<y[i]),i+1,t.first,t.second));
        else upd(r,F(tx|(c>x[i]),ty|(c<y[i]),i+1));
    }
    return g[tx][ty][i][now][l] = r;
}

4、更加复杂的例子

牛客挑战赛40C-小V和字符串

首先考虑 $F$是怎么计算的 $.$

不难发现我们有一种按位置扫描的暴力 $:$记录当前的答案和当前两个串 $A,B$的 $1$的个数之差 $.$

然后按这个数位 $DP,$再记个两个 $0/1$表示字典序即可 $.$

$DP$部分代码 $:$

int visc[N][2][2][N<<1]; int dpcnt[N][2][2][N<<1];
inline int Cnt(int pos,int bs,int ab,int sta){
    if (pos == n+1) return sta == 1002 ? 1 : 0;
    if (visc[pos][bs][ab][sta]) return dpcnt[pos][bs][ab][sta]; visc[pos][bs][ab][sta] = 1;
    int tot = 0;
    for (int va = 0; va <= 1; ++va) for (int vb = 0; vb <= 1; ++vb){
        if (va>vb&&!ab) continue;
        if (vb>S[pos]-'0' && !bs) continue;
        upd(tot,Cnt(pos+1,bs||(vb<(S[pos]-'0')),ab||(va<vb),sta + va - vb));
    }
    return dpcnt[pos][bs][ab][sta] = tot;
}

int viss[N][2][2][N<<1]; int dps[N][2][2][N<<1];
inline int Sum(int pos,int bs,int ab,int sta){
    if (pos == n+1) return 0;
    if (viss[pos][bs][ab][sta]) return dps[pos][bs][ab][sta]; viss[pos][bs][ab][sta] = 1;
    int tot = 0,vv;
    for (int va = 0; va <= 1; ++va) for (int vb = 0; vb <= 1; ++vb){
        if (va>vb&&!ab) continue;
        if (vb>S[pos]-'0' && !bs) continue;
        upd(tot,Sum(pos+1,bs|(vb<S[pos]-'0'),ab|(va<vb),sta + va - vb));
        vv = 0;
        if (va == vb) vv = 0;
        else if (sta == 1002){
            if (va == vb) vv = 0; else vv = P-pos;
        }
        else if (sta < 1002){
            if (va) vv = pos; else vv = P-pos;
        }
        else if (sta > 1002){
            if (va) vv = P-pos; else vv = pos;
        }
        upd(tot,(LL)vv*Cnt(pos+1,bs|(vb<S[pos]-'0'),ab|(va<vb),sta + va - vb)%P);
    }
    return dps[pos][bs][ab][sta] = tot;
}

一些例题

P1836 数页码

P4999 烦人的数学作业

[CQOI2016]手机号码

[SCOI2009] windy 数