zoj3707 斐波那契素数

题意：求第k个斐波那契素数对应的斐波那契数列及其之后的值%x==0的最小的斐波那契数值，(%m),
数据范围：(1<=k<=10^6),(3<=X<=100),(10<=M<=10^6)
解题思路：
斐波那契数列的特性：( 0,1,1,2,3,5,8... )
1、gcd(fib(n),fib(m))=fib(gcd(n,m))
2、如果fib(k)能被x整除，则fib(k*i)都可以被x整除。
3、f(0)+f(1)+f(2)+…+f(n)=f(n+2)-1
4、f(1)+f(3)+f(5)+…+f(2n-1)=f(2n)
5、f(2)+f(4)+f(6)+…+f(2n) =f(2n+1)-1
6、[f(0)]^2+[f(1)]^2+…+[f(n)]^2=f(n)·f(n+1)
7、f(0)-f(1)+f(2)-…+(-1)^n·f(n)=(-1)^n·[f(n+1)-f(n)]+1
8、f(m+n)=f(m-1)·f(n-1)+f(m)·f(n)
9、[f(n)]^2=(-1)^(n-1)+f(n-1)·f(n+1)
10、f(2n-1)=[f(n)]^2-[f(n-2)]^2
11、3f(n)=f(n+2)+f(n-2)
12、f(2n-2m-2)[f(2n)+f(2n+2)]=f(2m+2)+f(4n-2m) [ n〉m≥-1,且n≥1]
还有一个结论：
计算(a/b)%c  其中b能整除a
如果b与c互素，则(a/b)%c=a*b^(phi(c)-1)%c
如果b与c不互素，则(a/b)%c=(a%bc)/b
对于b与c互素和不互素都有(a/b)%c=(a%bc)/b成立


斐波那契素数的定义：( 2,3,5,13,89... )
对于斐波那契数列的第i个元素s[i]，如果gcd(s[i],s[j])==1 (1<=j<i)，那么s[i]是一个斐波那契素数。


斐波那契素数的特性(可根据斐波那契数列求出斐波那契素数)：
对于斐波那契数列(0,1,1,2,3,5,8...)，下标从0开始，下标为3,4,5,7,11,13...的数就是斐波那契素数。
即除了前两个外，斐波那契数列的素数下标对应的就是斐波那契素数，而如果下标不是素数，那么对应的也不是斐波那契素数。


所以这一题就先筛出斐波那契素数下标(3,4,5,7,11,13...)，然后就可以知道第k个斐波那契素数对应的斐波那契数列的下标，

根据矩阵快速幂计算出斐波那契数值，判断是否%x==0，如果不是就继续找下一个，直到找到为止。最后a/x%m=(a%(x*m))/x 

参考资料：http://blog.csdn.net/acdreamers/article/details/8989772

#include <stdio.h>
#include <iostream>
#include <string.h>
#include <math.h>
using namespace std;
typedef long long LL;
const LL N=16000000;
LL p[N];
bool prime[N];
LL k=1;
void isprime()
{
    LL i,j;
    p[0]=1;
    memset(prime,true,sizeof(prime));
    for(i=2;i<N;i++)
    {
        if(prime[i])
        {
            p[k++]=i;
            for(j=i+i;j<N;j+=i)
            {
                prime[j]=false;
            }
        }
    }
    p[1]=3;p[2]=4;
}
typedef struct
{
    LL m[2][2];
}Matrix;
Matrix per={1,0,0,1};//单位矩阵
Matrix a={1,1,1,0};
Matrix multi(Matrix a,Matrix b,LL MOD)
{
    Matrix c;
    LL i,j;
    for(i=0;i<2;i++)
    {
        for(j=0;j<2;j++)
        {
            c.m[i][j]=0;
            for(k=0;k<2;k++)
            {
                c.m[i][j]+=a.m[i][k]*b.m[k][j];
                c.m[i][j]%=MOD;
            }
        }
    }
    return c;
}
Matrix matrix_mod(LL k,LL MOD)//矩阵快速幂求 p矩阵的k次方
{
    Matrix p=a,ans=per;
    while(k)
    {
        if(k&1)
        {
            ans=multi(ans,p,MOD);
            k--;
        }
        k>>=1;
        p=multi(p,p,MOD);
    }
    return ans;
}
int main()
{
    LL K,X,M,t,i,ret,r;
    isprime();
    scanf("%lld",&t);
    while(t--)
    {
        scanf("%lld%lld%lld",&K,&X,&M);
        Matrix ans;

        for(i=p[K];;i++)
        {
            ans=matrix_mod(i-1,X);//求出的ans.m[0][0]就是fib[i]
            if((ans.m[0][0]%X==0))
            {
                r=i;
                break;
            }
        }
        ans=matrix_mod(r-1,M*X);
        ret=ans.m[0][0]/X;
        printf("%lld\n",ret);
    }
}
/*
input:
1
1 3 10
output:
1
*/