【bzoj2432】【NOI2011】兔农

题目描述

农夫栋栋近年收入不景气，正在他发愁如何能多赚点钱时，他听到隔壁的小 朋友在讨论兔子繁殖的问题。 问题是这样的：第一个月初有一对刚出生的小兔子，经过两个月长大后，这 对兔子从第三个月开始，每个月初生一对小兔子。新出生的小兔子生长两个月后 又能每个月生出一对小兔子。问第 n 个月有多少只兔子？ 聪明的你可能已经发现，第 n 个月的兔子数正好是第 n 个 Fibonacci(斐波那 契)数。栋栋不懂什么是 Fibonacci 数，但他也发现了规律：第 i+2 个月的兔子数 等于第 i 个月的兔子数加上第 i+1 个月的兔子数。前几个月的兔子数依次为： 1 1 2 3 5 8 13 21 34 … 栋栋发现越到后面兔子数增长的越快，期待养兔子一定能赚大钱，于是栋栋 在第一个月初买了一对小兔子开始饲养。 每天，栋栋都要给兔子们喂食，兔子们吃食时非常特别，总是每 k 对兔子围 成一圈，最后剩下的不足 k 对的围成一圈，由于兔子特别害怕孤独，从第三个月 开始，如果吃食时围成某一个圈的只有一对兔子，这对兔子就会很快死掉。 我们假设死去的总是刚出生的兔子，那么每个月的兔子数仍然是可以计算的。 例如，当 k=7 时，前几个月的兔子数依次为： 1 1 2 3 5 7 12 19 31 49 80 … 给定 n，你能帮助栋栋计算第 n 个月他有多少对兔子么？由于答案可能非常 大，你只需要告诉栋栋第 n 个月的兔子对数除 p 的余数即可。

输入

输入一行，包含三个正整数 n, k, p。

输出

输出一行，包含一个整数，表示栋栋第 n 个月的兔子对数除 p 的余数。

样例输入

6 7 100

样例输出

7

题解：

　　矩阵快速幂+......万恶的分类讨论。

　　%%%%http://blog.csdn.net/u011265346/article/details/46331419。

　　

#include<cstdio>
#include<map>
#include<iostream>
#include<cstring>
using namespace std;
map<long long ,int >mp;
typedef long long ll;
const int N=(int ) 1e6+50;
inline ll powmod(ll a,ll b,ll p){
    ll ans=1;
    while(b){
        if(b&1) ans=a*ans%p;
        a=a*a%p;
        b>>=1;
    }return ans;
}
int vis[12*N];
ll n,k,p,phi_k;
ll fib[12*N];
ll step[N];
int cnt,from;
bool circle;
inline ll gcd(ll a,ll b){
    return b==0?a:gcd(b,a%b);
}
inline ll phi(ll x){
    ll ans=1;
    for(ll i=2;i*i<=x;i++)
        if(x%i==0){
            ans*=i-1;
            x/=i;
            while(x%i==0)
                x/=i,ans*=i;
        }
    return ans*(x==1?1:x-1);
}
inline void init(){
    phi_k=phi(k);
    fib[1]=fib[2]=1;
    for(int i=3;i<=6*k;i++){
        fib[i]=(fib[i-1]+fib[i-2])%k;
        if(!vis[fib[i]])
            vis[fib[i]]=i;
    }
    for(ll i=1,j;;){
        mp[i]=++cnt;
        ll t=gcd(i,k);
        if(t>1) break;
        else{
            j=powmod(i,phi_k-1,k);
            if(!vis[j]){
                break;
            }
            else{
                i=i*fib[vis[j]-1]%k;
                step[cnt]=(ll)vis[j];
                if(mp.count(i)){
                    circle=true;
                    from=mp[i];break;
                }
            }
        }
    }
    step[1]-=2;
}
struct Matrix{
    ll a[4][4];
    Matrix(){memset(a,0,sizeof(a));}
    void e(){
        a[1][2]=a[2][1]=a[2][2]=a[3][3]=1;
    }
    void f(){
        a[1][1]=a[2][2]=a[3][3]=1;a[3][2]=-1;
    }
    friend Matrix operator *(Matrix x,Matrix y){
        Matrix c;
        for(int i=1;i<=3;i++)
            for(int j=1;j<=3;j++){
                for(int k=1;k<=3;k++)
                    (c.a[i][j]+=x.a[i][k]*y.a[k][j])%=p;
                (c.a[i][j]+=p)%=p;
            }
        return c;
    }
    friend Matrix operator ^(Matrix x,ll b){
        Matrix ans;
        for(int i=1;i<=3;i++) ans.a[i][i]=1;
        while(b){
            if(b&1) ans=ans*x;
            b>>=1;
            x=x*x;
        }return ans;
    }
    void print(){
        for(int i=1;i<=3;i++){
            for(int j=1;j<=3;j++)printf("%lld ",a[i][j]);puts("");
        }
    }
}a,b;
ll ans;
inline void solve(){
    if(circle){
        n-=2;
        a.e(),b.f();
        Matrix now;
        now.a[1][1]=now.a[1][2]=now.a[1][3]=1;
        int i;
        for(i=1;i<from&&n>=step[i];n-=step[i],i++)
            now=now*(a^step[i])*b;
        if(i<from) { 
            now=now*(a^n);
            ans=now.a[1][2];
            return ;
        }
        else{
            ll all_cic=0;
            for(i=from;i<=cnt;i++)
                all_cic+=step[i];
            ll cic=n/all_cic;
            n-=cic*all_cic;
            Matrix c;
            for(i=1;i<=3;i++)   c.a[i][i]=1;
            for(i=from;i<=cnt;i++)
                c=c*(a^step[i])*b;
            now=now*(c^cic);
            for(i=from;n>=step[i];n-=step[i],i++)
                now=now*(a^step[i])*b;
            now=now*(a^n);
            ans=now.a[1][2];return;
        }
    }
    else{
        
        n-=2;
        a.e(),b.f();
        Matrix now;
        now.a[1][1]=now.a[1][2]=now.a[1][3]=1;
        int i;
        for(i=1;step[i]&&n>=step[i];n-=step[i],i++){
            now=now*(a^step[i])*b;
        }
        now=now*(a^n);ans=now.a[1][2];return ;
    }
}
int main(){
    scanf("%lld%lld%lld",&n,&k,&p);
    if(n==1){
        printf("%lld\n",1%p);
        return 0;
    }
    init();
    solve();
    printf("%lld\n",ans);
}