HDU4565 So Easy!

  1 /*
  2  HDU4565 So Easy!
  3  http://acm.hdu.edu.cn/showproblem.php?pid=4565
  4  数论 快速幂 矩阵快速幂
  5  题意：求[（a+sqrt(b))^n ]%m [ ]代表向上取证
  6  联想到斐波那契数列，所以第一感觉是矩阵快速幂
  7  后来发现根本不用这么麻烦，我们认为a+b*sqrt(c)中的a为实部
  8  b为虚部，就能直接用二元乘法模拟矩阵快速幂了，因此这两种方法
  9  是等价的。于是乎，我们就解决了精度问题。
 10  然而即使我们得出了结果中的a+b*sqrt(c)也不能这么算
 11  因为[b*sqrt(c)]%m不和[b%m*sqrt(c)]%m同余，因此只能另辟蹊径。
 12  注意到题目中的(a-1)^2< b < a^2 <=> 0<a-sqrt(b)<1
 13  所以 0<(a+sqrt(b))^n<1 而 (a+sqrt(b))^n与(a-sqrt(b))^n是公轭的
 14  所以其和只剩下了实部的二倍。因此只需求出实部，将其乘以2就是答案。
 15  */
 16 #include <cstdio>
 17 #include <algorithm>
 18 #include <cstring>
 19 #include <cmath>
 20 #include <vector>
 21 #include <queue>
 22 #include <iostream>
 23 #include <map>
 24 #include <set>
 25 #define test
 26 using namespace std;
 27 const int Nmax=1005;
 28 //const long long mod=2147483647LL;
 29 //double eps=1e-8;
 30 long long a,b,n,m;
 31 struct Num
 32 {
 33     long long a;
 34     long long b;
 35     long long c;
 36     Num(long long _a,long long _b,long long _c)
 37     {
 38         a=_a;
 39         b=_b;
 40         c=_c;
 41     }
 42     friend Num operator * (Num x,Num y)
 43     {
 44         long long a1=x.a%m,b1=x.b%m,a2=y.a%m,b2=y.b%m,c=x.c;
 45         long long a=(a1*a2%m+((c*b1*b2)%m))%m;
 46         long long b=((a1*b2)%m+((b1*a2)%m))%m;
 47         return Num(a,b,c);
 48     }
 49     friend Num qpow(Num base,long long n)
 50     {
 51         Num ans(1LL,0LL,base.c);
 52         //ans.print();
 53         while(n>0LL)
 54         {
 55             if(n&1LL)
 56                 ans=ans*base;
 57             base=base*base;
 58             n>>=1;
 59             //ans.print();
 60         }
 61         //printf("%lld %lld %lld\n",ans.a,ans.b,ans.c);
 62         return ans;
 63     }
 64     void print()
 65     {
 66         printf("a:%lld b:%lld c:%lld\n",a,b,c);
 67     }
 68 };
 69
 70 void work()
 71 {
 72     Num base(a,1LL,b);
 73     Num ans=qpow(base,n);
 74     //ans.print();
 75     long long a=ans.a,b=ans.b,c=ans.c;
 76     //long long res=(long long)ceil(a+b*sqrt(c));//不能这么写，因为整数和小数直接不能模
 77     //while(a<=0)
 78         //a+=m;
 79     //while(b<=0)
 80         //b+=m;
 81     //ans.print();
 82     //long long res=(long long)ceil(a+b*sqrt(c));
 83     //printf("res:%lld\n",res);
 84     //res=res%m;
 85     //printf("%lld %lld %ll\n",a,b);
 86     //ans.print();
 87     //if(2LL*a!=ceil(a+b*sqrt(c)))
 88         //printf("NO\n");
 89     //res=(2LL*a)%m;
 90     //printf("ans:%lld myans:%lld\n",(2LL*a)%m,(long long)ceil(a+b*sqrt(c))%m);
 91     long long res=2LL*a%m;
 92     printf("%lld\n",res);
 93
 94 }
 95 int main()
 96 {
 97     #ifdef test
 98     //freopen("hdu4565.in","r",stdin);
 99     #endif
100     while(scanf("%lld%lld%lld%lld",&a,&b,&n,&m)==4)
101         work();
102     return 0;
103 }