原题大意

给出一个大数$N$其形式为：

其中 $c_i\in{0,1,\cdots,9} \ , \ 1 \le l_i \le 10^9 \, \ 1 \le k \le 10^5$

求 $\lfloor \frac{N}{M} \rfloor \pmod {10007} $

思路

我们考虑将向下取整符号去掉，

对于 $\text{ans}:=\lfloor \frac{n}{m} \rfloor \pmod {M} $ 可以化成：

我们记 $ \text{inv}(m) $ 为 $m$ 在模 $M$ 意义下的逆元 ，并令：

又由 $(*)$ 式可知：

这里 $G(l_i,M)$ 表示 $c_i\cdot(10^0+10^1+\cdots+10^{l_i-1})\pmod M$
可以使用几何级数的递推公式求出：

同理：

然后将 $r_m$ 和 $r_M$ 改成递推的形式（从高位到低位），就可以在 $O(k)$ 的时间复杂度内求出答案。

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
const ll M=10007;

ll qpow(ll a, ll b, ll m)
{
    ll res=1;
    a%=m;
    while(b)
    {
        if(b&1)res=res*a%m;
        a=a*a%m;
        b>>=1;
    }
    return res;
}

ll inv(ll x)
{
    return qpow(x,M-2,M);
}

ll sum(ll l, ll m)
{
	if(l==0)return 0;
	if(l==1)return 1;
	if(~l&1)return sum(l>>1, m)*(qpow(10, l>>1, m)+1)%m;
	return (sum(l-1,m)*10+1)%m;
}

int main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);

    ll k,m; cin >> k >> m;
    vector<pll> A(k);
    for(int i=k-1; i>=0; i--)
        cin >> A[i].first >> A[i].second;

    ll a=0, rm=0, p=1, pm=1;
    
    // 这里是从低位到高位递推的
    for(auto [c,l]:A)
    {
    	ll ca=c*sum(l,M)%M;
    	a=(a+ca*p)%M;
    	p=p*qpow(10,l,M)%M;
    	
    	ll cm=c*sum(l,m)%m;
    	rm=(rm+cm*pm)%m;
    	pm=pm*qpow(10,l,m)%m;
    }
    
    ll ans=(a-(rm%M)+M)%M*inv(m)%M;
    cout<<ans;

    return 0;
}

扩展

M不是质数

如果M不是质数，但 $gcd(m,M)≠1$ 那么逆元 $ \text{inv}(m) $ 可以用扩展欧几里得算法求出；
但 $gcd(m,M)≠1$ 时，逆元根本不存在，这是 $(**)$ 式可以变形成：

一个关于q的线性同余方程组，可以先提取公因子 $g=gcd(m,M)$ 归一化，

此时 $gcd(\frac{m}{g} , \frac{M}{g})=1$ ，就可以用扩展欧几里得算法求出特解，

在加上 $k$ 个 $\frac{M}{g}$ ，就可以得到通解：

但是 $k$ 是由 $n$ 唯一确定的，只通过 $r_m \ , \ r_M $ 而不计算 n/m 时是无法直接确定k的，应该换个思路。

显然这里 $gcd(n,M)>1$ 对于周期为 $M$ 的 $r_M$ 以及周期为 $m$ 的 $r_m$ 来说，将周期扩大到其最小公倍数 $mM$ 。
这样当n加上m时 $\lfloor n/m \rfloor$ 就增加1，继续增加到 $nM$ 后此时商就变成了M，在模M意义下又回到了0，当然是符合向下取整后取余的。
具体地，我们设n对 $mM$ 的商为Q余数为R，那么有：

我们要算的是 $\lfloor \frac{n}{m} \rfloor \pmod M$。将 $n$ 代入，得：

这样完美避开了逆元和k无法确定的问题。

m和l变大

如果 $1 \le m,l \le 10^{18}$ ，最暴力的方法就是直接使用大数模拟完成，例如以下程序：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
#include <bits/stdc++.h>
using namespace std;

using i128 = __int128;
using ll = long long;

i128 mul(i128 a, i128 b, i128 m)
{
    i128 res=0;
    while(b)
    {
        if(b&1) res=(res+a)%m;
        a=(a+a)%m;
        b>>=1;
    }
    return res;
}

pair<i128,i128> calc(i128 l,i128 m)
{
    if(l==0) return {1%m,0};
    if(l==1) return {10%m,1%m};
    auto [p,s]=calc(l>>1,m);
    i128 p2=mul(p,p,m), s2=(mul(s,p,m)+s)%m;
    if(~l&1) return {p2,s2};
    else return {mul(p2,10,m),(mul(s2,10,m)+1)%m};
}

int main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    int k;ll m,M;
    cin >> k >> m >> M;

    i128 mM=(i128)m*M, cur=0;
    for(int i=0;i<k;i++)
    {
        ll c,l; cin >> c >> l;
        auto [p,g]=calc(l,mM);
        i128 add=mul((i128)c,g,mM);
        cur=(mul(cur,p,mM)+add)%mM;
    }

    i128 ans=(cur/m)%M;
    cout<<(ll)ans;
}

但是很不幸，由于防止过程中乘法溢出，这里使用了模数乘法，导致最后的时间复杂度变成 $O(k \log l \cdot \log (m\cdot M))$ ，显然会T。

如果预处理 $10^i$ 模 $m\cdot M$ 的各项/和，但由于范围过大防止溢出，还是只能通过模数乘法来计算，时间复杂度仍然不变。

拆模

引理：

今有一长度为 $l$ 的十进制数 $n$ ，

那么其 $\lfloor n/m \rfloor \mod M$ 的结果可以表示为：

我们记数字的前i位（高i位）为 $ni:=\sum{j=l-i-1}^{l-1}d_j10^j$ ，那么有：

令 $n0=d{l-1}$ 不难发现有递推公式：

这样就可以递推处理商了，而且很好地避免了逆元等情况。

为了处理乘法范围过大，同时复杂度又不能接受，由上面的引理我们可以想到维护 $q\ ,\ r$ 模M，通过之前的递推公式，我们知道从高位到低位读入时有：

那么其商和余数的递推公式可以写成：

这里使用预处理 $10^i$ 分别模M和m的个项和合，然后把长度 $l$ 用二进制拆成若干个 $2^i$ 的和，这样就可以在 $O(k \log l)$ 的时间复杂度内求出答案。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
#include <bits/stdc++.h>
using namespace std;

using i128 = __int128;
using ll = long long;

const int LOG=60;
vector<ll> p10m(LOG), p10M(LOG);
vector<ll> s10m(LOG), s10M(LOG);

void init(ll m, ll M)
{
    p10m[0]=10%m,p10M[0]=10%M;
    s10m[0]=1%m,s10M[0]=1%M;

    for(int i=1;i<LOG;i++)
    {
        p10m[i]=(i128)p10m[i-1]*p10m[i-1]%m;
        p10M[i]=(i128)p10M[i-1]*p10M[i-1]%M;
        s10m[i]=((i128)s10m[i-1]*p10m[i-1]+s10m[i-1])%m;
        s10M[i]=((i128)s10M[i-1]*p10M[i-1]+s10M[i-1])%M;
    }
}

pair<i128,i128> calc(i128 l,i128 mod,
                     const vector<ll>& p10,
                     const vector<ll>& s10)
{
    i128 p=1,s=0;
    for(int i=0;l;i++,l>>=1)
        if(l&1)
        {
            p=p*p10[i]%mod;
            s=(s*p10[i]+s10[i])%mod;
        }
    return {p,s};
}

int main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    int k;ll m,M;
    cin >> k >> m >> M;
    init(m,M);  
    i128 r=0, q=0;
    for(int i=0;i<k;i++)
    {
        ll c,l; 
        cin >> c >> l;
        auto [p,s]=calc(l,m,p10m,s10m);
        auto [P,S]=calc(l,M,p10M,s10M);  
        i128 t=(i128)p*r+(i128)c*s;
        r=t%m;
        q=((i128)P*q+t/m) % M;
    }
    cout << (ll)q;
}

但是这样写WA了，原因是对于 $q{i+1}$ 的向下取整的分子应该是除了 $r_i$ 剩下都不应该取模，以上程序计算 $10^{l_i} r_i+c_i\sum{j=0}^{l_i-1}10^j$ 时，对于10的次幂以及和都做了取模处理，这显然是错误的。

如果使用模数为 $mM$ ，同样地，中间过程计算的积肯定会溢出，此时又用懒乘法，导致时间复杂度又变回去了。

下面给出正确至少可以ac的代码 虽然是pypy版本 ，直接完整地算出其分子：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
import sys

def sol():
    it=iter(map(int, sys.stdin.read().split()))
    try:
        K=next(it); m=next(it); M=next(it)
    except StopIteration: return

    mM = m * M
    P, S = [10%mM], [1%mM]
    for i in range(62):
        S.append((S[i]*P[i]+S[i])%mM)
        P.append((P[i]*P[i])%mM)

    r, q = 0, 0
    for _ in range(K):
        c, l = next(it), next(it)
        p, s, i = 1, 0, 0
        l = l
        while l:
            if l&1:
                s=(s*P[i]+S[i])%mM
                p=(p*P[i])%mM
            l>>=1
            i+=1
        t = r*p+c*s
        q = (q*(p%M)+(t%mM//m))%M # 这里利用之前的结论扩大模数mM
        r = t%m
        
    sys.stdout.write(str(q))

if __name__ == "__main__":
    sol()

得益于pypy的JIT，这个代码既可以简洁地保持高精度，又可以在2.5s内跑完。

拆彻底