ホームに戻る
TopCoder の傾向と対策４（Ｃ＋＋編：数え上げ）
０、はじめに
Div2 の hard は数え上げの問題がでます。
単に数え上げという視点で考えてみます。
１、1000000007 で割った余りを答えよ問題対策
TopCoder では整数の６４ビット制限を考慮してくれる。
よって、long long より大きな数値を考えなくても良い。
この数値を超える答えが出る場合には、
答えを 1000000007 で割った余りを答えよ、とされる。
足し算においては足し算を行うごとに % 1000000007 する。
引き算のときは先に引かれる数に 1000000007 を足しておく。
そこから引く数を引いて % 1000000007 する。
こうすると結果がマイナスにならず結果がおかしくならない。
掛け算する場合は足し算と同様である。
A * B をするとき、
A、B それぞれを 1000000007 で割ったときの商と余りを
商A、余A、商B、余Bとする。
A * B = (商A * 1000000007 + 余A) * (商B * 1000000007 + 余B)
= (商A * 商B * 1000000007 + (商A * 余B + 商B * 余A))
* 1000000007 + 余A * 余B
よって、
(A * B) % 1000000007 = (A % 1000000007) * (B % 1000000007)
ゆえに、掛け算のときもそのたびに % 1000000007 すればよい。
和、差、積においては 1000000007 のような素数以外の数字でも成り立つ。
割り算の場合は複雑である。
まず、余りについての標記の説明をする。
７を４で割ると余りは３である。
これを数式で書くといかのようになる。
7 ≡ 3 (mod 4)
さて、フェルマーの小定理は以下のようである。
a^(p-1) ≡ 1 (mod p) // ただし、p は素数で、a は p で割りきれない
m/n (m は n で割り切れる)の計算をする場合、
フェルマーの小定理より次のように書ける。
n^(p - 2) * n ≡ 1 (mod p)
このとき n の逆元は 1 / n^(p - 2) である。
よって、p = 1000000007 のとき、
(m / n) % 1000000007 = m * n^(1000000007 - 2) % 1000000007
簡単にいうと割り切れる m / n を 1000000007 で割った余りは、
(m * n^1000000005) % 1000000007 とできる。
これは 1000000007 が素数であるからできることである。
素数でなければ計算式の中で割り算を使うことができない。
さて、n^1000000005 のところを考えてみる。
単純に 1000000005 回掛け算をすると時間が足りない。
このとき次のようにすると圧倒的に計算量を縮めることができる。
1000000005 は２進数で 111011100110101100101000000101 と書ける。
よって、
n 、 n^2 、 n^4 、 n^8 ... を先に計算しておき、
そのうちの 1000000005 の2進数を下の桁から 1 であるものを見ていき、
n * n^4 * n^512 * n^2048 * ...
と 1 であるもののみ掛け算すると３０回ほどのループで計算できる。
コードで書くと次のようになる。
long long bp[40];
long long binpow(long long a, long long p){
bp[0] = 1;
bp[1] = a;
for(int i = 2; i << 40; ++i){
bp[i] = (bp[i - 1] * bp[i - 1]) % 1000000007;
}
int d = 1;
long long res = 1;
while(p){
if(p & 1){
res = (res * bp[d]) % 1000000007;
}
d++;
p = p >> 1;
}
return res;
}
次に、a^b を 1000000007 で割った余りを考える。
b が 1000000005 のような値であれば上記までの説明で計算できる。
しかし、b が100!のような値であった場合には使えない。
単純に((a%1000000007)^(100!%1000000007))%1000000007は間違いである。
b がとても大きな値である場合もフェルマーの小定理を元に考えることができる。
例えば、a^2000000000 を考える。
a^(p-1) ≡ 1 (mod p)
であるから、
a^2000000000 = a^(1000000006 + 999999994)
= a^(1000000007-1)*a^999999994 ≡ 1 (mod 1000000007)*a^999999994
つまり、
a^bを1000000007で割った余りは、
((a%1000000007)^(b%1000000006))%1000000007
とできる。
もちろんフェルマーの小定理が元なので100000000７が素数であること。
a は p で割りきれないことが必要である。
また、1000000006 は素数で無いので、b の中で割り算ができないことに注意する。
まとめ、
・和と積はそのつど p で割った余りを使えば良い。（p は素数でなくても良い）
・差はマイナスにならないように先に p を足しておく。（p は素数でなくても良い）
・商は割り切れる場合のみ、割る数の p-2 乗を掛ける。（p は素数であること！）
・累乗の指数部分は、p-1 で割った余りを使う。（基数がpで割れない！p は素数であること！）
Java の BigInteger を使うのもルール違反では無い。
また Python は特別なことをせずに大きな値を計算できる。
よって、時には BigInteger や Pyhton を使う戦略もアリである。
ただ、値が大きくなるたび計算が重たくなる。
やはり途中で適度に % 1000000007 をとったほうが良い。
もちろん最後に % 1000000007 するのを忘れずに。
２、場合の数
場合の数は k^n や n! や nPr や nCr を組み合わせる問題。
該当する場合を計数する問題は良く出る。
計数する数はほとんどの場合が巨大な数値になる。
よって、1000000007 で割った余りを返す場合の対策が必要である。
問題
昇順の整数の数列がある。同じ値が続いても良い。
Xの部分は値が不明である。
Xに自由に整数を埋めていいとするとき、
ありえる数列の可能性は何通りか？
答えを1000000007で割った余りを求めよ。
1,3,3,X,45,87,X,X,X,2017,X,X,91021503
考え方
例えば最初のXは3から45の可能性があり43通りの可能性がある。
次の3つのXは87から2017の間に3つの値が入る。
これは1931のうちから重複組み合わせで3つの値を選ぶことになる。
よって、nCrの(1931+3-1)C3で計算できる。
次の2つのXは(91019487+2-1)C2で計算できる。
よって、43*(1932)C3*(91019488)C2通りが答えになる。
途中の計算では1000000007の余りを考慮した計算が必要。
３、行列累乗
問題
f(0) = 1、f(1) = 1、f(x) = f(x - 1) + f(x - 2) のとき、
f(1000000000) を 1000000007 で割った余りを求めよ。
解を出すには 1000000000 回のループを必要としそうだが、
以下の方法を使うとループを 9 回にまで縮めることができる。
#define MOD 1000000007
vector<vector<long long> > mul(vector<vector<long long> > a, vector<vector<long long> > b){
vector<vector<long long> > v(a.size(), vector<long long>(b[0].size(), 0));
for(long long i = 0; i < a.size(); i++){
for(long long k = 0; k < b.size(); k++){
for(long long j = 0; j < b[0].size(); j++){
v[i][j] = (v[i][j] + (a[i][k] * b[k][j]) % MOD) % MOD;
}
}
}
return v;
}
vector<vector<long long> > binpow(vector<vector<long long> > a, long long n){
vector<vector<long long> > v(a.size(), vector<long long>(a.size(), 0));
for(long long i = 0; i < a.size(); i++){
v[i][i] = 1;
}
while(n > 0){
if(n & 1){
v = mul(v, a);
}
a = mul(a, a);
n >>= 1;
}
return v;
}
long long n = 1000000000;
vector<vector<long long> > a(2, vector<long long>(2, 0));
a[0][0] = 1; a[0][1] = 1;
a[1][0] = 1; a[1][1] = 0;
a = binpow(a, n + 1);
std::cout << a[1][0] << std::endl;
この方法を使えば
f(0) = 1、f(1) = 2、f(2) = 3、
f(x) = f(x - 1) + 2 * f(x - 2) + 3 * f(x - 3) + 4
といったさらに複雑な計算も瞬時にできる。
４、包除原理
包除原理は事象の重複を取り除くために有効な考え方である。
例えば、ＡとＢという２つの事象が起こりうる場合、
少なくともＡかＢのどちらかが起こる確率を求めよというとき。
Ａの起きる確率とＢの起きる確率を足し算するのはときに間違いになる。
なぜならＡとＢが同時に起きる場合が考慮されていない。
もしＡとＢが同時に起きる可能性があるのであれば、
Ａの起きる確率＋Ｂの起きる確率−ＡとＢが同時に起きる確率
が少なくともＡかＢのどちらかが起こる確率になる。
これは事象がＡＢＣである場合
Ａ＋Ｂ＋Ｃ−ＡＢ−ＡＣ−ＢＣ＋ＡＢＣ
になる。
事象がＡＢＣＤである場合、
Ａ＋Ｂ＋Ｃ＋Ｄ−ＡＢ−ＡＣ−ＡＤーＢＣ−ＢＤ−ＣＤ
＋ＡＢＣ＋ＡＢＤ＋ＡＣＤ＋ＢＣＤ−ＡＢＣＤ
になる。
これは単に満たす事象が奇数であれば足し算し、
満たす事象が偶数であれば引き算をすればいいと考えられる。
メビウス関数
n=60 のとき
メビウス関数はmapで次のようになる。
1->1,2*3->1,2*5->1,3*5->1
2->-1,3->-1,5->-1,2*3*5->-1
それ以外は0
実行時間はsqrt(n)である。
素因数分解したときの要素が奇数か偶数かの判定が可能。
map<int,int> moebius(int n){
map<int,int> ret;
vector<int> primes;
for(int i = 2; i*i <= n; i++){
if(n%i == 0){
primes.push_back(i);
while(n%i == 0)n /= i;
}
}
if(n != 1)primes.push_back(n);
int m = primes.size();
for(int i = 0; i < (1<<m); i++){
int mu = 1;
int d = 1;
for(int j = 0; j < m; j++){
if(i>>j & 1){
mu *= -1;
d *= primes[j];
}
}
ret[d] = mu;
}
return ret;
}
