WARush

SRMの結果とか、解けた問題のコードを書いていきます

SRM566 Div1 Medium " PenguinEmperor"

問題

http://community.topcoder.com/stat?c=problem_statement&pm=12338

パーシーは皇帝ペンギンになりたいと思っている。
まず、彼は自分にその資格があるのか試すため、長いたびに出る必要がある。

ペンギン皇国にはいくつかの都市がある。
全ての都市はとても高い山を取り囲むように、円状に位置していた。
都市は時計回りに0からnumCities-1の番号がつけられている。

パーシーは0番の都市に住んでおり、そこが彼のたびの始まりの場所である。
旅の最初の日、彼は0番目の都市の隣の都市へ移動する。
二日目は出発する都市から2つ離れた都市へ移動する。
旅のk日目は、現在いる都市からkだけ離れた都市へ移動する。
各日で、パーシーは山の周りを時計回りに移動するか、
反時計回りに移動するかを選ぶ事が出来る。

あなたは都市の数を表すnumCitiesが与えられる。
また、パーシーが旅を始めてからの日数daysPassedが与えられる。
daysPassed日目に0番目の都市に戻ってくるような旅程のパターン数を返せ。
2つの旅程において、任意のk日目で、違う都市に旅をしていれば、
それは違う旅程とする。

制約

2 <= numCities <= 350
1 <= dayPassed <= 10^18


考えた事

こ、これは蟻本で習った累乗行列!
練習のチャンス!

numCities日ごとの移動は実質同じなので、
numCities日後のパターン数の変化っぷりを行列にして、
それをdaysPassed / numCitiesだけ掛けてやる。
これはものすごくでかい数になるが、
繰り返し2乗法でlog化する事ができる。

余った日数は350未満になるので、
普通に計算してあげればよい。

2時間後

はぁ・・はぁ・・やっと実装できた。
よしSampleテストだ!

いつまでたってもいつまでたってもテストが終わらない。

行列の掛け算は350^3
それをたぶんlog10^18ぐらいやるんだろうな
あれ?350^3の時点で5000万だから無理だろこれ・・

行列の持ち方を工夫する

掛けるとnumCities日後に変化する行列は循環行列というものらしく、
一番上の列だけ持っておけば表現できる。
さらに、同じ循環行列同士の掛け算の結果はまた循環行列となり、
計算量は列数をNとすると、O(N^2)で出来る。

これを利用する事で全体の計算量を、都市の数をC、日数をDとすると、
O(C^2 * logD)にすることができる。


ソースコード

typedef vector< long long > Vec;

const int MOD = 1000000007;

class PenguinEmperor {
public:

    // 循環行列 a * b
    Vec mul_mod( const Vec& a, const Vec& b ){
        int n = a.size();
        Vec r( n, 0 );
        for( int i = 0; i < n; i++ ){
            for( int j = 0; j < n; j++ ){
                r[i] = ( a[j] * b[(i + j) % n] % MOD + r[i] ) % MOD;
            }
        }
        return r;
    }

    // 循環行列 a ^ n
    Vec pow_mod( Vec a, long long n ){
        Vec r( a.size(), 0 );
        r[0] = 1;
        while( n > 0 ){
            if( n & 1 ) r = mul_mod( a, r );
            a = mul_mod( a, a );
            n >>= 1;
        }
        return r;
    }

    // パターン数をd日後にする(都市数C)
    Vec pass( Vec a, int d, int C ){
        Vec temp( C, 0 );
        for( int i = 1; i <= d; i++ ){
            for( int j = 0; j < C; j++ ){
                int c1 = (j + i) % C;
                int c2 = (j - i + C) % C;
                if( c1 == c2 ){
                    temp[j] = a[c1];
                }else{
                    temp[j] = (a[c1] + a[c2]) % MOD;
                }
            }
            a = temp;
        }
        return a;
    }

    int countJourneys(int C, long long D ){
        
        long long a = D / C;
        int b = D % C;

        Vec m( C, 0 );
        m[0] = 1;

        // C日後のパターン数の循環行列を取得
        m = pass( m, C, C );

        // D / C日後のパターン数を計算
        m = pow_mod( m, a );

        // 循環行列は
        // a(0), a(N-1), a(N-2), ... a(2), a(1)
        // となっているためスワップ
        for( int i = 1; i < C - i; i++ ){
            swap( m[i], m[C-i] );
        }

        // 余った日数を計算
        m = pass( m, b, C );    

        return (int)m[0];
    }
};