快速幂（算法）的原理
快速幂算法
- 快速幂
- 数学原理
- 算法实现
- OJ题展示
- 不用高精度计算
- 二进制指数的高精度计算
- 数学题
- 等差数列和等比数列
- 计数原理
快速幂
求 ( a b ) % n (a^b)\%n (ab)%n的结果（即 a a a的 b b b次方，再除以 n n n得到的余数）。
利用程序求解时，可以利用循环语句，将 b b b个 a a a相乘，并在每次乘法运算（要防止溢出）后取模 n n n，即可求出结果，时间复杂度 O ( b ) O(b) O(b)。
即使是这样简单的问题，当 b b b的取值很大时，程序的效率也不高，比如 b = 1 0 9 b = 10^9 b=109，计算机要在1s内解决一个这样的问题都比较困难。如果要求解很多个，显然会超时。因此，需要思考算法的优化问题。
利用分治思想进行优化：
若 b b b是偶数，则 a b = a b 2 × a b 2 a^b = a^{\frac{b}{2}}\times a^{\frac{b}{2}} ab=a2b×a2b；若 b b b是奇数， a b = a b 2 × a b 2 × a a^b = a^{\frac{b}{2}}\times a^{\frac{b}{2}}\times a ab=a2b×a2b×a（这里 b 2 \frac{b}{2} 2b取下整）。即我们只要求出 a b 2 a^{\frac{b}{2}} a2b，然后再通过一次乘法运算，即可求出 a b a^b ab。继续将 b 2 \frac{b}{2} 2b分解成 b 4 \frac{b}{4} 4b， b 8 \frac{b}{8} 8b……最终 b b b会变成1，并且只需要 [ l o g 2 b ] [log_2b] [log2b]次就可以将 b b b变成1。我们称这种快速计算幂运算的方法为快速幂算法。
数学原理
在做这块的题的时候被极致的数值给震撼（ex）了。所以十分有必要补一下相关的数学知识来应对极致的数值。毕竟不可能每个题的每个步骤都用高精度计算代替。
《深入浅出程序设计》进阶篇的部分参考（侵权删，为方便复习删减部分）：
定义1：对于整数 a a a和 b b b，满足 b > 0 b>0 b>0，则存在唯一的整数 q q q和 r r r，满足 a = b q + r a = bq + r a=bq+r，其中 0 ≤ r < b 0 \leq r < b 0≤r<b。其中称 q q q为商、 r r r为余数。余数用 a m o d b a\mod b amodb或者 ( a % b ) (a \% b) (a%b)表示。
带模算式的运算：
(1)
(
a
+
b
)
m
o
d
M
=
(
a
m
o
d
M
+
b
m
o
d
M
)
m
o
d
M
(a + b) \mod M = (a \mod M + b \mod M) \mod M
(a+b)modM=(amodM+bmodM)modM
(2)
(
a
−
b
)
m
o
d
M
=
(
a
m
o
d
M
−
b
m
o
d
M
)
m
o
d
M
(a - b) \mod M = (a \mod M - b \mod M) \mod M
(a−b)modM=(amodM−bmodM)modM
即使 a ≥ b a \geq b a≥b，也不能保证 a m o d M ≥ b m o d M a \mod M \geq b \mod M amodM≥bmodM，假如12个苹果分给5个小朋友剩下2个，而分15个苹果就不剩了。
在c++中，当 a % M − b % M a \% M - b \% M a%M−b%M小于0时， ( a % M − b (a \% M - b % M) \% M (a%M−b的结果是负数，那么就可以写成 ( a % M − b % M + M ) % M (a \% M - b \% M + M) \% M (a%M−b%M+M)%M以保证结果是非负数。
(3) ( a × b ) m o d M = ( ( a m o d M ) × ( b m o d M ) ) m o d M (a \times b) \mod M = ((a \mod M) \times (b \mod M)) \mod M (a×b)modM=((amodM)×(bmodM))modM
类似乘法结合律，可带入定义证明。除法见乘法逆元。
定义2：若两数 x x x和 y y y除以 b b b的余数相等，则称 x x x和 y y y模 b b b同余，记作 x ≡ y ( m o d b ) x \equiv y \pmod b x≡y(modb)。
依据这种分类方法，可以得到同余的以下几个最基础的性质。
- 反身性： x ≡ x ( m o d M ) x \equiv x \pmod M x≡x(modM)。
- 对称性：若 x ≡ y ( m o d M ) x \equiv y \pmod M x≡y(modM)，则 y ≡ x ( m o d M ) y \equiv x \pmod M y≡x(modM)。
- 传递性：若 x ≡ y ( m o d M ) x \equiv y \pmod M x≡y(modM)， y ≡ z ( m o d M ) y \equiv z \pmod M y≡z(modM)，则 x ≡ z ( m o d M ) x \equiv z \pmod M x≡z(modM)。
根据上述模运算的计算方法，同余还有以下性质：
4. 同加性：若
x
≡
y
(
m
o
d
M
)
x \equiv y \pmod M
x≡y(modM)，则
x
+
z
≡
y
+
z
(
m
o
d
M
)
x + z \equiv y + z \pmod M
x+z≡y+z(modM)。
5. 同乘性：若
x
≡
y
(
m
o
d
M
)
x \equiv y \pmod M
x≡y(modM)，则
x
×
z
≡
y
×
z
(
m
o
d
M
)
x\times z \equiv y\times z \pmod M
x×z≡y×z(modM)。
6. 同幂性：若
x
≡
y
(
m
o
d
M
)
x \equiv y \pmod M
x≡y(modM)，则
x
n
≡
y
n
(
m
o
d
M
)
x^n \equiv y^n \pmod M
xn≡yn(modM)。
同余还有以下这个非常重要的性质：
7.
x
≡
y
(
m
o
d
b
)
⇔
b
∣
(
x
−
y
)
x \equiv y \pmod b \Leftrightarrow b | (x - y)
x≡y(modb)⇔b∣(x−y)。
以分苹果的例子解释这个性质。
12
12
12个苹果分给
5
5
5人，最后剩下
2
2
2个；
22
22
22个苹果分给5
人，依然剩下
2
2
2个。也就是
22
≡
12
(
m
o
d
5
)
22 \equiv 12 \pmod 5
22≡12(mod5)。可以发现，
22
22
22个苹果比12个苹果多出来的那些苹
果，一定平分到所有
5
5
5个小朋友手上，每个小朋友多分到
2
2
2个苹果，一共多了
5
×
2
=
10
5\times2 = 10
5×2=10个苹果；如
果每个小朋友多分到
x
x
x个苹果，那么总共就多了
5
x
5x
5x个苹果，并且
5
5
5一定是
5
x
5x
5x的一个约数，也可
表示为
5
∣
(
22
−
12
)
5 | (22 - 12)
5∣(22−12)。类似地，也可证明右式可以推导出左式。
模运算的其他乘法性质（忘了从哪里抄来的）：
- 若 a ≡ b ( m o d m ) a \equiv b \pmod{m} a≡b(modm)且 c ≡ d ( m o d m ) c \equiv d \pmod{m} c≡d(modm)，那么 a c ≡ b d ( m o d m ) ac \equiv bd \pmod{m} ac≡bd(modm)。
- 分配律： ( a + b ) × c m o d m = ( ( a × c ) m o d m + ( b × c ) m o d m ) m o d m (a + b)\times c\bmod m=((a\times c)\bmod m+(b\times c)\bmod m)\bmod m (a+b)×cmodm=((a×c)modm+(b×c)modm)modm。
根据结合律可以了解到，快速幂算法是合理的，对中间结果取模不会影响最终结果。
算法实现
快速幂算法的非递归实现：
int quickPow(int a, int b, int n) {//int可能存在负数的情况
int ret = 1;
while (b) {
if (b % 2 == 1)//b的中间值是奇数的情况
ret = ret * a % n;//*作为乘法的优先级高于%
a = a * a % n;
b = b / 2;
}
return ret;
}
比如 a 13 a^{13} a13，因为 1 3 ( 10 ) = 110 1 ( 2 ) 13_{(10)}=1101_{(2)} 13(10)=1101(2)，所以 a 13 = a 8 × a 4 × a a^{13}=a^8\times a^4\times a a13=a8×a4×a。
可以看到， a 13 a^{13} a13可以按照 13 13 13的二进制数进行拆分：
根据公式 ( a × b ) m o d M = ( ( a m o d M ) × ( b m o d M ) ) m o d M (a \times b) \mod M = ((a \mod M) \times (b \mod M)) \mod M (a×b)modM=((amodM)×(bmodM))modM，即使每个 a a a的若干整数次幂都对 n n n求模， ( a 13 ) % n (a^{13})\%n (a13)%n的结果依旧不会有影响。
这个特性同样可以应用到实际应用中：比如指数有200位，此时就可以先将指数转换成二进制，再进行计算。
如果运算过程中出现负数，可尝试将int换成long long或unsigned long long。
long long的存储范围是
[
−
2
63
,
2
63
−
1
]
[-2^{63},2^{63} - 1]
[−263,263−1]，即
[
−
9223372036854775808
,
9223372036854775807
]
[-9223372036854775808,9223372036854775807]
[−9223372036854775808,9223372036854775807]。
而unsinged long long的存储范围是
[
0
,
2
64
]
[0,2^{64}]
[0,264]，即
[
0
,
18446744073709551615
]
[0,18446744073709551615]
[0,18446744073709551615]。
如果用了unsigned long long和long long也溢出的话，只能用高精度计算题。
快速幂算法的递归实现（一般用不上）：
int quickPow(int a, int b, int n) {
if (b == 1)return a;
if (b % 2 == 0) {//b是偶数
int t = quickPow(a, b / 2, n);
return t * t % n;
}
else {//b是奇数
int t= quickPow(a, b / 2, n);
t = t * t % n;
t = t * a % n;
return t;
}
}
一般更建议使用非递归模式，可以减少函数栈帧的压力。
但无论是哪一种，都要注意用模运算的性质来降低数值。
快速幂算法将复杂运算拆分成能用已知符号表示的运算，配合高精度计算，能处理更大数量级的幂运算和模运算。
OJ题展示
不用高精度计算
题目链接：1616：A 的 B 次方
这题的题意就是求
a
b
(
m
o
d
m
)
a^b\pmod m
ab(modm)。但数据量可能太大，所以每一步都要用long long的同时，还要利用模运算的乘法性质对结果进行限制。
AC参考程序：
#ifndef _CRT_SECURE_NO_WARNINGS
#define _CRT_SECURE_NO_WARNINGS 1
#endif
#include<iostream>
#include<string>
#include<algorithm>
#include<vector>
using namespace std;
typedef long long ll;
ll quickPow(ll a, ll b, ll m) {
ll res = 1LL;//让计算机将1看做long long型数据
a %= m;
while (b) {
if (b % (2LL))
res = res * a % m;
a = a * a % m;
b /= 2LL;
}
return res;
}
void ac() {
ll a, b, m;
cin >> a >> b >> m;
cout << quickPow(a, b, m);
}
int main() {
ac();
return 0;
}
二进制指数的高精度计算
题目链接：1234：2011
这个OJ题要求的 ( 201 1 b ) % 10000 (2011^b)\%10000 (2011b)%10000的指数部分 b b b有200位，c语言没有任何内置类型能存储这么庞大的数据，所以可以用高精度计算将指数转换成二进制。
之后就是快速幂算法：
#ifndef _CRT_SECURE_NO_WARNINGS
#define _CRT_SECURE_NO_WARNINGS 1
#endif
#include<iostream>
#include<vector>
#include<algorithm>
#include<string>
using namespace std;
void div2(string& st) {
int num = 0;
string ans;
for (auto& a : st) {
num = num * 10 + a - '0';
ans = ans + (char)((num / 2) + '0');
num %= 2;
}
while (ans[0] == '0' && ans.size() > 1)
ans.erase(0, 1);
st = ans;
}
int quickPow(int a, string b, int n) {
int res = 1;
for (size_t i = b.size() - 1; i != -1; i--) {
if (b[i] - '0')
res = res * a % n;
a = a * a % n;
}
return res;
}
void ac() {
int k;
string num, bnum;//binary num，二进制
cin >> k;
while (k--) {
cin >> num;
while (num != "0") {//高精度求二进制
if ((num[num.size() - 1] - '0') % 2)
bnum = '1' + bnum;//栈的思想应用
else
bnum = '0' + bnum;
div2(num);
}
cout << quickPow(2011, bnum, 10000) << endl;//快速幂算法
bnum.clear();
}
}
int main() {
ac();
return 0;
}
数学题
快速幂算法更多的还是用来解数学题。数学的难度取决于出题人的心情。
等差数列和等比数列
题目链接：1615：【例 1】序列的第 k 个数
高中数学题知识点复习（公式参考）：
-
等差数列的公差： a n + 1 − a n = d a_{n+1}-a_n=d an+1−an=d， d d d为常数。
-
等差中项： a m + n 2 = a m + a n 2 a_{\frac{m+n}{2}}=\frac{a_m+a_n}{2} a2m+n=2am+an。
-
等差数列 { a n } \{a_n\} {an}的通项公式： a n = a 1 + ( n − 1 ) d = a m + ( n − m ) d a_n=a_1+(n-1)d=a_m+(n-m)d an=a1+(n−1)d=am+(n−m)d， n > m ≥ 1 n>m\geq 1 n>m≥1。
-
等差数列的前 n n n项和公式： S n = ( a 1 + a n ) × n 2 S_n=\frac{(a_1+a_n)\times n}{2} Sn=2(a1+an)×n，通项公式可代入。
-
等比数列的公比： a n + 1 a n = q \frac{a_{n+1}}{a_n}=q anan+1=q， q ≠ 0 q\neq0 q=0。
-
等比中项： G = a n ⋅ a m G=\sqrt{a_n\cdot a_m} G=an⋅am。
-
等比数列 { a n } \{a_n\} {an}的通项公式： a n = a 1 q n − 1 = a m q n − m a_n=a_1q^{n-1}=a_mq_{n-m} an=a1qn−1=amqn−m， n > m ≥ 1 n>m\geq 1 n>m≥1。
-
等比数列的前 n n n项和公式： S n = a 1 ( 1 − q n ) 1 − q = a 1 − a n q 1 − q S_n=\frac{a_1(1-q^n)}{1-q}=\frac{a_1-a_nq}{1-q} Sn=1−qa1(1−qn)=1−qa1−anq， q ≠ 1 q\neq1 q=1。
-
等差数列乘等比数列组成的新数列 { a n ⋅ b n } = ( A n + B ) q n − 1 \{a_n\cdot b_n\}=(An+B)q^{n-1} {an⋅bn}=(An+B)qn−1的前 n n n项和公式：
S n = ( k ⋅ n + m ) ⋅ q n − m S_n=(k\cdot n+m)\cdot q^n-m Sn=(k⋅n+m)⋅qn−m， k = a q − 1 ， m = b − k q − 1 k=\frac{a}{q-1}，m=\frac{b-k}{q-1} k=q−1a，m=q−1b−k，证明用错位相减法，这里省略。
用到了再补充。
参考程序：
#ifndef _CRT_SECURE_NO_WARNINGS
#define _CRT_SECURE_NO_WARNINGS 1
#endif
#include<iostream>
#include<string>
#include<algorithm>
#include<vector>
using namespace std;
/*
http://ybt.ssoier.cn:8088/problem_show.php?pid=1615
这题纯数值。即使是快速幂算法，也要通过模运算的乘法性质
对每一步操作进行极致的限制。
*/
typedef unsigned long long ull;
const ull pp = 200907;
ull quickPow(ull a, ull b, ull n) {
ull res = 1;
a %= n;
while (b) {
if (b % 2)
res = res * a % n;
a = a * a % n;
b /= 2;
}
return res;
}
void ac() {
int T;
cin >> T;
while (T--) {
ull a, b, c;
ull k;
cin >> a >> b >> c >> k;
if (b - a == c - b)
cout << ((a % pp + ((k - 1) * (c - b)) % pp) % pp) << endl;
else
cout << ((a % pp) * quickPow(c / b, k - 1, pp)) % pp << endl;
}
}
int main() {
ac();
return 0;
}
计数原理
题目链接：1618：越狱
根据乘法原理，总的方案数为
m
×
m
×
.
.
.
×
m
=
m
n
m\times m\times ...\times m=m^n
m×m×...×m=mn。
想要没人越狱，则满足条件的方案数为
m
×
(
m
−
1
)
×
(
m
−
1
)
=
m
×
(
m
−
1
)
n
−
1
m\times (m-1)\times (m-1)=m\times (m-1)^{n-1}
m×(m−1)×(m−1)=m×(m−1)n−1。
有人越狱的方案数为 m n − ( m × ( m − 1 ) n − 1 ) m^n - (m\times (m-1)^{n-1}) mn−(m×(m−1)n−1) 。
对这个数学式求模得到结果：
(
m
n
−
(
m
×
(
m
−
1
)
n
−
1
)
)
%
p
=
(
(
m
n
%
p
−
(
(
m
%
p
)
×
(
(
m
−
1
)
n
−
1
%
p
)
)
%
p
+
x
×
p
)
%
p
(m^n - (m\times (m-1)^{n-1}))\%p=((m^n\%p-((m\%p)\times ((m-1)^{n-1} \%p))\%p + x\times p)\%p
(mn−(m×(m−1)n−1))%p=((mn%p−((m%p)×((m−1)n−1%p))%p+x×p)%p（推导请移步其他博主）。
x
∗
p
x*p
x∗p是为了防止表达式出现负数，一般情况
x
x
x取1即可。
参考程序：
#ifndef _CRT_SECURE_NO_WARNINGS
#define _CRT_SECURE_NO_WARNINGS 1
#endif
#include<iostream>
#include<vector>
using namespace std;
typedef unsigned long long ull;
const ull p=100003ull;
ull quickPow(ull a,ull b,ull n){
ull ans=1ull;
while(b){
if(b%2ull)
ans=ans*a%n;
a=a*a%n;
b/=2ull;
}
return ans;
}
void ac1(){
ull n,m;
cin>>m>>n;
cout<<( quickPow(m,n,p)- ( (m%p)*quickPow(m-1,n-1,p) )%p +p)%p;
}
int main() {
ac1();
return 0;
}