4-快速幂-矩阵快速幂

快速幂-矩阵快速幂

快速幂就是快速算底数的n次幂。其时间复杂度为 O(log₂N)， 与朴素的O(N)相比效率有了极大的提高

(0)同余定理

(a×b)%c=(a%c×b%c)%c

(a+b)%c=(a%c + b%c)%c

(a-b)%c=(a%c - b%c + c)%c

快速幂-矩阵快速幂

(1)快速幂

1. 法一

$a^{2k}=(a^k{)}^2$ 幂次为偶数情况下

$a^{2k+1}=(a^k{)}^2\ast a$ 幂次为奇数情况下

代码以递归形式实现$a^n$

long long mod=1e9+7;
long long quick_pow1(long long a,long long n)
{
    //a^n
    if (n<0)
    {
        cerr << "error quick pow" << endl;
        return -1;
    }
    if (a==0)
        return 0;
    if (n==0)
        return 1;
    a%=mod;
    if ((n&1)==1)
        return a*quick_pow1(a,n-1)%mod;
    else
    {
        long long half=quick_pow1(a,n/2);
        return half*half%mod;
    }
}

1. 法二

求$a^k$则将幂次k按照2的幂次展开，即方法：先将k转换为二进制数，把二进制数首先写成加权系数展开式，计算展开式中的每个加数，$a^k$即展开成：以a为底，幂次为2的n次幂的数，相乘，前提是k>=0

例子：求$a^{29}$，则29D=11101B=$2^0$+$2^2$+$2^3$+$2^4$=1+4+8+16=29D
即$a^{29}$=$a^1$ $a^4$ $a^8$ * $a^{16}$

代码以迭代形式实现$a^n$

long long mod=1e9+7;
long long quick_pow2(long long a,long long n)
{
    if (n<0)
    {
        cerr << "error quick pow" << endl;
        return -1;
    }
    if (a==0)
        return 0;
    a%=mod;
    long long ans=1;
    while (n!=0)
    {
        if ((n&1)==1)
        {
            ans*=a;
            ans%=mod;
        }
        a*=a;
        a%=mod;
        n>>=1;
    }
    return ans;
}

(3)矩阵快速幂

例子：求斐波那契数列

$$\{\begin{array}{c}f(i)=f(i-1)+f(i-2)\\ f(0)=f(1)=1\end{array}$$

常规解法是直接dp模拟一遍，时间复杂度O(n)，空间复杂度用滚动数组可做到O(1)

解：显然$\left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]$ $\left[\begin{array}{c}f0\\ f1\end{array}\right]$ = $\left[\begin{array}{c}f1\\ f2\end{array}\right]$
这里构建了一个矩阵A，使得 A $\left[\begin{array}{c}f0\\ f1\end{array}\right]$ = $\left[\begin{array}{c}f1\\ f2\end{array}\right]$
又可以推导得$\left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]$ $\left[\begin{array}{c}f1\\ f2\end{array}\right]$ = $\left[\begin{array}{c}f2\\ f3\end{array}\right]$
即最终公式1: A $\left[\begin{array}{c}f(i-2)\\ f(i-1)\end{array}\right]$ = $\left[\begin{array}{c}f(i-1)\\ f(i)\end{array}\right]$，A = $\left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]$
且矩阵无交换律，但是有结合律，即A $\left[\begin{array}{c}f0\\ f1\end{array}\right]$ = $\left[\begin{array}{c}f1\\ f2\end{array}\right]$ 且 A $\left[\begin{array}{c}f1\\ f2\end{array}\right]$ = $\left[\begin{array}{c}f2\\ f3\end{array}\right]$
可以推导得 A $\left[\begin{array}{c}f1\\ f2\end{array}\right]$ = A (A $\left[\begin{array}{c}f0\\ f1\end{array}\right]$) = $\left[\begin{array}{c}f2\\ f3\end{array}\right]$ = $A^2$ $\left[\begin{array}{c}f0\\ f1\end{array}\right]$
即最终通项公式2: $A^n$ * $\left[\begin{array}{c}f0\\ f1\end{array}\right]$ = $\left[\begin{array}{c}f(n)\\ f(n+1)\end{array}\right]$ ， 这里斐波那契额数列A = $\left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]$ , 将要求的f(k)代入公式中的f(n)求即可

所以转化成新问题：求$A^n$，此时A是一个矩阵。此时只需结合快速幂的方法二即可。

tips：为了设计方便，通常将通项公式中的所有矩阵都拓展到A的维度
即此题求斐波那契数列，通项公式优化为$(\left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]{)}^n$ * $\left[\begin{array}{cc}f(0)& 0\\ f(1)& 0\end{array}\right]$ = $\left[\begin{array}{cc}f(n)& 0\\ f(n+1)& 0\end{array}\right]$

求斐波那契数列，用矩阵快速幂方法

/*
* f(0)=f(1)=1
* f(i)=f(i-1)+f(i-2) 求f(n)，通过矩阵快速幂
* 则矩阵通项公式
A=[0 1]   即 A^n * [f(0)] = [f(n)  ]
  [1 1]            [f(1)] = [f(n+1)]
*/
//矩阵A的维度是M*M，这里是2*2
const int M=2;
//矩阵A的类
class Ma
{
public:
    int a[M][M];
    Ma()
    {
        memset(a,0,sizeof(a));
    }
    void become_unit_matrix()
    {
        //设置成单位矩阵
        for (int i=0;i<M;++i)
            for (int j=0;j<M;++j)
                a[i][j]= i==j ? 1 : 0;
    }
    void show()
    {
        cout << "----" << endl;
        cout << "Matrix is " << M << '*' << M << endl;
        for (int i=0;i<M;++i)
        {
            for (int j=0;j<M;++j)
                cout << a[i][j] << ' ';
            cout << endl;
        }
        cout << "----" << endl;
    }
    Ma operator*(const Ma &B) const
    {
        //重载两个矩阵相乘
        Ma ans;
        for (int i=0;i<M;++i)
            for (int j=0;j<M;++j)
                for (int k=0;k<M;++k)
                    ans.a[i][j]+=(this->a[i][k])*(B.a[k][j]);
        return ans;
    }
    Ma operator^(int n) const
    {
        //重载幂次运算，并结合快速幂
        Ma ans;
        //单位矩阵乘任何矩阵还是任何矩阵本身
        ans.become_unit_matrix();
        Ma A=*this;
        while(n!=0)
        {
            if ((n&1)==1)
                ans=ans*A;
            A=A*A;
            n>>=1;
        }
        return ans;
    }
};

int get_fibo(int n)
{
    // 0 1 2 3 4 5  6  7  8
    // 1 1 2 3 5 8 13 21 34
    Ma A;
    A.a[0][0]=0;
    A.a[0][1]=A.a[1][0]=A.a[1][1]=1;
    //A.show();
    Ma F;
    F.a[0][0]=F.a[1][0]=1;
    //F.show();
    Ma ans=(A^n)*F;
    //ans.show();
    return ans.a[0][0];
}