CGOを使ったGoパッケージをRaspberryPy用にCrossCompileしたい
と、思ったのですが、Macでやろうと思うとなんだか面倒だったのでクロスコンパイル用のコンテナを作った。
かなりニッチだと思いますが、せっかく作ったので整理して下記に置いた。
github.com
使い方
事前にbuildしたいソース一式をgithubに配置して置く必要があります。
使い方は簡単で、まずコンテナをビルドし、
docker build . -t gorasp:build
そしてバイナリ出力先をマウントしつつ実行。(下記ではgorasp_builder直下のbuild)
docker run -it -v "$(pwd)"/build:/go/build gorasp:build
あとは、github repositryのURLを聞かれるので入力し、CGOを使っている場合はその旨入力すればクロスコンパイルされる。
余談
この作業をするに併せて、久々にラズパイのOSを書き込んだんですが昨今はimagerなんてものがあるんですね。
マジで楽すぎて感動した。文明の進歩に感謝...。
bit DP:TSP(ABC180E)
atcoder.jp
TSP!
Traveling Salesman Problem
巡回セールスマン問題(じゅんかいセールスマンもんだい、英: traveling salesman problem、TSP)は、都市の集合と各2都市間の移動コスト(たとえば距離)が与えられたとき、全ての都市をちょうど一度ずつ巡り出発地に戻る巡回路のうちで総移動コストが最小のものを求める(セールスマンが所定の複数の都市を1回だけ巡回する場合の最短経路を求める)組合せ最適化問題である。
全探索を行えば勿論最適解は求まるが、あり得る訪問経路の数は通りになるのですぐに組合せ爆発する。(例えばなら)
wikipediaにあるようにTSPはいわゆるNP困難な問題で、多項式時間で解けるアルゴリズムは見つかっていないので現実的な時間で最適解を求められるかどうかは都市の数()に依存する。
AtCoderの時間制約的にはが10以上程度になると既に全探索では厳しくなってくるが、程度ならばDPによる高速化で十分に解ける場合が多い。
ちなみに全探索版で解ける問題はABC183Cで出てましたね。
atcoder.jp
どう解くか
というわけでDPによる高速化をしていく。
重要なポイントは、現状態から終端までの最適な訪問順は「これまでどの都市を訪問してきたか」ということと「今いる都市」という情報のみに依存し、「これまで訪問した都市の訪問順序」には依らないということ。つまり残りの未訪問の都市集合と、スタート地点が同じなら、最適な経路は当然一意に定まる。
※このエントリの図が非常に分かりやすく参考になりました。
dalekspritner.hatenablog.com
よってあり得る状態は、「これまでどの都市を訪問してきたか」という状態と、「今いる都市」であるで表現でき、を状態-に至るまでの最小コストとするとを求めれば良いということになる。
遷移は、から今いる都市を除いた状態(現状態の前状態としてあり得る状態)をと表すと、は,()の中から最小のものを選べば良いので下記のようになる。
bitによる状態Sの表現
上記までで遷移は分かったので、あとは,を配列のindexで取り扱えるような値で表せればok。ここでがなんだったかを振り返ると、ある都市を訪問したか否かを表現できれば良いのでbitで取り扱えば良さそう。
例えば、都市1,3,4を訪問した場合だと下記のように表せば良い。
これで各状態を整数値で取り扱えるようになったので、前述の遷移を実装すれば終了。
計算量
あり得る状態数はであり、各状態に対して遷移を集約するので計算量はになる。
実装例
ここまでのロジックを忠実に実装すると下記のようになるが、実際最内ループの条件分岐は無くても動く。(初期値をINFにしていることでどちらにしても更新されないので。)
#include <algorithm> #include <cmath> #include <cstring> #include <deque> #include <iomanip> #include <iostream> #include <map> #include <queue> #include <set> #include <stack> #include <string> #include <vector> #define rep(i, n) for (int i = 0; i < (n); i++) #define reps(i, n, s) for (int i = (s); i < (n); i++) using ll = long long; using namespace std; constexpr long long INF = 1LL << 60; constexpr int MAX_N = 17; ll dp[(1 << MAX_N)][MAX_N]; int main() { cin.tie(0); ios::sync_with_stdio(false); int n; cin >> n; vector<ll> x(n), y(n), z(n); rep(i, n) cin >> x[i] >> y[i] >> z[i]; // calc cost vector<vector<ll>> cost(n, vector<ll>(n, 0)); rep(i, n) { reps(j, n, i + 1) { ll c = abs(x[i] - x[j]) + abs(y[i] - y[j]); cost[i][j] = c + max(0LL, z[j] - z[i]); cost[j][i] = c + max(0LL, z[i] - z[j]); } } // TSP fill(dp[0], dp[(1 << MAX_N)], INF); dp[0][0] = 0; reps(S, 1 << n, 1) { rep(v, n) { if ((S >> v) & 1) { rep(j, n) { ll S_hat = S - (1 << v); if (S_hat >> j & 1 || S_hat == 0) { dp[S][v] = min(dp[S][v], dp[S_hat][j] + cost[j][v]); } } } } } cout << dp[(1 << n) - 1][0] << endl; return 0; }
bit DP楽しい。
LongestCommonSubsequence(EDP-F)
atcoder.jp
有名なDPの問題なので知ってはいたが、ややこしそうなので後回しにしていた。
先日買ったけんちょんさん(@drken1215)の本を進めていたら章末問題に載っていたので、もりもり鍛えるべく重い腰を上げて解いてみた。
アルゴリズムとデータ構造をげっつしたのでもりもり鍛える pic.twitter.com/lmv1eW06JN
— KLEE (@kashee337) October 1, 2020
自分的には結構ややこしくて遷移の理解に苦しんだ。
文字列系アルゴリズムへの苦手意識が拭えない。
最長共通部分列問題(LCS Problem)
二つの異なる系列が与えられ、その共通部分列を答える問題。
ここでいう系列としては文字列などの記号列を取り扱う場合が多い。
また部分列は、順序は保持する必要はあるが隣接関係は崩れても良いというあたりがややこしい。
詳しくはwikiを参照。
最長共通部分列問題 - Wikipedia
どう解くか
dpの基本である部分問題に分割してそれを順に解いていくということをする。(当たり前..)
具体例として、EducationalDPコンテストの入力である二つの文字列,を考える。
まず、,の,文字目をそれぞれ,で表すとする。
そしてを,の最長共通部分列長さと定義する。
そうすると、,文字目に注目し①ならば、その文字を除いた前の文字列,に今の一致した分の長さを加えれば良いので、
となる。
次に②の場合を考えると、dp配列の値は,のLCSと,のLCSのうちより長いほうの値がそのまま現状態の長さになり、遷移は下記のようになる。
この遷移を適用し、dp配列を埋めていくと下記表のようになる。
ここで矢印は更新に適用された遷移であり、""は空文字を表している。
この表の上では斜めの遷移が①(注目文字の追加)の遷移、縦横の遷移が②(注目文字の削除)となっている。
これは,を一致させるにはを削除する必要があることを考えると分かりやすい。
実装例
遷移がここまで書ければあとは実装するだけだが、上記までの処理ではLCSの長さしか分からない。
復元が求められる場合は、表の右下から見ていって丸が着いた部分に遭遇したら斜め左上に遷移し、そうで無い部分のときは縦横のうち自信と同じ値を持つ方向に遷移するというポリシーで表の左上まで辿っていけばLCSそのものが得られる。
この時、縦横のうち両方とも同じ値の場合にどちらを優先するかで復元されるLCSは異なる物になる場合がある(LCSは一般に複数パターンあり得る)。
復元も含めて実装した例は下記のようになる。
#include <iostream> #include <algorithm> #include <string> #include <vector> #define rep(i, n) for (int i = 0; i < (n); i++) using namespace std; int main() { string s,t; cin>>s>>t; int n=s.size(); int m=t.size(); vector<vector<int>>dp(n+1,vector<int>(m+1,0)); //各状態でのLCS長さをメモ rep(i,n){ rep(j,m){ if(s[i]==t[j]){//同じ文字なら dp[i+1][j+1]=dp[i][j]+1; } else{ dp[i+1][j+1]=max(dp[i][j+1],dp[i+1][j]); } } } //後ろから辿っていきLCSを復元 int i=n,j=m; string ans=""; while(i>0&&j>0){ if(s[i-1]==t[j-1]){ ans+=s[i-1]; i--;j--; } else{ if(dp[i][j]==dp[i-1][j]){//横方向優先 i--; } else{ j--; } } } //後ろ向きにメモしたので逆順に直す reverse(ans.begin(),ans.end()); cout<<ans<<endl; return 0; }
ややこしい...
DP~imos法による高速化(ABC179D)
素朴に配るDPを実装するとになりTLEとなる(なってしまった....)。
例えば下記の入力だと遷移数がマス目の長さと同じオーダになるので通らない。
200000 1 1 200000
こういう時は累積和でdpを高速化するのが定番らしい。
解説動画が非常に分かりやすいが、貰うdpで解説されていたので自分が本番の時に書こうとしていた配るdpで実装する場合どうすれば良いのか調べてACしたのでメモする。
■解説動画
www.youtube.com
結論から書くとimos法で累積和をシミュレートしながら更新していけば良い。
imos法による配るdpの高速化
imos法自体については考案者の方のページが非常に分かりやすい。
imoz.jp
自分の理解だと、状態が遷移する境界値だけメモしておきシミュレートすることで区間系累積和の高速化などに応用できる手法という理解。
区間数を,マス目の長さをとすると、境界値のメモは境界を生成する区間数に比例するとなり、シミュレートは頭から足してくだけなのでとなる。
よって全体としてはで累積和が計算できる。
今回のD問題ではマス目の頭から見ていき、各位置毎に自分の位置から右側へ区間の境界値を配りシミュレートしながら更新していけば良い。
よってマス目の長さ分だけ遷移区間K個分の境界値のメモを行うので、全体のオーダはとなり間に合う。
実装例
全体のコードを載せるとロジックが見にくいのでmain部分のみ。全コードはここ。
int main() { int n, k; cin >> n >> k; vector<pair<int, int>> v(k); rep(i, k) { cin >> v[i].first >> v[i].second; } vector<mint> dp(n + 1, 0); dp[0] = 1; reps(i, n, 1) { //現在マスの値を更新 dp[i] += dp[i - 1]; //現在マスの値を配る for (auto p : v) { //境界インデックスを計算 int l = i + p.first; int r = i + p.second + 1; if (l > n) continue;//左境界の領域外判定 dp[l] += dp[i]; if(r > n+1) continue;//右境界の領域外判定 dp[r] -= dp[i]; } } cout << dp[n] << endl; return 0; }
累積和用の配列とdp用の配列を共通にして実装したが、バグを仕込みかねないので分けた方がいいのかもしれない。
包除原理(ABC172E)
包除原理を応用する問題が出題されて本番中は手も足も出なかった。
いろんな解説を見てACしたので理解をまとめる。
atcoder.jp
包除原理とは
ものすごくざっくり書くと、複数の集合から構成される全体集合を重複なく求めるための考え方。基本的には全体を足して、重複する部分に関しては別途足し引きして辻褄を合わせる。
詳細は下記を参照。
mathtrain.jp
どう応用したか
今回の例では、ユニークな数字で構成される2つの数列の内、任意のでとなるパターンを数え上げる。
このような数列を直接求めるのは難しいので、全体の数列の組み合わせを全通り求めてそこから条件を満たさないものを除いていくことを考える。
全体集合
全体の数列を表す集合は、m個の中からn個を2数列分並べればよいので次の式で求まる。
条件を満たさない集合
今回の条件を満たさない集合とは、任意のでとなるが1つ以上ある集合のことなのでこれを数え上げていく。
ここで、直接「ちょうど箇所ので重複する集合」を求められれば良いが、これも難しい。「少なくとも箇所ので重複する集合」なら数えやすそうだが、これを全て足してしまうと重複する部分が発生する。例えば「少なくとも箇所ので重複する集合」には「少なくとも箇所ので重複する集合」も含まれている。
ここで、包除原理を使って重複を削除し「箇所ので値が重複する集合,」を数え上げればよい。
上記を計算するために、「少なくとも箇所ので重複する集合」について考える。
これを数え上げるのは簡単で、次のように考える。
- まず、長さの数列の中で重複させる位置を箇所選ぶ ->
- 次に、その箇所に重複させる数字を個の中から選んで並べる ->
- 最後に、残ったの位置に当てはめる数字を個の中から選んで並べる(数列2つ分を選ぶため二乗となることに注意) ->
よって重複がのときの集合は下記のようになる
ここまでくれば、を求めるにはを包除原理に従って足し合わせればよい。
よって今求めたい集合は全体からこれを引けばよいので、
となる。ここで、のときは、
となるので、上式を整理すると結局、
を求めればよいことになる。
実装例
#include <bits/stdc++.h> #define rep(i, n) for (int i = 0; i < (n); i++) #define reps(i, n, s) for (int i = (s); i < (n); i++) using ll = long long; using namespace std; constexpr long long MAX = 5100000; constexpr long long INF = 1LL << 60; constexpr int MOD = 1000000007; ll fac[MAX], finv[MAX], inv[MAX]; void COMinit() { fac[0] = fac[1] = 1; finv[0] = finv[1] = 1; inv[1] = 1; for (int i = 2; i < MAX; i++) { fac[i] = fac[i - 1] * i % MOD; inv[i] = MOD - inv[MOD % i] * (MOD / i) % MOD; finv[i] = finv[i - 1] * inv[i] % MOD; } } ll nCk(int n, int k) { if (n < k) return 0; if (n < 0 || k < 0) return 0; return fac[n] * (finv[k] * finv[n - k] % MOD) % MOD; } ll nPk(ll n, ll k) { if (n == 0 || k == 0) return 1; return fac[n] * finv[n - k] % MOD; } class mint { long long x; public: mint(long long x = 0) : x((x % MOD + MOD) % MOD) {} mint operator-() const { return mint(-x); } mint& operator+=(const mint& a) { if ((x += a.x) >= MOD) x -= MOD; return *this; } mint& operator-=(const mint& a) { if ((x += MOD - a.x) >= MOD) x -= MOD; return *this; } mint& operator*=(const mint& a) { (x *= a.x) %= MOD; return *this; } mint operator+(const mint& a) const { mint res(*this); return res += a; } mint operator-(const mint& a) const { mint res(*this); return res -= a; } mint operator*(const mint& a) const { mint res(*this); return res *= a; } mint pow(ll t) const { if (!t) return 1; mint a = pow(t >> 1); a *= a; if (t & 1) a *= *this; return a; } // for prime MOD mint inv() const { return pow(MOD - 2); } mint& operator/=(const mint& a) { return (*this) *= a.inv(); } mint operator/(const mint& a) const { mint res(*this); return res /= a; } friend ostream& operator<<(ostream& os, const mint& m) { os << m.x; return os; } }; int main() { cin.tie(0); ios::sync_with_stdio(false); ll n, m; cin >> n >> m; COMinit(); // nCk*mPk*(m-kPn-k)^2 mint ans = 0; rep(k, n + 1) { mint now = nPk(m - k, n - k); now *= now; now *= nPk(m, k); now *= nCk(n, k); if (k % 2 != 0) now = -now; ans += now; } cout << ans << endl; return 0; }
参考
以下、参考にさせて頂いたサイトです。
www.youtube.com
qiita.com
回文検出(Manacher)
何を解きたいか
ご存じのように回文とは,のような逆順でも変化しない文字列のことです。このとき、ある文字列の部分文字列で作られる回文の最大長さや数を問う問題を解きたい場合を考えます。ナイーブに解こうとすると二重のfor文で全ての部分文字列を生成し、をチェックすればよいですが、となるので短い文字列にしか適用できません。
この問題を回文の特性を上手に利用してで解く方法があり、それが今回説明するManacherのアルゴリズムです。
どう解くか
上述したように回文の特性を上手いこと利用します。 説明を簡単にするために奇数長の回文のみを検出対象とします。偶数への拡張は別途最後にコード例を示します。
今回は、具体例としてという文字列で考えます。
この文字列を頭から見ていって、各文字を中心に両端に何文字同じ文字を持つかをメモしていきます。これは、メモされた値をとすると、その文字を中心にの長さの回文があるという情報をメモしたことになります。
具体的に図で見ていくと次のような形で処理が進んでいきます。
◆STEP1
1文字目はで、右側にしか文字が存在せず回文の長さは自分自身だけのため、1をメモします。
◆STEP2
2文字目はで、両側にが存在するため、回文の長さは自分自身+1として2をメモします。
◆STEP3
3文字目はで、両側には同じ文字は存在しないため、1をメモします。
◆STEP4
4文字目はで、両側は矢印のの部分まで同じ文字が続くので、4をメモします。
この先STEP5,6と進めていきたいところですが、文字列をよく見るとSTEP4でを中心にという回文が構成されていることが分かっているため、STEP5,6の結果はSTEP2,3と同じになり探索不要なことが分かります。
◆STEP5,6
このように、ある位置で長さの回文が検出されたときまでの計算結果を用いることで探索を効率化しようというのが基本アイデアです。
ただここで、利用できる過去の計算結果の範囲には気を付ける必要があります。
例えばSTEP7ではわかりやすくそういったことが起きます。
◆STEP7
このSTEPのメモ値は3となり、着目文字の対称にある左端ののメモ値とは一致しません。
ある文字を中心に一文字ずつ両側に同じ文字があるかどうかを探索していくとき、現在着目してる文字を含む最長の回文(STEP7では)より内側までしか対称性が保証されません。そのため、最長の回文より外の範囲へ探索範囲が伸びた時点でそこから先は未知の世界のため、真面目に探索する必要があります。
実装例
これらを実装したコード例を下記に示します。
#include <bits/stdc++.h> using namespace std; int main() { string s = "abacababa"; int n = s.size(); vector<int> rad(n); int c = 0, r = 0; while (c < n) { // cを中心に同じ文字がどこまで連続するか while (0 <= c - r && c + r < n && s[c - r] == s[c + r]) r++; rad[c] = r; //回文の長さに応じて利用可能な範囲を確認しつつメモ int k = 1; while (0 <= c - k && k + rad[c - k] < r) { rad[c + k] = rad[c - k]; k++; } //すでに計算が終わった分だけ中心と探索半径をずらす c += k; r -= k; } int Lmax = *max_element(rad.begin(), rad.end()) * 2 - 1; cout << Lmax << endl; return 0; }
偶数長への拡張
偶数長回文への拡張は簡単で、元の文字列の各文字の間にを挟んで同じアルゴリズムで検出すればよいです。そうするとのインデックスには偶数長の回文の長さが、元の文字のインデックスには奇数長の回文の情報がメモされます。
ただし、最後に分の長さを差し引きしなければならないことに注意が必要です。
#include <bits/stdc++.h> using namespace std; int main() { string _s = "abacabaaa"; int n = _s.size(); //'#'の追加 string s(2 * n + 1, '#'); for (int i = 0; i < n; i++) s[2 * i + 1] = _s[i]; n = 2 * n + 1; vector<int> rad(n); int c = 0, r = 0; while (c < n) { // cを中心に同じ文字がどこまで連続するか while (0 <= c - r && c + r < n && s[c - r] == s[c + r]) r++; rad[c] = r; //回文の長さに応じて利用可能な範囲を確認しつつメモ int k = 1; while (0 <= c - k && k + rad[c - k] < r) { rad[c + k] = rad[c - k]; k++; } //すでに計算が終わった分だけ中心と探索半径をずらす c += k; r -= k; } //'#'分の補正 for (int i = 0; i < n; i++) { if (i % 2 == 0) rad[i] = (rad[i] - 1) / 2; else rad[i] /= 2; } int Lmax = *max_element(rad.begin(), rad.end()) * 2 - 1; cout << Lmax << endl; return 0; }
参考
以下、参考にさせて頂いたサイトです。
snuke.hatenablog.com