HDU 3978 斐波那契循环节 ( 不会写啊!!)

传送门

题意:给出f(f(f...f(n)...)) 总共嵌套k次。问最后模p的值是多少。

首先应该明白的是这个题有循环节的。一个数模N的循环节就是这个数分解成素因子乘积的形式p1^a1*p2^a2*p3^a3...后,斐波那契模pi^ai的循环节的最大公约数。

那么一个素数的k次幂的循环节=斐波那契模上这个素数的循环节乘上p^(k-1)。

而一个素数p的循环节 如果p>5并且是5的二次剩余,那么循环节就是(p-1)的因子,否则就是2*(p+1)的因子。所以2 3 5 的时候需要特判一下。

知道这些就能求每一次嵌套的循环节了,通过矩阵连乘即可得出答案。

和这道题差不多 https://blog.csdn.net/henucm/article/details/102681157

#include <iostream>
#include <string.h>
#include <algorithm>
#include <stdio.h>
#include <math.h>
 
using namespace std;
typedef long long LL;
 
const int M=2;
 
struct Matrix
{
    LL m[M][M];
};
 
Matrix per= {1,0,0,1};
 
Matrix multi(Matrix a,Matrix b,LL MOD)
{
    Matrix c;
    int i,j,k;
    for(i=0; i<M; i++)
    {
        for(j=0; j<M; j++)
        {
            c.m[i][j]=0;
            for(k=0; k<M; k++)
            {
                c.m[i][j]+=a.m[i][k]*b.m[k][j]%MOD;
            }
            c.m[i][j]%=MOD;
        }
    }
    return c;
}
 
Matrix power(Matrix a,LL k,LL MOD)
{
    Matrix ans=per,p=a;
    while(k)
    {
        if(k&1)
        {
            ans=multi(ans,p,MOD);
            k--;
        }
        k>>=1;
        p=multi(p,p,MOD);
    }
    return ans;
}
 
LL gcd(LL a,LL b)
{
    return b? gcd(b,a%b):a;
}
 
LL quick_mod(LL a,LL b,LL m)
{
    LL ans=1;
    a%=m;
    while(b)
    {
        if(b&1)
        {
            ans=ans*a%m;
            b--;
        }
        b>>=1;
        a=a*a%m;
    }
    return ans;
}
 
//勒让德符号
int legendre(int a,int p)
{
    if(quick_mod(a,(p-1)>>1,p)==1) return 1;
    else                           return -1;
}
 
const int N=1000005;
const int NN=50005;
 
bool prime[N];
int p[N];
int num[NN],pri[NN];
int num1[NN],pri1[NN];
int arr[NN];
int loop[N];
int k,cnt,c;
 
void isprime()
{
    k=0;
    int i,j;
    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;
            }
        }
    }
}
 
void find(int n,int pri[],int num[])
{
    cnt=0;
    int t=(int)sqrt(1.0*n);
    for(int i=0; p[i]<=t; i++)
    {
        if(n%p[i]==0)
        {
            int a=0;
            pri[cnt]=p[i];
            while(n%p[i]==0)
            {
                a++;
                n/=p[i];
            }
            num[cnt]=a;
            cnt++;
        }
    }
    if(n>1)
    {
        pri[cnt]=n;
        num[cnt]=1;
        cnt++;
    }
}
 
void dfs(int dept,int product=1)
{
    if(dept==cnt)
    {
        arr[c++]=product;
        return;
    }
    for(int i=0; i<=num1[dept]; i++)
    {
        dfs(dept+1,product);
        product*=pri1[dept];
    }
}
 
int find_loop(int n)
{
    find(n,pri,num);
    int cnt1=cnt;
    LL ans=1;
    for(int i=0; i<cnt1; i++)
    {
        c=0;
        int record=1;
        if(pri[i]==2)
            record=3;
        else if(pri[i]==3)
            record=8;
        else if(pri[i]==5)
            record=20;
        else
        {
            if(legendre(5,pri[i])==1)
                find(pri[i]-1,pri1,num1);
            else
                find(2*(pri[i]+1),pri1,num1);
            dfs(0,1);
            sort(arr,arr+c);
            for(int k=0; k<c; k++)
            {
                Matrix A;
                A.m[0][0]=1;
                A.m[0][1]=1;
                A.m[1][0]=1;
                A.m[1][1]=0;
                Matrix a=power(A,arr[k]-1,pri[i]);
                int x=(a.m[0][0]+a.m[0][1])%pri[i];
                int y=(a.m[1][0]+a.m[1][1])%pri[i];
                if(x==1&&y==0)
                {
                    record=arr[k];
                    break;
                }
            }
        }
        for(int k=1; k<num[i]; k++)
            record*=pri[i];
        ans=ans/gcd(ans,record)*record;
    }
    return ans;
}
 
void Solve(int p,int k)
{
    loop[0]=p;
    for(int i=1; i<=k; i++)
        loop[i]=find_loop(loop[i-1]);
}
 
int work(int n,int k,int p)
{
    int t=n;
    LL ret,MOD;
    Matrix ans;
    Matrix A;
    A.m[0][0]=1;
    A.m[0][1]=1;
    A.m[1][0]=1;
    A.m[1][1]=0;
    Solve(p,k);
    for(int i=k; i>=0; i--)
    {
        MOD=loop[i];
        ans=power(A,t,MOD);
        ret=(ans.m[1][0]+ans.m[1][1])%MOD;
        t=ret;
    }
    return ret;
}
 
int main()
{
    isprime();
    int T,n,k,p,tt=1;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d%d",&n,&k,&p);
        printf("Case #%d: %d\n",tt++,work(n,k,p));
    }
    return 0;
}

 

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值