题解 P5130 【纯粹的弹幕地狱】

· · 题解

随手水~

首先要求出每一次子弹命中的概率。

设子弹起始位置为(x_1,y_1)，目标位置为(x_2,y_2)，不难发现子弹能造成伤害当且仅当\gcd(|x_1-x_2|,|y_1-y_2|)=1。

所以考虑枚举i=|x_1-x_2|和j=|y_1-y_2|。注意到i<j和i>j的情况实际对称，所以只考虑其中一边，最后乘2即可。为了规避i=1时特殊的情况，这里的i从2开始计算，而特殊情况会稍后计算。

确定了i和j之后，可以推出共有4(n-i+1)(n-j+1)种可能的位置（乘以4是因为目标可能在起始点的右上、右下、左下、左上四种方向），所以得到：

2\times 4\sum_{i=2}^n\sum_{j=0}^{i-1}(n-i+1)(n-j+1)[\gcd(i,j)=1]

这个式子肯定可以莫比乌斯反演无疑了；但是这里提供一种其他做法。

2\times 4\sum_{i=2}^n\sum_{j=0}^{i-1}(n-i+1)(n-j+1)[\gcd(i,j)=1] =8\sum_{i=2}^n(n-i+1)\sum_{j=0}^{i-1}(n-j+1)[\gcd(i,j)=1] =8\sum_{i=2}^n(n-i+1)\left[(n+1)\varphi(i)-\sum_{j=0}^{i-1}j[\gcd(i,j)=1]\right]

后面那个求和到底是什么东西？

注意到一个性质：对于i>1，\gcd(j,i)=1 \Longleftrightarrow \gcd(i-j,i)=1，也就是说所有i以内与i互质的自然数可以配对出现，每一对的和为i。考虑到一共有\varphi(i)个与i互质的数，也就是一共会出现\frac{\varphi(i)}{2}对，那么有：

\sum_{j=0}^{i-1}j[\gcd(i,j)=1]=\frac{i\varphi(i)}{2}

于是所求式为：

8\sum_{i=2}^n(n-i+1)\left[(n+1)\varphi(i)-\frac{i\varphi(i)}{2}\right]

直接O(n)求就可以了。

还剩下i=1的情况没有考虑。这时候有且仅有j=0和j=1是合法的，分别有4n^2和4n(n+1)种位置，直接加上即可。

上面的和是合法的方案数，除以总共的方案数就是概率，最后乘以\frac{ak(k+1)(2k+1)}{6}+\frac{bk(k+1)}{2}+ck就得到了总伤害的期望。

代码：

#include <iostream>
#include <cstdlib>
#include <cstdio>
#include <cstring>
#define mod 19260817

using namespace std;
const int MAXN = 10000005;
typedef long long ll;
ll n,k,a,b,c;
int prime[750005],tot = 0;
ll phi[MAXN];
bool vis[MAXN];
inline void make_table(int n)
{
    phi[1] = 1;
    for(int i = 2; i<=n; i++)
    {
        if(!vis[i])
        {
            prime[++tot] = i;
            phi[i] = i-1;
        }
        for(int j = 1; j<=tot&&i*prime[j]<=n; j++)
        {
            vis[i*prime[j]] = true;
            if(i%prime[j]==0)
            {
                phi[i*prime[j]] = phi[i]*prime[j];
                break;
            }
            else phi[i*prime[j]] = phi[i]*(prime[j]-1);
        }
    }
} 
inline ll qpow(ll a, ll b)
{
    ll res = 1;
    while(b)
    {
        if(b&1) res = res*a%mod;
        a = a*a%mod;
        b >>= 1;
    }
    return res;
}
inline ll inv(ll x)
{
    return qpow(x,mod-2);
}

int main()
{
    cin >> n >> a >> b >> c >> k;
    make_table(n+5);
    ll cnt = 0;
    ll inv2 = inv(2);
    for(int i = 2; i<=n; i++)
        cnt = (cnt+((n+1)*phi[i]%mod-(phi[i]*i%mod*inv2)%mod+mod)%mod*(n-i+1)%mod*8%mod)%mod;
    cnt = (cnt+n*n%mod*4%mod+4*n%mod*(n+1)%mod)%mod;
    ll totnum = (n+1)*(n+1)%mod*((n+1)*(n+1)%mod)%mod;
    ll prob = cnt*inv(totnum)%mod;
    ll dmg = (a*(k*(k+1)%mod*(2*k+1)%mod*inv(6)%mod)%mod+
             b*(k*(k+1)%mod*inv(2)%mod)%mod+c*k%mod)%mod;
    cout << prob*dmg%mod << endl;
    return 0;
}

这个方法能否拓展到低于线性的时间尚未可知，我暂时没有好的想法。不过莫比乌斯反演倒是肯定行。