HDU4565 So Easy!

/*
 HDU4565 So Easy!
 http://acm.hdu.edu.cn/showproblem.php?pid=4565
 数论 快速幂 矩阵快速幂
 题意：求[（a+sqrt(b))^n ]%m [ ]代表向上取证
 联想到斐波那契数列，所以第一感觉是矩阵快速幂
 后来发现根本不用这么麻烦，我们认为a+b*sqrt(c)中的a为实部
 b为虚部，就能直接用二元乘法模拟矩阵快速幂了，因此这两种方法
 是等价的。于是乎，我们就解决了精度问题。
 然而即使我们得出了结果中的a+b*sqrt(c)也不能这么算
 因为[b*sqrt(c)]%m不和[b%m*sqrt(c)]%m同余，因此只能另辟蹊径。
 注意到题目中的(a-1)^2< b < a^2 <=> 0<a-sqrt(b)<1
 所以 0<(a+sqrt(b))^n<1 而 (a+sqrt(b))^n与(a-sqrt(b))^n是公轭的
 所以其和只剩下了实部的二倍。因此只需求出实部，将其乘以2就是答案。
 */
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
#include <queue>
#include <iostream>
#include <map>
#include <set>
#define test
using namespace std;
const int Nmax=1005;
//const long long mod=2147483647LL;
//double eps=1e-8;
long long a,b,n,m;
struct Num
{
    long long a;
    long long b;
    long long c;
    Num(long long _a,long long _b,long long _c)
    {
        a=_a;
        b=_b;
        c=_c;
    }
    friend Num operator * (Num x,Num y)
    {
        long long a1=x.a%m,b1=x.b%m,a2=y.a%m,b2=y.b%m,c=x.c;
        long long a=(a1*a2%m+((c*b1*b2)%m))%m;
        long long b=((a1*b2)%m+((b1*a2)%m))%m;
        return Num(a,b,c);
    }
    friend Num qpow(Num base,long long n)
    {
        Num ans(1LL,0LL,base.c);
        //ans.print();
        while(n>0LL)
        {
            if(n&1LL)
                ans=ans*base;
            base=base*base;
            n>>=1;
            //ans.print();
        }
        //printf("%lld %lld %lld
",ans.a,ans.b,ans.c);
        return ans;
    }
    void print()
    {
        printf("a:%lld b:%lld c:%lld
",a,b,c);
    }
};

void work()
{
    Num base(a,1LL,b);
    Num ans=qpow(base,n);
    //ans.print();
    long long a=ans.a,b=ans.b,c=ans.c;
    //long long res=(long long)ceil(a+b*sqrt(c));//不能这么写，因为整数和小数直接不能模
    //while(a<=0)
        //a+=m;
    //while(b<=0)
        //b+=m;
    //ans.print();
    //long long res=(long long)ceil(a+b*sqrt(c));
    //printf("res:%lld
",res);
    //res=res%m;
    //printf("%lld %lld %ll
",a,b);
    //ans.print();
    //if(2LL*a!=ceil(a+b*sqrt(c)))
        //printf("NO
");
    //res=(2LL*a)%m;
    //printf("ans:%lld myans:%lld
",(2LL*a)%m,(long long)ceil(a+b*sqrt(c))%m);
    long long res=2LL*a%m;
    printf("%lld
",res);

}
int main()
{
    #ifdef test
    //freopen("hdu4565.in","r",stdin);
    #endif
    while(scanf("%lld%lld%lld%lld",&a,&b,&n,&m)==4)
        work();
    return 0;
}