「刷题」一个人的数论

一个人的数论,这题也昨晚好久了,是另外一道神题。

#include<cstdio>
#include<iostream>
using namespace std;
typedef long long ll;
const ll mod=1e9+7,maxn=1005,maxd=105;
ll d,w,p[maxn],pk[maxn],z[maxn],a[maxd][maxd],b[maxd],x[maxd];
ll qw(ll a,ll b)
{
    ll ans=1;
    for(;b;b>>=1,a=a*a%mod) if(b&1) ans=ans*a%mod;
    return ans;
}
ll abs(ll a)
{
    return a>0?a:-a;
}
void Guass()
{
    for(int i=1;i<=d+1;i++)
    {
        for(int j=1;j<=d+1;j++) a[i][j]=qw(i,j);
        for(int j=1;j<=i;j++) b[i]=(b[i]+qw(j,d))%mod;
    }
    for(int i=1;i<=d+1;i++)
    {
        int q=i;
        for(int j=i+1;j<=d+1;j++)
            if(abs(a[j][i])>abs(a[q][i]))
                q=j;
        for(int j=1;j<=d+1;j++) swap(a[i][j],a[q][j]);
        swap(b[i],b[q]);
        if(a[i][i]==0) continue;
        ll inv=qw(a[i][i],mod-2);
        for(int j=1;j<=d+1;j++)
            if(i!=j)
            {
                ll r=a[j][i]*inv%mod;
                for(int k=1;k<=d+1;k++) a[j][k]=((a[j][k]-a[i][k]*r%mod)%mod+mod)%mod;
                b[j]=((b[j]-1LL*b[i]*r%mod)%mod+mod)%mod;
            }
    }
    for(int i=d+1;i>=1;i--)
    {
        for(int j=i+1;j<=d+1;j++) b[i]=((b[i]-x[j]*a[i][j]%mod)%mod+mod)%mod;
        ll inv=qw(a[i][i],mod-2);
        x[i]=b[i]*inv%mod;
    }
}
ll donit()
{
    ll ans=0;
    for(int i=1;i<=d+1;i++)
    {
        ll res=x[i];
        for(int j=1;j<=w;j++) 
        {
            ll k1=pk[j]*i%(mod-1),k2=(k1+d-i)%(mod-1);
            ll t=((qw(p[j],k1)-qw(p[j],k2))%mod+mod)%mod;
            res=res*t%mod;
        }
        ans=(ans+res)%mod;
    }
    return ans;
}
int main()
{
    scanf("%lld%lld",&d,&w);
    Guass();
    for(int i=1;i<=w;i++) scanf("%lld%lld",&p[i],&pk[i]);
    printf("%lld
",donit());
    return 0;
}
View Code
原文地址:https://www.cnblogs.com/Lrefrain/p/11222774.html