题目描述

B 君有两个好朋友,他们叫宁宁和冉冉。

有一天,冉冉遇到了一个有趣的题目:输入 b;d;n,求$\lfloor \left ( \frac{b+\sqrt{d}}{2} \right ) ^n \rfloor \mathrm{mod} \ p$

其中p=7528443412579576937。

输入格式

一行三个整数 b;d;n。

输出格式

一行一个数表示模 7528443412579576937 之后的结果。

样例输入

1 5 9

样例输出

76

数据范围

其中 $0<b^2 \le d<(b+1)^2 \le 10^{18}$,并且 b mod 2=1,d mod 4=1

菜鸡的咆哮

不难发现$ \frac{b+\sqrt{d}}{2} ​$与一元二次方程的解$ \frac{-b+\sqrt{\Delta}}{2a} ​$非常相似

我们阔以小心翼翼的写出另外一个解$ \frac{b-\sqrt{d}}{2} $

令$ x_1=\frac{b+\sqrt{d}}{2},x_2=\frac{b+\sqrt{d}}{2}​$

$x_1+x_2=\frac{b+\sqrt{d}}{2}+\frac{b-\sqrt{d}}{2} =b​$

$x_1x_2=\frac{b+\sqrt{d}}{2} \times \frac{b-\sqrt{d}}{2} =\frac{b^2-d}{4}​$

有了$x_1,x_2​$之间的关系,我们尝试对$( \frac{b+\sqrt{d}}{2}) ^n+ (\frac{b-\sqrt{d}}{2}) ^n​$降次:

$x_1^n+x_2^n=(x_1+x_2)(x_1^{n-1}+x_2^{n-1})-x_1x_2(x_1^{n-2}+x_2^{n-2})​$

设$f[i]=(\frac{b+\sqrt d}{2})^i+(\frac{b-\sqrt d}{2})^i​$且一定是整数

可以得到递推式$f[i]=bf[i-1]+\frac{d-b^2}{4}f[i-2]​$且$f[0]=2,f[1]=b​$

然后我们就可以用矩阵快速幂愉快地递推了

但是当我发现那个p竟然大得出奇QAQ

于是我们使用快速乘防止溢出

但要注意在+p是也存在溢出的风险所以要判断一下,小于零才加

不要忘了在最后减去$ (\frac{b-\sqrt{d}}{2}) ^n$但是啷个减哎???

因为$b^2 \le d<(b+1)^2\Rightarrow b \le \sqrt d<b+1$

因为是向下取整,当n为偶数是不影响答案,

但当n为基数且$b\ne\sqrt d$时答案要减一。

代码

#include<bits/stdc++.h>
#define size 5
#define mod 7528443412579576937
#define ll long long
using namespace std;
ll b,d,n;
ll Add(ll x,ll y){
    return x+y>=mod||x+y<0?x+y-mod:x+y;
}
ll check(ll x,ll y){
    return x-y<0?x-y+mod:x-y;
}
ll qkmul(ll x,ll y){
    ll d=(ll)(x*(long double)y/mod+0.5);
    return check(x*y,d*mod);
}
struct matrix{
    int row,col;
    ll v[size][size];
    matrix(){memset(v,0,sizeof(v));}
    matrix operator *(const matrix &b)const{
        matrix c;
        c.row=row;c.col=b.col;
        for(int i=1;i<=c.row;i++)
        for(int j=1;j<=c.col;j++)
        for(int k=1;k<=col;k++){
            c.v[i][j]=Add(c.v[i][j],qkmul(v[i][k],b.v[k][j]));
        }
        return c;
    }
};
ll qkpow(matrix f,matrix base){
        while(n){
            if(n&1) f=f*base;
            base=base*base;
            n>>=1;
        }
        return f.v[1][1];
}
matrix f,base;
void init(){
    f.row=1;f.col=2;
    f.v[1][1]=b;f.v[1][2]=2;
    base.row=base.col=2;
    base.v[1][1]=b;base.v[2][1]=(d-b*b)/4;
    base.v[1][2]=1;
}
int main(){
    cin>>b>>d>>n;
    int xx=n;
    init();
    if(n==0){
        cout<<1<<endl;
        return 0;
    }
    n-=1;
    ll ans=qkpow(f,base);
    if(xx%2==0&&b*b!=d) ans-=1;
    if(ans<0) ans+=mod;
    cout<<ans<<endl;
    return 0;
}