地下的太阳

题目描述

结论

显然对于e数组只需要知道最后一行，很容易得到结论：
我们需要知道前m-1行所有数的和，设为t，那么$e[m][i]=n*i^m+t$
如何计算t呢？
对于计算前i行的和，假设已经知道前i-1行的和为t，现在计算t’。
显然$t'=t*n+\sum_{j=1}^nj^i$
涉及到求自然数幂和，可以用第一类斯特林数的方法预处理，详见我blog中的“计算自然数幂和的各种方法”（我懒得去贴链接了）
于是有了最后一行，每次就可以得到A和B。

计算答案

$\sum_{i=1}^s\sum_{j=1}^sA[i]*B[j]*(|i-j|+1)^2$

其中
$(|i-j|+1)^2$
$=(i-j)^2+1+2|i-j|$
$=i^2+j^2-2ij+1+2|i-j|$
对于前四项分别计算，以+1为例
$\sum_{i=1}^s\sum_{j=1}^sA[i]*B[j]$
$=\sum_{i=1}^sA[i]*\sum_{j=1}^sB[j]$
这显然容易计算，另外三项同理
观察最后一项
讨论i>=j的情况，是这样的
$\sum_{i=1}^s\sum_{j=1}^iA[i]*B[j]*2(i-j)$
拆括号，变成$2*i-2*j$
以$2*i$为例
$=\sum_{i=1}^sA[i]*2i*\sum_{j=1}^iB[j]$
那我只要知道B的前缀和一样容易计算
$2*j$同理
$i<j$同理

#include<cstdio>
#include<algorithm>
#include<cmath>
#define fo(i,a,b) for(i=a;i<=b;i++)
using namespace std;
typedef long long ll;
const int maxm=5000+10,maxs=30000+10,mo=1004535809;
ll f[maxm],su[maxm][maxm],a[maxs],b[maxs],aa[maxs],bb[maxs],sa[maxs],sb[maxs],ssa[maxs],ssb[maxs];
ll i,j,k,s,t,m,ca,ans;
ll n,l,r;
ll quicksortmi(ll x,ll y){
    if (!y) return 1;
    ll t=quicksortmi(x,y/2);
    t=(ll)t*t%mo;
    if (y%2) t=(ll)t*x%mo;
    return t;
}
void work(ll n,ll m){
    f[0]=n%mo;
    f[1]=(ll)(n%mo)*((n+1)%mo)%mo*quicksortmi(2,mo-2)%mo;
    ll i,t;
    ll j;
    fo(i,2,m-1){
        f[i]=1;
        fo(j,n+1-(ll)i,n+1)
            if (j%((ll)i+1)==0) f[i]=(ll)f[i]*((j/((ll)i+1))%mo)%mo;else f[i]=(ll)f[i]*(j%mo)%mo;
        fo(j,0,i-1){
            t=(ll)su[i][j]*f[j]%mo;
            if (((ll)i+j)%2==1) t=mo-t;
            (t+=mo)%=mo;
            f[i]=(f[i]-t)%mo;
            (f[i]+=mo)%=mo;
        }
    }
}
ll F(ll x){
    return ((ll)quicksortmi((x%mo),m)*(n%mo)%mo+t)%mo;
}
int main(){
    scanf("%lld%lld",&n,&m);
    su[0][0]=1;
    fo(i,1,m) su[i][0]=0,su[i][i]=1;
    fo(i,1,m)
        fo(j,1,i-1)
            su[i][j]=(su[i-1][j-1]+(ll)(i-1)*su[i-1][j]%mo)%mo;
    work(n,m);
    t=f[1];
    fo(i,2,m-1){
        t=(ll)t*(n%mo)%mo;
        t=(t+(ll)(n%mo)*f[i]%mo)%mo;
    }
    //t=746187567;
    scanf("%lld",&ca);
    l=1;r=n;
    while (ca--){
        scanf("%lld",&s);
        fo(i,1,s) a[i]=F(l+(ll)i-1);
        l+=s;
        fo(i,1,s) b[i]=F(r-(ll)s+(ll)i);
        r-=s;
        ans=0;
        fo(i,1,s) sa[i]=(sa[i-1]+a[i])%mo;
        fo(i,1,s) sb[i]=(sb[i-1]+b[i])%mo;//sum a and b
        fo(i,1,s) ans=(ans+(ll)a[i]*i%mo*i%mo*sb[s]%mo)%mo;//calc a[i]*i*i*b[j]
        fo(i,1,s) ans=(ans+(ll)b[i]*i%mo*i%mo*sa[s]%mo)%mo;//calc a[i]*b[j]*j*j
        ans=(ans+(ll)sa[s]*sb[s]%mo)%mo;//calc a[i]*b[j]
        //i>=j
        fo(i,1,s) ans=(ans+(ll)a[i]*2%mo*i%mo*sb[i]%mo)%mo;//calc a[i]*2*i*b[j](j<=i)
        fo(i,1,s) bb[i]=(ll)b[i]*(-2)%mo*i%mo;
        fo(i,1,s) ssb[i]=(ssb[i-1]+bb[i])%mo;//sum b[i]*(-2)*i
        fo(i,1,s) ans=(ans+(ll)a[i]*ssb[i]%mo)%mo;//calc a[i]*b[j]*2*j(j<=i)
        //i<=j
        fo(i,1,s) ans=(ans+(ll)b[i]*2%mo*i%mo*sa[i]%mo)%mo;//calc a[i]*b[j]*2*j(j>=i)
        fo(i,1,s) aa[i]=(ll)a[i]*(-2)%mo*i%mo;
        fo(i,1,s) ssa[i]=(ssa[i-1]+aa[i])%mo;//sum a[i]*(-2)*i
        fo(i,1,s) ans=(ans+(ll)b[i]*ssa[i]%mo)%mo;//calc a[i]*2*i*b[j](j>=i)
        fo(i,1,s) bb[i]=(ll)b[i]*i%mo;
        fo(i,1,s) ssb[i]=(ssb[i-1]+bb[i])%mo;//sum b[i]*i
        fo(i,1,s) ans=(ans-(ll)2*a[i]%mo*i%mo*ssb[s]%mo)%mo;//calc -2*a[i]*i*b[j]*j
        (ans+=mo)%=mo;
        printf("%lld\n",ans);
    }
}