10/18/2018 快速幂和矩阵快速幂

快速幂：

    通常来说，假如让我们计算x^n的值的时候，我们会想到x * x * x * .... * x，这种方法的效率是O(n)。但是有更高效的方法：

1. 如果当前的指数是偶数，我们就把指数拆成两半，得到两个相同的数，然后把这两个相同的数相乘，可以得到原来的数。

2. 如果当前的指数是奇数，我们把当前的指数拆成两半（floor(n / 2))，得到两个相同的数，然后把这两个数相乘之后再乘一个base，可以得到原来的数

  上述的方法里运用到的思想是Divide and Conquer，这种方法叫做快速幂，效率为O(log(n))。

---

    public static long quickPower(int n, int p) {
        long ans = 1;
        while (p != 0) {
            if (p % 2 == 1) {
                ans *= n;
            }
            n *= n;
            p >>= 1;
        }
        //can also return the result of n^(-p) by using 1.0 / ans
        return ans;
    }

---

矩阵快速幂：

   基本思想和快速幂一样，也是运用Divide and Conquer的思想，假设有矩阵A和幂n。若幂是偶数，则将幂分成两半，得到两个相同的矩阵（A^(n / 2))，最终把这个两个矩阵相乘得到最终答案。若幂是奇数，则将幂分成两半（floor(n / 2))，得到两个相同的矩阵，最终再乘一次A，得到最终答案。

首先假设我们有两个矩阵A和B，A和B可以相乘仅当A的列数等于B的行数，所以一个假如一个矩阵想要可以自己与自己相乘，它必须是方阵（n * n)。我们先写出矩阵相乘的代码：

    public static long[][] matrixMultiply(long[][] m1, long[][] m2) {
        long[][] ans = new long[m1.length][m2[0].length];

        for (int i = 0; i < ans.length; i++) {
            for (int j = 0; j < ans[0].length; j++) {
                for (int k = 0; k < m2.length; k++) {
                    ans[i][j] += m1[i][k] * m2[k][j];
                }
            }
        }
        return ans;
    }

之后我们再将“矩阵相乘”和“快速幂”运用在一起，得到下面的代码：

    public static long[][] quickMatrixMultiply(long[][] m1, int p) {
        long[][] ans = new long[m1.length][m1.length];
        while (p != 0) {
        //Why is ans and m1
            if (p % 2 == 1) {
                ans = matrixMultiply(ans, m1);
            }
            m1 = matrixMultiply(m1, m1);
            p >>= 1;
        }
        return ans;
    }

最终矩阵快速幂的效率是O(n^3(log(n)))。

---

应用：

快速幂的一个常见应用是用来计算Fibonacci Numbers。我们知道Fibonacci Numbers的递推式是

假如把它看成一个 1 * 2 的矩阵，即，则我们能构造矩阵 T =  *  = ，即

                                                          

这样，我们就可以用矩阵快速幂来计算Fibonacci了。

Sidenote: what matters is the second column in the matrix

---

Example:

GNY Regional 2015 Immortal Porpoises https://open.kattis.com/problems/porpoises 

import java.util.*;
public class ImmortalPorpoises {
	static long b = 1000000000;
	public static long[][] matrixMultiplication(long[][] m1, long[][] m2) {
		long[][] ans = new long[m1.length][m2[0].length];
		for (int i = 0; i < ans.length; i++) {
			for (int j = 0; j < ans[0].length; j++) {
				for (int k = 0; k < m2.length; k++) {
					//(a + b) % m = ((a % m) + (b % m)) % m 
					ans[i][j] = (ans[i][j] % b + m1[i][k] * m2[k][j] % b) % b;
				}
			}
		}
		return ans;
	}
	public static long[][] quickMatrixMultiplication(long[][] m, long n) {
		//Notice this is the identity matrix
		long[][] ans = new long[m.length][m.length];
		for (int i = 0; i < ans.length; i++) {
			ans[i][i] = 1;
		}
		while (n != 0) {
			if (n % 2 == 1) 
				ans = matrixMultiplication(ans, m);
			m = matrixMultiplication(m, m);
			//equivalent to n /= 2
			n >>= 1;
		}
		return ans;
	}
	public static long getFibo(long n) {
		if (n == 1 || n == 2) {
			return 1;
		}
		long res = 0;
		long[][] start = {{0,1}, {1,1}};
		//Actually, why is this n - 2? A: You can have T(n) = T(n - 1) * A^m (0 <= m) and you can 
		//observe that when n = 3 you only need one transformation. Thus, m = n - 2.
		long[][] temp = quickMatrixMultiplication(start, n - 2);
		res += (temp[0][1] + temp[1][1]) % b;
		return res;
	}
	public static void main(String args[]) {
		Scanner sc = new Scanner(System.in);
		int kase = sc.nextInt();
		for (int i = 0; i < kase; i++) {
			long index = sc.nextLong();
			long power = sc.nextLong();
			System.out.println(index + " " + getFibo(power));
		}
		sc.close();
	}
}