zoukankan      html  css  js  c++  java
  • 51nod1229-序列求和V2【数学,拉格朗日插值】

    正题

    题目链接:http://www.51nod.com/Challenge/Problem.html#problemId=1229


    题目大意

    给出(n,k,r)

    [sum_{i=1}^ni^kr^i ]

    (1leq Tleq 20,1leq n,rleq 10^{18},1leq kleq 2000)


    解题思路

    如此明显的式子直接开推

    [S_k=sum_{i=1}^ni^kr^i,rS_k=sum_{i=2}^{n+1}(i-1)^kr^i ]

    [(r-1)S_k=n^kr^{n+1}-r+sum_{i=2}^nleft((i-1)^k-i^k ight)r^i ]

    二项式展开((i-1)^k)

    [(r-1)S_k=n^kr^{n+1}-r+sum_{i=2}^nsum_{j=0}^{k-1}(-1)^{k-j}inom{k}{j}r^i ]

    然后把(j)提到前面去

    [(r-1)S_k=n^kr^{n+1}-r+sum_{j=0}^{k-1}(-1)^{i-j}inom{k}{j}sum_{i=2}^nr^i ]

    [Rightarrow S_k=frac{n^kr^{n+1}-r+sum_{j=0}^{k-1}(-1)^{k-j}inom{k}{j}(S_j-r)}{r-1} ]

    这样(S_k)就可以(O(k^2))递推了。

    当然当(r=1)的时候,不能再使用这个公式,此时(sum_{i=1}^ni^k)是很经典的问题,直接拉格朗日插值插出一个(k+1)次多项式即可。

    此题到这里就圆满结束了,但是以直觉判断上面那个式子可以卷积,拆开组合数然后疯狂跳步一下就是

    [(r-1)frac{S_k-r}{k!}=n^kr_{n+1}-r-(r-1)frac{r}{k!}+sum_{j=0}^{k-1}frac{(-1)^{k-j}}{(k-j)!}frac{S_j-r}{j!} ]

    [H(x)=sum_{i=0}^infty (n^ir_{n+1}-r-(r-1)frac{r}{i!})x^i ]

    [G(x)=sum_{i=0}^infty frac{S_i-r}{i!},F(x)=sum_{i=1}^{infty}frac{(-1)^i}{i!} ]

    那么有

    [(r-1)G(x)=H(x)+F(x)G(x) ]

    [Rightarrow G(x)=frac{H(x)}{r-1-F(x)} ]

    然后多项式求逆即可,时间复杂度(O(klog k)),虽然这题的模数不能用,但是可以顺便解决掉序列求和V5

    但是最近写的多项式求逆有点多,咕了,所以上面的式子如果有错我也没办法(((、


    code

    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #define ll long long
    using namespace std;
    const ll N=2100,P=1e9+7;
    ll T,n,k,r,inv[N],fac[N],pre[N],suf[N],s[N];
    ll power(ll x,ll b){
    	ll ans=1;b%=P-1;x%=P;
    	while(b){
    		if(b&1)ans=ans*x%P;
    		x=x*x%P;b>>=1;
    	}
    	return ans;
    }
    ll C(ll n,ll m)
    {return fac[n]*inv[m]%P*inv[n-m]%P;}
    signed main()
    {
    	scanf("%lld",&T);
    	inv[0]=fac[0]=inv[1]=1;
    	for(ll i=2;i<N;i++)inv[i]=P-inv[P%i]*(P/i)%P;
    	for(ll i=1;i<N;i++)fac[i]=fac[i-1]*i%P,inv[i]=inv[i-1]*inv[i]%P;
    	while(T--){
    		scanf("%lld%lld%lld",&n,&k,&r);r%=P;
    		if(r==1){
    			ll ans=0;k+=2;n%=P;
    			pre[0]=suf[k+1]=1;
    			for(ll i=1;i<=k;i++)pre[i]=pre[i-1]*(n-i)%P;
    			for(ll i=k;i>=1;i--)suf[i]=suf[i+1]*(n-i)%P;
    			for(ll i=1,p=0;i<=k;i++){
    				p=(p+power(i,k-2))%P;
    				(ans+=p*pre[i-1]%P*suf[i+1]%P*inv[i-1]%P*(((k-i)&1)?P-inv[k-i]:inv[k-i])%P)%=P;
    			}
    			printf("%lld
    ",(ans+P)%P);
    		}
    		else{
    			ll z=power(r,n+1),invr=power(r-1,P-2);
    			s[0]=(z-r+P)*invr%P;n%=P;
    			for(ll i=1,pw=n;i<=k;i++,pw=pw*n%P){
    				s[i]=(z*pw-r+P)%P;s[i-1]=(s[i-1]-r+P)%P;
    				for(ll j=0;j<i;j++)
    					(s[i]+=(((i-j)&1)?(P-1):(1))*s[j]%P*C(i,j)%P)%=P;
    				s[i]=s[i]*invr%P;
    			}
    			printf("%lld
    ",(s[k]+P)%P);
    		}
    	}
    	return 0;
    }
    
  • 相关阅读:
    Spring Boot之@ImportResource、配置类、占位符表达式
    Spring Boot之测试类报错问题
    Spring Boot之@EnableAutoConfiguration源码解读
    Spring Boot之第一个程序及执行原理。
    eclipse中git使用中的冲突解决
    python画国旗
    第十六周进度
    个人课程总结
    人月神话之阅读笔记03
    第十五周进度
  • 原文地址:https://www.cnblogs.com/QuantAsk/p/15339896.html
Copyright © 2011-2022 走看看