GYM100526I Interesting Integers
题目大意
定义一种 \(Gabonacci\) 数列:
\[\begin{array}{c}
G_1=a\\
G_2=b\\
G_i=G_{i-1}+G_{i-2}
\end{array}
\]
给定一个正整数 \(n\) ,求最小的 \(a\) , \(b\) 使得 \(n\) 是该数列上的一个数。 \(2\le n\le 10^9\)
题解
看到这道题,首先的想法是可以先与处理出类似于矩阵快速幂一样的东西,使得我们可以在 \(O(1)\) 的时间内算出对于一种 \(a\) 与 \(b\) 在 \(x\) 位之后的数字。因为 \(Fibonacc\) 数列的第 \(50\) 位早已超过了 \(10^9\) ,所以我们大约处理 \(50\) 位左右即可。
作者这里直接利用了矩阵乘法,预处理出前 \(40\) 位的矩阵。而这里的 \(Gabonacci\) 数列和 \(Fibonacc\) 数列的状态转移矩阵是完全一模一样的,应该没人不知道吧,给大家看一下:
\[\begin{bmatrix}
1&1 \\
1&0
\end{bmatrix}
\]
而初始矩阵就易得为:
\[\begin{bmatrix}
b&a
\end{bmatrix}
\]
即:
\[\begin{bmatrix}
G_{i-1}&G_{i-2}
\end{bmatrix}
*
\begin{bmatrix}
1&1 \\
1&0
\end{bmatrix}
=
\begin{bmatrix}
G_i&G_{i-1}
\end{bmatrix}
\]
然后我们考虑每一个预处理好的矩阵,我们如何通过状态转移矩阵和最后的值来求出初始的矩阵呢?我们可以回顾矩阵乘法,发现对于某一个转移矩阵:
\[\begin{bmatrix}
x&y \\
y&z
\end{bmatrix}
\]
\[G_i=a*x+b*y
\]
这不就是求一个不定方程的特殊解吗?我们可以利用扩展欧几里得,计算出这个解。
代码如下:
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int MAXN=45;
struct Matrix
{
int n,m,h[3][3];
Matrix()
{
n=m=0;
memset(h,0,sizeof(h));
}
void Set1(const int x)
{
n=m=x;
for(int i=1;i<=x;++i)
h[i][i]=1;
}
};
Matrix operator * (const Matrix a,const Matrix b)
{
Matrix ans;
ans.n=a.n;
ans.m=b.m;
for(int i=1;i<=a.n;++i)
{
for(int j=1;j<=b.m;++j)
{
for(int k=1;k<=a.m;++k)
ans.h[i][j]+=a.h[i][k]*b.h[k][j];
}
}
return ans;
}
int t,n;
Matrix stp[MAXN],tmp,ans;
int l,r,mid,tag;
int gcd(int a,int b)
{
if(b==0)
return a;
return gcd(b,a%b);
}
int lcm(int a,int b)
{
return a*b/gcd(a,b);
}
void exgcd(int a,int b,int &x,int &y,int z)
{
if(b==0)
{
x=z;
y=0;
return ;
}
exgcd(b,a%b,x,y,z);
int tmp=x;
x=y;
y=tmp-a/b*y;
return ;
}
signed main()
{
stp[1].n=stp[1].m=2;
stp[1].h[1][1]=stp[1].h[1][2]=stp[1].h[2][1]=1;
for(int i=2;i<=MAXN-5;++i)
stp[i]=stp[i-1]*stp[1];
tmp.n=1;
tmp.m=2;
cin>>t;
while(t--)
{
cin>>n;
for(int i=MAXN-5;i>=1;--i)
{
int a=stp[i].h[1][1],b=stp[i].h[2][1];
int x,y,z=gcd(a,b);
if(n%z!=0)
continue;
exgcd(a,b,x,y,z);
x*=n/z;
y*=n/z;
int addx=lcm(a,b)/a;
int addy=lcm(a,b)/b;
int k=abs(x-y)/(addx+addy);
if(x<y)
{
if(x+addx*k<y-addy*k)
++k;
x+=addx*k;
y-=addy*k;
}
else
{
if(x-addx*k<y+addy*k)
--k;
x-=addx*k;
y+=addy*k;
}
if(x<=0||y<=0)
continue;
printf("%lld %lld\n",y,x);
break;
}
}
return 0;
}