garakutagoya

興味、関心のあったこと、そして私の気持ちなどを残していきたい

Typical DP Contest_E_数

Atcoder DP 動的計画法 桁DP

コンテストの名前の通りDPの問題.どのようにDPテーブルを更新するかが少し特徴的なのでブログにまとめておく.新規性は1ミリもありません.

atcoder.jp

問題の説明

N以下の正整数であり,十進表記したときの各桁の数の和がDの倍数であるものの個数を[tex: mod 1,000,000,007 ]で求めよ.

条件

  • 1≦N≦10^10000
  • 1≦D≦100

 

考えた事

見ての通りNがとんでもなく大きい.もちろんN毎に判定してcountするようなO(N)解法が全く持って間に合わないのは火を見るよりも明らかだろう.

この問題に限らず,一般に調べるべきもののサイズが非常に大きいもの*1であるとき,string型に直して各桁について調べるアプローチが典型である.今回もそのようなアプローチを取る.つまり,Nを最大10000桁の数とみて,その各々の桁の情報を使って問題を解く方針だ.

今回ならどのように見て行けばいいだろうか.もちろん素直に各桁の数字を決めていって全部試す方法では意味がない.どうすべきか.

例えばN=1235, D=4であったとして,一々これを1,2,3,4......1235と試すのはどう考えても無駄が大きい.これを変えたいのだ.ここでDPの考えが出てくる.先にDPテーブルをどう作ればいいかを見てもらおう.

dp_{(確定桁数),(N未満確定かどうか),(問題の条件式)}

として,DPテーブルに上の条件で満たしている数の個数を書き入れていくのだ

基本はこの形である.数を上の桁から決定していくように眺めるのだ.この場合だと千の位が0か1のどちらかで決定する.

イメージとしては配るDPを想像してもらえばいい.上の桁(例えば0???)を決定していく事で,下の桁の情報が増えていく形だ.

千の位が0である時,どうやってもNを下回る事が確定しているので,残りの百の位,十の位,一の位で4の倍数となるものは千の位が0となるものから来ることが分かる.残りの桁の決定でも同様である.

千の位が0であるとき,もちろんN未満確定なので,

dp_{2,1,(問題の条件式)'}       +=     dp_{1,1,(問題の条件式)'}

となる.明らかに桁の情報が増えたとしても元を辿れば前の桁の情報からきている.

今回であれば,各桁の合計がDで割り切れるかどうかなので,問題の条件式は各桁の合計jとして

dp_{1,1,j\bmod D}

となる.上の例であるならば千の位は0として進めているので,j=0となる.N=1235としているので1桁目は0,1のいずれしか存在しない.つまり,

dp_{1,1,j}=1(j=0, 1)

dp_{1,1,j}=(otherwise)

の形である.

そして,上の桁が0確定であるものを,dp_{2,1,(0+0)\bmod D},dp_{2,1,(0+1)\bmod D},dp_{2,1,(0+2)\bmod D},...dp_{2,1,(0+9)\bmod D}に配っていくのだ.

ここで0~9は何かというと二桁目,つまり百の位の数字であり,00~09の段階で条件を満たしているものは何個かというように見ている.

今回であると

dp_{2,1,0} =3

dp_{2,1,1} =2

dp_{2,1,2} =2

dp_{2,1,3} =2

となる.見て分かるようにさっきまで雑然としていたものが徐々にまとまってきている.

まだ終わりではない.今度は1桁目,つまり千の位が1のものを考える.千の位が1であるとき,N未満確定かどうかは分からないため,

dp_{1,0,1}=1となる.今度はこれを配っていく.

千の位が1であり,Nを超えないように百の位を決めるとき,百の位の数として考えられる候補は0,1,2であり,このうち0,1は明らかにN未満確定する.

(1+0)%D=1,(1+1)\bmod D=2である事に注意すると,

dp_{2,1,1} += dp_{1,0,1}

dp_{2,1,2} += dp_{1,0,1}

のように配られ,結果,

dp_{2,1,0} =3

dp_{2,1,1} =3

dp_{2,1,2} =3

dp_{2,1,3} =2

となる.最後に,千の位が1,百の位が2であるものを考えると,これはN未満確定ではないため,

dp_{2,0,(1+2)\bmod D}+=dp_{1,0,1}

と配られる.つまり,

dp_{2,0,3}=1

となる.以上をまとめると,

dp_{2,1,0} =3

dp_{2,1,1} =3

dp_{2,1,2} =3

dp_{2,1,3} =2

dp_{2,0,0}=0

dp_{2,0,1}=0

dp_{2,0,2}=0

dp_{2,0,3}=1

となる.これはN以下となるよう上2桁を決める時の候補が00~09,10,11,12の計12個である事と一致している.

ポイントとなるのは桁数の合計を記録するのではなく,桁数の合計をmod Dで取った数を記録している事にある.こうする事で配列の幅をDとし,情報の整理と計算量削減に成功している.

乱暴な言い方かもしれないが,D=4であるとき,022?と026?と035?と116?が持つ情報量は全く同じであるから,1括りにまとめておくのが楽だというような解釈をすればいいだろう.これは楽であるだけでなく,こう1括りにまとめてしまう事で,計算量の幅をO(10^{Nの桁数})

から大きく落とし込める点でも有益である.

この纏めるというところが桁DPの真髄かもしれない.まず,前の桁を決定し問題の条件式に応じてまとめ,その次の桁にその情報(個数)を配り,またまとめ,また配りを繰り返していく.

今回であればまとめる作業は[tex:*2\bmod  D]であり,これはただ剰余を求めているだけであるからO(1)で行える.本来は0~9で一々見て,各々の桁の和について一つ一つ考えるため,計算量はO(N)であったところから比べると大きな進歩であろう.これを踏まえてコードに書き起こすと以下の通り.

ll n = N.size();
    vector<vector<vector<mint>>> dp(
        n + 1, vector<vector<mint>>(
                   2, vector<mint>(D, 0)));  // dp[桁数][N未満確定か][j % D]
    dp.at(0).at(0).at(0) = 1;  // 0桁目は0とみなすため,dp[0][0][0]=1
    rep(i, n) {
        rep(j, D) {
            rep(k, 10) {
                dp.at(i + 1).at(1).at[(j + k) % D] += dp.at(i).at(1).at(j);  
//上から見てiより前の桁でN未満確定から考える.
            }
            ll keta = N.at(i) - '0';  // Nの上から見てi桁目
            rep(m, keta) {
                dp.at(i + 1).at(1).at[(j + m) % D] += dp.at(i).at(0).at(j);
            }
            dp.at(i + 1).at(0).at[(j + keta) % D] += dp.at(i).at(0).at(j);
        }
    }

 

初期条件として0桁目で0であるとき,Nの0桁目は0であるとみなして,dp_{0,0,0}=1としている.

後はコードにある通り,上から桁を決め,Dの剰余の回数配り,新しい桁としてあり得るのは0~9であるとしてループを回す.

計算量はO(n×D×10)=O(10^8)となり,制限時間に間に合う.

答えとしては,dp_{n,1,0},dp_{n,0,0}を合計したものから1を引いたものとなる.

dp_{n,1,0}はn桁目*3まで決めてN未満確定のもの(0000~1234)で,各桁の合計をDで割った余りが0となるものの合計であり,dp_{n,0,0}はn桁目まで決めて,N未満確定でないもの(これは実質N=1235のみを指す)の中で,各桁の合計をDで割った余りが0となるものの合計である.*4

1を引く理由としては初期条件の部分でdp_{n,0,0}=1としており,それが巡り巡って,dp_{n,1,0}の中に0000の個数1個分を余分に含む形となっているからである.問題文が「N以下の非負整数で～」となっていた場合,この処理は行わなくていい.

いかがだっただろうか.正直,自分のプログラミングの力が足りなく,例に大きく頼った未熟な説明となってしまった.これが今の自分の限界である.

とはいえ,今の自分の力で精一杯,桁DPのアルゴリズムを説明できていると思う.これを読んだ誰かの理解の一助となる事を願っております.今回も最後までお読みいただきありがとうございました.

プログラム

#include <bits/stdc++.h>
using namespace std;
#include <atcoder/all>
using namespace atcoder;
using ll = long long;
using mint = modint1000000007;
#define rep(i, n) for (ll i = (0); i < (ll)(n); i++)
int main() {
    ll D;
    cin >> D;
    string N;
    cin >> N;
    ll n = N.size();
    vector<vector<vector<mint>>> dp(
        n + 1, vector<vector<mint>>(
                   2, vector<mint>(D, 0)));  // dp[桁数][N未満確定か][j % D]
    dp.at(0).at(0).at(0) = 1;  // 0桁目は0とみなすため,dp[0][0][0]=1
    rep(i, n) {
        rep(j, D) {
            rep(k, 10) {
                dp.at(i + 1).at(1).at[(j + k) % D] += dp.at(i).at(1).at(
                    j);  //上から見てiより前の桁でN未満確定から考える.
            }
            ll keta = N.at(i) - '0';  // Nの上から見てi桁目
            rep(m, keta) {
                dp.at(i + 1).at(1).at[(j + m) % D] += dp.at(i).at(0).at(j);
            }
            dp.at(i + 1).at(0).at[(j + keta) % D] += dp.at(i).at(0).at(j);
        }
    }
    mint sum = dp.at(n).at(0).at(0) + dp.at(n).at(1).at(0) - 1;
    cout << sum.val() << endl;
}

*1:多倍長整数というそうです

*2:今までの桁の合計\bmod D)+(新しい桁の数字

*3:もちろん,全ての桁を決定した事になる

*4:ここでは実質的に0か1となる

AtCoder Grand Contest 007_B_Construct Sequences

Atcoder 青色

青色の問題.久々にアルゴリズムもへったくれもない問題です.数学好きなのでこういった問題の方が好きです.

atcoder.jp

問題の説明

1からNを並べた順列{p_i}(1≦i≦N)が与えられる.以下の条件を満たす数列{a_N}および,{b_N}を出力せよ.

条件

  • 1≦a_i,b_i≦10^9
  • 正整数列{a_i}は狭義単調増加
  • 正整数列{b_i}は狭義単調減少
  • a_{p_i}+b_{p_i}は狭義単調増加
  • 1≦N≦200,000

 

考えた事

問題文にある通り,答え方としてはp_iに合わせてa_i,b_iを作りだし,答えればよい.

まず考えるべきは全探索である.極端な話であるが*1,この問題の入力としてあり得る数列p_i全てに合わせa_iとb_iを先に全て作りだせば答えとなる.

この指針には2つ問題がある.まず第一に計算量である.当然のことながら,上の方法を行う場合,計算量はO(200000!)となりどう考えても無理がある.

次に解法としての適当さである.そもそもこの問題は答えを1つ列挙する問題である.裏を返せば全部試すなと読めるのではなかろうか.答えが複数種類ある事を暗に示唆されている問題で愚鈍に全て試すのは明らかに効率がよくない.

ではどうすべきか.こういった問題において効果的なのは条件を徐々に厳しくする形で問題を見つめる事である.全ての条件を一気に考えても分からなくなってにゃにゃ！？となっておしまいである.我々は人間だ.カワイイカワイイにゃんこと違って物事を整理して考えねばいけない種族だ.がんばろう.

この場合だと「a_{p_i}+b_{p_i}は狭義単調増加」はめんどうそうである.ひとまず置いて置こう.他から考える.とはいえ他の条件は狭義単調増加,狭義単調減少な数列を作るだけなのですぐできるだろう.なぜわざわざすぐできる事を分けて先に書いたか.ここが今回の大きなポイントである.

つまり今回の問題はp_iによらず,先にa_i,b_iの数列を作っておき,p_iに合わせてそれを微調整するように考えればいいのである.

a_i+b_iがどの順で大きくなるかはp_iが入力されるまで分からない.ならば先にどんな数列を作っておくとよいだろうか.自分なら先に作っておくべき数列はiによらずa_i+b_iが同じ大きさになる数列と考える.

そして後からp_iの順にそれぞれの数列の項に+iずつしていけば「a_{p_i}+b_{p_i}は狭義単調増加」の条件も満たす.ここで注意しなければならないのは,+iによる変化量である.これを行うと一番小さいa_{p_1}とa_{p_{20000}}とで19999の差が出来る.「{a_i}は狭義単調増加」という条件を満たしたままこれを行うには,先にa_iとa_{i+1}とで20,000の差をつけておけばいい.

つまり,a_1=20000,a_2=40000,a_3=60000......といった具合だ.そしてこれによって先に作るべきb_iの詳細も分かる.「iによらずa_i+b_iが同じ大きさになる数列」を作ればいいのだからb_i-b_{i+1}=20000となるように作ればいい.

vll a(N);
    vll b(N);
    rep(i, N) { a.at(i) = 20000 * (i + 1); }
    rep(i, N) { b.at(i) = 20000 * (N - i); }

このようにして作ったa_iとb_iが以上である.これによって「a_{p_i}+b_{p_i}は狭義単調増加」以外の全てを満たし,「iによらずa_i+b_iが同じ大きさになる数列」を作成できた.

後は,

rep(i, N) {
        a.at(p.at(i) - 1) += i;
        b.at(p.at(i) - 1) += i;
    }

として,p_iの順に1ずつ差をつけていけば,元々の条件を満たしながら「a_{p_i}+b_{p_i}は狭義単調増加」の条件をみたす.

そしてこれは高々,(2×10^4)^2=4×10^8であり,a_i,b_i≦10^9の制約を満たしている.以上をもって数列a_i,b_iを作成する事に成功した.

こういった問題は典型解がなく,難しいが非常に楽しいものである.コツというか自分の考えとしては条件を緩める,特にp_iのような変数絡みの条件を後回しにすると,解放のアプローチが得やすいだろう.ここまでお読みいただきありがとうございました.

プログラム

#include <bits/stdc++.h>
using namespace std;
using ll = long long;
using vll = vector<long long>;
#define reps(i, a, n) for (ll i = (a); i < (ll)(n); i++)
#define rep(i, n) for (ll i = (0); i < (ll)(n); i++)
int main() {
    ll N;
    cin >> N;
    vll p(N);
    rep(i, N) { cin >> p.at(i); }
    vll a(N);
    vll b(N);
    rep(i, N) { a.at(i) = 20000 * (i + 1); }
    rep(i, N) { b.at(i) = 20000 * (N - i); }
    rep(i, N) {
        a.at(p.at(i) - 1) += i;
        b.at(p.at(i) - 1) += i;
    }
    rep(i, N) { cout << a.at(i) << ' '; }
    cout << endl;
    rep(i, N) { cout << b.at(i) << ' '; }
    cout << endl;
}

*1:極端とつけている時点で察しはつくかもしれない

第１２回日本情報オリンピック 予選（過去問)_D_ 暑い日々 (Hot days)

Atcoder DP 動的計画法

DPの問題.前日の情報をどうするか,そこがカギとなります.

atcoder.jp

問題の説明

温度条件がT度の下でD日服を選んだ場合の派手さの差の絶対値の合計の最大値を求めよ.

条件

  • 1≦D≦200
  • 1≦T≦60

 

考えた事

派手さの差の絶対値の合計を考える.つまり前の情報を遷移させる必然性を伺わせる事から,この問題の解法としてその場その場で最適な服を選び続けるGreedyかDPが思いつくだろう.

事実としてこの問題はGreedyに近い発想でも解く事が出来る*1.

keidaroo.hatenablog.com

この記事にあるように,温度条件を満たしている服の中で1日毎に派手さが最大と最小のものを取り,前の日に選んだ服との派手さの差が大きくなる方を順々に選び続ければいいとなる.

上の方法が成立する理由は,差の絶対値を足し続ける手法であるため,基本端っこと端っこ(最大と最小)を選び続ける方が差の絶対値が大きくなるからである.この場合の計算量はO(D)であり,非常に高速である.

ただ,残念ながら自分はそれが思いつかなかった.ポンコツやね.ではどうするか.愚直にDPをしていこう.また手順に則ってやっていく.

ステップ1.全探索出来るかどうか

全探索する場合の計算量の見積もりをしよう.とりあえず服を全部選べると仮定した場合が都合がよい*2.

この場合,N枚の服から1枚を選ぶことをD日繰り返すため,計算量はO(N^D)となる.

指数時間,Dの上限は200,当然全探索では間に合いません.

ステップ2.状態の遷移を考える

この問題の場合,前日に選んだ服と今日選んだ服の派手さの差の絶対値をそこまでの派手さの差の絶対値の合計に足し合わせる事で遷移するだろう.

ステップ3.変数を見つけ,配列の表を書く

遷移する変数は何になるだろうか.自分は最初,日数,前日に選んだ服,当日の服を変数にするといいと考えた.三次元DPである.この場合,服を選べない場合と同じ服しか選べない場合を差別化するために初期値は-INF=-10^{18}で設定した.

注意したいのは1日目で,この時にはいずれにせよ派手さの差の絶対値の合計値ゼロである.しかし,温度条件を満たしている服と満たしていない服が存在するため,前者には0,後者はそのまま-INFとした.

もちろんこの発想は間違っていない.ただ図に起こしてみるとどうであろうか.

これは入力例1を図に起こしたものである.iが前日に選んだもの(存在しない場合は0),jがその日に選んだものである.

見ての通り遷移がややこしい.なぜか.j→iとなるからである.ある日に選んだものが次の日に影響するため,このようにswapする.だからややこしいのである.

また,dp_{2,3,2}のような後に生きない無駄なデータがある.そしてどのデータが無駄か必要か判別するのもややこしい.使うデータを判別するならばdp_{2,2}列の中の最大値を使えばよい(前々日に選んだ服はどうでもいいため).しかし毎回毎回それを調べるのも億劫である.

そして自分は何がしたかったのか0以上のデータを判別してからこれをやろうとしたため,ややこしすぎて実装出来なくなってしまった.さてどうしようか.

ステップ4.漸化式の作成

大事なのは変数である.上の中でいらない変数はなんであろうか.答えは昨日選んだ服である.なぜいらないか.理由は簡単である,もし服を選べない場合,そこは-INFとなっているため,昨日の服の選択の如何に関わらず,遷移する中での派手さの合計の最大値を求めさえすれば,問題は起きないという事だ.

何が言いたいかというと上の図のi(0≦i≦4)を全て圧縮して以下の図のようにすればいい.

iを全て圧縮して二次元DPとし,dp_{2,2}にはdp_{1,j}に全てabs(C_j-C_2)を足し,その最大値だけを記録していけばよいとなる.ポイントとなるのは初期値であり,初期値が-INFであるからこそ,前日に選べた服と選べなかった服を差別化できる.一回の派手さの差の絶対値を足したとて,-∞は-∞のまま.一々「前の日に使えた服は何だっけ～、ううううう.」などと悩まなくてよいのだ.

したがって変数は2つでよい.一応もう少ししっかり書くと,A_j≦T_d≦B_jの時,

dp_{d,j}=max(dp_{d, j},dp_{d-1, i} + abs(C_i-C_j))

となるだろう.iについては服の数だけループを回して全て調べあげればよい.

もちろん条件を満たしていない場合,dp_{d,j}=dp_{d, j}となるため,一々書く必要もない.

今回の場合,「最高気温が天気予報に従ったときに着るのに適した服が,D日間のどの日に対しても1つ以上存在することが保証されている.」とあるため,必ず1つは条件を満たす.d行全てが-INFで終わる事はあり得ない事が保証されているため,このようにして考えていい.まとめると,

dp_{d,j}=max(dp_{d, j},dp_{d-1, i} + abs(C_i-C_j))(A_j≦T_d≦B_j)

dp_{d,j}=dp_{d,j}=-INF(otherwise)

となる.初期条件としては

dp_{1,j}=0(A_j≦T_0≦B_j)

dp_{1,j}=-INF(otherwise)

とするとよいだろう.

ステップ5.漸化式を元に配列を埋める

ここまで書けば実装は難しくない.

vvll dp(D + 1, vll(N, -INF));

としてDP配列を作る.とりあえず-INFで置いて置き,

rep(i, N) {
        if (A.at(i) <= T.at(0) && B.at(i) >= T.at(0)) {
            dp.at(1).at(i) = 0;
        }
    }

として0を埋めればよい.

後は3重ループを回し,配列を埋めればよい.

reps(i, 2, D + 1) {
        rep(j, N) {
            if (A.at(j) <= T.at(i - 1) && B.at(j) >= T.at(i - 1)) {
                rep(k, N) {
                    CHMAX(dp.at(i).at(j),
                          dp.at(i - 1).at(k) + abs(C.at(j) - C.at(k)));
                }
            }
        }
    }

漸化式の部分はこのようになる.計算量はO(D×N×N)=O(D×N^2)となり,高々O(10^6)となるため,十分間に合う.

ステップ6.求めたい配列の値を見て答えを導く

dp_{D,j}をjを動かして,最大値を出力すればよい.遷移がややこしくタフな問題であった.お疲れ様でした.

プログラム

#include <bits/stdc++.h>
using namespace std;
#include <atcoder/all>
using namespace atcoder;
using ll = long long;
using ld = long double;
using P = pair<int, int>;
using Graph = vector<vector<int>>;
using vi = vector<int>;
using vl = vector<long>;
using vll = vector<long long>;
using vvi = vector<vi>;
using vvl = vector<vl>;
using vvll = vector<vll>;
using vvvll = vector<vvll>;
using vs = vector<string>;
using pii = pair<long long, long long>;
using mint = modint1000000007;
const long double EPS = 1e-10;
const long long INF = 1e18;
const long double PI = acos(-1.0L);
#define reps(i, a, n) for (ll i = (a); i < (ll)(n); i++)
#define rep(i, n) for (ll i = (0); i < (ll)(n); i++)
#define rrep(i, n) for (ll i = (1); i < (ll)(n + 1); i++)
#define repd(i, n) for (ll i = n - 1; i >= 0; i--)
#define rrepd(i, n) for (ll i = n; i >= 1; i--)
#define ALL(n) begin(n), end(n)
#define IN(a, x, b) (a <= x && x < b)
#define INIT                          \
    std::ios::sync_with_stdio(false); \
    std::cin.tie(0);
template <class T>
inline T CHMAX(T& a, const T b) {
    return a = (a < b) ? b : a;
}
template <class T>
inline T CHMIN(T& a, const T b) {
    return a = (a > b) ? b : a;
}

int main() {
    ll D, N;
    cin >> D >> N;
    vll T(D);
    rep(i, D) { cin >> T.at(i); }
    vll A(N);
    vll B(N);
    vll C(N);
    rep(i, N) { cin >> A.at(i) >> B.at(i) >> C.at(i); }
    vvll dp(D + 1, vll(N, -INF));
    rep(i, N) {
        if (A.at(i) <= T.at(0) && B.at(i) >= T.at(0)) {
            dp.at(1).at(i) = 0;
        }
    }

    reps(i, 2, D + 1) {
        rep(j, N) {
            if (A.at(j) <= T.at(i - 1) && B.at(j) >= T.at(i - 1)) {
                rep(k, N) {
                    CHMAX(dp.at(i).at(j),
                          dp.at(i - 1).at(k) + abs(C.at(j) - C.at(k)));
                }
            }
        }
    }
    ll max = -INF;
    rep(i, N) { CHMAX(max, dp.at(D).at(i)); }
    cout << max << endl;
}

*1:後述するがこれでも前の情報を記録する必要があり,厳密にはDPである.ただ,計算の工夫によって大きく変わる事から紹介する.

*2:こういう場合,どのようなパターンが入力されるか分からないので計算量的に最悪を考えておくといい

AtCoder Beginner Contest 015_D_高橋くんの苦悩

Atcoder 水色 DP 動的計画法

DPの問題.ナップサック問題+α.スクリーンショットの重要度ってなんやねん.

atcoder.jp

問題の説明

枚数がK枚,幅がWに収まる条件の下で画像を選ぶときの重要度の合計の最大値を求めよ.

条件

  • 1≦W≦10,000
  • 1≦N≦50
  • 1≦K≦N

 

考えた事

枚数制限がない場合,これは典型的なナップサック問題と呼ばれるものとなる.

ja.wikipedia.org

(競プロでは)ナップサック問題はDPで解く事が一般的である*1.

garakutagoya.hatenablog.com

上の記事に記載の通り,DPの基本は変数の捜索にある.ナップサック問題でのDPによる使い方をかいつまんで話すと,荷物の添字をiとし,重量の総和をjとすると,「i番目までの荷物で重量j以下となるように荷物を選んだときの価値の最大値」としてdp_{i,j}テーブルを作ればよい.

今回であれば幅と重要度がそれに対応する.幅の制限を満たしながら重要度を最大とするように考えれば上のナップサック問題と何一つ変わらない.読者の皆様には簡単でしょう.

今回のポイントは枚数制限の要素が追加されている事だ.ナップサック問題でいうところの制限がもう1つ追加されている点に難しさがある.そこに注意してDPを考えてみよう.

ステップ1.全探索出来るかどうか

全探索する場合の計算量の見積もりをしよう.幅の制限は一先ず置いておき,0~K個選び,それ以外を選ばないくらいで考えると楽でいい.

\sum_{k=0}^K {}_nC_k=2^K

であるから計算量はO(2^K)となりお察しの通り宇宙時間である.無念.

ステップ2.状態の遷移を考える

この問題の場合,前までの状態に写真を足せるかどうかで見ればよいだろう.その場合,幅,重要度はもちろんのこと写真の枚数が1枚増える事に注目するのが今回の大きなポイントである.オタクの部屋と違って普通の人は部屋に物を増やしたがらない.

ステップ3.変数を見つけ,配列の表を書く

遷移する変数を考えてみよう.この場合鍵となるのは小さい範囲で状況の変化を考えてみる事だ.もしスクリーンショットを選んだとしよう,その場合に変化するものは枚数,幅の合計,重要度の合計,そして画像を選んだという事実である.枚数,幅,選択した画像が分かっていれば重要度の合計を求める事が出来る.したがって上の3つを変数として考えてみよう.

イメージとしてはf(枚数,幅の合計,使用可能な画像の種類)=(重要度の合計)と思ってもらえれば差し支えない.

枚数を0≦k≦K,使用可能な画像の添字を0≦i≦N,幅の合計を0≦j≦Wとして表を書いてみる.この場合,表を一度に書くと3次元を記述する事になってしまうため,kの値で分けて書くとよい*2.

f:id:hikki-dominion:20220416134833j:plain

入力例1を用いて作成した表を以上に載せた.遷移している部分を赤と青で記しているが,3次元であるためどこが遷移しているのかまだ難しいだろう.私自身もこれだけでは60がどこから来たものか分からず考えてしまった.

次のステップで遷移について考えてみよう.

ステップ4.漸化式の作成

ここからは漸化式の作成に入る.キーとなるのは状態の変化の仕方が画像を選んだ/選んでいないで二分される事だ.

まずはi枚目の画像を選べない場合を仮定しよう.この場合,j\lt A_iが条件となる.選べない場合,枚数も,幅の合計も一切変わらない.変わる部分といえばi枚目を使わないという事はi-1で考えても変わらない,iをi-1にしてもよい点だろう.

つまり,

  • k→k
  • i→i-1
  • j→j

としてよいのである.したがって,

dp_{k,i,j}=dp_{k,i-1,j}(j\lt A_i)

となる.

次にi枚目の画像を選べる場合を考える.もし選ばないならばk,i,jの遷移は選べない場合と全く変わらない.

選ぶ場合,使用可能な枚数は1枚減り,幅もA_i減る.i番目を使う事が分かっていれば,足し算の線形性と同じ画像は2度使わない制約から,後はi-1の場合で考えてよい.そして重要度はB_i上がる.そしてDPには選ぶ選ばないいずれかの重要度の合計の最大値を記載する事を繰り返せばいい*3.

そのため,

  • k→k-1
  • i→i-1
  • j→j-A_i

dp_{k,i,j}=max(dp_{k,i-1,j}, dp_{k-1,i-1,j-A_i}+B_i)(j≧A_i)

となる.まとめると,

dp_{k,i,j}=dp_{k,i-1,j}(j\lt A_i)

dp_{k,i,j}=max(dp_{k,i-1,j}, dp_{k-1,i-1,j-A_i}+B_i)(j≧A_i)

となる.初期条件としては

dp_{k,i,j}=0(k=0∨i=0∨j=0)

とすると分かりやすいだろう.

ステップ5.漸化式を元に配列を埋める

実装をするにあたってどうしたらいいだろう.ここまで読んでいただけた方ならお分かりの通り,3次元配列を用いる.

vvvll dp(K + 1, vvll(N + 1, vll(W + 1, 0)));

としてDP配列を作る.初期条件はdp_{k,i,j}=0(k=0∨i=0∨j=0)だけで十分であるが,(何もない)=0としても齟齬は生じないためそうしている.距離コストの最小値を求める場合,この部分がINFになる事もあるため注意が必要だ.

後は3重ループを回し,配列を埋めればよい.

if (j < ab.at(i - 1).at(0)) {
                    dp.at(k).at(i).at(j) = dp.at(k).at(i - 1).at(j);
                } else {
                    dp.at(k).at(i).at(j) =
                        max(dp.at(k).at(i - 1).at(j),
                            dp.at(k - 1).at(i - 1).at(j - ab.at(i - 1).at(0)) +
                                ab.at(i - 1).at(1));

漸化式の部分はこのようになる.計算量はO(K×N×W)となり,高々O(10^7)となるため,やや危ないが間に合う.

ステップ6.求めたい配列の値を見て答えを導く

dp_{K,N,W}を見て,それを出力すればいい.なかなかタフな問題であった.

プログラム

#include <bits/stdc++.h>
using namespace std;
#include <atcoder/all>
using namespace atcoder;
using ll = long long;
using ld = long double;
using P = pair<int, int>;
using Graph = vector<vector<int>>;
using vi = vector<int>;
using vl = vector<long>;
using vll = vector<long long>;
using vvi = vector<vi>;
using vvl = vector<vl>;
using vvll = vector<vll>;
using vvvll = vector<vvll>;
using vs = vector<string>;
using pii = pair<long long, long long>;
using mint = modint1000000007;
const long double EPS = 1e-10;
const long long INF = 1e18;
const long double PI = acos(-1.0L);
#define reps(i, a, n) for (ll i = (a); i < (ll)(n); i++)
#define rep(i, n) for (ll i = (0); i < (ll)(n); i++)
#define rrep(i, n) for (ll i = (1); i < (ll)(n + 1); i++)
#define repd(i, n) for (ll i = n - 1; i >= 0; i--)
#define rrepd(i, n) for (ll i = n; i >= 1; i--)
#define ALL(n) begin(n), end(n)
#define IN(a, x, b) (a <= x && x < b)
#define INIT                          \
    std::ios::sync_with_stdio(false); \
    std::cin.tie(0);
template <class T>
inline T CHMAX(T& a, const T b) {
    return a = (a < b) ? b : a;
}
template <class T>
inline T CHMIN(T& a, const T b) {
    return a = (a > b) ? b : a;
}

int main() {
    ll W;
    cin >> W;
    ll N, K;
    cin >> N >> K;
    vvll ab(N, vll(2));
    rep(i, N) {
        rep(j, 2) { cin >> ab.at(i).at(j); }
    }
    vvvll dp(K + 1, vvll(N + 1, vll(W + 1, 0)));
    rrep(k, K) {
        rrep(i, N) {
            rrep(j, W) {
                if (j < ab.at(i - 1).at(0)) {
                    dp.at(k).at(i).at(j) = dp.at(k).at(i - 1).at(j);
                } else {
                    dp.at(k).at(i).at(j) =
                        max(dp.at(k).at(i - 1).at(j),
                            dp.at(k - 1).at(i - 1).at(j - ab.at(i - 1).at(0)) +
                                ab.at(i - 1).at(1));
                }
            }
        }
    }
    cout << dp.at(K).at(N).at(W) << endl;
}

*1:他の解き方について知りたい方は

第2回：ナップサック問題を色々な方法で解いてみた【ブレインパッドの数理最適化ブログ】 - Platinum Data Blog by BrainPad

を参照されたい.

*2:平面毎にスライスした後,高さで積分し,立体の体積を求める問題を思い浮かべてもらうと分かりやすいだろう.大学受験ではしばしば出てくる

*3:これも線形性による

AtCoder_Beginner_Contest_240_C Jumping Takahashi

Atcoder 茶色 DP 動的計画法

いよいよDP.最初は優しめの問題から.

atcoder.jp

問題の説明

高橋君がN回a_iまたはb_i段のジャンプを行います.N回のジャンプの後にXにいる可能性があるか答えなさい.

条件

  • 1≦N≦100
  • 1≦a_i\lt b_i≦100
  • 1≦X≦10,000
  • 入力は全て整数.

 

考えた事

まず,DPとは何をする手法なのか.端的に言えば値を記憶する配列を作り,それを使い回すものであろう.この問題においてそれは出来るのだろうか.

例えば1回目にa_1段上がり,2回目にb_2段上がったものと,1回目にa_1段上がり,2回目にa_2段上がったものでは当然現在位置が異なる.Greedyと違って今回はX段まで上る最適な方法は存在しない(入力によって方法を変えなければならない)ため,X段に辿り着けるかどうかは全てを試さなければ分からない.

ところがこれを全探索で実装しようとするとTLEエラーとなってしまう.全てのiについてa_iとb_iを選んだ場合を考えると計算量はO(2^N)乗となり,N=100が上限では計算が間に合わない.

なぜ全探索ではダメなのか.無駄が多いのである.最後だけa_Nとb_Nで異なるが,それ以外は全く同じものを選んだ場合を考えてもらいたい.全探索ではこれを一々計算し直すため,結果として非常に遅い手法となってしまう.

これを解消するのがDPである.配列にN-1回目までに行ける部分を記憶させる事で,N回目に行ける部分はN-1回目までに行ける部分にa_Nまたはb_Nを足せば求められる.DPがメモ化再帰と呼ばれる所以である.*1裏を返せばN-1回の状態を記憶させる事に意味のない問題の場合,DPは不適当な解法だといえる.

f:id:hikki-dominion:20220411182404j:plain

状態の遷移を上の図*2に表した.i回目に座標jにいる場合は1を,そうでない場合は0を記している.0回目(登り始める前)は0にいるため,0にのみ1がつき,そうでない部分は0が記される.

1回目のジャンプを行うと0"+"a_1または0"+"b_1のジャンプを行う.ここで重要なのはi回目のジャンプで行ける先を求める場合,i-1回までの結果にそのままa_iまたはb_iを足せば求められるという点だ.これが可能であるからこそ,これを帰納的にN回繰り返してN回目に行ける場所が求められる.

少し前の繰り返しになるが,「繰り返し」が出来るかどうかがDPの使用可否を決める.どれを変数とすれば状態を遷移させられるか,繰り返しを作る事が出来るか,もし出来るならばどのような遷移となるか,これを考える事が非常に重要である.

競技プログラミングにおいては条件のところで変数のヒントが与えられている事も多い.今回においてパラメータ(変数)はi(1≦i≦N)とj(1≦j≦100×100=10,000)であり,N,X,a_i,b_iは定数である.問題文にパラメータとしてiは既に記されており,後は自分でjを置けばDPの遷移を考えられる.この知見は大事にしたい.

さあ実装だ.DPの実装は上の配列を埋める事に終始する.配列を埋めるにはどうすればいいのか.ここで出てくるのが漸化式である.

i-1回目までの遷移とi回目の遷移を結び付ける式があれば配列を埋めて行ける.DPの最終工程はこれである.

今回であればdp_{i-1,j}=1ならばdp_{i,j+a_i}=dp_{i,j+b_i}=1となる.それ以外は0となる.つまり,

dp_{i,j+a_i}=dp_{i,j+b_i}=1(dp_{i-1,j}=1)

dp_{i,j}=1=0(otherwise.)

となる.まとめてコードにすると,

vvll dp(N + 1, vll(10001, 0));
    dp.at(0).at(0) = 1;
    rrep(i, N) {
        rep(j, 10001) {
            if (dp.at(i - 1).at(j)) {
                dp.at(i).at(j + a.at(i - 1)) = 1;
                dp.at(i).at(j + b.at(i - 1)) = 1;
            }
        }
    }

となる.プログラム上は0からスタートしているため,iの添え字が1つずれているのには注意したい.0を配列のデフォルトにし,iとjを動かして,上の条件を満たすものだけに1をつけていく.

計算量に注目するとこれはO(N×10,000)であり,全探索のO(2^N)と比べ,大幅に短縮されている事が分かるだろう.

今回は制限に余裕があるため行っていないが,例えばj+a_{i-1}≦X,j+b_{i-1}≦Xをif文に追加し,X以上となった場合は考えない*3ように実装すれば計算量はO(N×X)となり更に短縮できる.

DPの流れをまとめると

  1. 全探索出来るかどうかを考える.(計算量が大きな基準となる)
  2. 状態が遷移するか,前の状態から今の状態を帰納的に導き出せるかを考える.(前の状態に足せるかどうかを見ると良さげ.)
  3. 変数を見つけ,小さな値の範囲でいいから配列を書いてみる.
  4. それを基に漸化式を作る.
  5. 漸化式を元に変数を動かし配列を全て埋めていく.
  6. 求めたい値に対応する配列を見て,答えを導き出す.

となるだろう.最後はもちろん問題によって異なるが今回であればdp_{N,X}である.これが1か0かがジャンプで辿り着けるかどうかを決める為,それを参照して出力すればいい.

プログラム

#include <bits/stdc++.h>
using namespace std;
#include <atcoder/all>
using namespace atcoder;
using ll = long long;
using ld = long double;
using P = pair<int, int>;
using Graph = vector<vector<int>>;
using vi = vector<int>;
using vl = vector<long>;
using vll = vector<long long>;
using vvi = vector<vi>;
using vvl = vector<vl>;
using vvll = vector<vll>;
using vs = vector<string>;
using pii = pair<long long, long long>;
using mint = modint1000000007;
const long double EPS = 1e-10;
const long long INF = 1e18;
const long double PI = acos(-1.0L);
#define reps(i, a, n) for (ll i = (a); i < (ll)(n); i++)
#define rep(i, n) for (ll i = (0); i < (ll)(n); i++)
#define rrep(i, n) for (ll i = (1); i < (ll)(n + 1); i++)
#define repd(i, n) for (ll i = n - 1; i >= 0; i--)
#define rrepd(i, n) for (ll i = n; i >= 1; i--)
#define ALL(n) begin(n), end(n)
#define IN(a, x, b) (a <= x && x < b)
#define INIT                          \
    std::ios::sync_with_stdio(false); \
    std::cin.tie(0);
template <class T>
inline T CHMAX(T& a, const T b) {
    return a = (a < b) ? b : a;
}
template <class T>
inline T CHMIN(T& a, const T b) {
    return a = (a > b) ? b : a;
}

int main() {
    ll N, X;
    cin >> N >> X;
    vll a(N);
    vll b(N);
    rep(i, N) {
        ll aa, bb;
        cin >> aa >> bb;
        a.at(i) = aa;
        b.at(i) = bb;
    }
    vvll dp(N + 1, vll(10001, 0));
    dp.at(0).at(0) = 1;
    rrep(i, N) {
        rep(j, 10001) {
            if (dp.at(i - 1).at(j)) {
                dp.at(i).at(j + a.at(i - 1)) = 1;
                dp.at(i).at(j + b.at(i - 1)) = 1;
            }
        }
    }

    if (dp.at(N).at(X)) {
        cout << "Yes" << endl;
    } else {
        cout << "No" << endl;
    }
}

*1:正確にはメモ化+分割統治法である.

*2:この図の事をDPテーブルと言う人もいる

*3:これを枝刈りという

AtCoder Beginner Contest 091_C_2D Plane 2N Points

Atcoder Greedy 水色

Greedyの問題.どうマッチングさせると仲良しさんが作られるか.

atcoder.jp

問題の説明

赤い点の外側にある青い点を結ぶことで作る事の出来るペアの個数の最大値を求めよ.

ただしX=(a, b),Y=(c, d)でa\lt cかつb\lt dのとき,XはYの外側にあるという.

条件

  • 1≦N≦100
  • 赤い点と青い点のx座標,y座標は全て相異なる.

 

考えた事

最初考えた事は,赤い点からどのように青い点を選べばペアを多く作れるかであった.そういった中で,なるべくx座標の差を小さくしたかったため,赤い点のx座標の小さい順から考え,ペアを作る事の出来る中でx座標が最小の青い点を選ぶことを繰り返すものであった.

しかしこれでは上手くいかない.文章からも明らかであるが上の状態ではy座標の情報が一切反映されていない.我々が目指している事は闇雲に外側の点を選ぶのではなく,外側は外側でもなるべく内側によっている点を選ぶことである.もっというとなるべくペア間の距離を短く結ぶことにある.

なぜなら距離が離れれば離れるほどその点はより外側にあるという事になり,ペアを作る効率がよくない.どれだけ(原点から)外側にあるかの度合いを測る尺度として距離があり,その距離の長短とペアが対応している事は興味深い.

f:id:hikki-dominion:20220405205304j:plain

イメージにおこすとこのようになる.上のxy座標が最初に考えたペアの作り方である.見て分かるようにどう考えても効率的でない.

翻って下はどうであろうか.ペア間の距離が短くなっているのもさることながら,何より作成可能なペアの個数が増えている.おそらくこれを目指す事が最善手であると考えた.

ではこのようにするためにはどうすべきであろうか.自分は青い点に注目した.x座標の小さい青い点から(外側から),ペアを作成可能な赤い点を探していく手法だ.

そういった中でどれとペアを結ぶべきか.やはりここでも距離が出てくる.距離の計算は\sqrt{(c-a)^2+(d-b)^2}で求める事が出来る.そして,今はx座標が小さい方から見ているため,これ以上x座標の距離を縮めようとしてもあまり距離は縮まらない.むしろ注目している青い点に対し,赤い点のy座標がギリギリまで近い点を選ぶことでその距離を最小化出来る.このようにしてGreedyを行うと上手く行くのである.

「あえてはじめは距離が離れているところを選ぶのはありじゃないか.」,と思うかもしれない.しかし,今考えているのはx座標が小さい順である事に注意したい.最初にあえて自分に近い点を選ばず,遠くの点を選んだとしよう.そうした時,次の青い点は[tex;x]座標が離れているぶんもっと遠くなる.そうなってしまった場合,その行為に何の得もない.

また,それによって次のペアの距離が縮まったと仮定しても,その場合はそもそもはじめの青い点は次のペアの赤い点とペアを結べないか,どちらの赤い点も入れ替えてペアを作る事が可能かのどちらかでしかないため,「はじめにあえて距離が離れているものを選ぶ」方針ではGreedyに対し,ペアの数を減らす事さえあっても増やす事はありえない.

加えてGreedyの基本姿勢が「その場での最善」を優先し,それを繰り返す思考であるから「その場での最善」を優先せず,かといって後の方で取り返せる保証のない方針は得策ではない.*1

考えなければならないのは一度ペアに使った点はもう使えないという事.今回は一度ペアになった赤い点のx座標を大きくし,彼方に飛ばす事でそれに対応した.

if (k != -1) {
            count++;
            ab.at(k).at(0) = 2000;
        }

青い点はループで全探索するため,ペアの成否を考えなくてよい.扱いやすい.

プログラム

#include <bits/stdc++.h>
using namespace std;
int main() {
    int N;
    cin >> N;
    vector<vector<int>> ab(N, vector<int>(2));
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < 2; j++) {
            cin >> ab.at(i).at(j);
        }
    }
    vector<vector<int>> cd(N, vector<int>(2));
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < 2; j++) {
            cin >> cd.at(i).at(j);
        }
    }
    sort(cd.begin(), cd.end());
    int count = 0;
    for (int i = 0; i < N; i++) {
        int y = -1;
        int k = -1;
        for (int j = 0; j < N; j++) {
            if (cd.at(i).at(0) > ab.at(j).at(0) &&
                cd.at(i).at(1) > ab.at(j).at(1) && ab.at(j).at(1) > y) {
                y = ab.at(j).at(1);
                k = j;
            }
        }
        if (k != -1) {
            count++;
            ab.at(k).at(0) = 2000;
        }
    }
    cout << count << endl;
}

*1:Atcoder社の公式解説ではこの部分がより丁寧に記述されている.

AtCoder Regular Contest 006_C_積み重ね

Atcoder Greedy 緑色

Greedyはブログにまとめてなかったなと思ったので.ゴミ屋敷の住人必見！？

atcoder.jp

問題の説明

高橋君がN個の箱を積み上げる時に出来る箱の山の最小個数を求めなさい.箱にはその重さw_i(1≦i≦N)以下の箱を重ねる事ができる.

条件

  • 1≦N≦500
  • 1≦w_i≦100,000

 

考えた事

Greedyの基本は「その場での最善」を繰り返していく事にあります.これを意識して考えるとこの問題は分かりやすいかもしれないです.

では,この問題における「その場での最善」とは何か.結論からいうと今重ねられる山の中で,1番軽いものに重ね,重ねられない場合は山を作る事を繰り返す事となる(方針A).

ここでポイントとなるのはなるべく重ねる事である.例えば8,1,7のように積み上げるとして,8の上に1を載せるのはちょっともったいなく感じるかもしれない(方針B).しかし,それは気のせいでしかない.

f:id:hikki-dominion:20220405154151j:plain

 

上(方針A)と下(方針B)の比較を見て分かる通り,1の上に8を載せようが載せまいが出来る状況は全く同じなのである.

また,

f:id:hikki-dominion:20220405154327j:plain

とした場合,上の方針Aでは山が1つなのに対し,下の方針Bでは山が2つ出来てしまう.

一応下の山の方には8が残っているというメリットがあるものの,例えば次に8以下の数字が来た場合,上ならば新しい山を作り,下ならばその8の上に重ねる事で全く同じ状況*1となり,8より大きな数が来た場合は重ねられないため,どちらも山が1つずつ増えてしまい,どちらにせよ8を取っておく意味が全くない事が分かるだろう.

以上から方針Bは却下され,方針Aが最善と分かる.*2大事なのは方針Bが上振れても方針Aと変わらず,順番次第では山一つ増やしてしまう事だ.これが分かれば後は山のてっぺんの数字を記録する配列を作り,てっぺんが変わったら書き換え,新たに山が出来たらcountを増やし,てっぺんの数字を記録するプログラムを書けば容易に解けるだろう.

プログラム

#include <bits/stdc++.h>
using namespace std;
int main() {
    int N;
    cin >> N;
    vector<int> w(N);
    for (int i = 0; i < N; i++) {
        int we;
        cin >> we;
        we--;
        w.at(i) = we;
    }
    vector<int> y(100000, 0);
    int count = 0;
    for (int i = 0; i < N; i++) {
        bool flag = false;
        for (int j = w.at(i); j < 100000; j++) {
            if (y.at(j) == 1) {
                y.at(j) = 0;
                y.at(w.at(i)) = 1;
                flag = true;
                break;
            }
        }
        if (flag == false) {
            count++;
            y.at(w.at(i)) = 1;
        }
    }
    cout << count << endl;
}

*1:8,1,7の順に来た場合を連想した方は聡明である.

*2:詳しい証明は背理法と帰納法により与えられる.背理法の肝は上に記している