zoukankan      html  css  js  c++  java
  • luoguP4000 斐波那契数列

    题目链接

    luoguP4000 斐波那契数列

    题解

    根据这个东西
    https://www.cnblogs.com/sssy/p/9418732.html
    我们可以找出%p意义下的循环节
    然后就可以做了
    人傻,自带,大,常数

    代码

    #include<bits/stdc++.h> 
    using namespace std; 
    #define LL long long 
    const LL maxn = 1000007; 
    LL dp[maxn * 10]; 
    LL prime[maxn],s = 0; 
    bool vis[maxn];
    void init_prime() { 
        for(LL i = 2;i < maxn;i ++) { 
            if(!vis[i]) prime[s ++] = i; 
            for(LL j = 0;j < s && i * prime[j] < maxn;j ++) { 
                vis[i * prime[j]] = 1; 
                if(i%prime[j] == 0) break;
            } 
        } 
    } 
    LL pow_mod(LL a1,LL b1){
        LL ret = 1;
        for(;b1;b1 >>= 1,a1 *= a1) 
            if(b1 & 1) ret = ret * a1; 
        return ret; 
    } 
    LL pow_mod2(LL a,LL b,LL p) { 
        LL ret = 1; 
        for(;b;b >>= 1,a = a * a % p)  
            if(b & 1) ret = ret * a % p; 
        return ret; 
    } 
    LL gcd(LL a,LL b) { return b ? gcd(b,a % b) : a; } 
    bool f(LL n,LL p) { 
        if(pow_mod2(n,(p - 1) >> 1,p) != 1) return false; 
        return true;
    } 
    struct matrix{
        LL x1,x2,x3,x4; 
    }; 
    matrix matrix_a,matrix_b;
    matrix M2(matrix aa,matrix bb,LL mod){
        matrix tmp;
        tmp.x1 = (aa.x1 * bb.x1 % mod + aa.x2 * bb.x3 % mod) % mod; 
        tmp.x2 = (aa.x1 * bb.x2 % mod + aa.x2 * bb.x4 % mod) % mod; 
        tmp.x3 = (aa.x3 * bb.x1 % mod + aa.x4 * bb.x3 % mod) % mod; 
        tmp.x4 = (aa.x3 * bb.x2 % mod + aa.x4 * bb.x4 % mod) % mod; 
        return tmp; 
    }
    matrix M(LL n,LL mod){
        matrix a,b;
        a=matrix_a;b=matrix_b;
        for(;n;n >>= 1,a = M2(a,a,mod))  
            if(n & 1) b = M2(b,a,mod); 
        return b; 
    }  
    LL fac[100][2],l,x,fs[1000]; 
    void dfs(LL count,LL step) { 
        if(step == l) {  
            fs[x ++] = count; 
            return ; 
        } 
        LL sum = 1; 
        for(LL i = 0;i < fac[step][1];++ i) { 
            sum *= fac[step][0]; 
            dfs(count * sum,step + 1);
        }
        dfs(count,step + 1); 
    } 
    LL solve2(LL p){
        if(p < 1e6 && dp[p])   return dp[p]; 
        bool ok = f(5,p);//判断5是否为p的二次剩余 
        LL t; 
        if(ok) t = p - 1; 
        else t = 2 * p + 2; 
        l = 0;
        for(LL i = 0;i < s;i ++){
            if(prime[i] > t / prime[i])  break;  
            if(t % prime[i] == 0) { 
                LL count = 0; 
                fac[l][0] = prime[i];
                while(t % prime[i] == 0) count ++ , t /= prime[i]; 
                fac[l ++][1] = count;
            }
        }
        if(t > 1) fac[l][0]=t, fac[l++][1]=1; 
        x = 0;
        dfs(1 , 0); //求t的因子 
        sort(fs,fs + x);  
        for(LL i = 0;i < x;i ++) { 
            matrix m1 = M(fs[i],p); 
            if(m1.x1 == m1.x4 && m1.x1 == 1 && m1.x2 == m1.x3 && m1.x2 == 0) { 
                if(p < 1e6) dp[p] = fs[i]; 
                return fs[i]; 
            }  
        } 
    } 
    LL get_M(LL n){
        LL ans = 1,cnt; 
        for(LL i = 0;i < s;i ++) { 
            if(prime[i] > n / prime[i]) break; 
            if(n % prime[i] == 0) { 
                 LL count = 0;
                while(n % prime[i] == 0)  count ++,n /= prime[i]; 
                cnt = pow_mod(prime[i],count - 1); 
                cnt *= solve2(prime[i]); 
                ans = (ans / gcd(ans , cnt)) * cnt;
            } 
        }  
        if(n > 1){ 
            cnt = 1; 
            cnt *= solve2(n); 
            ans = ans / gcd(ans,cnt) * cnt;  
        }  
        return ans; 
    } 
    char Ss[30000007]; 
    struct bign { 
        int z[30000007],l; 
        void init() { 
     		memset(z,0,sizeof(z)); 
            scanf("%s",Ss + 1); 
            l = strlen(Ss + 1); 
            for(int i = 1;i <= l;i ++) 
                z[i] = Ss[l - i + 1] - '0'; 
        } 
        LL operator % (const long long & a) const { 
            LL b = 0; 
            for (int i = l;i >= 1;i --) 
                b = (b * 10 + z[i]) % a; 
            return b; 
        } 
    }z;  
    LL m1; 
    struct Matrix { 
        LL a[3][3]; 
        Matrix () { memset(a,0,sizeof a); }   	
        Matrix  operator * (const Matrix & p) const { 
            Matrix ret; 
            for(int i = 0;i <= 1;++ i) 
                for(int j = 0;j <= 1;++ j) 
                    for(int k = 0;k <= 1;++ k) 
                        ret.a[i][j] = (ret.a[i][j] + ( a[i][k] * p.a[k][j] ) % m1 ) % m1; 
            return ret; 
        } 
    }; 
    LL solve(LL x) { 
        Matrix p,q; 
        p.a[0][0] = 1; p.a[0][1] = 1; p.a[1][0] = 1; 
        q.a[0][1] = 1;  
        for(;x;x >>= 1,p = p * p) 
            if(x & 1) q = q * p; 
        return q.a[0][0]; 
    } 
    main() { 
        LL t,n;
        init_prime(); 
        matrix_a.x1 = matrix_a.x2 = matrix_a.x3 = 1;   
        matrix_a.x4 = 0; 
        matrix_b.x1 = matrix_b.x4 = 1;  
        matrix_b.x2 = matrix_b.x3 = 0;
        dp[2] = 3;dp[3] = 8;dp[5] = 20; 
        z.init(); scanf("%lld",&n); m1 = n; 
        n = get_M(n); 
        //printf("%lld
    ",n % m1); 
        n = z % n; 
        printf("%lld
    ",solve(n)); 
        return 0; 
    } 
    
  • 相关阅读:
    datanode报错Problem connecting to server
    使用命令查看hdfs的状态
    Access denied for user root. Superuser privilege is requ
    ElasticSearch默认的分页参数 size
    SparkStreaming Kafka 维护offset
    【容错篇】Spark Streaming的还原药水——Checkpoint
    251 Android 线性与相对布局简介
    250 Android Studio使用指南 总结
    249 如何解决项目导入产生的中文乱码问题
    248 gradle更新问题
  • 原文地址:https://www.cnblogs.com/sssy/p/9420285.html
Copyright © 2011-2022 走看看