BAEKJOON 12728: n제곱 계산

 

2008년 Google code jam round A에 나왔던 문제이다. 이 문제에서는 매우 정밀한 수의 연산을 요구하기 때문에 double이나 float, Decimal 등을 이용한 연산이 옳은 결과를 내지 못하도록 설계되어 있다.

 

Solution

먼저, $(3+\sqrt{5}{)}^n$에 대해서 생각해 보자. $\sqrt{5}$가 살짝 거슬린다. 이를 해결하기 위해 $(3-\sqrt{5}{)}^n$을 더해준다. 그러면 $\sqrt{5}$가 존재하는 모든 항을 없앨 수 있다. 조금 더 자세히 설명하자면,

$(3+\sqrt{5}{)}^n+(3-\sqrt{5}{)}^n=\sum _{i=0}^n(3{)}^{n-i}(\sqrt{5}{)}^i+\sum _{i=0}^n(3{)}^{n-i}(\sqrt{5}{)}^i(-1{)}^i$이고, 이는 $\sum _{i=0}^{\lfloor \frac{n}{2}\rfloor }(3{)}^{n-2i}(5{)}^i$로 표현할 수 있다.

 

이제, $(3+\sqrt{5}{)}^n$을 정리하였을 때를 생각해 보자. $\sqrt{5}$가 붙은 부분과 그렇지 않은 부분으로 나뉠 것이다. 그렇다면 유리수 부분을 $a_n$, 무리수 부분을 $b_n\sqrt{5}$라 하자.

 

이 포스팅에서 처음 도출한 식인 $\sum _{i=0}^{\lfloor \frac{n}{2}\rfloor }(3{)}^{n-2i}(5{)}^i$으로 돌아가자. 이는 $2a_n$과 같다. 그리고 $(3-\sqrt{5}{)}^n$는 항상 (0,1]의 범위에 들어간다. 그러므로, 우리가 구하고자 하는 $(3-\sqrt{5}{)}^n$의 값은 $2a_n-1$과 같다. 이를 이용하여 답을 구해줄 수 있다.

 

이제 $a_n$을 구해보자.

$(3+\sqrt{5}{)}^{n+1}=(3+\sqrt{5})(a_n+b_n\sqrt{5})=(3a_n+5b_n)+(a_n+3b_n)\sqrt{5}$이다.

그러므로, $a_{n+1}=3a_n+5b_n$, $b_{n+1}=a_n+3b_n$ 가 된다.

 

이를 행렬로 표현하여 나타냄으로서 빠르게 $a_n$을 구할 수 있다. 위 점화식을 행렬로 나타내면

$\left(\begin{array}{c}{a}_n\\ b_n\end{array}\right)={\left(\begin{array}{cc}3& 5\\ 1& 3\end{array}\right)}^n\left(\begin{array}{c}{a}_0\\ b_0\end{array}\right)$ 이다. 여기서 $(3+\sqrt{5}{)}^0=1$ 이므로, $a_0=1$, $b_0=0$ 이다.

 

행렬의 거듭제곱은 분할 정복을 통해 $O(\log (n))$에 수행하여 문제를 해결할 수 있다.

 

정답 코드

#include <bits/stdc++.h>
using namespace std;
using ll=long long;
const int mod=1000;
void solve(){
    int n;cin>>n;
    int A[2][2]={{3,5},{1,3}},B[2][2]={{3,5},{1,3}};n--;
    while(n){
        int a,b,c,d;
        if(n&1){
            a=(B[0][0]*A[0][0]+B[0][1]*A[1][0])%mod;
            b=(B[0][0]*A[0][1]+B[0][1]*A[1][1])%mod;
            c=(B[1][0]*A[0][0]+B[1][1]*A[1][0])%mod;
            d=(B[1][0]*A[0][1]+B[1][1]*A[1][1])%mod;
            B[0][0]=a;B[0][1]=b;B[1][0]=c;B[1][1]=d;
        }
        a=(A[0][0]*A[0][0]+A[0][1]*A[1][0])%mod;
        b=(A[0][0]*A[0][1]+A[0][1]*A[1][1])%mod;
        c=(A[1][0]*A[0][0]+A[1][1]*A[1][0])%mod;
        d=(A[1][0]*A[0][1]+A[1][1]*A[1][1])%mod;
        A[0][0]=a;A[0][1]=b;A[1][0]=c;A[1][1]=d;
        n>>=1;
    }
    ll x=(B[0][0]*2+999)%mod;
    if(x<10)cout<<"00";
    else if(x<100)cout<<"0";
    cout<<x;
}
int main(){
    ios_base::sync_with_stdio(false);cin.tie(NULL);
    int T;cin>>T;
    for(int i=1;i<=T;i++){
        cout<<"Case #"<<i<<": ";
        solve();
        cout<<"\n";
    }
    return 0;
}