ComentOJ-C72-4033 「火鼠的皮衣 -不焦躁的内心-」
给定 $_math_inline$n, a, b, p$math_inline_$ （p不一定是质数），求 $_math_inline$\sum\limits_{i=0}^{\lfloor\frac{n}{2}\rfloor} a^i b^{n-2i} {n\choose 2i}$math_inline_$ 对 $_math_inline$p$math_inline_$ 取模的结果，T组询问。
首先对求和公式进行转化；
$_math_inline$ \begin{equation} \begin{aligned} ans &= \sum\limits_{i=0}^{\lfloor\frac{n}{2}\rfloor} a^i b^{n-2i} {n\choose 2i} \\ &= \sum\limits_{i=0}^{n}(\sqrt{a})^ib^{n-i}{n \choose i}\left[i\%2==0\right] \end{aligned} \end{equation} $math_inline_$考虑把 $_math_inline$\left[i\%2==0\right]$math_inline_$ 改写成 $_math_inline$\frac{1^i + {(-1)}^i}{2}$math_inline_$
$_math_inline$ \begin{equation} \begin{aligned} ans &= \sum\limits_{i=0}^{n}(\sqrt{a})^ib^{n-i}{n \choose i}\left[i\%2==0\right] \\ &= \sum\limits_{i=0}^{n}(\sqrt{a})^ib^{n-i}{n \choose i}\frac{1^i + {(-1)}^i}{2} \\ &= \frac{1}{2}\left(\sum\limits_{i=0}^{n}(\sqrt{a})^ib^{n-i}{n \choose i} + \sum\limits_{i=0}^{n}(-\sqrt{a})^ib^{n-i}{n \choose i}\right) \end{aligned} \end{equation} $math_inline_$根据二项式定理 $_math_inline${(a+b)}^n=\sum\limits^n_{k=0}{n \choose k}x^{n-k}y^k=\sum\limits^n_{k=0}{n \choose k}x^ky^{n-k}$math_inline_$
$_math_inline$ \begin{equation} \begin{aligned} ans &= \frac{1}{2}\left(\left(b+\sqrt{a}\right)^n+\left(b-\sqrt{a}\right)^n\right) \end{aligned} \end{equation} $math_inline_$因为 a,b 固定，可以设 Fx 为当 n=x 时的答案，那么 Fx 一定是满足一个二项递推式，即 $_math_inline$Fx=A\times F_{x-1}+B\times F_{x-2}$math_inline_$ 。因为 a 不一定有二次剩余，所以不能用快速幂来做，p也不一定是质数，所以不能用 BM 算法求出 A,B（有一个求逆元的步骤）。不过可以通过手算 $_math_inline$F_0,F_1,F_2,F_3$math_inline_$ 然后解一个二元一次方程组来算出 A,B。
$_math_inline$ \begin{aligned} F_0 &= 1 \\ F_1 &= b \\ F_2 &= b^2+a \\ F_3 &= b^3+3ab \end{aligned} $math_inline_$或者
考虑特征根方程： $_math_inline$x^2=Ax+B$math_inline_$ ，即 $_math_inline$x^2-Ax-B=0$math_inline_$ ，方程的两个解理应为答案式子中的两个特征根 $_math_inline$b+\sqrt{a},b-\sqrt{a}$math_inline_$ 。根据韦达定理，可知 A 为两根之和 $_math_inline$(b+\sqrt{a})+(b-\sqrt{a})=2b$math_inline_$ ，-B 为两根之积，即 $_math_inline$B=-(b+\sqrt{a})(b-\sqrt{a})=a-b^2$math_inline_$ 。那么求出 $_math_inline$F_0=1, F_1=b$math_inline_$ 之后用矩阵乘法就可以在 $_math_inline$O(log(n))$math_inline_$ 的时间内求出。
矩阵为
$_math_inline$ \left ( \begin{matrix} 1 & b \end{matrix} \right ) \left ( \begin{matrix} 0 & a-b^2 \\ 1 & 2b \end{matrix} \right ) $math_inline_$或
$_math_inline$ \left ( \begin{matrix} b & 1 \end{matrix} \right ) \left ( \begin{matrix} 2b & 1 \\ a-b^2 & 0 \end{matrix} \right ) $math_inline_$Code 1
1/* ***********************************************
2Author : Akvicor
3Created Time : Mon Oct 28 14:54:35 2019
4File Name : D.cpp
5************************************************ */
6
7#include <bits/stdc++.h>
8
9#define FAST_IO ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
10#define endl '\n'
11#define ASB using namespace std; typedef long long ll; namespace AkvicorS {
12#define ASE } int main() { return AkvicorS::sol(); }
13
14ASB
15
16int t;
17ll n, a, b, p;
18
19inline ll mul(ll x, ll y){ return (__int128)x*y%p; }
20
21struct Mat{
22 int r, c;
23 ll num[2][2];
24 Mat(int _r=2, int _c=2){r=_r;c=_c;}
25 Mat operator * (const Mat &obj) const{
26 Mat res(r, obj.c);
27 for(int i = 0; i < res.r; ++i){
28 for(int j = 0; j < res.c; ++j){
29 res.num[i][j] = 0;
30 for(int k = 0; k < c; ++k){
31 res.num[i][j] = (res.num[i][j] + mul(num[i][k], obj.num[k][j]))%p;
32 }
33 }
34 }
35 return res;
36 }
37}f, g;
38
39int sol(){
40 FAST_IO;
41
42 cin >> t;
43 while(t--){
44 cin >> n >> a >> b >> p;
45 a %= p; b %= p;
46 f.r = 1;
47 f.c = g.r = g.c = 2;
48
49 f.num[0][0] = 1; f.num[0][1] = b;
50
51 g.num[0][0] = 0; g.num[0][1] = (a+p-mul(b, b)) % p;
52 g.num[1][0] = 1; g.num[1][1] = 2*b%p;
53 for(; n; n >>= 1, g = g*g) if(n&1) f = f*g;
54 cout << f.num[0][0] << endl;
55 }
56
57 return 0;
58}
59
60ASE
Code 2
1/* ***********************************************
2Author : Akvicor
3Created Time : Mon Oct 28 14:54:35 2019
4File Name : D.cpp
5************************************************ */
6
7#include <bits/stdc++.h>
8
9#define FAST_IO ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
10#define endl '\n'
11#define ASB using namespace std; typedef long long ll; namespace AkvicorS {
12#define ASE } int main() { return AkvicorS::sol(); }
13
14ASB
15
16int t;
17ll n, a, b, p;
18typedef pair<ll, ll> pll;
19
20#define fi first
21#define se second
22
23ll mul(ll x, ll y){ return (__int128)x*y%p; }
24
25pll operator * (pll a, pll b){
26 return make_pair( ( mul(a.fi, b.fi) + mul( mul(a.se, b.se), ::AkvicorS::a) ) % p, ( mul(a.se, b.fi) + mul(a.fi, b.se) ) % p);
27}
28
29pll qpow(pll x, ll y){
30 pll ans = make_pair(1, 0);
31 for(; y; y>>=1, x=x*x) if(y&1) ans = ans * x;
32 return ans;
33}
34
35int sol(){
36 FAST_IO;
37
38 cin >> t;
39 while(t--){
40 cin >> n >> a >> b >> p;
41 cout << qpow(make_pair(b, 1), n).fi << endl;
42 }
43
44 return 0;
45}
46
47ASE
除另有声明外，本博客文章均采用 知识共享 (Creative Commons) 署名 4.0 国际许可协议 进行许可。转载请注明原作者与文章出处。