[BOJ 28365] 점화식과 주기

꽤나 대수학적 지식을 필요로 하는 문제이다. 
문제에서 주어진 $S, T$를 그대로 활용하자. 그럼 다음이 성립한다. 
\[\begin{bmatrix} x_{n+1} \\ x_{n} \end{bmatrix} = \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}^{n} \begin{bmatrix} x_{1} \\ x_{0} \end{bmatrix}\]
그럼
\[ \begin{bmatrix} x_{S+T+1} \\ x_{S+T} \end{bmatrix}- \begin{bmatrix} x_{S+1} \\ x_{S} \end{bmatrix}= \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}^{S}\left( \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}^{T} \begin{bmatrix} x_{1} \\ x_{0} \end{bmatrix}- \begin{bmatrix} x_{1} \\ x_{0} \end{bmatrix}\right)\]
따라서 우리가 찾아야 하는 $S, T$는
\[\begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}^{S}\left( \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}^{T} \begin{bmatrix} x_{1} \\ x_{0} \end{bmatrix}- \begin{bmatrix} x_{1} \\ x_{0} \end{bmatrix}\right) \equiv \begin{bmatrix} 0 \\ 0 \end{bmatrix}\pmod{p}\]
이다. 우리는 $b\not\equiv 0\pmod{p}$ iff $ \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}$ is invertible 임을 안다. 
먼저 $b\not\equiv 0\pmod{p}$라고 가정해보자. 또, 알아야 할 사실은 $|GL_{2}(\mathcal{F}_{p})|=(p^2-1)(p^2-p)$이라는 것이다. 따라서 \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix}는 order가 유한하다는 것을 알 수 있다. 
만약 위 행렬의 서로 다른 eigenvalue가 2개 존재하면 위 행렬은 대각화 가능하고 
\[ \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix} \simeq \begin{bmatrix} r & 0 \\ 0 & s \end{bmatrix}\]
이며 $r^{p-1} \equiv s^{p-1} \equiv 1\pmod{p}$ 이기 때문에 위 행렬의 order는 $p-1$을 나눈다. 
행렬이 eigenvalue를 갖지만 중복된 값으로 1개 가진다고 하자. 그럼 Jordan form에 의해서
\[ \begin{bmatrix} a & b \\ 1 & 0 \end{bmatrix} \simeq \begin{bmatrix} r & 1 \\ 0 & r \end{bmatrix}\]
가 된다. 이 행렬의 $n$제곱을 구해보면
\[ \begin{bmatrix} r & 1 \\ 0 & r \end{bmatrix}^{n} = \begin{bmatrix} a^{n} & na^{n-1} \\ 0 & a^{n} \end{bmatrix}\]
만약 $n=p(p-1)$이라면 오른쪽의 행렬은 단위 행렬이 된다. 
즉, 행렬이 eigenvalue를 가지면 그 행렬의 order는 $p(p-1)$을 나눈다. 
 
eigenvalue를 가지지 않는다고 하자. 그럼 이 행렬의 특성 다항식은 $\mathcal{F}_{p}$에서 irreducible이다. 특성 다항식의 차수가 2이므로 $\mathcal{F}_{p^{2}}$에서 이 다항식은 완전히 split 된다. 따라서 $GL_{2}(\mathcal{F}_{p^{2}})$에서 eigenvalue를 distinct하게 가지게 된다. 즉, 이 행렬의 order는 $p^{2}-1$을 나눈다. 
 
우리에게 주어진 행렬이 eigenvalue를 갖는지 확인하기 위해서는 $x^{2}-ax-b\equiv 0\pmod{p}$의 solution이 존재하는지 확인하면 된다. 행렬의 order가 무엇을 나누는지 확인했다면 각 수 $p(p-1),p^{2}-1$의 약수를 하나씩 골라 행렬을 거듭제곱 하면서 order을 찾아내도록 하자. 약수를 찾는 과정에서 pollard rho 알고리즘을 사용하자. 게다가, 주어진 조건에서 $p=O(10^9)$이기 때문에 $p(p-1),p^{2}-1=O(10^{18})$로 충분히 int64 범위 내에 존재한다. 
 
그럼 이제 $b\equiv 0\pmod{p}$ 인 경우를 고려하자. 이 경우 $x_{n}$은 단지 등비수열이 된다. 
$x_{1} \equiv 0\pmod{p}$이라고 하자. 만약 $x_{0}\equiv 0\pmod{p}$이면 $S=0,T=1$이고 그렇지 않으면 $S=1,T=1$이다. 
$x_{1}\not\equiv 0\pmod{p}$라고 하자. 만약 $a\equiv 0\pmod{p}$이면 $S=2,T=1$이다. 그렇지 않으면 $a^{T}x_{1}\equiv x_{1}\pmod{p}$인 최소의 $T$를 찾는다. 즉, $a\pmod{p}$의 order를 찾으면 된다. $a$의 order는 반드시 $p-1$을 나누기에 위와 같은 방법을 사용하면 된다. 찾은 order를 $k$라고 하면 정답은 $S=1, T=k$ 이다. 
 
만약 $x_{0} \equiv x_{1}\equiv 0\pmod{p}$가 된다면 정답은 항상 $S=0, T=1$이라는 것을 유의하자. 
코드는 아래에 있다. 
 

# 코드 보기

#include <bits/stdc++.h>
 
typedef long long ll;
ll mod;
 
ll fpow(ll n,ll k,ll mod){
    ll s=1;
    for(;k;k>>=1){
        if(k&1) s=(__int128)s*n%mod;
        n=(__int128)n*n%mod;
    }
    return s;
}
 
ll miller(ll n,ll a){
    ll t=n-1;
    while(t%2==0) if(fpow(a,t/=2,n)==n-1) return 1;
    t=fpow(a,t,n);
    return t==n-1 || t==1;
}
 
ll prime(ll n){
    if(n==1) return 0;
    int m=12,i;
    ll a[12]={2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
    for(i=0;i<m;i++){
        if(a[i]==n) return 1;
        if(!miller(n,a[i])) return 0;
    }
    return 1;
}
 
ll gcd(ll a,ll b){
    return b?gcd(b,a%b):a;
}
 
ll rf(ll x,ll n){
    return ((__int128)x*x+13)%n;
}
 
ll pollard_rho(ll n){
    ll x=2,y=2,t=1,c=13;
    if(prime(n)) return n;
    if(n%2==0) return 2;
    while(t==1 || t==n){
        x=rf(x,n);
        y=rf(rf(y,n),n);
        t=gcd(x>y?x-y:y-x,n);
        if(t==n) x=y=rand()%n;
    }
    return prime(t)? t:pollard_rho(t);
}
 
ll solve(ll s[][2], ll n){
    if(n == 1) return 0;
    ll p = 0;
    ll i,j;
    ll r;
    while(n-1){
        r=pollard_rho(n);
        n/=r;
        s[p++][0]=r;
        s[p-1][1]=1;
    }
    for(i=0;i<p;i++){
        for(j=0;j<p-i-1;j++){
            if(s[j][0]>s[j+1][0]){
                ll t=s[j][0];
                s[j][0]=s[j+1][0];
                s[j+1][0]=t;
            }
        }
    }
    ll t[1001][2], pt = 0;
    t[pt][0] = s[0][0];
    t[pt][1] = 1;
    for(i = 1;i < p;i++){
        if(s[i][0] == s[i - 1][0]) t[pt][1]++;
        else{
            pt++;
            t[pt][0] = s[i][0];
            t[pt][1] = 1;
        }
    }
    pt++;
    for(i = 0;i < pt;i++){
        s[i][0] = t[i][0];
        s[i][1] = t[i][1];
    }
    return pt;
}
 
void mat_mul(ll A[2][2], ll B[2][2]){
    ll C[2][2];
    C[0][0] = ((A[0][0] * B[0][0]) % mod + (A[0][1] * B[1][0]) % mod) % mod;
    C[0][1] = ((A[0][0] * B[0][1]) % mod + (A[0][1] * B[1][1]) % mod) % mod;
    C[1][0] = ((A[1][0] * B[0][0]) % mod + (A[1][1] * B[1][0]) % mod) % mod;
    C[1][1] = ((A[1][0] * B[0][1]) % mod + (A[1][1] * B[1][1]) % mod) % mod;
    A[0][0] = C[0][0]; A[0][1] = C[0][1]; A[1][0] = C[1][0]; A[1][1] = C[1][1];
}
 
void mat_fpow(ll A[][2], ll k){
    ll S[2][2] = {1, 0, 0, 1};
    for(;k;k >>= 1){
        if(k & 1) mat_mul(S, A);
        mat_mul(A, A);
    }
    A[0][0] = S[0][0]; A[0][1] = S[0][1]; A[1][0] = S[1][0]; A[1][1] = S[1][1];
}
 
ll min;
ll a, b;
 
void f(ll s[][2], ll now, ll end, ll val){
    if(val >= min) return;
    if(now == end){
        ll A[2][2] = {a, b, 1, 0};
        mat_fpow(A, val);
        if(A[0][0] == 1 && A[0][1] == 0 && A[1][0] == 0 && A[1][1] == 1) min = val;
        return;
    }
    ll t = 1;
    for(int i = 0;i <= s[now][1];i++){
        f(s, now + 1, end, val * t);
        t *= s[now][0];
    }
}
 
void g(ll s[][2], ll now, ll end, ll val){
    if(val >= min) return;
    if(now == end){
        if(fpow(a, val, mod) == 1) min = val;
        return;
    }
    ll t = 1;
    for(int i = 0;i <= s[now][1];i++){
        g(s, now + 1, end, val * t);
        t *= s[now][0];
    }
}
 
int main(){
    ll i, j, len;
    ll s[1001][2], t[1001][2], ls, lt;
    ll x0, x1;
    scanf("%lld %lld %lld %lld %lld", &mod, &a, &b, &x0, &x1);
    if(b % mod == 0) {
        if(x1 % mod == 0){
            if(x0 % mod == 0) printf("0 1");
            else printf("1 1");
            return 0;
        }
        if(a % mod == 0){
            printf("2 1");
            return 0;
        }
        
        min = mod - 1;
        len = solve(s, mod - 1);
        //for(i = 0;i < len;i++) printf("%lld %lld\n",s[i][0],s[i][1]);
        g(s, 0, len, 1);
        printf("1 %lld", min);
        return 0;
    }
    min = mod * mod;
    ll inv4 = fpow(4, mod - 2, mod);
    ll D = (a * a % mod * inv4 % mod + b) % mod;
    ll solvab = fpow(D, (mod - 1) / 2, mod);
    if(solvab != mod - 1){
        len = solve(s, (mod - 1) / 2);
        s[len][0] = mod; s[len][1] = 1;
        s[len + 1][0] = 2; s[len + 1][1] = 1;
        len += 2;
        //for(i = 0;i < len;i++) printf("%lld %lld\n",s[i][0],s[i][1]);
    }
    else{
        ll twos = 0;
        ll A = mod - 1, B = mod + 1;
        while(A % 2 == 0){
            A /= 2;
            twos++;
        }
        while(B % 2 == 0){
            B /= 2;
            twos++;
        }
        ls = solve(s, A);
        lt = solve(t, B);
        len = ls + lt;
        for(i = ls;i < len;i++){
            s[i][0] = t[i - ls][0];
            s[i][1] = t[i - ls][1];
        }
        s[len][0] = 2; s[len][1] = twos;
        len++;
        //for(i = 0;i < len;i++) printf("%lld %lld\n",s[i][0],s[i][1]);
    }
    f(s, 0, len, 1);
    if(x1 % mod == 0 && x0 % mod == 0) min = 1;
    printf("0 %lld", min);
}

# 코드 닫기