632
技術社區[雲棲]
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;
}
最後更新：2017-04-03 16:48:47