NOD 1147 还原连分数

题意：给出非完全平方数，问你第k项连分数还原后的分数，对1e9+7取余。

对于连分数，我们可以表示为：

对于无理数，ai一定是无穷数列，反之，对于有理数，ai一定是有穷数列。

对于上式中的p与q，有递推式：

而对于sqrt(n)来说，ai中的首项为一个单独的整数，除了它后面的都会循环。

对于这个问题，当然是先求出数列ai，求这个直接倒啊倒。

关键是有了这个ai数列，如何求到第k位的结果。这个当然有循环节，然后在一个循环节内是可以先处理出结果。

然后再计算有多少个这样的循环节，这部分就可以快速幂，剩下的部分就直接multi暴力吧。

此过程中还有一个重要的点，就是：

怎么处理的？

在本题由于sqrt(n)的整数部分没有存进a数组，我们知道p1=a1，p2=a1a2+1，所以我们把上述矩阵表示为：

所以这样完全就解决问题了，如果是q，最后面那个矩阵应该是0，1

#include <iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define M 1000000007
using namespace std;
const double eps=1e-8;
const int MAX=2;
typedef struct
{
    long long m[MAX][MAX];
} Matrix;
Matrix I= {1,0,
           0,1
          };
Matrix yl= {0,1,
            1,0
           };
Matrix data,ans,P,Q;
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;
}
long long gcd(long long a,long long b)
{
    return b?gcd(b,a%b):a;
}
int getxun(long long n,long long *a)
{
    double k=sqrt(n*1.0);
    int cnt=0;
    long long q=1,p=(long long)k,tmp=1;
    double first=k-p;
    while(1)
    {
        tmp=q;
        q=n-p*p;
        long long G=gcd(tmp,q);
        q/=G;
        tmp/=G;
        k=(sqrt(1.0*n)+p)*tmp/q;
        a[cnt++]=(long long)k;
        p=abs(p-a[cnt-1]*q);
        if(fabs(k-a[cnt-1]-first)<eps) break;
    }
    return cnt;
}
int main()
{
    int n,k,f;
    long long s[10000];
    long long p,q;
    while(~scanf("%d%d",&n,&k))
    {
        P=I,ans=I,Q=I;
        f=getxun(n,s);
        int m=k/f,yu=k%f;
        if(m)
        {
            for(int i=f-1; i>=0; i--)
            {
                data=yl;
                data.m[0][0]=s[i];
                P=matrixmul(P,data);
            }
            Q=quickpow(m);
        }
        for(int i=yu-1; i>=0; i--)
        {
            data=yl;
            data.m[0][0]=s[i];
            ans=matrixmul(ans,data);
        }
        ans=matrixmul(ans,Q);
        p=ans.m[0][0];
        q=ans.m[0][1];
        long long ans1=(((long long)sqrt(n)*p)%M+q)%M;
        cout<<ans1<<"/"<<p<<endl;
    }
    return 0;
}