HDU3221-Brute—force Algorithm（欧拉函数）

Brute-force Algorithm

                                                                        Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)
                                                                                                   Total Submission(s): 2762    Accepted Submission(s): 737

Problem Description

Professor Brute is not good at algorithm design. Once he was asked to solve a path finding problem. He worked on it for several days and finally came up with the following algorithm:

Any fool but Brute knows that the function “funny” will be called too many times. Brute wants to investigate the number of times the function will be called, but he is too lazy to do it.

Now your task is to calculate how many times the function “funny” will be called, for the given a, b and n. Because the answer may be too large, you should output the answer module by P.

 

Input

There are multiple test cases. The first line of the input contains an integer T, meaning the number of the test cases.

For each test cases, there are four integers a, b, P and n in a single line.
You can assume that 1≤n≤1000000000, 1≤P≤1000000, 0≤a, b<1000000.

 

Output

For each test case, output the answer with case number in a single line.

 

Sample Input

  

   3
3 4 10 3
4 5 13 5
3 2 19 100
  

 

Sample Output

  

   Case #1: 2
Case #2: 11
Case #3: 12
  

 

Source

2009 Asia Shanghai Regional Contest Host by DHU

题意：f[n]=f[n-1]*f[n-2]，f[1]=a，f[2]=b，求出f[n]，并且答案对mod取余。
解题思路：可以发现a的指数和b的指数均类似于斐波那契数列。用矩阵的快速幂可以求出第n项a和b的指数分别是多少。但是这个指数会非常大，存不下来，需要对一个数去模。可以通过公式:A^n%mod=A^( n%Phi[mod] + Phi[mod] )%mod(n>=Phi[mod])。Phi[mod]表示不大于mod的数中与mod互质的数的个数，可以用欧拉函数来求：Phi[mod]=mod*(1-1/q1)*(1-1/q2)*(1-1/q3)*....*(1-1/qk)（q1,q2,q3...qk表示mod的素因子）。

#include <cstdio>
#include <cmath>
#include <algorithm>
#include <cstring>
#include <iostream>

using namespace std;

#define ll long long

ll a,b,mod,Phi;
int visit[1000090],k,prime[1000900];

struct Matric
{
    ll a[2][2];
    Matric()
    {
        memset(a,0,sizeof a);
    }
}x;

void init()
{
    memset(visit,0,sizeof visit);
    for(int i=2;i<=2000;i++)
    {
        if(!visit[i])
        {
            for(int j=i*2;j<=1000000;j+=i)
                visit[j]=1;
        }
    }
    k=0;
    for(int i=2;i<=1000000;i++)
        if(!visit[i])
            prime[k++]=i;
}

Matric mypow(Matric p,Matric q)
{
    Matric sum;
    for(int i=0;i<=1;i++)
    {
        for(int j=0;j<=1;j++)
        {
            for(int k=0;k<=1;k++)
            {
                sum.a[i][j]+=p.a[i][k]*q.a[k][j];
                if(sum.a[i][j]>Phi)
                    sum.a[i][j]=sum.a[i][j]%Phi+Phi;
            }
         }
     }
     return sum;
}

Matric getsum(Matric p,ll n)
{
     Matric sum;
     sum.a[0][0]=1;
     sum.a[0][1]=2;
     while(n)
     {
         if(n&1) sum=mypow(sum,p);
         p=mypow(p,p);
         n/=2;
     }
     return sum;
}

ll qpow(ll a,ll n)
{
    ll ans;
    ans=1;
    while(n)
    {
        if(n&1)
        {
            ans=ans*a;
            ans%=mod;
        }
        a=a*a;
        a%=mod;
        n/=2;
    }
    return ans;
 }

int main()
{
    init();
    int t,cas=0;
    ll n;
    scanf("%d",&t);
    while(t--)
    {
        scanf("%lld %lld %lld %lld",&a,&b,&mod,&n);
        printf("Case #%d: ",++cas);
        if(n==1) {printf("%lld\n",a%mod);continue;}
        else if(n==2) {printf("%lld\n",b%mod); continue;}
        else if(n==3) {printf("%lld\n",a*b%mod);continue;}
        else if(n==4) {printf("%lld\n",(a*b%mod)*b%mod);continue;}
        if(mod==1) {printf("0\n");continue;}
        Phi=1;
        ll num=mod;
        for(int i=0;i<k;i++)
        {
            if(prime[i]>mod) break;
            if(mod%prime[i]==0) {Phi*=(prime[i]-1);num/=prime[i];}
        }
        Phi*=num;
        x.a[0][1]=x.a[1][0]=x.a[1][1]=1;
        Matric ans=getsum(x,n-4);
        ll num1=ans.a[0][0];
        ll num2=ans.a[0][1];
        ll aa=qpow(a,num1);
        ll bb=qpow(b,num2);
        printf("%lld\n",aa*bb%mod);
    }
    return 0;
}