Codeforces 739D Recover a functional graph
変な問題だけど解法はそこそこきれい。
パス長とサイクル長のペアたちが満たすべき条件は、
- 各サイクル長に対し、パス長0の点の数はサイクル長の倍数(サイクルをいくつつくるかの条件)
- 各サイクル長に対し、パス長xの点が存在しないならパス長x+1の点も存在しない
である。この条件を満たすときの構成は容易。
まず、サイクル長が不定でパス長が定まっているような頂点であって、サイクル長もパス長も定まっているどの頂点のパス長よりもパス長が長いものたちは、全部同じサイクルにくっつけるのがよいことがわかる。このときのサイクル長を全探索する。
すると、条件は、「まだパス長やサイクル長が定まっていない頂点たちから、パス長xでサイクル長yの点が何個欲しい」という形で書ける。これを列挙し、パス長やサイクル長を定めていない頂点から適切に辺をはって最大流を流す。欲しい頂点と使える頂点が全部対応づいていればあとは適切に構成すればOK。
パス長もサイクル長も定まっている頂点は何をしてもOKなのでわざわざ辺をはらずに足りないときにとってきて使うことにすると辺の数が減り、計算量は全体でO(N^3)となる。
最初はもうちょっとgreedyでいけると思ってもう少し小さくフローを使おうとしていたが、嘘なことに気付き、そのあと全体をフローに突っ込んでしまえばいいことに気付いた。迷走である。ICPCでやったら破門。
パス長0,正,不定、サイクル長0,不定の6通りを全部違う変数で持とうとしたら実装が大変なことになったのでまとめて書き直した。これでも大変だけどだいぶましになってると思う。
本番でAC1人とかでもまあ時間かければ解けるっぽいけど、時間をかけないと解けないのでダメ。4時間くらいかかった。嘘に走ったり実装方針を立てるのに失敗せずにスパッと最短距離で実装できる能力が欲しいなぁ。
#include<stdio.h> #include<vector> #include<algorithm> #include<string> #include<iostream> #include<stdlib.h> using namespace std; int get(string s) { if (s == "?")return 400; int t = 0; for (int i = 0; i < s.size(); i++)t = t * 10 + s[i] - '0'; return t; } typedef pair<int, int>pii; vector<int>dat[500][500]; vector<int>fr; vector<int>zdat[500][500]; vector<int>zfr; int ans[300]; #define SIZE 1000 vector<int>pat[SIZE]; vector<int>cap[SIZE]; vector<int>rev[SIZE]; void adde(int a, int b, int c) { pat[a].push_back(b); pat[b].push_back(a); cap[a].push_back(c); cap[b].push_back(0); rev[a].push_back(pat[b].size() - 1); rev[b].push_back(pat[a].size() - 1); } int frv[SIZE], fre[SIZE]; bool flag[SIZE]; void dfs(int node) { if (flag[node])return; flag[node] = true; for (int i = 0; i<pat[node].size(); i++) { if (cap[node][i]>0 && (!flag[pat[node][i]])) { frv[pat[node][i]] = node; fre[pat[node][i]] = i; dfs(pat[node][i]); } } } int nod; int maxflow(int st, int go) { int ret = 0; for (;;) { fill(flag, flag + nod, false); dfs(st); if (!flag[go])return ret; int mini = 1000000000; int now = go; for (;;) { if (now == st)break; int t = fre[now]; now = frv[now]; mini = min(mini, cap[now][t]); } ret += mini; //flow++; now = go; for (;;) { if (now == st)break; int t = fre[now]; int nx = now; now = frv[now]; cap[now][t] -= mini; cap[nx][rev[now][t]] += mini; } } } int main() { int num; scanf("%d", &num); int mmxx = 0; for (int i = 0; i < num; i++) { string za, zb; cin >> za >> zb; int a = get(za), b = get(zb); if (a < 400 || b < 400)zdat[b][a].push_back(i); else zfr.push_back(i); if (b < 400 && a < 400)mmxx = max(mmxx, a); } vector<pii>adv; for (int i = mmxx + 1; i <= num; i++) { for (int j = 0; j < zdat[400][i].size(); i++) { adv.push_back(make_pair(i, zdat[400][i][j])); } zdat[400][i].clear(); } //for (int i = 0; i < adv.size(); i++)printf("....%d %d\n", adv[i].first,adv[i].second); for (int p = 1; p <= num; p++) { for (int i = 0; i < 500; i++)for (int j = 0; j < 500;j++)dat[i][j] = zdat[i][j]; fr = zfr; int sum = 0; for (int i = 0; i < 500; i++)sum += dat[p][i].size(); for (int i = 0; i < adv.size(); i++)dat[p][adv[i].first].push_back(adv[i].second); vector<pii>nee; for (int i = 1; i <= num; i++) { int s = 0; for (int j = 0; j < 500; j++)s += dat[i][j].size(); if (s == 0)continue; int maxi = 0; for (int j = 1; j <= num; j++)if (dat[i][j].size() > 0)maxi = j; for (int j = 1; j <= maxi; j++)if (dat[i][j].size() == 0)nee.push_back(make_pair(i, j)); int x = max(1, int(dat[i][0].size() + i - 1) / i)*i; for (int j = 0; j < x - dat[i][0].size(); j++)nee.push_back(make_pair(i, 0)); } if (nee.size() > num)continue; for (int i = 0; i < SIZE; i++) { pat[i].clear(); rev[i].clear(); cap[i].clear(); } nod = SIZE; for (int i = 0; i < nee.size(); i++) { adde(i, nee[i].first + 330, 1); adde(i, nee[i].second + 660, 1); } for (int i = 0; i < nee.size(); i++)adde(998, i, 1); for (int i = 0; i <= num; i++)adde(i + 330, 999, dat[i][400].size()); for (int i = 0; i <= num; i++)adde(i + 660, 999, dat[400][i].size()); int ret = maxflow(998, 999); if (ret + fr.size() < nee.size())continue; //printf("%d:\n", p); for (int i = 0; i < nee.size(); i++) { //printf(" %d %d\n", nee[i].first, nee[i].second); bool f = false; for (int j = 0; j < pat[i].size(); j++) { if (pat[i][j] != 998 && cap[i][j] == 0) { if (pat[i][j] < 660) { dat[nee[i].first][nee[i].second].push_back(dat[nee[i].first][400][dat[nee[i].first][400].size() - 1]); dat[nee[i].first][400].pop_back(); } else { dat[nee[i].first][nee[i].second].push_back(dat[400][nee[i].second][dat[400][nee[i].second].size() - 1]); dat[400][nee[i].second].pop_back(); } f = true; break; } } if (!f) { dat[nee[i].first][nee[i].second].push_back(fr[fr.size() - 1]); fr.pop_back(); } } for (int i = 0; i <= num; i++) { for (int j = 0; j < dat[i][400].size(); j++)dat[i][1].push_back(dat[i][400][j]); if (i != 0)for (int j = 0; j < dat[400][i].size(); j++)dat[p][i].push_back(dat[400][i][j]); else for (int j = 0; j < dat[400][i].size(); j++)dat[1][i].push_back(dat[400][i][j]); } for (int i = 0; i < fr.size(); i++)dat[1][0].push_back(fr[i]); //for (int i = 0; i <= num; i++)for (int j = 0; j <= num; j++)for (int k = 0; k < dat[i][j].size();k++)printf("%d %d %d\n", i, j, dat[i][j][k]); for (int i = 1; i <= num; i++) { for (int j = 0; j < dat[i][0].size(); j++) { ans[dat[i][0][j]] = dat[i][0][(j + 1) % i + (j / i)*i]; } for (int j = 1; j <= num; j++) { for (int k = 0; k < dat[i][j].size(); k++) { if (dat[i][j - 1].size() == 0) { goto l01; } //printf(" %d %d\n", dat[i][j][k], dat[i][j - 1][0]); ans[dat[i][j][k]] = dat[i][j - 1][0]; } } } for (int i = 0; i < num; i++)printf("%d ", ans[i] + 1); printf("\n"); return 0; l01:; } printf("-1\n"); }
Codeforces 750H New Year and Snowy Grid
新年初AC。GoodByeのボス。今年(去年)のは骨がない。
s-t間に2本の辺素パスが取れるという条件は、Mengerの定理より1点を取り除くとs-tが非連結になるという条件と同値なので、1点白いマスを除いてs-tを非連結にできるか判定すればいい。これは、追加する黒い点の周りの連結成分たちだけを見て8近傍でつなげていって、右上と左下からつづく連結成分が十分近くにあるかをうまいこと前処理とかして判定してやればOK。
簡単な割に解かれていない。TLが無駄に長いの何なんだろう。コードは長いけど実装はつらくない。
fflushを忘れてILE(初めてみた)、左上と右下の連結成分を眺めるリストに入れるのを忘れてWA。
#include<stdio.h> #include<vector> #include<algorithm> #include<set> using namespace std; typedef pair<int,int>pii; #define SIZE 2000000 class unionfind { public: int par[SIZE]; int ran[SIZE]; void init() { for(int i=0;i<SIZE;i++) { par[i]=i; ran[i]=0; } } int find(int a) { if(a==par[a])return a; else return par[a]=find(par[a]); } void unite(int a,int b) { a=find(a); b=find(b); if(a==b)return; if(ran[a]>ran[b]) { par[b]=a; } else { par[a]=b; } if(ran[a]==ran[b])ran[b]++; } }; unionfind uf; int map[1006][1006]; int dat[1006][1006]; void dfs(int x,int y,int c) { if(dat[x][y]!=0)return; dat[x][y]=c; for(int i=-1;i<=1;i++)for(int j=-1;j<=1;j++)if(map[x+i][y+j])dfs(x+i,y+j,c); } int main() { int mx,my,query; scanf("%d%d%d",&mx,&my,&query); for(int i=5;i<=my+3;i++)map[2][i]=1; for(int i=2;i<=my;i++)map[mx+3][i]=1; for(int i=5;i<=mx+3;i++)map[i][2]=1; for(int i=2;i<=mx;i++)map[i][my+3]=1; for(int i=3;i<=mx+2;i++) { for(int j=3;j<=my+2;j++) { char z; scanf(" %c",&z); if(z=='#')map[i][j]=1; } } int pt=1; for(int i=0;i<=mx+5;i++) { for(int j=0;j<=my+5;j++) { if(map[i][j]!=0&&dat[i][j]==0)dfs(i,j,pt++); } } set<pii>se; for(int i=0;i<=mx+5;i++) { for(int j=0;j<=my+5;j++) { if(map[i][j]==0)continue; for(int x=-2;x<=2;x++) { for(int y=-2;y<=2;y++) { if(map[i+x][j+y]!=0&&dat[i][j]!=dat[i+x][j+y]) { se.insert(make_pair(dat[i][j],dat[i+x][j+y])); } } } } } /*for(int i=0;i<=mx+6;i++) { for(int j=0;j<=my+6;j++) { printf("%c",'A'+dat[i][j]); } printf("\n"); }*/ int st=dat[2][my+3],go=dat[mx+3][2]; uf.init(); for(int p=0;p<query;p++) { int n; scanf("%d",&n); vector<pii>vec; for(int i=0;i<n;i++) { int za,zb; scanf("%d%d",&za,&zb); za+=2,zb+=2; dat[za][zb]=pt++; map[za][zb]=1; vec.push_back(make_pair(za,zb)); } vector<int>lis; for(int i=0;i<n;i++) { int nx=vec[i].first,ny=vec[i].second; for(int x=-1;x<=1;x++) { for(int y=-1;y<=1;y++) { if(map[nx+x][ny+y]) { uf.unite(dat[nx][ny],dat[nx+x][ny+y]); lis.push_back(dat[nx+x][ny+y]); } } } } lis.push_back(st); lis.push_back(go); if(uf.find(st)==uf.find(go)) { goto l01; } for(int i=0;i<n;i++) { int nx=vec[i].first,ny=vec[i].second; if(uf.find(dat[nx][ny])!=uf.find(st)&&uf.find(dat[nx][ny])!=uf.find(go))continue; for(int x=-2;x<=2;x++) { for(int y=-2;y<=2;y++) { if(map[nx+x][ny+y]) { if((uf.find(dat[nx][ny])==uf.find(st)&&uf.find(dat[nx+x][ny+y])==uf.find(go))||(uf.find(dat[nx][ny])==uf.find(go)&&uf.find(dat[nx+x][ny+y])==uf.find(st))) { goto l01; } } } } } for(int i=0;i<lis.size();i++) { for(int j=0;j<lis.size();j++) { if((uf.find(lis[i])==uf.find(st)&&uf.find(lis[j])==uf.find(go))||uf.find(lis[i])==uf.find(go)&&uf.find(lis[j])==uf.find(st)) { if(se.find(make_pair(lis[i],lis[j]))!=se.end())goto l01; } } } printf("YES\n"); if(false) { l01:; printf("NO\n"); } for(int i=0;i<lis.size();i++) { uf.par[lis[i]]=lis[i]; uf.ran[lis[i]]=0; } for(int i=0;i<vec.size();i++) { int x=vec[i].first,y=vec[i].second; uf.par[dat[x][y]]=dat[x][y]; uf.ran[dat[x][y]]=0; dat[x][y]=0; map[x][y]=0; } fflush(stdout); } }
ICPC アジア地区予選つくば大会 参加記
参加記です。Hziwara, hogloidとCxiv-Dxivで参加していました。
- 一日目
家を出る。まだTXにすら乗っていない時間にもうつくばについている人々が大量にいて集合時間を間違えていないか不安になる。調べたら大丈夫そうなので北千住でラーメンを食べる。
https://tabelog.com/tokyo/A1324/A132402/13161350/
食べログで評価が結構高かったので行った。北千住から徒歩5分。
煮干ラーメンなんだけど、煮干特有のきつい感じがなくて、でもしっかり煮干のうまみがある。たまに柚子が口に入って、これがまたスープに合う。おいしかったです。
なんだかんだでつくばに着く。Tシャツを着る前にKLab魔剤を飲みに行ったら係の人に「Tシャツ着ろ」って言われる。何を思ったのか、「But I want to drink it.」とか言ってそのまま飲む(意味不明)。
practiceは暇なのでおやつを食べたりほかのチームにちょっかいを出しに行ったりする。Javaチャレンジがないらしく結構早く終わった。
そのあと懇親会とかあって、みさわさんとフローの話をしたりきゅうりの服で遊んだり害悪しているうちに終わり、ホテルに着く。まともなホテルで何より。
ARCはFだけ考えて、わかったけど実装する時間はなかった。寝る。翌日ちゃんと起きるために目覚ましを5回くらいかける。
- 二日目
起床に成功する。なんと朝食を食べる。この時点で勝利である。(そんなことはない。)
なんだかんだでコンテストが始まる。英字キーボードに慣れていないのでGから読む。正確にいうと、
hogloid「最初のほうDEGwerがやってもしょうがないだろうから真ん中辺から読んで」
ぼく「はいじゃあFくらいから読みます」
コンテスト開始
Fの問題文が長いことに気付いたぼく「Gから読みます」
最悪である。
Gは読んですぐ解法がわかり若干不安になる。制約500000を50000だと思い込んで動的構築segtreeで二分探索とか言う。詰めずに後ろを読み進める。←あまりよくない
最初のほうが通されている。少しパソコンが開いたすきに書き始めるが、途中でさすがに量が多いので詰め直したほうがいいという結論に至り、交代して真面目に詰める。setでいいことに気付く。
再びパソコンに触れ、細かいことを口に出しながら実装するがサンプルが合わない。印刷するとバグに気付くが出すとREとか言われる。結局印刷デバッグを何回かやったら通った。3回もWAを出して、クソだなぁみたいな気分になる。
後ろのほうは読み進めた時とかに若干方針を考えたものがあったので、Iの雰囲気をふじわらさんに伝えてそれっぽいものを書いてもらう(証明とかは全然していない)。四角形のほう(のちょっとした一般化)しか渡しておらず、そりゃ通るわけがない。Jはざっと見たときにどうせ山登ればいけるだろう(幾何なので僕は実装しないけど)と思っていたので、まあ幾何やって山登りするだけなので時間できたらやっておいてくださいみたいなことを言う。
その間にhogloidと2人でHを考える。無向辺からなる木を縮約してDAGかどうか見るのでinfinite判定ができるというのを思いつき、解けそうな見た目の問題になったのでhogloidに投げると、hogloidが残りの部分を思いついてhogloidが実装をはじめる。量は多いけどまあ通せるやろみたいな気分でいる。
ふじわらさんとIがバグなのかどうかわからんなぁとか言っているうちに、縦横が互いに素に近い場合そもそも四角形にする必要がないやんということ(当たり前である)に気付き、それがわかるとオーダー解析の方針も一瞬で見える。三角形の場合と四角形の場合を両方使うとどちらもsqrtくらいのオーダーになるという話だった。見切り発車からの方針転換が起こっており危険なのだが、幸いちょっと変えるだけで済むらしい。ふじわらさんに伝える。
なんやかんやでHIが通り、JKが残ったのでふじわらさんがJをやる。Kをひたすら考える。全く分からないのであきらめて3並列でJのデバッグをしようとか思ったくらいの時刻にJが通る。のこりは暇である。凍結前トップで、終了直前にSJTUのFHの提出が両方増えており、ペナルティ的にもこれが両ACでなければ1位だという気分になる。
あんぱんがおいしかった。
解説は長く暇だった。聞きたかったKは、知識ゲーで解けるわけなしだった。SJTUはなんでこんなものを知っているんだ。
結果発表である。SJTUのどちらかが間に合っていないのを祈っていたし終了直前に提出が増えているので通っていないと思っていたが、両方通っていた。終了9分前と6分前にACしたらしく、負けを認めた。まあ選抜的には中国は関係ないので選抜としては実質1位(これはコンテスト終了時にはわかっていた)なので喜ぶ。
あとは企業ブースで問題を解き続けておなかがすいたのでそこらへんに生えているペンネを食べて、電車に長々と乗っておしまい。
去年と違って中難度が大量に並んでいたので得意セットだった。終了15分前くらいにできることがなくなるという一番いい形だった。得意分野が違うので3人いれば中難度は誰かは解けて、チーム戦は最高という気分になった。なにせ結局1問しか実装してないからね。一人ハマってもその間他が進み続けるからなんだかんだでカバーできるというのもチーム戦ならではの安心要素で、そういうふるまい方まで含めてやっぱり面白いなぁと。HとかIとか問題として面白かったし。特にHは大好き(しかも実装すらせずに済んで最高)。
まあ個人のふるまいとしてはGで無駄にペナルティを使ったので実装面ではあまりよくなかったし反省要素はあるけど、結局完数単独の2位だったし、解法は結構出したのでそこそこ貢献できたと思っています。まあなんだかんだで結果が良かったし満足しています。
WFがんばるぞい。
AOJ 2347 Sunny Graph
組合せ最適化で最大マッチングのところを読んで実装したので、解いた。
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2347
問題概要: グラフがある。頂点1を含む長さ3以上の奇サイクルをとり、残りの頂点で完全マッチングを作れるか判定せよ。
解法:
まずEdmondsの最大マッチングアルゴリズムを走らせる。奇サイクルは1点+完全マッチングに分解できるので、得られたマッチングのサイズが(頂点数-1)/2でなければNoを出力。
そうでない場合、このグラフのGallai-Edmonds分解を考える。Gallai-Edmonds分解とは、頂点集合を「ある完全マッチングが存在し、それに含まれない点」「その隣接点」「それ以外」の3つに分ける分割のこと。
条件より、頂点1はGallai-Edmonds分解で最初の集合に入る。Gallai-Edmonds分解において、最初の集合(に頂点を制限したグラフ)の全ての連結成分は、どの1点を除いても完全マッチングを持つという性質を持つ。そのようなグラフは全ての耳の長さが奇数の耳分解を持つので、頂点数が3以上なら奇サイクルとマッチングへの分解を持つ。
よって、Gallai-Edmonds分解において頂点1が最初の集合に入り、かつその連結成分のサイズが3以上なら、(残りの成分は2番目の集合の点とうまく結び付ければ)、条件を満たす。逆に、サイズが1なら、(その点は2番目の集合の点としか結ばれておらず、2番目の集合の点は最初の集合の点以外と(Tutte条件より)結ばれないため)、条件を満たさない。
Gallai-Edmonds分解は、単にEdmondsのアルゴリズムでの交互木の根からの距離の偶奇を求めれば求まるので、単に頂点1が(根からの距離偶数の)あるサイズ3以上の花に含まれるかを判定すればよい。
計算量は最大マッチングパート以外は微小で、O(N^3)となる。
#include<stdio.h> #include<vector> #include<algorithm> using namespace std; int match[222]; int par[222]; int root[222]; vector<int>pat[222]; int state[222];//0: root 1: odd 2: even 3: out of forest bool flag; bool vis[222]; bool shr[222]; int num;//scanfで直接入れる void dfs(int node) { if (vis[node])return; vis[node] = true; if (flag)return; for (int i = 0; i < pat[node].size(); i++) { if (flag)return; if (state[pat[node][i]] == 3) { par[pat[node][i]] = node; state[pat[node][i]] = 1; state[match[pat[node][i]]] = 2; dfs(match[pat[node][i]]); } else if (state[pat[node][i]] % 2 == 0) { if (root[node] == root[pat[node][i]])continue; vector<int>p1, p2; int v1 = node, v2 = pat[node][i]; for (;;) { p1.push_back(v1); if (state[v1] == 0)break; p1.push_back(match[v1]); v1 = par[match[v1]]; } for (;;) { p2.push_back(v2); if (state[v2] == 0)break; p2.push_back(match[v2]); v2 = par[match[v2]]; } if (v1 != v2) { flag = true; reverse(p1.begin(), p1.end()); reverse(p2.begin(), p2.end()); for (int i = 0; i < p1.size() - 1; i += 2) { match[p1[i]] = p1[i + 1]; match[p1[i + 1]] = p1[i]; } for (int i = 0; i < p2.size() - 1; i += 2) { match[p2[i]] = p2[i + 1]; match[p2[i + 1]] = p2[i]; } match[node] = pat[node][i]; match[pat[node][i]] = node; return; } else { reverse(p1.begin(), p1.end()); reverse(p2.begin(), p2.end()); int r = -1; for (int i = 0; i < min(p1.size(), p2.size()); i++) { if (p1[i] == p2[i]) { r = root[p1[i]]; } } int x1, x2; for (int i = 0; i < p1.size(); i++)if (root[p1[i]] == r)x1 = i; for (int i = 0; i < p2.size(); i++)if (root[p2[i]] == r)x2 = i; for (int i = x1 + 2; i < p1.size() - 1; i += 2) { par[p1[i]] = p1[i + 1]; } for (int i = x2 + 2; i < p2.size() - 1; i += 2) { par[p2[i]] = p2[i + 1]; } if (root[node] != r)par[node] = pat[node][i]; if (root[pat[node][i]] != r)par[pat[node][i]] = node; vector<int>nxt; fill(shr, shr + num, false); for (int i = x1; i < p1.size(); i++)shr[root[p1[i]]] = true; for (int i = x2; i < p2.size(); i++)shr[root[p2[i]]] = true; for (int i = 0; i < num; i++) { if (shr[root[i]]) { root[i] = r; if (state[i] == 1) { state[i] = 2; nxt.push_back(i); } } } for (int i = 0; i < nxt.size(); i++)dfs(nxt[i]); } } } } void maxmatching()//patに双方向に辺を入れて呼ぶとmatch[i]に最大マッチングでのマッチング相手が入る { for (int i = 0; i < num; i++)match[i] = i; for (;;) { for (int i = 0; i < num; i++) { if (match[i] == i)state[i] = 0; else state[i] = 3; root[i] = par[i] = i; vis[i] = false; } flag = false; for (int i = 0; i < num; i++) { if (state[i] == 0)dfs(i); } if (!flag)break; } } int main() { int way; scanf("%d%d", &num, &way); for (int i = 0; i < way; i++) { int za, zb; scanf("%d%d", &za, &zb); za--, zb--; pat[za].push_back(zb); pat[zb].push_back(za); } maxmatching(); int m = 0; for (int i = 0; i < num; i++)if (match[i] == i)m++; int r = root[0]; int cnt = 0; for (int i = 0; i < num; i++)if (root[i] == r)cnt++; if (m == 1 && cnt >= 3)printf("Yes\n"); else printf("No\n"); }
Codeforces 497E Subsequences Return
- 問題概要
k進法で0からn-1までの各桁の和mod kという数列があるので、異なる部分列の個数を求めよ。
k<=30,n<=10^18
- 解法
珍しくデータ構造系でない。良問。
まず、部分列の個数は、各文字について、現れる次の位置に配るDPをやると重複なく数えられる。これでO(nk)のDPができる。
これを行列累乗にするが、列が変なので、log_k^n回に分けて行列を累乗していく。
行列のi行j列を、iではじまってその区間の後はjをとるような列の数 としておく。iではじまってその中でおわる部分列の個数も別に持っておく。
i回目には、0~k^i-1までの行列を求める。このとき、1個まえのランクの行列をかけるが、始まる数が1個ずれるので1個行と列をシフトさせてかける。
これを繰り返すと、各ステップでO(k^4)がかかり、O(k^4log_k^n)で解ける。
こういう実装の細部を詰めるのが異常に苦手でできたコードからは想像もつかない量の時間を食うのやめたい。精進します。
#include<stdio.h> #include<vector> #include<algorithm> using namespace std; typedef long long ll; ll mod = 1000000007; typedef vector<vector<ll> >vvi; typedef pair<vvi, vector<ll> >pvv; int num; pvv merge(vvi a, vvi b, vector<ll>ta, vector<ll>tb, int stp) { vvi ret; ret.resize(num); for (int i = 0; i < num; i++)ret[i].resize(num); vector<ll>tr = ta; for (int i = 0; i < num; i++) { for (int j = 0; j < num; j++) { for (int k = 0; k < num; k++) { ret[i][k] += a[i][j] * b[(j + num - stp) % num][(k + num - stp) % num]; ret[i][k] %= mod; } } } for (int i = 0; i < num; i++) { for (int j = 0; j < num; j++) { tr[i] += a[i][j] * tb[(j + num - stp) % num]; tr[i] %= mod; } } return make_pair(ret, tr); } vvi dat[100]; vector<ll>td[100]; int main() { ll len; scanf("%I64d%d", &len, &num); dat[0].resize(num); for (int i = 0; i < num; i++)dat[0][i].resize(num); for (int i = 0; i < num; i++)dat[0][0][i] = dat[0][i][i] = 1; td[0].resize(num); td[0][0] = 1; vvi e; e.resize(num); for (int i = 0; i<num; i++)e[i].resize(num); for (int i = 0; i<num; i++)e[i][i] = 1; vector<ll>te; te.resize(num); ll now = 1; for (int p = 1;; p++) { if (double(now)*double(num)>(2e18))break; now *= num; vvi z = e; vector<ll>t = te; for (int i = 0; i < num; i++) { pvv x = merge(z, dat[p - 1], t, td[p - 1], i); z = x.first; t = x.second; } dat[p] = z; td[p] = t; } /*for (int i = 0; i < num; i++) { for (int j = 0; j < num; j++) { printf("%I64d ", dat[2][i][j]); } printf("\n"); } printf("\n"); for (int i = 0; i < num; i++)printf("%I64d ", td[2][i]); printf("\n");*/ vector<int>v; for (;;) { if (len == 0)break; v.push_back(len%num); len /= num; } vvi n = e; vector<ll>tn = te; for (int i = int(v.size()) - 1; i >= 0; i--) { for (int j = 0; j < v[i]; j++) { int d = j; for (int k = int(v.size()) - 1; k >= i + 1; k--)d += v[k]; d %= num; pvv z = merge(n, dat[i], tn, td[i], d); n = z.first; tn = z.second; } } ll ans = 1; for (int i = 0; i < num; i++)ans += tn[i]; printf("%I64d\n", ans%mod); }
Codeforces 687D Dividing Kingdom II
- 問題概要
頂点数1000以下の単純無向グラフがあり、辺には番号がついている。辺の番号のある区間の辺のみ使って二部グラフをつくるとき、使えない辺のコストの最大値の最小値を1000回くらい求めよ。
- 解法
10^9を実装しても通るらしい。許さない。
区間を平方分割し、各区間についてUFとかで答えを求める。
区間をマージするとき、必要なのはその区間の二部グラフにならないように辺を上からとっていったような最大全域森の辺だけ(他はあっても最大値の最小値には影響しない)なので、segtreeにこの区間たちをのせると合計O(N^2)でsegtreeが構築できて、平方分割の真ん中の区間はsegtreeを見ることでさぼる。端は愚直にやる。クエリ当たりO(Nlog N)となる。
これを頑張って1時間くらいでデバッグまで終わらせたのは成長だなぁと思いながら順位表を見たら異常に解かれていて何事かと思ったら10^9が通るんですね...ゆるさん。
#include<stdio.h> #include<vector> #include<algorithm> using namespace std; typedef pair<int, int>pii; typedef pair<pii, int>pi3; typedef pair<vector<pi3>, int>pvi; #define SIZE 2000 class unionfind { public: int par[SIZE]; int ran[SIZE]; int ren[SIZE]; void init() { for (int i = 0; i<SIZE; i++) { par[i] = i; ran[i] = 0; ren[i] = 1; } } int find(int a) { if (a == par[a])return a; else return par[a] = find(par[a]); } void unite(int a, int b) { a = find(a); b = find(b); if (a == b)return; if (ran[a]>ran[b]) { par[b] = a; ren[a] += ren[b]; } else { par[a] = b; ren[b] += ren[a]; } if (ran[a] == ran[b])ran[b]++; } }; unionfind uf; vector<pi3>merge(vector<pi3>va, vector<pi3>vb) { vector<pi3>vec; int pt = 0; for (int i = 0; i < va.size(); i++) { for (;;) { if (pt == vb.size())break; if (vb[pt].second > va[i].second)vec.push_back(vb[pt++]); else break; } vec.push_back(va[i]); } for (int i = pt; i < vb.size(); i++)vec.push_back(vb[i]); return vec; } pvi get(vector<pi3>vec) { uf.init(); vector<pi3>ret; for (int i = 0; i < vec.size(); i++) { int s = vec[i].first.first, t = vec[i].first.second; if (uf.find(s * 2) != uf.find(t * 2) && uf.find(s * 2) != uf.find(t * 2 + 1)) { uf.unite(s * 2, t * 2 + 1); uf.unite(s * 2 + 1, t * 2); ret.push_back(vec[i]); } if (uf.find(s * 2) == uf.find(t * 2)) { return make_pair(ret, vec[i].second); } } return make_pair(ret, -1); } pi3 dat[1000000]; vector<pi3>sortv(int lb, int ub) { vector<pii>zv; for (int j = lb; j <= ub; j++) { zv.push_back(make_pair(dat[j].second, j)); } sort(zv.begin(), zv.end()); reverse(zv.begin(), zv.end()); vector<pi3>z; for (int i = 0; i < zv.size(); i++) { z.push_back(dat[zv[i].second]); } return z; } #define SIZE 1024 class segtree { public: pvi seg[SIZE * 2]; void init(vector<pvi>vec) { seg[0].second = -1; for (int i = 0; i < vec.size(); i++) { seg[SIZE + i] = vec[i]; } for (int i = vec.size(); i < SIZE; i++)seg[SIZE + i].second = -1; for (int i = SIZE - 1; i >= 1; i--) { pvi t = get(merge(seg[i * 2].first, seg[i * 2 + 1].first)); seg[i].first = t.first; seg[i].second = max(t.second, max(seg[i * 2].second, seg[i * 2 + 1].second)); } } pvi calc(int beg, int end, int node, int lb, int ub) { if (ub < beg || end < lb)return seg[0]; if (beg <= lb&&ub <= end)return seg[node]; pvi a = calc(beg, end, node * 2, lb, (lb + ub) / 2); pvi b = calc(beg, end, node * 2 + 1, (lb + ub) / 2 + 1, ub); pvi t = get(merge(a.first, b.first)); return make_pair(t.first, max(t.second, max(a.second, b.second))); } }; segtree tree; int main() { int num, way, query; scanf("%d%d%d", &num, &way, &query); for (int i = 0; i < way; i++) { int za, zb, zc; scanf("%d%d%d", &za, &zb, &zc); za--; zb--; dat[i] = make_pair(make_pair(za, zb), zc); } /*for (int i = 0; i < way; i++) { for (int j = i; j < way; j++) { printf("%d %d %d\n:::", i, j, get(sortv(i, j)).second); vector<pi3>v = get(sortv(i, j)).first; for (int k = 0; k < v.size(); k++)printf("%d %d ", v[k].first.first, v[k].first.second); printf("\n"); } }*/ vector<pvi>sgv; for (int i = 0; i < (way + num - 1) / num; i++) { sgv.push_back(get(sortv(num*i, min(way, num*(i + 1))))); } tree.init(sgv); for (int p = 0; p < query; p++) { int za, zb; scanf("%d%d", &za, &zb); za--; zb--; int lb = za / num, ub = zb / num; if (lb == ub) { pvi t = get(sortv(za, zb)); printf("%d\n", t.second); } else { pvi v1 = get(sortv(za, (lb + 1)*num - 1)); pvi v2 = tree.calc(lb + 1, ub - 1, 1, 0, SIZE - 1); pvi v3 = get(sortv(ub*num, zb)); pvi t = get(merge(v1.first,merge(v2.first,v3.first))); /*printf("%d %d %d %d\n", v1.second, v2.second, v3.second, t.second); for (int k = 0; k < v1.first.size(); k++)printf("%d %d ", v1.first[k].first.first, v1.first[k].first.second); printf("\n"); for (int k = 0; k < merge(v1.first, merge(v2.first, v3.first)).size(); k++)printf("%d %d ", merge(v1.first, merge(v2.first, v3.first))[k].first.first, merge(v1.first, merge(v2.first, v3.first))[k].first.second); printf("\n"); for (int k = 0; k < v3.first.size(); k++)printf("%d %d ", v3.first[k].first.first, v3.first[k].first.second); printf("\n"); for (int k = 0; k < t.first.size(); k++)printf("%d %d ", t.first[k].first.first, t.first[k].first.second); printf("\n");*/ printf("%d\n", max(max(t.second, v2.second), max(v1.second, v3.second))); } } }
Codeforces 587F Duff is Mad
- 問題概要
文字列がたくさんある。次のクエリを処理せよ。
クエリ: k個目文字列の中に出てくるr番目からl番目までの文字列のoccurrenceの合計を求めよ。
- 解法
安心安全の筋肉系文字列 & 平方分割。
長い文字列に関してはSAで答えを前計算(Aho-Corasickをつかったほうが良かった気がするけどライブラリがないのでSAで殴った)して累積和の形で持っておき、短い文字列に関してはクエリを先読みしてtrie木のEuler-tour列に平方分割系のデータ構造をのせて殴る(更新がO(N)回、質問がO(Nsqrt(N))回くるときはsegtreeやBITより平方分割のほうが速い)。
CFを埋めているせいで最近こういう筋肉系ばかりやっている。このブログが実装自慢ブログと化している。
こういう実装をバグらせなくなってきたのはとてもいいこと(今回もSAのライブラリ以外に200行以上本質部分のコードを書いたのにあまりバグが生えずびっくりした)だが、頭の整理と気力が追い付かないので一息に書ききれない。
SAのせいで発生する(Aho-Corasickだとたぶん大丈夫)MLEを避けるためにバケットサイズの比がいびつになった。これなら余計にlogつけても通ったかもしれない。
#include<stdio.h> #include<vector> #include<algorithm> #include<string> #include<iostream> using namespace std; vector<int>sais(vector<int>vec)//1-originにして末尾に0を加えてから呼ぶのを忘れないこと { if (vec.empty()) { vector<int>v; return v; } vector<int>dat; dat.resize(vec.size()); dat[vec.size() - 1] = 0; for (int i = vec.size() - 2; i >= 0; i--) { if (vec[i]>vec[i + 1])dat[i] = 1; else if (vec[i]<vec[i + 1])dat[i] = 0; else dat[i] = dat[i + 1]; } vector<vector<int> >sa; int maxi = 0; for (int i = 0; i<vec.size(); i++)maxi = max(maxi, vec[i]); sa.resize(maxi + 1); vector<int>vnum; vnum.resize(maxi + 1); for (int i = 0; i<vec.size(); i++)vnum[vec[i]]++; vector<bool>islms; islms.resize(vec.size()); fill(islms.begin(), islms.end(), false); for (int i = 0; i <= maxi; i++) { sa[i].resize(vnum[i]); fill(sa[i].begin(), sa[i].end(), -1); } vector<int>pt1, pt2; pt1.resize(maxi + 1); pt2.resize(maxi + 1); for (int i = 0; i <= maxi; i++) { pt2[i] = vnum[i] - 1; } for (int i = vec.size() - 1; i >= 1; i--) { if ((dat[i - 1] == 1 && dat[i] == 0) || i == vec.size() - 1) { sa[vec[i]][pt2[vec[i]]--] = i; islms[i] = true; } } for (int i = 0; i <= maxi; i++) { for (int j = 0; j<vnum[i]; j++) { if (sa[i][j]>0) { if (dat[sa[i][j] - 1] == 1) { sa[vec[sa[i][j] - 1]][pt1[vec[sa[i][j] - 1]]++] = sa[i][j] - 1; } } } } for (int i = 1; i <= maxi; i++) { for (int j = pt2[i] + 1; j<vnum[i]; j++)sa[i][j] = -1; pt2[i] = vnum[i] - 1; } for (int i = maxi; i >= 0; i--) { for (int j = vnum[i] - 1; j >= 0; j--) { if (sa[i][j]>0) { if (dat[sa[i][j] - 1] == 0) { sa[vec[sa[i][j] - 1]][pt2[vec[sa[i][j] - 1]]--] = sa[i][j] - 1; } } } } vector<int>d; d.resize(vec.size()); int cnt = 0; vector<int>bef; bef.push_back(0); bool fl = false; for (int i = 0; i <= maxi; i++) { for (int j = 0; j<vnum[i]; j++) { if (islms[sa[i][j]]) { vector<int>zk; zk.push_back(vec[sa[i][j]]); for (int k = sa[i][j] + 1; k<vec.size(); k++) { zk.push_back(vec[k]); if (islms[k])break; } if (zk != bef)cnt++; else if (vec[sa[i][j]] != 0)fl = true; d[sa[i][j]] = cnt; bef = zk; } } } vector<int>vt; for (int i = 0; i<vec.size(); i++)if (islms[i])vt.push_back(d[i]); vector<int>gv; vector<int>nv; if (fl) { gv = sais(vt); vector<int>v; for (int i = 0; i<vec.size(); i++) { if (islms[i])v.push_back(i); } for (int i = 0; i<gv.size(); i++) { nv.push_back(v[gv[i]]); } } else { gv = vt; nv.resize(gv.size()); int pt = 0; for (int i = 0; i<gv.size(); i++) { for (;;) { if (islms[pt]) { nv[gv[i]] = pt; pt++; break; } pt++; } } } for (int i = 0; i <= maxi; i++) { fill(sa[i].begin(), sa[i].end(), -1); pt1[i] = 0; pt2[i] = vnum[i] - 1; } for (int i = nv.size() - 1; i >= 0; i--) { sa[vec[nv[i]]][pt2[vec[nv[i]]]--] = nv[i]; } for (int i = 0; i <= maxi; i++) { for (int j = 0; j<vnum[i]; j++) { if (sa[i][j]>0) { if (dat[sa[i][j] - 1] == 1) { sa[vec[sa[i][j] - 1]][pt1[vec[sa[i][j] - 1]]++] = sa[i][j] - 1; } } } } for (int i = 1; i <= maxi; i++) { for (int j = pt2[i] + 1; j<vnum[i]; j++)sa[i][j] = -1; pt2[i] = vnum[i] - 1; } for (int i = maxi; i >= 0; i--) { for (int j = vnum[i] - 1; j >= 0; j--) { if (sa[i][j]>0) { if (dat[sa[i][j] - 1] == 0) { sa[vec[sa[i][j] - 1]][pt2[vec[sa[i][j] - 1]]--] = sa[i][j] - 1; } } } } vector<int>ret; for (int i = 0; i <= maxi; i++) { for (int j = 0; j<vnum[i]; j++) { ret.push_back(sa[i][j]); } } return ret; } vector<int>calclcp(vector<int>str, vector<int>sa)//lcp: SAのi-1番目とi番目のlcp { vector<int>rsa; rsa.resize(sa.size()); for (int i = 0; i < sa.size(); i++)rsa[sa[i]] = i; vector<int>lcp; lcp.resize(sa.size()); int now = 1; for (int i = 0; i < str.size() - 1; i++) { if (now != 0)now--; for (;;) { if (str[i + now] == str[sa[rsa[i] - 1] + now])now++; else { lcp[rsa[i]] = now; break; } } } return lcp; } #define SIZE 262144 class segtree { public: int seg[SIZE * 2]; void init() { for (int i = 1; i < SIZE * 2; i++)seg[i] = 1000000000; } void update(int a, int b) { a += SIZE; seg[a] = min(seg[a], b); for (;;) { a /= 2; if (a == 0)break; seg[a] = min(seg[a * 2], seg[a * 2 + 1]); } } int get(int beg, int end, int node, int lb, int ub) { if (ub < beg || end < lb)return 1000000000; if (beg <= lb&&ub <= end)return seg[node]; return min(get(beg, end, node * 2, lb, (lb + ub) / 2), get(beg, end, node * 2 + 1, (lb + ub) / 2 + 1, ub)); } }; segtree tree; vector<int>str, sa, lcp, rsa; void init() { sa = sais(str); lcp = calclcp(str, sa); tree.init(); rsa.resize(sa.size()); for (int i = 0; i < sa.size(); i++)rsa[sa[i]] = i; for (int i = 0; i < str.size(); i++)tree.update(i, lcp[i]); } int getlcp(int a, int b)//a文字目からとb文字目からのlcpの長さ { return tree.get(min(rsa[a], rsa[b]) + 1, max(rsa[a], rsa[b]), 1, 0, SIZE - 1); } #define B 1250 #define NB 81 typedef long long ll; ll rrui[NB][200010]; ll rans[NB][100010]; int stp[100000]; int toz[100000]; vector<int>kou[100000]; ll ans[100000]; class trie { public: int nex[100001][26]; vector<int>dat[100001]; int pt; void init() { pt = 1; fill(nex[0], nex[0] + 26, -1); } void adds(string s, int d) { int now = 0; for (int i = 0; i < s.size(); i++) { if (nex[now][s[i] - 'a'] == -1) { nex[now][s[i] - 'a'] = pt; fill(nex[pt], nex[pt] + 26, -1); pt++; } now = nex[now][s[i] - 'a']; } dat[now].push_back(d); } int getdest(string s) { int now = 0; for (int i = 0; i < s.size(); i++) { if (nex[now][s[i] - 'a'] == -1)break; now = nex[now][s[i] - 'a']; } return now; } vector<int>eul; int ord[100001]; int fin[100001]; void calceul(int node) { ord[node] = eul.size(); eul.push_back(node); for (int i = 0; i < 26; i++)if (nex[node][i] != -1)calceul(nex[node][i]); fin[node] = eul.size() - 1; } }; trie tr; class sqq { public: ll now[B*NB]; ll flag[NB]; void resolve(int b) { for (int i = 0; i < B; i++)now[B*b + i] += flag[b]; flag[b] = 0; } void add(int lb, int ub, int t) { int a = lb / B, b = ub / B; resolve(a); resolve(b); if (a == b) { for (int i = lb; i <= ub; i++)now[i] += t; } else { for (int i = lb; i < a*B + B; i++)now[i] += t; for (int i = a + 1; i <= b - 1; i++)flag[i] += t; for (int i = b*B; i <= ub; i++)now[i] += t; } } ll get(ll a) { return flag[a / B] + now[a]; } }; sqq bi; typedef pair<int, int>pii; vector<pii>que1[100000], que2[100000]; int main() { int num, query; scanf("%d%d", &num, &query); vector<string>vec; for (int i = 0; i < num; i++) { string s; cin >> s; vec.push_back(s); } for (int i = 0; i < num; i++) { stp[i] = str.size(); for (int j = 0; j < vec[i].size(); j++)str.push_back(vec[i][j] - 'a' + 1); str.push_back(27); } str.push_back(0); init(); tr.init(); for (int i = 0; i < vec.size(); i++) { tr.adds(vec[i], i); } int pt = 0; for (int i = 0; i < num; i++) { if (vec[i].size() >= B) { for (int j = stp[i]; j < stp[i] + vec[i].size(); j++) { rrui[pt][rsa[j] + 1]++; } for (int j = 1; j <= str.size(); j++)rrui[pt][j] += rrui[pt][j - 1]; toz[i] = pt; pt++; } else { for (int j = 0; j < vec[i].size(); j++) { string zs; for (int k = j; k < vec[i].size(); k++)zs.push_back(vec[i][k]); kou[i].push_back(tr.getdest(zs)); } } } for (int i = 0; i < num; i++) { int lb, ub; int beg = 0, end = rsa[stp[i]]; for (;;) { if (beg == end)break; int med = (beg + end) / 2; if (getlcp(sa[med], sa[rsa[stp[i]]]) < vec[i].size())beg = med + 1; else end = med; } lb = beg; beg = rsa[stp[i]], end = str.size() - 1; for (;;) { if (beg == end)break; int med = (beg + end + 1) / 2; if (getlcp(sa[med], sa[rsa[stp[i]]]) < vec[i].size())end = med - 1; else beg = med; } ub = beg; for (int j = 0; j < pt; j++) { rans[j][i + 1] += rrui[j][ub + 1] - rrui[j][lb]; } } for (int i = 0; i < pt; i++) { for (int j = 0; j < num; j++) { rans[i][j + 1] += rans[i][j]; } } for (int p = 0; p < query; p++) { int za, zb, zc; scanf("%d%d%d", &za, &zb, &zc); za--; zb--; zc--; if (vec[zc].size() >= B) { ans[p] = rans[toz[zc]][zb + 1] - rans[toz[zc]][za]; } else { que1[za].push_back(make_pair(zc, p)); que2[zb].push_back(make_pair(zc, p)); } } tr.calceul(0); for (int i = 0; i < num; i++) { for (int j = 0; j < que1[i].size(); j++) { for (int k = 0; k < kou[que1[i][j].first].size(); k++) { ans[que1[i][j].second] -= bi.get(tr.ord[kou[que1[i][j].first][k]]); } } if (vec[i].size() < B)bi.add(tr.ord[kou[i][0]], tr.fin[kou[i][0]], 1); for (int j = 0; j < que2[i].size(); j++) { for (int k = 0; k < kou[que2[i][j].first].size(); k++) { ans[que2[i][j].second] += bi.get(tr.ord[kou[que2[i][j].first][k]]); } } } for (int i = 0; i < query; i++) { printf("%I64d\n", ans[i]); } }