mcfx's blog

题解、Writeup、游记和碎碎念

BZOJ 2432: [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 的余数即可。

#include<bits/stdc++.h>

typedef long long ll;

#define fo0(i,n) for(int i=0,i##end=n;i<i##end;i++)
#define fo1(i,n) for(int i=1,i##end=n;i<=i##end;i++)
#define fo(i,a,b) for(int i=a,i##end=b;i<=i##end;i++)

void exgcd(int a,int b,int &x,int &y){if(b)exgcd(b,a%b,y,x),y-=(a/b)*x;else x=1,y=0;}

ll n;
int k,p,pos[1000000],f[6000010]={0,1,1},p2[1000000];

struct mat
{
    int a[3][3];
    inline void operator *=(mat x)
    {
        ll t[3][3];memset(t,0,sizeof(t));
        fo0(i,3)fo0(j,3)fo0(k,3)t[i][k]+=(ll)a[i][j]*x.a[j][k];
        fo0(i,3)fo0(j,3)a[i][j]=t[i][j]%p;
    }
    inline mat operator *(mat x)
    {
        mat r=*this;r*=x;return r;
    }
};

mat pow(mat a,ll p)
{
    mat r=(mat){1,0,0,0,1,0,0,0,1};
    for(;p;p>>=1,a*=a)if(p&1)r*=a;
    return r;
}

struct data
{
    int st,len;mat x;
}s[1000001];

inline mat get(ll n,int tt)
{
    mat r=(mat){1,0,0,0,1,0,0,0,1};
    if(!n)return r;
    for(;n>=s[tt].len;n-=s[tt++].len)r*=s[tt].x;
    return r*pow((mat){0,1,0,1,1,0,0,0,1},n);
}

int main()
{
    scanf("%lld%d%d",&n,&k,&p);
    for(int i=3;;i++)
    {
        f[i]=f[i-1]+f[i-2];
        if(f[i]>=k)f[i]-=k;
        if(!pos[f[i]])pos[f[i]]=i;
        if(f[i]==1&&f[i-1]==0)break;
    }
    int ns=1,sc=1;
    for(;!p2[ns];sc++)
    {
        s[sc].st=ns;
        p2[ns]=sc;
        int t1,t2;
        exgcd(k,ns,t1,t2);
        if(t2<0)t2+=k;
        t1=pos[t2];
        if((ll)ns*t2%k!=1||!t1)break;
        s[sc].len=t1;
        s[sc].x=pow((mat){0,1,0,1,1,0,0,0,1},t1-1)*(mat){0,1,0,1,1,0,0,-1,1};
        ns=(ll)ns*f[t1-1]%k;
    }
    mat bfm=(mat){1,0,0,0,1,0,0,0,1},lpm=(mat){1,0,0,0,1,0,0,0,1},ans=(mat){1,0,0,0,1,0,0,0,1};
    ll bf=0,lp=0;
    fo1(i,p2[ns]-1)bf+=s[i].len,bfm*=s[i].x;
    fo(i,p2[ns],sc-1)lp+=s[i].len,lpm*=s[i].x;
    if(!lp)lp=1,lpm=(mat){0,1,0,1,1,0,0,0,1};
    if(n<bf)
    {
        ans*=get(n,1);
    }
    else
    {
        ans*=bfm;n-=bf;
        ans*=pow(lpm,n/lp);
        ans*=get(n%lp,p2[ns]);
    }
    int a1=ans.a[0][1]+ans.a[2][1];
    while(a1>=p)a1-=p;while(a1<0)a1+=p;
    printf("%d\n",a1);
}

日期： 2016-12-20
标签： BZOJ 暴力矩阵快速幂
这是一篇旧文，原始文章及评论可在 https://oldblog.mcfx.us/archives/156/ 查看。

mcfx's blog

BZOJ 2432: [Noi2011]兔农

Input

Output

Sample Input

Sample Output

Solution

Code