[bzoj 4833]最小公倍佩尔数

## [ 传送门 ]()

Description 

Let ((1+sqrt2)^n=e(n)+f(n)cdotsqrt2) , both (e(n)) and (f(n)) are integers

Let (g(n)) be the gcd of (f(1),f(2),...,f(n))

given (n), (p), where (p) is a prime number

Calculate the value of

[ sum_{i=1}^{n}icdot g(i)     mod  p ]

(Tleq 210 ,sum nleq 3×10^6)

Solution 

[f(n)=2f(n-1)+f(n-2),f(0)=1,f(1)=1 ]

Similar to the (Fibonacci) sequence, we have

[ gcd(f(a),f(b))=f(gcd(a,b)) ]

It's hard to evaluate LCM directly,but we can get it by maximum inversion

[ lcm(S)=prod_{T⊆S,T≠∅}gcd(T)^{(−1)^{|T|−1}} ]

so we can find that

[ g(n)=prod_{T⊆S,T≠∅}f(gcd(T))^{(−1)^{|T|−1}} ]

The next step is the most important.

construct a function (h) ,which satisfies

[ f(n)=prod_{d|n}h(d) ]

we can calculate (h(1...n)) easily by (O(nlog n))

and

[ egin{equation}  egin{split}  g(n)&=prod_{d=1}^n h(d)^{∑_{T⊆S,T≠∅,d|gcd(T)}(−1)^{|T|+1}}\  &=prod_{d=1}^nh(d)  end{split}  end{equation} ]

Code 

#include<bits/stdc++.h>
#define ll long long
#define reg register
#define db double
using namespace std;
inline int read()
{
	int x=0,f=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x*f;
}
const int MN=1e6+5;
int f[MN],h[MN],n,P,ans;
int Mul(int x,int y){return 1ll*x*y%P;}
int Add(int x,int y){return (x+y)%P;}
int fpow(int x){int r=1,m=P-2;for(;m;m>>=1,x=Mul(x,x))if(m&1)r=Mul(r,x);return r;} 
int main()
{
	int T=read();
	while(T--)
	{
		n=read(),P=read();
		reg int i;
		f[0]=0;h[0]=h[1]=f[1]=1;
		for(i=2;i<=n;++i) h[i]=1,f[i]=Add(Mul(f[i-1],2),f[i-2]);
		for(i=2;i<=n;++i)
		{
			h[i]=Mul(f[i],fpow(h[i]));
			for(int j=i<<1;j<=n;j+=i) h[j]=Mul(h[j],h[i]);
		}
		for(ans=0,i=1;i<=n;++i) h[i]=Mul(h[i],h[i-1]),ans=Add(ans,Mul(h[i],i));
		printf("%d
",ans);
	}
	return 0;
}

---

Blog来自PaperCloud，未经允许，请勿转载，TKS！