hdu 3978 Evil teacher's Final Problem (斐波那契数列的循环节)

  1/* ***********************************************
  2Author :111qqz
  3Created Time :Tue 01 Nov 2016 08:31:22 PM CST
  4File Name :code/hdu/3978.cpp
  5************************************************ */
  6#include <cstdio>
  7#include <cstring>
  8#include <iostream>
  9#include <algorithm>
 10#include <vector>
 11#include <queue>
 12#include <set>
 13#include <map>
 14#include <string>
 15#include <cmath>
 16#include <cstdlib>
 17#include <ctime>
 18#define fst first
 19#define sec second
 20#define lson l,m,rt<<1
 21#define rson m+1,r,rt<<1|1
 22#define ms(a,x) memset(a,x,sizeof(a))
 23typedef long long LL;
 24#define pi pair < int ,int >
 25#define MP make_pair
 26using namespace std;
 27const double eps = 1E-8;
 28const int dx4[4]={1,0,0,-1};
 29const int dy4[4]={0,-1,1,0};
 30const int inf = 0x3f3f3f3f;
 31struct Mat
 32{
 33    LL mat[2][2];
 34    void clear()
 35    {
 36        ms(mat,0);
 37    }
 38}M,M1;
 39Mat mul (Mat a,Mat b,LL mod)
 40{
 41    Mat c;
 42    c.clear();
 43    for ( int i = 0 ; i < 2 ; i++)
 44        for ( int j = 0 ; j < 2 ; j++)
 45            for ( int k  = 0 ; k < 2 ; k++)
 46                c.mat[i][j] = (c.mat[i][j] + a.mat[i][k] * b.mat[k][j]%mod)%mod;
 47    return c;
 48}
 49Mat mat_ksm(Mat a,LL b,LL mod)
 50{
 51    Mat res;
 52    res.clear();
 53    for ( int i = 0 ; i < 2 ; i++) res.mat[i][i] = 1;
 54    while (b>0)
 55    {
 56        if (b&1) res = mul(res,a,mod);
 57        b = b >> 1LL;
 58        a = mul(a,a,mod);
 59    }
 60    return res;
 61}
 62LL gcd(LL a,LL b)
 63{
 64    return b?gcd(b,a%b):a;
 65}
 66const int N = 1E6+7;
 67bool prime[N];
 68int p[N];
 69void isprime() //一个普通的筛
 70{
 71    int cnt = 0 ;
 72    ms(prime,true);
 73    for ( int i = 2 ; i < N ; i++)
 74    {
 75        if (prime[i])
 76        {
 77            p[cnt++] = i ;
 78            for ( int j = i+i ; j < N ; j+=i) prime[j] = false;
 79        }
 80    }
 81}
 82LL ksm( LL a,LL b,LL mod)
 83{
 84   LL res = 1;
 85   while (b>0)
 86   {
 87       if (b&1) res = (res * a) % mod;
 88       b = b >> 1LL;
 89       a = a * a % mod;
 90   }
 91   return res;
 92}
 93LL legendre(LL a,LL p) //勒让德符号,判断二次剩余
 94{
 95    if (ksm(a,(p-1)>>1,p)==1) return 1;
 96    return -1;
 97}
 98LL pri[N],num[N];//分解质因数的底数和指数。
 99int cnt; //质因子的个数
100void solve(LL n,LL pri[],LL num[])
101{
102    cnt = 0 ;
103    for ( int  i = 0 ; p[i] * p[i] <= n ; i++)
104    {
105        if (n%p[i]==0)
106        {
107            int Num = 0 ;
108            pri[cnt] = p[i];
109            while (n%p[i]==0)
110            {
111                Num++;
112                n/=p[i];
113            }
114            num[cnt] = Num;
115            cnt++;
116        }
117    }
118    if (n>1)
119    {
120        pri[cnt] = n;
121        num[cnt] = 1;
122        cnt++;
123    }
124}
125LL fac[N];
126int cnt2; //n的因子的个数
127void get_fac(LL n)//得到n的所有因子
128{
129    cnt2 = 0 ;
130    for (int i =  1 ; i*i <= n ; i++)
131    {
132        if (n%i==0)
133        {
134            if (i*i!=n)
135            {
136                fac[cnt2++] = i ;
137                fac[cnt2++] = n/i;
138            }
139            else fac[cnt2++] = i;
140        }
141    }
142}
143LL find_loop(LL n)
144{
145    solve(n,pri,num);
146    LL ans = 1;
147    for ( int i = 0 ; i < cnt ; i++)
148    {
149        LL rec = 1;
150        if (pri[i]==2) rec = 3;
151        else if (pri[i]==3) rec = 8;
152        else if (pri[i]==5) rec = 20;
153        else
154        {
155            if (legendre(5,pri[i])==1)
156                get_fac(pri[i]-1);
157            else get_fac(2*pri[i]+2);
158            sort(fac,fac+cnt2);
159            for ( int j = 0 ; j < cnt2 ; j++) //挨个验证因子
160            {
161                Mat tmp = mat_ksm(M,fac[j],pri[i]); //下标从0开始，验证fac[j]为循环节，应该看fib[0]==fib[fac[j]]和fib[1]==fib[fac[j]+1]是否成立
162                tmp = mul(tmp,M1,pri[i]);
163                if (tmp.mat[0][0]==1&&tmp.mat[1][0]==1)
164                {
165                    rec = fac[j];
166                    break;
167                }
168            }

 1        }
 2        for ( int j = 1 ; j < num[i] ; j++)
 3            rec *=pri[i];
 4        ans = ans / gcd(ans,rec) * rec;
 5    }
 6    return ans;
 7}
 8void init()
 9{
10    M.clear();
11    M.mat[0][1] = M.mat[1][0] = M.mat[1][1] = 1;
12    M1.clear();
13    M1.mat[0][0] = M1.mat[1][0] = 1;
14}
15LL get_fib(LL n,LL mod)
16{
17    if (mod==1) return 0;
18    Mat ret;
19    ret = mat_ksm(M,n-1,mod);
20    ret = mul(ret,M1,mod);
21    return (ret.mat[1][0])%mod;
22}
23LL n,P;
24int k;feibonaqi
25map< LL,LL >mp;

 1LL  loop[N];
 2int main()
 3{
 4#ifndef  ONLINE_JUDGE 
 5    freopen("code/in.txt","r",stdin);
 6#endif
 7    isprime();
 8    int T;
 9    cin>>T;
10    init();
11    int cas = 0 ;
12    while (T--)
13    {
14	scanf("%lld%d%lld",&n,&k,&P);
15	printf("Case #%d: ",++cas);
16	if (n==0||n==1)
17	{
18	    printf("1\n");
19	    continue;
20	}
21	loop[0] = P;
22	LL cur = n;
23	for (int i = 1; i <= k ; i++) loop[i] = find_loop(loop[i-1]);
24	for (int i = k ; i >= 0 ; i--) cur = get_fib(cur,loop[i]);
25	printf("%lld\n",cur);
26    }
27#ifndef ONLINE_JUDGE  
28    fclose(stdin);
29#endif
30    return 0;
31}