zoj3759（待解决+算法木有问题+but需要java大数）

  1 //zoj3759
//http://www.bingfengsa.com/info/18001.html
//学习笔记(仅仅是填充学习博文中的知识点而以)
/*
【题目要求】：
找到：x(x+1)/2=ky^2      -x^2-x+2ky^2=0          (2x+1)^2-8k*y^2=1
      x^2=ky^2           x^2-ky^2=0;              特殊
      x(3x-1)/2=ky^2     3x^2-x-2ky^2=0         （6x-1)^2-24ky^2=1
      x(2x-1)=ky^2       2x^2-x-2ky^2=0          (4x-1)^2-8y^2=1
      x=ky^2             0x^2-x+ky^2=0;          特殊
右边是统一化成的z^2-ny^2=0的形式,z代换了x

【换元】：z^2-ny^2=1
【佩尔方程】：
显然这个方程有正整数解，则x，y才可能有整数解，
///////佩尔方程的知识/////////
一、下面的不定方程称为佩尔（Pell）方程：
x^2 - ny^2= 1  设n是正整数
如果d是完全平方数，那么只有平凡解(-1,0),(1,0)
如果不是，一定有无穷多组正整数解
二、递归解：
x[i+1]=x[1]*x[i]+n*y[1]*y[i]
y[i+1]=x[1]*y[i]+y[1]*x[i]
三、基本解(x1,y1)：
首先根据根号n的渐进连分数（h/k)表示，找出前几项，察看x=h,y=k带入后是否是一组解。
是，则得到了(x1,y1)
四、渐进连分数的生成（如截图）：
五、佩尔方程的解是指数增长的（很少）
六、佩尔方程最小解的生成
from : http://blog.csdn.net/accelerator_916852/article/details/20411727
bool pell( int D, int& x, int& y ) {//D就是上文的N
    int sqrtD = sqrt(D + 0.0);//这里应该处理精度才行
    if( sqrtD * sqrtD == D ) return false;
    int c = sqrtD, q = D - c * c, a = (c + sqrtD) / q;
    int step = 0;
    int X[] = { 1, sqrtD };
    int Y[] = { 0, 1 };
    while( true ) {
        X[step] = a * X[step^1] + X[step];
        Y[step] = a * Y[step^1] + Y[step];
        c = a * q - c;
        q = (D - c * c) / q;
        a = (c + sqrtD) / q;
        step ^= 1;
        if( c == sqrtD && q == 1 && step ) {
            x = X[0], y = Y[0];
            return true;
        }
    }
}
///////佩尔方程的知识/////////
【继续】 现在我们找到了z1，y1，能得到所有的z[i],y[i],现在我们通过这些解，找到最小的x、y正整数解就可以了
【特殊方程的解决】：(题目没要求计算)
 1、x^2-ky^2=0;k判断k是否是完全平方数即可
 2、0x^2-x+ky^2=0;直接解出x即可
 */

#include<iostream>
#include<string.h>
#include<stdio.h>
#include<stdlib.h>
#include<cmath>
#include<algorithm>
#include<queue>
#include<stack>
#include<set>
#include<map>
#define LL unsigned long long
#define maxn 510
#define INF 99999
using namespace std;

int k,kind;

bool pell( int D, LL& x, LL& y ) {//D就是上文的N
    LL sqrtD = (LL)sqrt(D + 0.0);//这里应该处理精度才行
    if( sqrtD * sqrtD == D || (sqrtD+1) * (sqrtD+1)==D ) return false;
    LL c = sqrtD, q = D - c * c, a = (c + sqrtD) / q;
    LL step = 0;
    LL X[] = { 1, sqrtD };
    LL Y[] = { 0, 1 };
    while( true ) {
        X[step] = a * X[step^1] + X[step];
        Y[step] = a * Y[step^1] + Y[step];
        c = a * q - c;
        q = (D - c * c) / q;
        a = (c + sqrtD) / q;
        step ^= 1;
        if( c == sqrtD && q == 1 && step ) {
            x = X[0], y = Y[0];
            return true;
        }
    }
}

void solve(int N,int a,int b)//写成x=(z+a)/b的形式
{
    LL z1,y1;
//    int zi,yi;
    LL ansx,ansy;
    bool ok=pell(N,z1,y1);
    if (!ok)
    {
        printf("Impossible!
");
        return ;
    }
    if ((z1+a)%b==0) {ansx=(z1+a)/b,ansy=y1;cout<<ansx<<" "<<ansy<<endl;}
    else printf("Impossible!
");
}
int main()
{
    while(cin>>kind)
    {
        if (kind==0) break;
        cin>>k;
        if (kind==3) solve(8*k,-1,2);
        else if (kind==5) solve(24*k,1,6);
        else if (kind==6) solve(8*k,1,4);
    }
    return 0;
}