HDU3292.No more tricks, Mr Nanguo(pell 方程理论 + 矩阵快速幂)

今天刚刚学完pell 方程理论记录一下了
证明略,给大家推荐两个证明比较好的博客和学习网站
oi-wike-Pell方程
这是一个很好的总结博客添加链接描述

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int N = 2;
const int mod = 8191;
struct Matrix{
    ll a[N][N];
    int n,m;
    Matrix(int b, int c){
        memset(a, 0,sizeof(a));
        n = b;
        m = c;
    }
    Matrix operator*(const Matrix &B)const{
        Matrix C(n,B.m);
        for(int i = 0; i < 2; i++){
            for(int j = 0; j < 2; j++){
                for(int k = 0; k < 2; k++){
                    C.a[i][j] = (C.a[i][j] + B.a[i][k]*a[k][j]%mod)%mod;
                }
            }
        }
        return C;
    }
};
const int maxn = 2e5 + 5;
ll a[maxn];
bool pell(ll &x,ll &y,ll d){
    ll m = (ll)sqrt((double)d);
    //当前的d不能为完全平方数 
    if(m*m == d)
        return false;
    //连分数数位
    int num = 0;
    double sq = sqrt(d);
    a[num++] = m;
    ll b = m;//当前整数部分
    ll c = 1;//连分数最终展开时的分母
    double temp;//连分数展开时的每一项
    do{
        c = (d-b*b)/c;
        temp = (sq+b)/c;
        a[num++] = (int)(floor(temp));
        b = a[num-1]*c-b;
    //当有一位等于整数两倍时结束
    }while(a[num-1] != 2*a[0]);
    ll p = 1,q = 0;
    //转化成p / q形式 
    for(int i = num - 2;i >= 0;i--){
        ll temp = p;
        p = q + p*a[i];
        q = temp;
    }
    if((num - 1) % 2){
    	//连分数长度为奇数时
        x=2*p*p+1;
        y=2*p*q;
    }
    else{//连分数长度为偶数时
        x=p;
        y=q;
    }
    return true;
}
signed main(){
	ll n,k;
	while(cin >> n >> k){
		ll x,y;
		if(pell(x, y, n)){
			Matrix base = Matrix(2, 2);
			Matrix ans = Matrix(2, 2);
			base.a[0][0] = x,base.a[0][1] = n*y%mod;
			base.a[1][0] = y,base.a[1][1] = x;
			ans.a[0][0] = x,ans.a[1][0] = y;
			k--;
			while(k){
				if(k & 1)ans = ans * base;
				base = base * base;
				k >>= 1;
			}
			cout << ans.a[0][0]%mod << '\n';
		}else{
			cout << "No answers can meet such conditions" << '\n';
		}
	}
} 
/*
2 999888
3 1000001
4 8373
*/