BZOJ 2813: 奇妙的Fibonacci

2813: 奇妙的Fibonacci

Time Limit: 20 Sec  Memory Limit: 512 MB
Submit: 497  Solved: 134
[Submit][Status][Discuss]

Description

Fibonacci数列是这样一个数列：

F1 = 1， F2 = 1， F3 = 2 . . .

Fi = Fi-1 + Fi-2 (当 i >= 3)

pty忽然对这个古老的数列产生了浓厚的兴趣，他想知道：对于某一个Fibonacci数Fi，

有多少个Fj能够整除Fi (i可以等于j)，他还想知道所有j的平方之和是多少。

Input

第一行一个整数Q,表示Q个询问。

第二行四个整数:Q₁, A, B, C

第i个询问Q_i = (Q_i-1 * A + B) mod C + 1(当i >= 2)

Output

Ai代表第i个询问有多少个F_j能够整除F_Qi。

Bi代表第i个询问所有j的平方之和。

输出包括两行：

第一行是所有的Ai之和。

第二行是所有的Bi之和。

由于答案过大，只需要输出除以1000000007得到的余数即可。

Sample Input

2
2 2 1 8

Sample Output

6
55

HINT

对于100%的数据保证：Q <= 3*10^6，C <= 10^7，A <= 10^7，B <= 10^7，1 <= Q1<= C

Source

分析：

打表可以发现答案就是约数个数以及约数平方和...

如果要严谨的证明就是这个定理：$gcd(f[i],f[j])=f[gcd(i,j)]$

证明如下：

定理1：$gcd(f[n],f[n+1])=1$

$gcd(f[n],f[n+1])$

$=gcd(f[n+1]-f[n],f[n])$

$=gcd(f[n],f[n-1])$

$……$

$=gcd(f[2],f[1])$

$=1$

定理2：$f[m+n]=f[m-1]f[n]+f[m]f[n+1]$

$f[m+n]$

$=f[n+m-1]+f[n+m-2]$

$=2*f[n+m-2]+f[n+m-3]$

$=a[x]*f[n+m-x]+b[x]*f[n+m-x-1]$

$=(a[x]+b[x])*f[n+m-x-1]+a[x]*f[n+m-x-2]$

如果$x=1$，则$a[x]=f[2],b[x]=f[1]$

如果$x=2$，则$a[x]=f[3],b[x]=f[2]$

如果$x=n$，则$a[x]=f[n+1],b[x]=f[n]$

所以$f[m+n]=f[m]*f[n+1]+f[n]*f[m-1]$

定理3：$gcd(f[m+n],f[n])=gcd(f[m],f[n])$

$gcd(f[m+n],f[n])$

$=gcd(f[m]*f[n+1]+f[n]*f[m-1],f[n])$

$=gcd(f[m]*f[n+1],f[n])$

$=gcd(f[n+1],f[n])*gcd(f[m],f[n])$

$=gcd(f[m],f[n])$

然后根据辗转相减法：

设$m=p _{1}*n+r _{1},n=p _{2}*r _{1}+r _{2},r_{1}=p _{3}*r _{2}+r _{3}......$

$gcd(m,n)=gcd(n,r _{1})=......=r _{x}$

$gcd(f[m],f[n])=gcd(f[n],f[r _{1}])=......f[r _{x}]=f[gcd(m,n)]$

然后约数和还有平方和线性筛的时候维护一下就好了...

因为f[2]=1，所以如果q是奇数约数和要+1，平方和要+4...

代码：

#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
//by NeighThorn
using namespace std;

const int maxn=10000000+5,mod=1e9+7;

int a,b,c,q,q1,tot,ans1,ans2,vis[maxn],pri[maxn],cnt[maxn],Min[maxn],sum[maxn],els[maxn];

inline void prework(void){
	cnt[1]=sum[1]=1;
	for(int i=2;i<=10000000;i++){
		if(!vis[i])
			pri[++tot]=i,Min[i]=els[i]=1,cnt[i]=2,sum[i]=1LL*i*i%mod+1;
		for(int j=1;j<=tot&&i*pri[j]<=10000000;j++){
			vis[i*pri[j]]=1;
			if(i%pri[j]==0){
				Min[i*pri[j]]=Min[i]+1;els[i*pri[j]]=els[i];
				cnt[i*pri[j]]=(cnt[i]/(Min[i]+1))*(Min[i*pri[j]]+1);
				sum[i*pri[j]]=(1LL*sum[i]*pri[j]%mod*pri[j]%mod+sum[els[i]])%mod;
				break;
			}
			cnt[i*pri[j]]=cnt[i]+1,Min[i*pri[j]]=1;
			els[i*pri[j]]=i;cnt[i*pri[j]]=cnt[i]<<1;
			sum[i*pri[j]]=(sum[i]+1LL*sum[i]*pri[j]%mod*pri[j]%mod)%mod;
		}
	}
}

signed main(void){
	ans1=ans2=0;prework();
	scanf("%d%d%d%d%d",&q,&q1,&a,&b,&c);a%=c,b%=c;
	for(int i=1;i<=q;i++){
		if(i>1)
			q1=(1LL*q1*a+b)%c+1;
		(ans1+=cnt[q1]+(q1&1))%mod,(ans2+=sum[q1]+(q1&1)*4)%=mod;
	}
	printf("%d\n%d\n",ans1,ans2);
	return 0;
}