HDU 3292 pell方程

来源：互联网发布：print 用法 python 编辑：程序博客网时间：2024/05/18 04:59

题意：求x^2+n*y^2=1 按x排序第k大的解。

题目就是求佩尔方程，用矩阵连乘的方法，当n为完全平方数的时候无解。

如果我们求出Pell方程的最小正整数解后，就可以根据递推式求出所有的解。

则根据上式我们可以构造矩阵，然后就可以快速幂了。

这样就可以求出第k大的解。

HDU3292题就要用到上面的矩阵方法求第k大的解。

#include <iostream>#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>using namespace std;#define M 8191const int MAX=2;typedef struct{    long long m[MAX][MAX];} Matrix;Matrix I= {1,0,           0,1          };Matrix P;Matrix matrixmul(Matrix a,Matrix b) //矩阵乘法{    int i,j,k;    Matrix c;    for (i = 0 ; i < MAX; i++)        for (j = 0; j < MAX; j++)        {            c.m[i][j] = 0;            for (k=0; k<MAX; k++)                c.m[i][j]+=((a.m[i][k]%(M))*(b.m[k][j]%(M)))%(M);            c.m[i][j] %=M;        }    return c;}Matrix quickpow(long long n){    Matrix m = P, b = I;    while (n >= 1)    {        if (n & 1)            b = matrixmul(b,m);        n = n >> 1;        m = matrixmul(m,m);    }    return b;}int main(){    long long n,k,x,y;    while(~scanf("%I64d%I64d",&n,&k))    {        if(n==4||n==9||n==16||n==25)        {            puts("No answers can meet such conditions");            continue;        }        for(y=1;; y++)        {            x=(long long)sqrt(double(y*y*n+1));            if(x*x==y*y*n+1)                break;        }        P.m[0][0]=P.m[1][1]=x%M;        P.m[0][1]=(n*y)%M,P.m[1][0]=y%M;        Matrix ans=quickpow(k-1);        printf("%I64d\n",((ans.m[0][0]*x)%M+(ans.m[0][1]*y)%M)%M);    }    return 0;}