目錄
快速冪
快速冪取模
矩陣快速冪
矩陣快速冪取模
[HDU1005練習](#Number Sequence)
快速冪
? 冪運算:\(x ^ n\)
? 根據其一般定義我們可以簡單實作其非負整數情況下的函式
定義法:
int Pow (int x, int n) {
int result = 1;
while(n--) {
result *= x;
}
return result;
}
? 不難看出此時演算法的時間復雜度是\(O(n)\),一旦n取較大數值,計算時間就會大大增加,極其容易出現超時的情況,
快速冪:
? 首先要在此列舉兩個前提:
-
計算機是通過二進制進行存盤資料的,對二進制資料可以直接進行操作,
-
\(2^n+2^n=2*2^n=2^{n + 1}\)
? 對于\(x ^ n\),其中n可以表示為m位的二進制形式,即\(n=n_12^0+n_22^1+n_32^3+\cdots+n_m2^{m-1}\)
? 那么\(x ^ n=x ^ {n_12^0+n_22^1+n_32^3+\cdots+n_m2^{m-1}}\)
? 即\(x^n=x ^ {n_12^0}x^{n_22^1}x^{n_32^3}\cdots x^{n_m2^{m-1}}\)
? 根據前提1,計算機可以直接對二進制格式存盤的資料進行讀取,那么我們就可以對\(n\)一個個讀取再對其進行累乘,
? 當取到第\(a\)位時,其乘數有通項公式:\(x^{n_a2^{a-1}}\)
? 通過標準庫math,用代碼實作:
int Pow (int x, int n) {
int result = 1;
int a = 1;
while(n) {
if(n & 1) result *= round( pow(x, 1 << (a - 1)) );//round的作用在于double轉int時防止丟失精度,對1進行位運算是一種計算2的n次冪的快速方法
n >>= 1;
a++;
}
return result;
}
? 但實際上這個呼叫了標準庫的代碼并沒有實作快速冪,因為仍然是采用pow()進行的運算
此處由 \(2^n+2^n=2\times2^n=2^{n + 1}\)
即\((x ^ {2 ^ {n}} )^2=x ^ {2\times 2 ^ {n}} =x ^ {2 ^ {n + 1}}\)
因此我們可以通過對前項進行二次冪運算得到后項
先得到首項\(f(1)=x^{2^{1-1}}=x\)
即令int f = x
具體實作:
int Pow (int x, int n) {
int result = 1;
int f = x;
while(n) {
if(n & 1) result *= f;
f *= f;
n >>= 1;
}
return result;
}
不難發現此時演算法的時間復雜度由其次數的二進制位數而決定,即\(O(m)\)也就是\(O(log_2n)\)
另外因為此演算法常常用于計算大數,所以int型別最好都換成long long型別,防止溢位,
快速冪取模
? 對\(x^n\)取模運算:\(x^n\mod p\),當n極大時,同樣其運算量也會增大
? 這就需要我們尋求一種快速解決的方法,此時可以運用我們上文所提到的快速冪演算法
? 引論1:\((np+a)\mod p=a\mod p\quad (n\in\mathbb{Z} )\)
? 證明如下:設\(f(n,p,a)=(np+a)\mod p\quad (n\in\mathbb{Z} )\)
? 則由定義有\((np+a)\mod p\\=\frac{np+d}{p}-([\frac{np+d}{p}+1]-1)\\ =d-p([\frac{d}{p}+1]-1)\)
? 顯而易見,\((np+a)\mod p=a\)與\(n\)無關
? 令\(n=0\)得\((np+a)\mod p=a\mod p\quad (n\in\mathbb{Z} )\\ Q.E.D\)
? 引論2:\((mn)\mod p=((m\mod p)(n\mod p))\mod p\)
? 證明如下:令\(\left\{ \begin{array}{c} m=(ip+c)& (i\in\mathbb{Z} )\\ n=(jp+d) & (j\in\mathbb{Z} )\end{array} \right.\)
? 則有\(\left\{ \begin{array}{c} m\mod p=c\\ n\mod p=d \end{array} \right.\)
? 原式\((m*n)%p\\=((ip+c)(jp+d))\mod p\\=(jip^2+jpc+ipd+cd)\mod p\\=(jip^2+jpc+ipd+cd)\mod p\\=((jip+jc+id)p+cd)\mod p\\因為(jip+jc+id)\in\mathbb{Z},由引論1得:\\=(cd)\mod p\\=((m\mod p)(n\mod p))\mod p\)
? 即\((mn)\mod p=((m\mod p)(n\mod p))\mod p\\ Q.E.D\)
?
? 因此對于\(x^n\mod p\)
? 可以寫成\((x ^ {n_12^0}x^{n_22^1}x^{n_32^3}\cdots x^{n_m2^{m-1}})\mod p\\ =((x ^ {n_12^0}\mod p)(x^{n_22^1}\mod p)(x^{n_32^3}\mod p)\cdots (x^{n_m2^{m-1}}\mod p))\mod p\)
? 并且由之前的推定\((x ^ {2 ^ {n}} )^2=x ^ {2\times 2 ^ {n}} =x ^ {2 ^ {n + 1}}\)
? 有\((x ^ {2 ^ {n}} \mod p)^2\mod p =(x ^ {2 ^ {n}})^2\mod p=x ^ {2 ^ {n + 1}}\mod p\)
? 代碼實作:
int Mod (int x, int n, int p) {
int result = 1;
int f = x % p;
while(n) {
if(n & 1) result = (result * f)%p;
f = (f * f)%p;
n >>= 1;
}
return result;
}
矩陣快速冪
? 當\(X\)為任意方塊矩陣,即\(X=\left| \begin{matrix} x_{11} & \cdots & x_{1n}\\ \vdots & \ddots & \vdots \\ x_{n1} & \cdots & x_{nn}\\ \end{matrix} \right|\)時
? 同理\(X^a\),有\(X^a=X ^ {a_12^0}X^{a_22^1}X^{a_32^3}\cdots X^{a_m2^{m-1}}\)
? 同樣的,任意方塊矩陣也適用于快速冪
? 代碼實作:
#include <iostream>
#include <cstring>
using namespace std;
#define R 2
struct Matrix{
int data[R][R];
};
Matrix multi(Matrix a,Matrix b,int rec){
Matrix result;
memset(&(result.data), 0, sizeof(result.data));
for(int i = 0; i < rec; i++)
for(int j = 0; j < rec; j++)
for(int k = 0; k < rec; k++)
result.data[i][j] += a.data[i][k] * b.data[k][j];
return result;
}
Matrix poww (Matrix x, int n) {
Matrix result;
memset(&(result.data), 0,sizeof( result.data));
for(int i = 0; i < R; i++) result.data[i][i] = 1;
Matrix f = x;
while(n) {
if(n & 1) result = multi(result, f, R);
f = multi(f, f, R);
n >>= 1;
}
return result;
}
void MatrixPrint(Matrix target) {
for(int i = 0; i < R; i++) {
for(int j = 0; j < R; j++) cout << target.data[i][j]<< " ";
cout << endl;
}
}
int main() {
Matrix result;
memset(&(result.data), 0,sizeof( result.data));
result.data[0][0] = 1;
result.data[0][1] = 2;
result.data[1][0] = 3;
result.data[1][1] = 4;
MatrixPrint(result);
cout << endl;
result = poww(result, 3);
MatrixPrint(result);
return 0;
}
輸出結果:
1 2
3 437 54
81 118
矩陣快速冪取模
? 運算定義:當\(X\)為任意方塊矩陣,即\(X=\left| \begin{matrix} x_{11} & \cdots & x_{1n}\\ \vdots & \ddots & \vdots \\ x_{n1} & \cdots & x_{nn}\\ \end{matrix} \right|\)時
? 則\(X\mod p=\left| \begin{matrix} x_{11} \mod p & \cdots & x_{1n}\mod p\\ \vdots & \ddots & \vdots \\ x_{n1}\mod p & \cdots & x_{nn}\mod p\\ \end{matrix} \right|\)
?
? 引論1:\((m+n)\mod p=((m\mod p)+(n\mod p))\mod p\)
? 證明如下:令\(\left\{ \begin{array}{c} m=(ip+c)& (i\in\mathbb{Z} )\\ n=(jp+d) & (j\in\mathbb{Z} )\end{array} \right.\)
? 則有\(\left\{ \begin{array}{c} m \mod p=c\\ n \mod p=d \end{array} \right.\)
? 原式\((m+n)\mod p\\=((ip+c)+(jp+d))\mod p\\=((i+j)p+c+d)\mod p\\因為(i+j)\in\mathbb{Z},由上文引論得:\\=(c+d)\mod p\\=((m\mod p)+(n\mod p))\mod p\)
? 即\((m+n)\mod p=((m\mod p)+(n\mod p))\mod p\\ Q.E.D\)
? 引論2:有方塊矩陣\(X=\left| \begin{matrix} x_{11} & \cdots & x_{1n}\\ \vdots & \ddots & \vdots \\ x_{nn} & \cdots & x_{nn}\\ \end{matrix} \right|\),\(Y=\left| \begin{matrix} y_{11} & \cdots & y_{1m}\\ \vdots & \ddots & \vdots \\ y_{n1} & \cdots & y_{nm}\\ \end{matrix} \right|\)
? 則有\(XY\mod p=((X\mod p)·(Y\mod p))\mod p\)
? 證明如下:
? 令\(X=\left| \begin{matrix} x_{11} & \cdots & x_{1n}\\ \vdots & \ddots & \vdots \\ x_{n1} & \cdots & x_{nn}\\ \end{matrix} \right|\),\(Y=\left| \begin{matrix} y_{11} & \cdots & y_{1n}\\ \vdots & \ddots & \vdots \\ y_{n1} & \cdots & y_{nn}\\ \end{matrix} \right|\)
? 則
\(X·Y\mod p=\left| \begin{matrix} (x_{11}y_{11} + \cdots +x_{1n}y_{n1})\mod p & \cdots & (x_{11}y_{1n} + \cdots +x_{1n}y_{nn})\mod p\\ \vdots & \ddots & \vdots \\ (x_{n1}y_{11} + \cdots +x_{nn}y_{n1})\mod p & \cdots & (x_{n1}y_{1n} + \cdots +x_{nn}y_{nn})\mod p\\ \end{matrix} \right|\)
\(=\left| \begin{matrix} ((x_{11}y_{11})\mod p + \cdots +(x_{1n}y_{n1})\mod p )\mod p & \cdots & ((x_{11}y_{1n})\mod p + \cdots +(x_{1n}y_{nn})\mod p\mod p\\ \vdots & \ddots & \vdots \\ ((x_{n1}y_{11})\mod p + \cdots +(x_{nn}y_{n1})\mod p)\mod p & \cdots & ((x_{n1}y_{1n})\mod p + \cdots +(x_{nn}y_{nn})\mod p)\mod p\\ \end{matrix} \right|\)
\(=\left| \begin{matrix} (((x_{11}\mod p)(y_{11}\mod p))\mod p + \cdots +((x_{1n}\mod p)(y_{n1}\mod p))\mod p )\mod p & \cdots & (((x_{11}\mod p)(y_{1n}\mod p))\mod p + \cdots +(((x_{1n}\mod p)(y_{nn}\mod p))\mod p)\mod p\\ \vdots & \ddots & \vdots \\ (((x_{n1}\mod p)(y_{11}\mod p))\mod p + \cdots +((x_{nn}\mod p)(y_{n1}\mod p))\mod p)\mod p & \cdots & (((x_{n1}\mod p)(y_{1n}\mod p))\mod p + \cdots +((x_{nn}\mod p)(y_{nn}\mod p))\mod p)\mod p\\ \end{matrix} \right|\)
\(=\left| \begin{matrix} ((x_{11}\mod p)(y_{11}\mod p) + \cdots +(x_{1n}\mod p)(y_{n1}\mod p))\mod p & \cdots & ((x_{11}\mod p)(y_{1n}\mod p) + \cdots +((x_{1n}\mod p)(y_{nn}\mod p))\mod p\\ \vdots & \ddots & \vdots \\ ((x_{n1}\mod p)(y_{11}\mod p) + \cdots +(x_{nn}\mod p)(y_{n1}\mod p))\mod p & \cdots & ((x_{n1}\mod p)(y_{1n}\mod p) + \cdots +(x_{nn}\mod p)(y_{nn}\mod p))\mod p\\ \end{matrix} \right|\)
\(=\left( \left| \begin{matrix} x_{11}\mod p & \cdots & x_{1n}\mod p\\ \vdots & \ddots & \vdots \\ x_{n1}\mod p & \cdots & x_{nn}\mod p\\ \end{matrix} \right|\left| \begin{matrix} y_{11}\mod p & \cdots & y_{1n}\mod p\\ \vdots & \ddots & \vdots \\ y_{n1}\mod p & \cdots & y_{nn}\mod p\\ \end{matrix} \right| \right)\mod p\)
\(=((X\mod p)·(Y\mod p))\mod p \\\text{即}XY\mod p=((X\mod p)·(Y\mod p))\mod p\\Q.E.D\)
說明任意矩陣也符合快速冪取模的演算法
具體實作:
#include <iostream>
#include <cstring>
using namespace std;
#define R 2
struct Matrix{
int data[R][R];
};
Matrix multi(Matrix a,Matrix b,int rec,int p){
Matrix result;
memset(&(result.data), 0,sizeof( result.data));
for(int i = 0; i < rec; i++)
for(int j = 0; j < rec; j++) {
for(int k = 0; k < rec; k++)
result.data[i][j] += a.data[i][k] * b.data[k][j];
result.data[i][j] %= p;
}
return result;
}
Matrix poww (Matrix x, int n, int p) {
Matrix result;
memset(&(result.data), 0,sizeof( result.data));
for(int i = 0; i < R; i++) result.data[i][i] = 1;
Matrix f = x;
while(n) {
if(n & 1) result = multi(result, f, R, p);
f = multi(f, f, R, p);
n >>= 1;
}
return result;
}
void MatrixPrint(Matrix target) {
for(int i = 0; i < R; i++) {
for(int j = 0; j < R; j++) cout << target.data[i][j]<< " ";
cout << endl;
}
}
int main() {
Matrix result;
memset(&(result.data), 0,sizeof( result.data));
result.data[0][0] = 1;
result.data[0][1] = 2;
result.data[1][0] = 3;
result.data[1][1] = 4;
MatrixPrint(result);
cout << endl;
result = poww(result, 3, 3);
MatrixPrint(result);
return 0;
}
輸出結果:
1 2
3 41 0
0 1
練習(HDU1005):
Number Sequence
Problem Description
A number sequence is defined as follows:
\(f(1) = 1, f(2) = 1, f(n) = (A * f(n - 1) + B * f(n - 2)) \mod 7.\)
Given A, B, and n, you are to calculate the value of f(n).
Input
The input consists of multiple test cases. Each test case contains 3 integers A, B and n on a single line (1 <= A, B <= 1000, 1 <= n <= 100,000,000). Three zeros signal the end of input and this test case is not to be processed.
Output
For each test case, print the value of f(n) on a single line.
Sample Input
1 1 3
1 2 10
0 0 0
Sample Output
2
5
Author
CHEN, Shunbao
Source
ZJCPC2004
? 由題意不難發現,題目輸入三個數A,B,n,要求我們求出廣義的斐波那契數列(generalized Fibonacci sequence)第n項的取模
? 即\(f(n) = (A * f(n - 1) + B * f(n - 2)) \mod 7 = (Aa_{n-1}+Ba_{n-2})= a_n \mod 7\)
? 因此通過前后項關系,可以構建矩陣\(\left| \begin{matrix} a_{n} & 0\\ a_{n-1} & 0\\ \end{matrix} \right|\mod 7=\left( \left| \begin{matrix} A & B\\ 1 & 0\\ \end{matrix} \right|\left| \begin{matrix} a_{n-1} & 0\\ a_{n-2} & 0\\ \end{matrix} \right| \right)\mod 7\)
\(=\left( \left| \begin{matrix} A & B\\ 1 & 0\\ \end{matrix} \right|^{n-2}\left| \begin{matrix} a_2 & 0\\ a_1 & 0\\ \end{matrix} \right| \right)\mod 7\)
\(=\left(\left( \left| \begin{matrix} A & B\\ 1 & 0\\ \end{matrix} \right|^{n-2}\mod 7\right)\left(\left| \begin{matrix} a_2 & 0\\ a_1 & 0\\ \end{matrix} \right|\mod7\right) \right)\mod 7\)
如此,就轉化為一個矩陣快速冪的問題
具體實作(謝謝提醒!原來我的solution忽略了\(n = 1\)和\(n = 2\)的情況,導致無限回圈):
#include <iostream>
#include <cstring>
using namespace std;
#define R 2
struct M{
int data[R][R];
};
M multi(M a,M b,int rec,int p){
M result;
memset(&(result.data), 0,sizeof( result.data));
for(int i = 0; i < rec; i++)
for(int j = 0; j < rec; j++) {
for(int k = 0; k < rec; k++)
result.data[i][j] += a.data[i][k] * b.data[k][j];
result.data[i][j] %= p;
}
return result;
}
M poww (M x, int n, int p) {
M result;
M f = x;
memset(&(result.data), 0,sizeof( result.data));
for(int i = 0; i < R; i++) result.data[i][i] = 1;
while(n) {
if(n & 1) result = multi(result, f, R, p);
f = multi(f, f, R, p);
n >>= 1;
}
return result;
}
M solve(int a, int b, int n, int p) {
M result;
M trans;
memset(&(result.data), 0,sizeof( result.data));
memset(&(trans.data), 0,sizeof( trans.data));
result.data[0][0] = 1;
result.data[1][0] = 1;
trans.data[0][0] = a;
trans.data[0][1] = b;
trans.data[1][0] = 1;
trans = poww(trans, n, p);
result = multi(trans, result, R, p);
return result;
}
int main() {
int a, b, n;
while(true) {
cin >> a >> b >> n;
if(!a && !b && !n) break;
if(n == 2 || n == 1 ){
cout << 1 << endl;
continue;
}
cout << solve(a, b, n - 2, 7).data[0][0] << endl;
}
return 0;
}
轉載請註明出處,本文鏈接:https://www.uj5u.com/houduan/83184.html
標籤:C++
下一篇:Spring——事務管理
