Solution：SP27519 Z124H - Zeros of the fundamental Fibonacci period

Argon_Cube · 2025-06-24 20:09:59 · 题解

请先阅读 Easy Version SP27484 Z124 - Zeros in Fibonacci period 的题解。

垃圾高精卡常题/tuu/tuu 下面看起来写的很少但是每一步卡常都要试好久，卡了两天终于过了。

和 SP27484 的做法是完全一样的，但是我们发现 p<10^{100}，现在我们需要高精了。

然后我们发现这个时候因为高精时间复杂度是 \mathrm O(T\log^3 p)，这能过？？？

首先我们要先整一个压位高精，这样常数就可以小到消掉一个 \log p。自己手写的十进制压 8 位高精应该就别想过了，我们需要使用二进制压 64 位的 GNU MP 高精库（参考 https://zhuanlan.zhihu.com/p/389419581 ）才有希望能卡过去。

然后我们直接把板子粘过来发现本地跑了 2.5\rm s。交上去 4.4\rm s 的时限 TLE 了？？？？

这个时候你发现直接高精是过不去 SP27484 的（这个题前两组数据和时限和 SP27484 一样），必须要特判 p<10^{18} 不用高精。交上去还是 TLE？？？

这个时候发现 SPOJ 是不开 O2 的，于是我们直接加一个 #pragma 手动开一下 O2，虽然还是过不去但是通过一些测试发现此时我们在 SPOJ 上能跑出来大约 8000 组数据（总共 10^4 组）。说不定再卡卡就过了呢！

pairiinf fib(intinf n,const intinf& p)
{
    if(n==0)
        return (pairiinf){0,1};
    if(n%2==1)
    {
        auto c=fib(n-1,p);
        return (pairiinf){(c.a+c.b)%p,c.a};
    }
    auto c=fib(n/2,p);
    return (pairiinf){c.a*(c.b*2+c.a)%p,(c.a*c.a+c.b*c.b)%p};
}

这个是我们一开始算 F_n\bmod p 的代码。其中 intinf 是一个封装好的 GMP 高精类。我们通过反复的删代码发现，本地的 2.5\rm s 里面竟然有几百毫秒在 n%2 n-1 n/2 上面？？这个时候我们阅读一下 gmp.h 的源码，发现 GMP 内部存储是按二进制位压 64 位从低到高存的。我们发现根本不需要对 n 操作，直接对 n 从高到低扫每一位并动态维护 a,b 就好了！于是我们就得到了：

pairiinf fib(const intinf& n,const intinf& p)
{
    intinf a=0,b=1;
    const auto& m=n.n[0];
    ford(i,m._mp_size-1,0)
        ford(j,63,0)
        {
            auto c=a*(a+b+b),d=a*a+b*b;
            a=c%p,b=d%p;
            if(m._mp_d[i]&1ull<<j)
                (b+=a)%=p,swap(a,b);
        }
    return {a%p,b%p};
}

然后还是过不去？？这个时候就开始玄学了，经过尝试发现乘法次数减少了好像耗时更长了，那我们减少一点取模，当 b+=a 时不取模，最后再把 a 取一下模，竟然真的变快了？

但是还是过不去，我们看看能不能减少一点赋值操作，然后经过许多次尝试发现下面的写法是本地唯一能卡进 2\rm s 的：

pairiinf fib(const intinf& n,const intinf& p)
{
    intinf a=0,b=1;
    const auto& m=n.n[0];
    ford(i,m._mp_size-1,0)
        ford(j,63,0)
        {
            auto d=a*a+b*b;
            (a*=(a+b+b))%=p,b=d%p;
            if(m._mp_d[i]&1ull<<j)
                b+=a,swap(a,b);
        }
    return {a%p,b};
}

竟然还是过不去？？那我把最后的 a 取模删了最后判 a==0||a==p 总行了吧？还是过不去？？？

无语了什么垃圾卡常题，这都不给我过吗？本地都能跑到 1.8\rm s 了没理由过不去吧？再交两发看看？过了？？？？？

#pragma GCC optimize(2)

#include <algorithm>
#include <iostream>
#include <cassert>
#include <cstring>
#include <vector>
#include <bitset>
#include <string>
#include <array>
#include <gmp.h>

typedef int i32;
typedef unsigned int u32;
typedef long long i64;
typedef unsigned long long u64;
typedef __int128_t i128;
typedef __uint128_t u128;

struct fastIO
{
    static const int BUFF_SZ = 1 << 18;
    char inbuf[BUFF_SZ], outbuf[BUFF_SZ];
    fastIO()
    {
        setvbuf(stdin, inbuf, _IOFBF, BUFF_SZ);
        setvbuf(stdout, outbuf, _IOFBF, BUFF_SZ);
    }
} IO;

u128 rd_from_str(const char *s)
{
    int len = strlen(s);
    u128 k = 0;
    for (int i = 0; i < len; ++i)
        k = (k << 1) + (k << 3) + (s[i] ^ 48);
    return k;
}

void wr(u128 x)
{
    if (x > 9)
        wr(x / 10);
    putchar(x % 10 + '0');
}

// area 2 : GMP Library
#define ONLINE_JUDGE

#ifdef ONLINE_JUDGE

#ifdef __cplusplus
extern "C"
{
#endif
    // SPOJ, codechef, atcoder, acwing, openjudge, ideone (64-bit), etc.
    const char *gmp_path = "/usr/lib/x86_64-linux-gnu/libgmp.so.10";
    void *__libc_dlopen_mode(const char *x, int y);
    void *__libc_dlsym(void *x, const char *y);
#define DLOPEN __libc_dlopen_mode
#define DLSYM __libc_dlsym

#ifdef __cplusplus
}
#endif

void *p = DLOPEN(gmp_path, 2);
void *get(const char *name) { return DLSYM(p, name); }
#define D(name) const auto my_##name = (decltype(name) *)get("__g" #name)
D(mpz_init);
D(mpz_init_set_ui);
D(mpz_init_set_si);
D(mpz_init_set);
D(mpz_init_set_str);
D(mpz_clear);
D(mpz_set);
D(mpz_set_ui);
D(mpz_set_si);
D(mpz_swap);
D(mpz_set_str);
D(mpz_get_str);
D(mpz_inp_str);
D(mpz_out_str);

D(mpz_neg);
D(mpz_abs);
D(mpz_add);
D(mpz_add_ui);
D(mpz_sub);
D(mpz_sub_ui);
D(mpz_mul);
D(mpz_mul_ui);
D(mpz_mul_si);
D(mpz_mul_2exp);
D(mpz_tdiv_r_ui);
D(mpz_tdiv_r);
D(mpz_tdiv_q_ui);
D(mpz_tdiv_q);
D(mpz_tdiv_qr);
D(mpz_tdiv_q_2exp);
D(mpz_divexact);
D(mpz_divexact_ui);
D(mpz_pow_ui);
D(mpz_sqrt);

D(mpz_xor);
D(mpz_and);
D(mpz_ior);
D(mpz_cmp);
D(mpz_sizeinbase);
D(mpz_popcount);
D(mpz_gcd);

D(mpz_probab_prime_p);
D(mpz_nextprime);

D(mpz_fac_ui);
D(mpz_root);

D(mpz_scan0);
D(mpz_scan1);
#undef D

#else

#define my_mpz_init mpz_init
#define my_mpz_init_set_ui mpz_init_set_ui
#define my_mpz_init_set_si mpz_init_set_si
#define my_mpz_init_set mpz_init_set
#define my_mpz_init_set_str mpz_init_set_str
#define my_mpz_clear mpz_clear
#define my_mpz_set mpz_set
#define my_mpz_set_ui mpz_set_ui
#define my_mpz_set_si mpz_set_si
#define my_mpz_swap mpz_swap
#define my_mpz_set_str mpz_set_str
#define my_mpz_get_str mpz_get_str
#define my_mpz_inp_str mpz_inp_str
#define my_mpz_out_str mpz_out_str
#define my_mpz_neg mpz_neg
#define my_mpz_abs mpz_abs
#define my_mpz_add mpz_add
#define my_mpz_add_ui mpz_add_ui
#define my_mpz_sub mpz_sub
#define my_mpz_sub_ui mpz_sub_ui
#define my_mpz_sub_si mpz_sub_si
#define my_mpz_mul mpz_mul
#define my_mpz_mul_ui mpz_mul_ui
#define my_mpz_mul_si mpz_mul_si
#define my_mpz_mul_2exp mpz_mul_2exp
#define my_mpz_tdiv_r_ui mpz_tdiv_r_ui
#define my_mpz_tdiv_r mpz_tdiv_r
#define my_mpz_tdiv_q_ui mpz_tdiv_q_ui
#define my_mpz_tdiv_q mpz_tdiv_q
#define my_mpz_tdiv_qr mpz_tdiv_qr
#define my_mpz_tdiv_q_2exp mpz_tdiv_q_2exp
#define my_mpz_divexact mpz_divexact
#define my_mpz_divexact_ui mpz_divexact_ui
#define my_mpz_pow_ui mpz_pow_ui
#define my_mpz_sqrt mpz_sqrt
#define my_mpz_xor mpz_xor
#define my_mpz_and mpz_and
#define my_mpz_ior mpz_ior
#define my_mpz_cmp mpz_cmp
#define my_mpz_sizeinbase mpz_sizeinbase
#define my_mpz_popcount mpz_popcount
#define my_mpz_gcd mpz_gcd
#define my_mpz_probab_prime_p mpz_probab_prime_p
#define my_mpz_nextprime mpz_nextprime
#define my_mpz_fac_ui mpz_fac_ui
#define my_mpz_root
#define my_mpz_scan0 mpz_scan0
#define my_mpz_scan1 mpz_scan1
#endif

// radix : 10
struct BigInteger_GMP
{
    mpz_t n;

    // 1.consturct and set

    BigInteger_GMP() { my_mpz_init(n); }
    BigInteger_GMP(long x) { my_mpz_init_set_si(n, x); }
    BigInteger_GMP(std::string s) { my_mpz_init_set_str(n, s.c_str(), 10); }
    BigInteger_GMP(const BigInteger_GMP &o) { my_mpz_init_set(n, o.n); }
    ~BigInteger_GMP() { my_mpz_clear(n); }

    void set(const BigInteger_GMP &o) { my_mpz_set(n, o.n); }
    void set_si(long x) { my_mpz_set_si(n, x); }

    // 2.input and output

    // input from fp, and set gmp
    void input(FILE *fp = stdin, int base = 10) { my_mpz_inp_str(n, fp, base); }
    // set gmp by str
    void set_str(char *str, int base = 10) { my_mpz_set_str(n, str, base); }
    // output gmp
    void print(FILE *fp = stdout) const { my_mpz_out_str(fp, 10, n); }
    void println(FILE *fp = stdout) const { my_mpz_out_str(fp, 10, n), fputc('\n', fp); }
    // get str from gmp
    void get_str(char *str, int base = 10) const { my_mpz_get_str(str, base, n); }

    u128 output_u128() const
    {
        char s[60];
        memset(s, 0, sizeof(s));
        get_str(s, 10);
        return rd_from_str(s);
    }

    // 3.arithmetics

    void operator=(const BigInteger_GMP &o) { my_mpz_set(n, o.n); }
    bool operator==(const BigInteger_GMP &o) const { return my_mpz_cmp(n, o.n) == 0; }
    bool operator!=(const BigInteger_GMP &o) const { return my_mpz_cmp(n, o.n) != 0; }
    bool operator<(const BigInteger_GMP &o) const { return my_mpz_cmp(n, o.n) < 0; }
    bool operator>(const BigInteger_GMP &o) const { return my_mpz_cmp(n, o.n) > 0; }
    bool operator<=(const BigInteger_GMP &o) const { return my_mpz_cmp(n, o.n) <= 0; }
    bool operator>=(const BigInteger_GMP &o) const { return my_mpz_cmp(n, o.n) >= 0; }

    friend BigInteger_GMP abs(const BigInteger_GMP &o)
    {
        BigInteger_GMP ret;
        my_mpz_abs(ret.n, o.n);
        return ret;
    }

    friend BigInteger_GMP neg(const BigInteger_GMP &o)
    {
        BigInteger_GMP ret;
        my_mpz_neg(ret.n, o.n);
        return ret;
    }

    BigInteger_GMP operator+(const BigInteger_GMP &o) const
    {
        BigInteger_GMP ret;
        my_mpz_add(ret.n, n, o.n);
        return ret;
    }

    BigInteger_GMP operator+(const unsigned long &o) const
    {
        BigInteger_GMP ret;
        my_mpz_add_ui(ret.n, n, o);
        return ret;
    }

    BigInteger_GMP operator-(const BigInteger_GMP &o) const
    {
        BigInteger_GMP ret;
        my_mpz_sub(ret.n, n, o.n);
        return ret;
    }

    BigInteger_GMP operator-(const unsigned long &o) const
    {
        BigInteger_GMP ret;
        my_mpz_sub_ui(ret.n, n, o);
        return ret;
    }

    BigInteger_GMP operator*(const long &o) const
    {
        BigInteger_GMP ret;
        my_mpz_mul_si(ret.n, n, o);
        return ret;
    }

    BigInteger_GMP operator*(const BigInteger_GMP &o) const
    {
        BigInteger_GMP ret;
        my_mpz_mul(ret.n, n, o.n);
        return ret;
    }

    BigInteger_GMP operator/(const unsigned long &o) const
    {
        BigInteger_GMP ret;
        my_mpz_tdiv_q_ui(ret.n, n, o);
        return ret;
    }

    BigInteger_GMP operator/(const BigInteger_GMP &o) const
    {
        BigInteger_GMP ret;
        my_mpz_tdiv_q(ret.n, n, o.n);
        return ret;
    }
    // mpz_divexact can only user when a is divided by b
    BigInteger_GMP divide_exact(const BigInteger_GMP& o) const
    {
        BigInteger_GMP ret;
        my_mpz_divexact(ret.n, n, o.n);
        return ret;
    }

    BigInteger_GMP operator%(const unsigned long &o) const
    {
        BigInteger_GMP ret;
        my_mpz_tdiv_r_ui(ret.n, n, o);
        return ret;
    }

    BigInteger_GMP operator%(const BigInteger_GMP &o) const
    {
        BigInteger_GMP ret;
        my_mpz_tdiv_r(ret.n, n, o.n);
        return ret;
    }

    void div_rem(const BigInteger_GMP &o, BigInteger_GMP &ret_div, BigInteger_GMP &ret_rem)
    {
        my_mpz_tdiv_qr(ret_div.n, ret_rem.n, n, o.n);
    }

    BigInteger_GMP operator<<(unsigned long x) const
    {
        BigInteger_GMP ret;
        my_mpz_mul_2exp(ret.n, n, x);
        return ret;
    }

    BigInteger_GMP operator>>(unsigned long x) const
    {
        BigInteger_GMP ret;
        my_mpz_tdiv_q_2exp(ret.n, n, x);
        return ret;
    }

    BigInteger_GMP operator^(unsigned long p) const
    {
        BigInteger_GMP ret;
        my_mpz_pow_ui(ret.n, n, p);
        return ret;
    }

    BigInteger_GMP &operator+=(const BigInteger_GMP &o) { return (*this) = (*this) + o, (*this); }
    BigInteger_GMP &operator-=(const BigInteger_GMP &o) { return (*this) = (*this) - o, (*this); }
    BigInteger_GMP &operator*=(const BigInteger_GMP &o) { return (*this) = (*this) * o, (*this); }
    BigInteger_GMP &operator/=(const BigInteger_GMP &o) { return (*this) = (*this) / o, (*this); }
    BigInteger_GMP &operator%=(const BigInteger_GMP &o) { return (*this) = (*this) % o, (*this); }
    BigInteger_GMP &operator<<=(unsigned long x) { return (*this) = (*this) << x, (*this); }
    BigInteger_GMP &operator>>=(unsigned long x) { return (*this) = (*this) >> x, (*this); }
    BigInteger_GMP &operator^=(unsigned long p) { return (*this) = (*this) ^ p, (*this); }

    size_t log(int k) const { return ((*this) == BigInteger_GMP(1)) ? 0 : my_mpz_sizeinbase(n, k); }

    BigInteger_GMP root(int k) const
    {
        BigInteger_GMP ret;
        my_mpz_root(ret.n, n, k);
        return ret;
    }
    friend BigInteger_GMP gcd(BigInteger_GMP a, BigInteger_GMP b)
    {
        BigInteger_GMP ret;
        my_mpz_gcd(ret.n, a.n, b.n);
        return ret;
    }

    bool is_prime() const { return my_mpz_probab_prime_p(n, 25); }
    BigInteger_GMP next_prime() const
    {
        BigInteger_GMP ret;
        my_mpz_nextprime(ret.n, n);
        return ret;
    }
    void factorial(unsigned long x) { my_mpz_fac_ui(n, x); }
    // scan from low to high
    unsigned long scan_first_r_zero(unsigned long start = 0) const { return my_mpz_scan0(n, start); }
    unsigned long scan_first_r_one(unsigned long start = 0) const { return my_mpz_scan1(n, start); }
};

#define rgall(arr)          (arr).begin(),(arr).end()
#define rgo1(arr,cnt)       (arr).begin()+1,(arr).begin()+1+(cnt)
#define rgcnt(arr,cnt)      (arr).begin(),(arr).begin()+(cnt)
#define rgany(arr,rgl,rgr)  (arr).begin()+(rgl),(arr).begin()+(rgr)
#define fori(i,a,b)         for(int i=(a);i<=(b);i++)
#define ford(i,a,b)         for(int i=(a);i>=(b);i--)
#define fori0(i,a,b)        for(int i=(a);i<(b);i++)
#define ford0(i,a,b)        for(int i=(a);i>(b);i--)
#define fr first
#define sc second

using namespace std;

typedef BigInteger_GMP intinf;

struct pairiinf
{
    intinf a,b;
};
struct pairi128
{
    __int128 a,b;
};

pairi128 fib128(long long n,const long long& p)
{
    if(!n)
        return (pairi128){0,1};
    if(n&1)
    {
        auto c=fib128(n-1,p);
        return (pairi128){(c.a+c.b)%p,c.a};
    }
    auto c=fib128(n>>1,p);
    return (pairi128){c.a*(c.b*2+c.a)%p,(c.a*c.a+c.b*c.b)%p};
}
pairiinf fib(const intinf& n,const intinf& p)
{
    intinf a=0,b=1;
    bool flag=false;
    const auto& m=n.n[0];
    ford(i,m._mp_size-1,0)
        ford(j,63,0)
        {
            if(flag)
            {
                auto d=a*a+b*b;
                (a*=(a+b+b))%=p,b=d%p;
            }
            if(m._mp_d[i]&1ull<<j)
            {
                flag=true;
                b+=a,swap(a,b);
            }
        }
    return {a,b};
}

int solve_i128(long long p)
{
    long long p0=((p+1)%5<3?p-1:(p+1)*2);
    if(p<6)
        return p-1;
    while(!(p0&1))
        p0=p0/2;
    auto c=fib128(p0,p);
    if(c.a==0)
        return 4;
    else if((c.b*2+c.a)%p==0)
        return 1;
    else
        return 2;
}

int main(int argc,char* argv[],char* envp[])
{
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    intinf maxp(string(97,'9')+"203"),pow18("1000000000000000000");//10**100-797
    cin>>argc;
    string ps;
    while(argc--)
    {
        cin>>ps;
        intinf p(ps),p0;
        if(p<pow18)
        {
            cout<<solve_i128(stoll(ps))<<'\n';
            continue;
        }
        p0=((p+1)%5<3?p-1:p+p+2);
        while(p0%2==0)
            p0>>=1;
        pairiinf c=fib(p0,p);
        if(c.a==0||c.a==p)
            cout<<4;
        else if(c.a*(c.b*2+c.a)%p==0)
            cout<<1;
        else
            cout<<2;
        cout<<'\n';
    }
    return 0;
}