SP8547 MAIN74 - Euclids algorithm revisited 题解

\(\text{Description}\)

记 \(S\) 为 \(a,b\quad(a\ge b)\) 用辗转相除法得到最大公约数时的迭代次数。

给定 \(S\)，求最小的满足条件的 \(a+b\pmod {10^9+7}\)

限制

\(0\le S\le10^{18}\)

\(\text{Solution}\)

初时看到这个题感觉无从下手，又因为 \(S\) 的范围很明显的是一个结论题。

所以我首先打了个暴力打表找找规律，然后就可以发现答案和斐波那契数列有着密不可分的关系

然后我就突然想起来一个性质： 斐波那契数列任意相邻的两项互质

证明如下：

记 \(f_k\) 表示斐波那契数列第 \(k\) 项的值，则

\[\gcd (f_k,f_{k-1}) = \gcd(f_{k-1},f_{k}-f_{k-1}) = \gcd(f_{k-1},f_{k-2})=\cdots=\gcd(f_2,f_1)=1 \]

又发现刚才证明的过程就是辗转相除法。

所以迭代的次数 \(S\) 对应的最小答案就是斐波那契数列的项数加一。

所以我们只要用矩阵乘法加速一下求斐波那契数列就可以了

时间复杂度 \(\Theta(T\log S)\)

但是这里要特判一下 \(S=0\) 和 \(S = 1\) 时的情况。

具体细节见代码。

\(\text{Code}\)

/*If you are full of hope,you will be invincible*/
#include <bits/stdc++.h>
#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/hash_policy.hpp>
#include<ext/pb_ds/priority_queue.hpp>
std::mt19937 hpy(time(nullptr)+(unsigned long long)(new char));
namespace STD{
	constexpr int HL=1<<20;
	char buf[HL],*t1=buf,*t2=buf;
	inline char getc(){return t1==t2&&(t2=(t1=buf)+fread(buf,1,HL,stdin),t1==t2)?EOF:*t1++;}
	template <typename Tp>
	inline void read(Tp &x) {
		x = 0;static int f = 1;char ch = getc();
		for(;!isdigit(ch);ch=getc()) if(ch=='-') f = 0;
		for(;isdigit(ch);ch=getc()) x=(x<<3)+(x<<1)+(ch^48);
		x=f?x:-x;
	}
	template <typename Tp,typename ...Rp> inline void read(Tp &x,Rp&... rst) {read(x),read(rst...);}
	char buff[HL],*T=buff;
	void flush(){fwrite(buff,1,T-buff,stdout);T=buff;}
	inline void putc(char ch){if(T==buff+HL)flush();*T++=ch;}
	template<typename Tp> inline void print(Tp x){
	    if(x<0)putc('-'),x=-x;
	    if(!x) return putc('0'),void();
	    static int st[20],tp;while(x)st[++tp]=x%10,x/=10;
	    while(tp)putc(st[tp]^48),--tp;
	}
	template <typename Tp> inline void print(const char &ch,Tp x) {print(x),putc(ch);}  
	template<typename T,typename ...Args> inline void print(const char &ch,T x,Args ...args) {print(ch,x);print(ch,args...);}
    template<typename Tp> inline Tp max(const Tp &a,const Tp &b){return a>b?a:b;}
    template<typename Tp> inline Tp min(const Tp &a,const Tp &b){return a<b?a:b;}
} using namespace STD;
using ll = long long;
using u64 = unsigned long long;
using u32 = unsigned int;
constexpr int mod = 1e9 + 7;
struct Martix {
    int a[3][3];
    Martix(){memset(a,0,sizeof a);}
    Martix(const int val){a[0][0]=a[1][1]=val;}
    friend Martix operator * (Martix a,Martix b) {
        Martix c;
        for(int k = 1;k <= 2;++k)
            for(int i = 1 ;i <= 2;++i) 
                for(int j = 1;j <= 2;++j)  
                    c.a[i][j] = (c.a[i][j]+1ll*a.a[i][k]*b.a[k][j]%mod)%mod;
        return c;
    }
};
inline Martix ksm(Martix a,u64 b) {
    Martix ans;
    ans.a[1][1] = ans.a[1][2] = 1;
    for(;b;b>>=1,a=a*a) if(b&1) ans=ans*a;
    return ans;
}
inline int Get(const u64 &num) {
    if(num <= 2) return 1;
    Martix Base,ans;
    Base.a[1][1] = Base.a[2][1] = Base.a[1][2] = 1;
    ans = ksm(Base,num-2);
    return ans.a[1][1] % mod;
}
inline int Solve(const u64 &n) {
    if(n == 0) return 0;
    if(n == 1) return 2;
    return Get(n+3);//由于要求和，所以是第 n+3 项
}
int main(int argc,const char **argv) {
    int T;read(T);
    while(T--) {
        u64 n;read(n);
        print('\n',Solve(n));
    }   
    return flush(),0;
}