結果

問題 No.42 貯金箱の溜息
ユーザー codershifthcodershifth
提出日時 2015-08-16 20:04:23
言語 C++11
(gcc 11.4.0)
結果
AC  
実行時間 47 ms / 5,000 ms
コード長 6,571 bytes
コンパイル時間 1,246 ms
コンパイル使用メモリ 167,524 KB
実行使用メモリ 6,944 KB
最終ジャッジ日時 2024-07-18 09:54:55
合計ジャッジ時間 1,760 ms
ジャッジサーバーID
(参考情報)
judge3 / judge1
このコードへのチャレンジ
(要ログイン)

テストケース

テストケース表示
入力 結果 実行時間
実行使用メモリ
testcase_00 AC 42 ms
6,816 KB
testcase_01 AC 47 ms
6,940 KB
testcase_02 AC 46 ms
6,944 KB
権限があれば一括ダウンロードができます

ソースコード

diff #

#include <bits/stdc++.h>

typedef long long ll;
typedef unsigned long long ull;

#define FOR(i,a,b) for(int (i)=(a);i<(int)(b);i++)
#define REP(i,n) FOR(i,0,n)
#define RANGE(vec) (vec).begin(),(vec).end()

using namespace std;

class ModInt {
    long long value_;

public:
    static long long Mod;
    static void set_mod(long long mod) { Mod= mod; }
    static long long  get_mod(void) { return Mod; }

    ModInt() : value_(0) {}
    ModInt(long long val) {
            if (val >= Mod)
                value_ = val % Mod;
            else if (val >= 0)
                value_ = val;
            else // val < 0
                value_ = (((-val)/Mod+1)*Mod+val) % Mod;
            assert(value_ < Mod);
    }
    ModInt(const ModInt &other) { value_ = other.value_; }
    ~ModInt() {}

#define OPT(OP)                                                         \
    ModInt operator OP (const ModInt &other) const {                    \
        return ModInt(value_ OP other.value_);                          \
    }                                                                   \
    template<class T>                                                   \
    ModInt &operator OP##=(const T &other) {                            \
        return (*this = (*this OP ModInt(other)));                      \
    }
    OPT(+) OPT(-) OPT(*)
#undef OPT

    bool operator==(const ModInt &other) const {
        return (value_ == other.value_);
    }
    bool operator!=(const ModInt &other) const {
        return !(*this == other);
    }
    ModInt &operator=(const ModInt &other) {
        value_ = other.value_;
        return *this;
    }
    // cast overload
    operator int() const { return value_; }
    operator long long() const { return value_; }
    long long value() const { return value_; }
    friend std::ostream& operator<<(std::ostream& os, const ModInt& m) {
        os<<m.value_;
        return os;
    }
};
#define OPT(OP,T)                                                        \
    ModInt operator OP (const ModInt &a, T b) { return a OP ModInt(b); } \
    ModInt operator OP (T a, const ModInt &b) { return ModInt(a) OP b; }
    OPT(*,int) OPT(+,int) OPT(-,int) OPT(*,long long) OPT(+,long long) OPT(-,long long)
#undef  OPT
long long ModInt::Mod = (long long)(1E+9)+9;

template<typename T>
T extgcd(T a, T b, T &x, T &y) {
    T d = a;
    if (b != 0)
    {
        d = extgcd(b, a%b, y, x);
        y -= (a/b)*x;
    }
    else
    {
        x=1; y=0;
    }
    return d;
}

ModInt inverse(const ModInt &a) {
        ll x,y;
        auto gcd = extgcd((ll)a, ModInt::Mod, x, y);
        if (gcd != 1)
            return 0;
        return x;
}
class SighOfSavingsBox {
public:
    void solve(void) {
            // M <= 10^18 なので M に対する DP ではダメ。どう高速化するか?
            //
            // f(x,m) := C[x] 以下の coin を使って m 円をちょうど払う方法の数
            //
            // とすると,
            //
            // f(x,m) = f(x-1,m) + f(x-1,m-C[x]) + f(x-1,m-2*C[x]) + ... + f(x-1,m-k*C[x]) (k*C[x] <= m) ...(*)
            // であり
            //
            // f(x,m-C[x]) = f(x-1,m-C[x]) + ... + f(x-1,m-k*C[x]) より
            // f(x,m) = f(x-1,m) + f(x,m-C[x]) と書ける
            //
            // このままだと m が大きすぎて計算できない。
            //
            // f(0,m) = 1 (0 次の多項式)であり p を未満の任意の整数とすると
            //
            // f(x,m*C[x]+p)     - f(x,(m-1)*C[x]+p) = f(x-1,m*C[x]+p)
            // f(x,(m-1)*C[x]+p) - f(x,(m-2)*C[x]+p) = f(x-1,(m-1)*C[x]+p)
            //    :
            // f(x,C[x]+p)       - f(x,p)          = f(x-1,C[x]+p)
            //
            // なので
            //
            // f(x,m*C[x]+p) = f(x,p) + f(x-1,m*C[x]+p) + ... + f(x-1,C[x]+p)
            //                          <---------------- m ---------------->
            //
            // g(x) を d 次多項式とすると
            // g(x-k) = x^d + G(x-k) と書けるので
            //
            // x                 x
            // ∑ g(y)  = x*x^d + ∑ G(x-k)  となり d+1 次式となる。
            // y=1               k=1
            //
            // よって f(x,m*C[x]+p) は再帰的に m についての x 次式となる
            //
            // ( f(x,m) が x 次式になるわけではないことに注意!!!
            //   なぜなら (*) における k は k=m/C[x] であってこの値はかなり大きい )
            //
            // この x 次多項式の係数ををラグランジェ補完によって決定してやればよい。
            //
            // この多項式 f(x,m*C[x]+p) は p に依存する。
            //
            const vector<int> C = {1,5,10,50,100,500};

            // 多項式を補完するためにサンプル点を準備
            // 多項式は p に依存するためサンプルは多めにとっておく
            int K = C.size();
            int x = K-1;
            int nSamples = C[x]*(K+1) + 1;
            vector<vector<ModInt>> dp(K+1, vector<ModInt>(nSamples));

            // 初期化
            // C[0]=1 で表現できるのは 1 通り
            fill(RANGE(dp[0]), 1);

            // O(nSamples*K)
            FOR(x,1,K)
            REP(m,nSamples)
            {
                dp[x][m] = dp[x-1][m];
                if (m >= C[x])
                    dp[x][m] += dp[x][m-C[x]];
            }


            int T;
            cin>>T;
            REP(_,T)
            {
                ll M;
                cin>>M;

                //
                // sample を使って f(x,m*C[x]+p) をラグランジェ補完して
                //  m = M/C[x]
                //  p = M%C[x]
                // を代入する
                ModInt res = 0;
                ll m = M/C[x];
                ll p = M%C[x];
                // O(K^2*log(K))
                REP(i, K)
                {
                    ModInt c = 1;
                    REP(j, K)
                    {
                        if (i == j)
                            continue;
                        c *= (m-j);
                        c *= inverse(i-j);
                    }
                    res += dp[x][i*C[x]+p]*c;
                }
                cout<<res<<endl;
            }
    }
};

#if 1
int main(int argc, char *argv[])
{
        ios::sync_with_stdio(false);
        auto obj = new SighOfSavingsBox();
        obj->solve();
        delete obj;
        return 0;
}
#endif
0