問題概要
1,5,10,50,100,500円玉硬貨がある。今$M[i]$円の買い物をしたい。
$T$回$M[i]$が与えられるので、$M[i]$円払うときの硬貨の組合せを$10^9+9$で割ったあまりで求めよ。
$T\leqq10^4$、$M[i]\leqq10^{18}$
解法
naiveなDPでは$O(M)$でこれでは計算が終わらない。
$dp(i):=i$円払う組合せとしてやればいいが、ここでDPを変形する。
$dp(x,i):=x$番目までの硬化を使用して$i$円払うときの組合せとして、
式は、$dp(x,i) = dp(x-1,i) + dp(x,i-c[x])$ となる。
実は$dp[x][a*c[x]+b]$は、$a$の$x$次式で表せる。(天下り的だけど)
したがって5次式であることがわかり、6点あれば多項式補間で解ける。
互いの倍数を考慮すると666で1周期があって、これを上手にまとめるのかなと思って頑張ったけどダメだった。
多項式補完できるのですね…
[証明]
$dp[x][a*c[x]+b]$は、$a$の$x$次式で表せる。これを示す。
$dp[0][any]=1$より、これは$0$次式
$dp[x-1][any]$が$x-1$次式のとき、$dp[x][any]$が$x$次式なら嬉しい。(帰納法が成り立つので！)
$dp[x][a*C[x]+b]$
$= dp[x][(a-1)*C[x]+b] + dp[x-1][a*C[x]+b]$
$dp[x][(a-1)*C[x]+b]$
$= dp[x][(a-2)*C[x]+b] + dp[x-1][(a-1)*C[x]+b]$
$dp[x][(a-2)*C[x]+b]$
$= dp[x][(a-3)*C[x]+b] + dp[x-1][(a-2)*C[x]+b]$
…
$dp[x][C[x]+b]$
$= dp[x][b] + dp[x-1][C[x]+b] $
となるので,
$ dp[x][a*C[x]+b]$
$= dp[x][b] + dp[x-1][a*C[x]+b] $
$+ dp[x-1][(a-1)*C[x]+b] $
$+ … + dp[x-1][C[x]+b]$ …(※)
ここで, 数学的帰納法の仮定より, dp[x-1][a*C[x-1]+b] という形の式は a の x-1 次式になります。また, C[x-1] は問題で与えられている数値上 C[x] の約数になっているので, dp[x-1][a*C[x]+b] という式は, dp[x-1][a*C[x-1]+b] と同じく a の x-1 次式になるはずです。
ここで (※) 式に戻ると, 右辺は 定数 + (a に関する x-1 次式 a 個の和) という形になっています。よって, 右辺は x 次式です。 以上より dp[x][a*C[x]+b] は x 次式になります。
(n付き二項間漸化式を考えるととても良い ref:mayokoさんの解説) したがって、seedとなるのは$b$の部分。 これをもとに、6点を使用して関数を補間すれば正しい。
[備考]
$M= a*500+b $
$ans(0)=dp[0*500+b], $
$ans(1)=dp[1*500+b], $
… ,
$ans(5)=dp[5*500+b]$
$ans(t) = ans(0)* \displaystyle \frac{(t-1)(t-2)(t-3)(t-4)(t-5)}{(0-1)(0-2)(0-3)(0-4)(0-5)}$
$+ans(1)* \displaystyle \frac{ (t-0)(t-2)(t-3)(t-4)(t-5)}{(1-0)(1-2)(1-3)(1-4)(1-5)}$
+ …,
と補完する。
ラグランジュ補間では低次の多項式補間が上手にできるが、高次の補間は上手にできない。その様な場合はスプライン補間を使うと良い。
硬貨数を$C$として、
計算量:$O(dp+TC^2)$
ソース
#include <bits/stdc++.h> using namespace std; using VS = vector<string>; using LL = long long; using VI = vector<int>; using VVI = vector<VI>; using VL = vector<LL>; using VVL = vector<VL>; #define FOR(i, s, e) for (LL(i) = (s); (i) < (e); (i)++) // 繰り返し二乗法 long long powmod(long long x, long long n, long long mod) { long long res = 1; while (n > 0) { if (n & 1) { res = (res*x) % mod; } x = (x*x) % mod; n >>= 1; } return res; } // フェルマーの小定理 a^(p-1)≡1 -> a^-1≡a^(p-2) long long fermat_inv_mod(long long a, long long p) { return powmod(a, p - 2, p); } int main() { cin.tie(0); ios_base::sync_with_stdio(false); int T; cin >> T; VI c({ 1,5,10,50,100,500 }); const LL mod = 1000000009; VL dp(3010, 0); dp[0] = 1; FOR(i, 0, 6) { FOR(m, 0, 3000) { if (m + c[i] <= 3000) (dp[m + c[i]] += dp[m]) %= mod; } } FOR(t, 0, T) { LL M; cin >> M; LL a = M / 500; LL b = M % 500; LL ans = 0; FOR(j, 0, 6) { LL L = dp[500 * j + b]; FOR(k, 0, 6) { if (j != k) { L *= (a - k) % mod; L %= mod; L *= fermat_inv_mod(j - k, mod); L %= mod; } } ans += L + mod; ans %= mod; } /* naive VL dp(M + 501, 0); dp[0] = 1; FOR(i, 0, 6) { FOR(m, 0, M) { (dp[m+c[i]] += dp[m])%=mod; } } debug(dp[M]);*/ cout << ans << "\n"; } return 0; }