HDU4565 So Easy! —— 共轭构造、二阶递推数列、矩阵快速幂

转载

mb5fdb1266ce6df 2018-02-07 14:33:00

文章标签 #include 浮点数 i++ 递推 .net 文章分类 后端开发

题目链接：https://vjudge.net/problem/HDU-4565

So Easy!

Time Limit: 2000/1000 MS (Java/Others) Memory Limit: 32768/32768 K (Java/Others)
Total Submission(s): 5525 Accepted Submission(s): 1841

Problem Description

　　A sequence S_n is defined as:

HDU4565 So Easy! —— 共轭构造、二阶递推数列、矩阵快速幂_i++

Where a, b, n, m are positive integers.┌x┐is the ceil of x. For example, ┌3.14┐=4. You are to calculate S_n.

　　You, a top coder, say: So easy!

Input

　　There are several test cases, each test case in one line contains four positive integers: a, b, n, m. Where 0< a, m < 2¹⁵, (a-1)²< b < a², 0 < b, n < 2³¹.The input will finish with the end of file.

Output

　　For each the case, output an integer S_n.

Sample Input

2 3 1 20132 3 2 2013 2 2 1 2013

Sample Output

414 4

Source

2013 ACM-ICPC长沙赛区全国邀请赛——题目重现

Recommend

zhoujiaqi2010

题解：

1.因为：0< a，(a-1)²< b < a²,。当a>=1时， a-1<根号b<a，那么 0<(a-根号b)<1；当0<a<1时， 1-a<根号b<a，那么那么 0<(a-根号b)<2*a-1<1。综上：0<(a-根号b)<1。所以0<(a-根号b)^n<1

2.这里假设b = 根号b，以方便描述。根据上述结论，那么可得：[(a+b)^n] = [(a+b)^n + (a-b)^n] = (a+b)^n + (a-b)^n 。

解释：

2.1 因为a-1<b<1，所以b必定是浮点数，那么(a+b)^n 也必定是浮点数，此时，再加上个大于0小于1的浮点数(a-b)^n，那么 [(a+b)^n + (a-b)^n] 有可能等于[(a+b)^n] ，也有可能等于[(a+b)^n] +1，这就要取决(a+b)^n的小数部分与(a-b)^n的小数部分之和是否大于1。

2.2 此时，就要将(a+b)^n+(a-b)^n展开进行分析。展开后可知，当b的指数为奇数时，正负抵消；当b的指数为偶数时，b^2k 为一个整数。综上：(a+b)^n+(a-b)^n为一个整数，即表明(a+b)^n的小数部分与(a-b)^n的小数部分之和刚好等于1，所以[(a+b)^n] = [(a+b)^n + (a-b)^n] = (a+b)^n + (a-b)^n 。

3.根据上述分析，问题转化为求：S[n] = (a+b)^n + (a-b)^n 。而这个式子可以看成是二阶齐次递推式的通项公式，可以根据通项公式反推回递推式，然后构造矩阵进行求解。具体如下：

HDU4565 So Easy! —— 共轭构造、二阶递推数列、矩阵快速幂_.net_02

HDU4565 So Easy! —— 共轭构造、二阶递推数列、矩阵快速幂_浮点数_03

以上来自：javascript:void(0)

代码如下：

HDU4565 So Easy! —— 共轭构造、二阶递推数列、矩阵快速幂_#include_04

1 #include <iostream>
 2 #include <cstdio>
 3 #include <cstring>
 4 #include <algorithm>
 5 #include <vector>
 6 #include <cmath>
 7 #include <queue>
 8 #include <stack>
 9 #include <map>
10 #include <string>
11 #include <set>
12 using namespace std;
13 typedef long long LL;
14 const int INF = 2e9;
15 const LL LNF = 9e18;
16 //const int MOD = 1e9+7;
17 const int MAXN = 1e6+100;
18 
19 int MOD;
20 const int Size = 2;
21 struct MA
22 {
23     LL mat[Size][Size];
24     void init()
25     {
26         for(int i = 0; i<Size; i++)
27         for(int j = 0; j<Size; j++)
28             mat[i][j] = (i==j);
29     }
30 };
31 
32 MA mul(MA x, MA y)
33 {
34     MA ret;
35     memset(ret.mat, 0, sizeof(ret.mat));
36     for(int i = 0; i<Size; i++)
37     for(int j = 0; j<Size; j++)
38     for(int k = 0; k<Size; k++)
39         ret.mat[i][j] += 1LL*x.mat[i][k]*y.mat[k][j]%MOD, ret.mat[i][j] = (ret.mat[i][j]%MOD+MOD)%MOD;
40     return ret;
41 }
42 
43 MA qpow(MA x, LL y)
44 {
45     MA s;
46     s.init();
47     while(y)
48     {
49         if(y&1) s = mul(s, x);
50         x = mul(x, x);
51         y >>= 1;
52     }
53     return s;
54 }
55 
56 
57 int main()
58 {
59     LL a, b, n, m, f[2];
60     while(scanf("%lld%lld%lld%lld", &a,&b,&n,&m)!=EOF)
61     {
62         MOD = (int)m;
63         f[0] = 2; f[1] = 2*a;
64         if(n<=1)
65         {
66             printf("%lld\n", f[n]%MOD);
67             continue;
68         }
69 
70         MA s;
71         memset(s.mat, 0, sizeof(s.mat));
72         s.mat[0][0] = 2*a; s.mat[0][1] = b-a*a;
73         s.mat[1][0] = 1; s.mat[1][1] = 0;
74 
75         s = qpow(s, n-1);
76         LL ans = (1LL*s.mat[0][0]*f[1]%MOD+1LL*s.mat[0][1]*f[0]%MOD+2*MOD)%MOD;
77         printf("%lld\n", ans);
78     }
79 }