zoukankan      html  css  js  c++  java
  • hdu 3221 Brute-force Algorithm(高速幂取模,矩阵高速幂求fib)

    http://acm.hdu.edu.cn/showproblem.php?pid=3221


    一晚上搞出来这么一道题。。Mark。


    给出这么一个程序。问funny函数调用了多少次。

    我们定义数组为所求:f[1] = a,f[2] = b, f[3] = f[2]*f[3]......f[n] = f[n-1]*f[n-2]。相应的值表示也可为a^1*b^0%p。a^0*b^1%p,a^1*b^1%p,.....a^fib[n-3]*b^fib[n-2]%p。即a,b的指数从n=3以后与fib数列一样。


    由于n非常大。fib[n]也想当大。

    a^fib[n]%p能够利用a^fib[n]%p = a^(fib[n]%phi[p]+phi[p])%p进行降幂,条件时fib[n]>=phi[p]。求fib[n]%phi[p]能够构造矩阵。利用矩阵高速幂求fib[n]%phi[p]。


    #include <stdio.h>
    #include <iostream>
    #include <map>
    #include <set>
    #include <list>
    #include <stack>
    #include <vector>
    #include <math.h>
    #include <string.h>
    #include <queue>
    #include <string>
    #include <stdlib.h>
    #include <algorithm>
    #define LL long long
    #define _LL __int64
    #define eps 1e-12
    #define PI acos(-1.0)
    #define C 240
    #define S 20
    using namespace std;
    const int maxn = 110;
    
    struct matrix
    {
        LL mat[3][3];
        void init()
        {
            memset(mat,0,sizeof(mat));
            for(int i = 0; i < 2; i++)
                mat[i][i] = 1;
        }
    } m;
    
    LL a,b,p,n,phi_p;
    LL fib[10000000];
    
    //phi[p]
    LL Eular(LL num)
    {
        LL res = num;
        for(int i = 2; i*i <= num; i++)
        {
            if(num%i == 0)
            {
                res -= res/i;
                while(num%i == 0)
                    num /= i;
            }
        }
        if(num > 1)
            res -= res/num;
        return res;
    }
    //矩阵相乘
    matrix mul_matrix(matrix x, matrix y)
    {
        matrix ans;
        memset(ans.mat,0,sizeof(ans.mat));
        for(int i = 0; i < 2; i++)
        {
            for(int k = 0; k < 2; k++)
            {
                if(x.mat[i][k] == 0) continue;
                for(int j = 0; j < 2; j++)
                {
                    ans.mat[i][j] = (ans.mat[i][j] + x.mat[i][k]*y.mat[k][j])%phi_p;
                }
            }
        }
        return ans;
    }
    //a^t%phi_p
    LL pow_matrix(LL t)
    {
        matrix a,b;
        a.mat[0][0] = a.mat[0][1] = a.mat[1][0] = 1;
        a.mat[1][1] = 0;
        b.init();
        while(t)
        {
            if(t&1)
                b = mul_matrix(a,b);
            a = mul_matrix(a,a);
            t >>= 1;
        }
        return b.mat[0][0];
    }
    //a^t%p
    LL pow(LL a, LL t)
    {
        LL res = 1;
        a %= p;
        while(t)
        {
            if(t&1)
                res = res*a%p;
            a = a*a%p;
            t >>= 1;
        }
        return res;
    }
    //a^fib[t]%p转化为a^(fib[t]%phi[p]+phi[p])%p,fib[t] >= phi[p]。
    LL solve(LL a, LL t)
    {
        fib[0] = 1;
        fib[1] = 1;
        int i;
        for(i = 2; i <= t; i++)
        {
            fib[i] = fib[i-1] + fib[i-2];
            if(fib[i] >= phi_p)
                break;
        }
        if(i <= t) //当满足条件fib[t] >= phi[p]时,进行降幂
        {
            LL c = pow_matrix(t) + phi_p;
            return pow(a,c);
        }
        else
            return pow(a,fib[t]);
    }
    
    int main()
    {
        int test;
        scanf("%d",&test);
        for(int item = 1; item <= test; item++)
        {
            scanf("%lld %lld %lld %lld",&a,&b,&p,&n);
            printf("Case #%d: ",item);
            if(n == 1)
            {
                printf("%lld
    ",a%p);
                continue;
            }
            if(n == 2)
            {
                printf("%lld
    ",b%p);
                continue;
            }
            if(n == 3)
            {
                printf("%lld
    ",a*b%p);
                continue;
            }
            if(p == 1)
            {
                printf("0
    ");
                continue;
            }
            phi_p = Eular(p);
            LL res = solve(a,n-3)*solve(b,n-2)%p;
            printf("%lld
    ",res);
        }
        return 0;
    }
    


  • 相关阅读:
    Samba 4.0 RC3 发布
    SymmetricDS 3.1.7 发布,数据同步和复制
    Express.js 3.0 发布,Node.js 的高性能封装
    GIFLIB 5.0.1 发布,C语言的GIF处理库
    jQuery UI 1.9.1 发布
    SVN Access Manager 0.5.5.14 发布 SVN 管理工具
    DynamicReports 3.0.3 发布 Java 报表工具
    HttpComponents HttpClient 4.2.2 GA 发布
    AppCan 2.0 正式发布,推移动应用云服务
    Ruby 2.0 的新功能已经冻结
  • 原文地址:https://www.cnblogs.com/yxwkf/p/5113611.html
Copyright © 2011-2022 走看看