結果
| 問題 | No.5025 Forest Fire |
| コンテスト | |
| ユーザー |
 Shki
|
| 提出日時 | 2026-10-20 07:17:18 |
| 言語 | C++23(gcc16) (gcc 16.2.0 + boost 1.92.0 + ACL) |
| 結果 |
WA
不安定
|
| 実行時間 | - |
| コード長 | 26,342 bytes |
| 記録 | |
| コンパイル時間 | 5,033 ms |
| コンパイル使用メモリ | 411,728 KB |
| 実行使用メモリ | 9,932 KB |
| スコア | 49,042,477 |
| 最終ジャッジ日時 | 2026-10-24 14:36:12 |
| 合計ジャッジ時間 | 129,872 ms |
|
ジャッジサーバーID （参考情報） |
judge2_0 / <nil> |
（要ログイン）
| ファイルパターン | 結果 |
|---|---|
| other | AC * 49 WA * 1 |
ソースコード
#include <bits/stdc++.h>
using namespace std;
// Heuristic solver for the N=50, K=100 stochastic forest-fire problem.
// C++17. All objective values below omit the constant divisor N*N.
#ifndef GREEDY_S
#define GREEDY_S 32
#endif
#ifndef GREEDY_END
#define GREEDY_END 1.20
#endif
using ll = long long;
using Clock = chrono::steady_clock;
struct RNG {
uint64_t x;
explicit RNG(uint64_t s=1234567890123456789ULL) : x(s) {}
inline uint64_t next() {
uint64_t z = (x += 0x9e3779b97f4a7c15ULL);
z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL;
z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL;
return z ^ (z >> 31);
}
inline int integer(int n) { return int(next() % uint64_t(n)); }
inline double real() { return (next() >> 11) * 0x1.0p-53; }
};
struct Solver {
int N, K, V, eligible=0;
vector<int> d, trees;
vector<array<int,4>> nbr;
vector<uint8_t> deg;
RNG rng;
Clock::time_point start;
double limit = 1.57;
static constexpr int MAXS=384;
vector<vector<uint8_t>> samples;
// Tarjan DFS scratch arrays.
vector<int> tin, low, sub, sepSum, sepSq, visitOrder;
int timer=0;
const vector<uint8_t> *curCut=nullptr, *curDry=nullptr;
vector<int> seen, compId, compSizes;
int stamp=0;
vector<int> queue_buf;
Solver(int n, int k, vector<int> dd) : N(n), K(k), V(n*n), d(move(dd)),
nbr(V), deg(V), rng(0xa920a6b592391fdeULL), start(Clock::now()),
tin(V), low(V), sub(V), sepSum(V), sepSq(V),
seen(V,0), compId(V,0), queue_buf(V) {
for(int v=0; v<V; ++v) {
if(d[v]>0) trees.push_back(v);
int i=v/N,j=v%N;
if(i>0) nbr[v][deg[v]++]=v-N;
if(i+1<N) nbr[v][deg[v]++]=v+N;
if(j>0) nbr[v][deg[v]++]=v-1;
if(j+1<N) nbr[v][deg[v]++]=v+1;
}
eligible=int(trees.size());
K=min(K,eligible);
// Unbiased random samples shared across all candidates (common random numbers).
samples.assign(MAXS, vector<uint8_t>(V));
for(int s=0; s<MAXS; ++s) {
for(int v=0; v<V; ++v)
samples[s][v] = (d[v]>0 && rng.integer(10)<d[v]);
}
}
inline double elapsed() const {
return chrono::duration<double>(Clock::now()-start).count();
}
inline bool timeLeft(double margin=0.0) const { return elapsed()+margin < limit; }
// Computes sum_C |C|^2 on one fixed sample, using a DFS flood fill.
inline ll evalOne(const vector<uint8_t> &cut, const vector<uint8_t> &dry) {
++stamp;
ll ans=0;
for(int v=0;v<V;v++) {
if(cut[v] || !dry[v] || seen[v]==stamp) continue;
int head=0,tail=1;
queue_buf[0]=v;
seen[v]=stamp;
while(head<tail) {
int u=queue_buf[head++];
for(int t=0;t<deg[u];t++) {
int w=nbr[u][t];
if(seen[w]!=stamp && !cut[w] && dry[w]) {
seen[w]=stamp;
queue_buf[tail++]=w;
}
}
}
ans += (ll)tail*tail;
}
return ans;
}
inline ll eval(const vector<uint8_t> &cut, int begin, int cnt) {
ll ans=0;
for(int s=begin;s<begin+cnt;s++) ans+=evalOne(cut,samples[s]);
return ans;
}
// A complete DFS computes the removal gain of every node in a component.
void dfsTarjan(int u, int p) {
tin[u]=low[u]=timer++;
sub[u]=1;
sepSum[u]=sepSq[u]=0;
visitOrder.push_back(u);
for(int k=0;k<deg[u];k++) {
int v=nbr[u][k];
if((*curCut)[v] || !(*curDry)[v]) continue;
if(tin[v]<0) {
dfsTarjan(v,u);
sub[u]+=sub[v];
low[u]=min(low[u],low[v]);
if(low[v]>=tin[u]) {
sepSum[u]+=sub[v];
sepSq[u]+=sub[v]*sub[v];
}
} else if(v!=p) low[u]=min(low[u],tin[v]);
}
}
// O(S*(V+edges)) all-vertices removal marginal gains.
void removalGains(const vector<uint8_t> &cut,int S,vector<ll>&gain) {
fill(gain.begin(),gain.end(),0LL);
for(int s=0;s<S;s++) {
curCut=&cut;
curDry=&samples[s];
fill(tin.begin(),tin.end(),-1);
timer=0;
visitOrder.clear();
for(int r=0;r<V;r++) {
if(cut[r] || !samples[s][r] || tin[r]>=0) continue;
int first=(int)visitOrder.size();
dfsTarjan(r,-1);
int sz=(int)visitOrder.size()-first;
int sq=sz*sz;
for(int k=first;k<(int)visitOrder.size();k++) {
int u=visitOrder[k];
int rest=sz-1-sepSum[u];
gain[u]+=(ll)sq-sepSq[u]-rest*rest;
}
}
}
}
inline int countCut(const vector<uint8_t>&a) {
int n=0; for(uint8_t c:a) n+=c; return n;
}
// Greedily complete an existing (possibly empty) solution.
// Has a deterministic per-run time cap; safely returns a valid partial solution.
void greedyFill(vector<uint8_t>&cut, int S, double endTime) {
vector<ll> gain(V);
int selected=countCut(cut), iterations=0;
double blockStart=elapsed();
int currentS=S;
while(selected<K && elapsed()<endTime) {
removalGains(cut,currentS,gain);
int best=-1;
ll bestG=-1;
for(int v:trees) if(!cut[v] && gain[v]>bestG) {
best=v; bestG=gain[v];
}
if(best<0) break;
cut[best]=1;
++selected; ++iterations;
if((iterations&3)==0 && selected<K) {
double now=elapsed();
double stepCost=(now-blockStart)/4.;
double budget=endTime-now;
if(stepCost*(K-selected)>0.85*budget && currentS>12) {
int proposed=(int)(currentS*(0.83*budget)/
max(1e-9,stepCost*(K-selected)));
currentS=max(12,min(currentS,proposed));
}
blockStart=now;
}
}
}
// A full left-right or top-bottom blocking path; 8-connected blocked
// path separates 4-connected open cells. DP permits |delta row|<=1.
vector<int> barrier(bool horizontal,int target,int variation) {
const int INF=1000000000;
vector<array<int,50>> prev(N);
vector<int> cost(N*N,INF);
auto at=[&](int along,int across)->int {
return horizontal ? across*N+along : along*N+across;
};
for(int across=0;across<N;across++) {
int v=at(0,across);
int dist=abs(across-target);
cost[across]=(d[v] ? 1000+3*d[v] : 0) + dist*(variation?12:25);
prev[0][across]=-1;
}
for(int col=1;col<N;col++) {
for(int row=0;row<N;row++) {
int best=INF,par=-1;
for(int dr=-1;dr<=1;dr++) {
int pr=row+dr;
if(pr<0 || pr>=N) continue;
int val=cost[(col-1)*N+pr];
if(val<best) best=val,par=pr;
}
int v=at(col,row);
cost[col*N+row]=best+(d[v] ? 1000+3*d[v] : 0) +
abs(row-target)*(variation?12:25);
prev[col][row]=par;
}
}
int bestrow=0,best=INF;
for(int row=0;row<N;row++) {
int val=cost[(N-1)*N+row]+abs(row-target)*30;
if(val<best) best=val,bestrow=row;
}
vector<int> path;
int row=bestrow;
for(int col=N-1;col>=0;col--) {
path.push_back(at(col,row));
row=prev[col][row];
}
return path;
}
vector<uint8_t> pathCut(const vector<int>&a,const vector<int>&b={}) {
vector<uint8_t> result(V,0);
for(int v:a) if(d[v])result[v]=1;
for(int v:b) if(d[v])result[v]=1;
return result;
}
void addHelpfulFill(vector<uint8_t>&cut,int cntSample) {
// Short greedy fill is useful especially where water saves barrier budget.
if(countCut(cut)>=K) return;
vector<ll> gain(V);
for(int iter=countCut(cut);iter<K;iter++) {
removalGains(cut,cntSample,gain);
int best=-1;ll bestVal=-1;
for(int v:trees) if(!cut[v] && gain[v]>bestVal)
best=v,bestVal=gain[v];
if(best<0) break;
cut[best]=1;
}
}
// Cost of restoring a previously felled tree. Components and sizes are
// computed once per sample; each cut vertex then needs <=4 neighbor probes.
void restorationCosts(const vector<uint8_t>&cut,int S,vector<ll>&cost) {
fill(cost.begin(),cost.end(),0LL);
vector<int> cutList=collectCut(cut);
for(int s=0;s<S;s++) {
const auto& dry=samples[s];
fill(compId.begin(),compId.end(),-1);
compSizes.clear();
for(int v=0;v<V;v++) {
if(cut[v] || !dry[v] || compId[v]>=0) continue;
int cid=(int)compSizes.size();
int head=0,tail=1;
queue_buf[0]=v;compId[v]=cid;
while(head<tail) {
int u=queue_buf[head++];
for(int k=0;k<deg[u];k++) {
int w=nbr[u][k];
if(!cut[w] && dry[w] && compId[w]<0) {
compId[w]=cid;
queue_buf[tail++]=w;
}
}
}
compSizes.push_back(tail);
}
for(int v:cutList)if(dry[v]) {
int ids[4],cnt=0;
ll joined=1,old=0;
for(int k=0;k<deg[v];k++) {
int w=nbr[v][k],c=compId[w];
if(c<0)continue;
bool seenBefore=false;
for(int i=0;i<cnt;i++)if(ids[i]==c)seenBefore=true;
if(seenBefore)continue;
ids[cnt++]=c;
ll z=compSizes[c];
joined+=z;old+=z*z;
}
cost[v]+=joined*joined-old;
}
}
}
// SA solution lists allow fast proposal updates without mutating originals.
struct Sol {
vector<uint8_t> cut;
ll value=LLONG_MAX;
};
vector<int> collectCut(const vector<uint8_t>&cut) {
vector<int> v;
for(int u:trees)if(cut[u])v.push_back(u);
return v;
}
// Score on stronger, independent validation samples to reduce overfitting.
ll validation(const vector<uint8_t>&cut) {
return eval(cut,128,192);
}
void normalize(vector<uint8_t>&cut) {
int have=countCut(cut);
if(have>K) {
vector<int> used=collectCut(cut);
while(have>K) {
int v=used[rng.integer((int)used.size())];
if(cut[v]) {cut[v]=0;--have;}
}
}
if(have<K) {
// Add trees that cannot hurt the objective; score against all samples.
vector<ll> g(V);
removalGains(cut,24,g);
vector<pair<ll,int>> cand;
for(int v:trees) if(!cut[v]) cand.emplace_back(g[v],v);
sort(cand.begin(),cand.end(),greater<pair<ll,int>>());
for(auto [value,v]:cand) {
if(have==K)break;
cut[v]=1;have++;
}
}
}
// Pick a swap size 1..7 from a time-dependent smooth distribution.
// Larger changes are common at the beginning; 1-point refinements dominate
// later, but every allowed swap size retains nonzero probability.
int sampleSwapSize(double progress) {
int maxk=min({7,K,eligible-K});
if(maxk<=0) return 0;
double lambda=0.95+1.40*progress;
double weights[8], sum=0;
for(int k=1;k<=maxk;k++) {
weights[k]=0.97*exp(-lambda*(k-1))+0.03;
sum+=weights[k];
}
double roll=rng.real()*sum;
for(int k=1;k<=maxk;k++) {
roll-=weights[k];
if(roll<0) return k;
}
return maxk;
}
// Exactly k distinct old cuts are replaced by k distinct previously
// uncut trees. Local/clustered and priority-guided choices are mixed.
bool mutateSwapRange(vector<uint8_t>&x, const vector<uint8_t>&base,
const vector<int>&cutVertices, int k,
const vector<int>&lowRestore,
const vector<int>&highGain) {
if(k<1 || (int)cutVertices.size()<k || eligible-(int)cutVertices.size()<k)
return false;
x=base;
int anchor=cutVertices[rng.integer((int)cutVertices.size())];
bool cluster=(k>=2 && rng.integer(100)<47);
bool guide=(!lowRestore.empty() && rng.integer(100)<(k==1?67:45));
int clusterRadius=2+(k+1)/2;
vector<int> removed; removed.reserve(k);
for(int i=0;i<k;i++) {
int a=-1;
if(cluster) {
for(int tries=0; tries<28;tries++) {
int v=cutVertices[rng.integer((int)cutVertices.size())];
if(!x[v])continue;
if(abs(v/N-anchor/N)+abs(v%N-anchor%N)<=clusterRadius) {
a=v;break;
}
}
}
if(a<0 && guide) {
int cap=min((int)lowRestore.size(),max(35,k*8));
for(int tries=0;tries<16;tries++) {
int v=lowRestore[rng.integer(cap)];
if(x[v]) {a=v;break;}
}
}
if(a<0) {
for(int tries=0;tries<20;tries++) {
int v=cutVertices[rng.integer((int)cutVertices.size())];
if(x[v]) {a=v;break;}
}
}
if(a<0) for(int v:cutVertices) if(x[v]) {a=v;break;}
if(a<0)return false;
x[a]=0;
removed.push_back(a);
}
for(int i=0;i<k;i++) {
int a=removed[i],b=-1;
// Preserve the strong priority-guided one-point move from the
// original algorithm; larger moves preferentially stay local.
if(guide && k==1 && !highGain.empty()) {
int cap=min((int)highGain.size(),80);
for(int tries=0;tries<16;tries++) {
int v=highGain[rng.integer(cap)];
if(!base[v] && !x[v]) {b=v;break;}
}
}
// Keep edits geographically coherent, especially for k>1.
if(b<0 && rng.integer(100)<(k>=2?75:57)) {
int radius=rng.integer(100)<75?2:6;
int r0=a/N,c0=a%N;
for(int tries=0;tries<28;tries++) {
int r=r0+rng.integer(2*radius+1)-radius;
int c=c0+rng.integer(2*radius+1)-radius;
if(r<0||r>=N||c<0||c>=N)continue;
int v=r*N+c;
if(d[v]>0 && !base[v] && !x[v]) {b=v;break;}
}
}
if(b<0 && !highGain.empty() && rng.integer(100)<65) {
int cap=min((int)highGain.size(),max(80,k*16));
for(int tries=0;tries<20;tries++) {
int v=highGain[rng.integer(cap)];
if(!base[v] && !x[v]) {b=v;break;}
}
}
if(b<0) {
for(int tries=0;tries<30;tries++) {
int v=trees[rng.integer(eligible)];
if(!base[v] && !x[v]) {b=v;break;}
}
}
if(b<0) for(int v:trees) if(!base[v] && !x[v]) {b=v;break;}
if(b<0)return false;
x[b]=1;
}
return true;
}
// Shift structured separators; always preserves the number of cuts.
bool mutate(vector<uint8_t> &x, const vector<uint8_t> &base, int type,
const vector<int>&cutVertices) {
x=base;
if(cutVertices.empty()) return false;
if(type==2) { // shift all cuts on one row/column to an adjacent row/column
bool horizontal=rng.integer(2);
int line=-1;
// Choose one of the most populated lines.
vector<pair<int,int>> ls;
for(int z=0;z<N;z++) {
int n=0;
for(int w=0;w<N;w++) {
int v=horizontal?z*N+w:w*N+z;
n+=base[v];
}
if(n>=5)ls.push_back({n,z});
}
if(ls.empty())return false;
sort(ls.begin(),ls.end(),greater<pair<int,int>>());
line=ls[rng.integer(min((int)ls.size(),4))].second;
int off=rng.integer(2)*2-1;
if(line+off<0||line+off>=N)return false;
vector<int> removed,added;
for(int w=0;w<N;w++) {
int a=horizontal?line*N+w:w*N+line;
int b=horizontal?(line+off)*N+w:w*N+line+off;
if(!base[a] || a==b)continue;
if(d[b] && !base[b]) {removed.push_back(a);added.push_back(b);}
}
if(removed.size()<3)return false;
for(int v:removed)x[v]=0;
for(int v:added)x[v]=1;
return x!=base;
}
// Local segment cluster move: shift a connected block of cut sites.
int seed=cutVertices[rng.integer((int)cutVertices.size())];
int move=rng.integer(4);
int dr[4]={-1,1,0,0}, dc[4]={0,0,-1,1};
int arow=seed/N,acol=seed%N;
int rr=arow+dr[move],cc=acol+dc[move];
if(rr<0||rr>=N||cc<0||cc>=N)return false;
vector<int> group{seed};
vector<uint8_t> mark(V,0);mark[seed]=1;
int len=3+rng.integer(10);
for(int idx=0;idx<(int)group.size() && (int)group.size()<len;idx++) {
int u=group[idx];
int dir0=rng.integer(8);
for(int j=0;j<8;j++) {
int q=(j+dir0)%8;
int r=u/N + (q/3)-1;
int c=u%N + (q%3)-1;
if(r<0||r>=N||c<0||c>=N)continue;
int v=r*N+c;
if(v==u||!base[v]||mark[v])continue;
mark[v]=1;group.push_back(v);
if((int)group.size()==len)break;
}
}
vector<int> from,to;
for(int u:group) {
int r=u/N+dr[move],c=u%N+dc[move];
if(r<0||r>=N||c<0||c>=N)continue;
int v=r*N+c;
if(d[v] && !base[v]) {from.push_back(u);to.push_back(v);}
}
if(from.size()<2)return false;
for(int u:from)x[u]=0;
for(int u:to)x[u]=1;
return x!=base;
}
vector<uint8_t> run() {
if(K==0)return vector<uint8_t>(V,0);
vector<uint8_t> best(V,0);
ll bestVal=LLONG_MAX;
const int VALID_START=128, VALID_SIZE=192;
auto consider=[&](const vector<uint8_t>&x) {
ll v=eval(x,VALID_START,VALID_SIZE);
if(v<bestVal) {bestVal=v;best=x;}
};
// Always create a valid global separator first (strong on dense forests).
{
vector<uint8_t> middle(V,0);
for(int v:trees) if(v/N==N/2 || v%N==N/2)middle[v]=1;
normalize(middle);
consider(middle);
}
// Build a full, independent stochastic greedy solution before the
// global candidates. This avoids skipping greedy on slow hardware.
{
vector<uint8_t> cut(V,0);
greedyFill(cut, GREEDY_S, min(0.96,elapsed()+0.90));
#ifdef PROFILE
cerr << "greedy selected " << countCut(cut) << " at " << elapsed() << "\n";
#endif
if(countCut(cut)<K)normalize(cut);
consider(cut);
}
#ifdef PROFILE
cerr << "baseline " << elapsed() << " best " << bestVal/192.0 << "\n";
#endif
// Shortlist global separators. Cheap screening avoids wasting the
// 2-second budget on weak or redundant geometric candidates.
if(elapsed()<1.02) {
vector<vector<int>> horizontal,vertical;
for(int t: {13,17,20,23,25,27,30,33,36}) {
horizontal.push_back(barrier(true,t,0));
vertical.push_back(barrier(false,t,0));
}
vector<pair<ll,vector<uint8_t>>> top;
for(int h=0;h<(int)horizontal.size();h++) {
for(int v=0;v<(int)vertical.size();v++) {
if(elapsed()>1.25)break;
vector<uint8_t> cut=pathCut(horizontal[h],vertical[v]);
if(countCut(cut)>K)continue;
top.push_back({eval(cut,0,6),move(cut)});
}
}
for(int a=0;a<(int)horizontal.size();a++)
for(int b=a+1;b<(int)horizontal.size();b++) {
if(elapsed()>1.25)break;
auto cut=pathCut(horizontal[a],horizontal[b]);
if(countCut(cut)<=K)
top.push_back({eval(cut,0,6),move(cut)});
}
for(int a=0;a<(int)vertical.size();a++)
for(int b=a+1;b<(int)vertical.size();b++) {
if(elapsed()>1.25)break;
auto cut=pathCut(vertical[a],vertical[b]);
if(countCut(cut)<=K)
top.push_back({eval(cut,0,6),move(cut)});
}
sort(top.begin(),top.end(),[](auto&a,auto&b){return a.first<b.first;});
if((int)top.size()>12)top.resize(12);
for(auto &item:top) {
if(elapsed()>1.33)break;
item.first=eval(item.second,0,48);
}
sort(top.begin(),top.end(),[](auto&a,auto&b){return a.first<b.first;});
for(int i=0;i<min(4,(int)top.size());i++) {
if(elapsed()>1.36)break;
auto cut=top[i].second;
if(countCut(cut)<K) normalize(cut);
consider(cut);
}
}
#ifdef PROFILE
cerr << "structured " << elapsed() << " best " << bestVal/192.0 << "\n";
#endif
// Ensure a valid fallback solution even on slower systems.
if(countCut(best)<K)normalize(best);
// Anneal on a medium-resolution objective, and validate record candidates
// using independent samples. A deliberate safety margin is left to output.
if(elapsed()>1.42)return best;
vector<uint8_t> cur=best;
const int SA_S=64;
ll curVal=eval(cur,0,SA_S);
vector<uint8_t> candidate;
const double annealStart=elapsed(), annealEnd=limit-0.07;
if(annealEnd>annealStart+0.03) {
double duration=annealEnd-annealStart;
double totalScale=double(curVal)/SA_S;
double T0=max(1.0,totalScale*0.006);
double T1=max(0.2,totalScale*0.000003);
vector<int> used=collectCut(cur);
int iter=0;
int lastRefresh=-1000, acceptCount=0;
ll saBest=curVal;
vector<int> lowestRestore,highestGain;
vector<ll> restore(V), gain(V);
auto refreshRanking=[&]() {
const int S=24;
restorationCosts(cur,S,restore);
removalGains(cur,S,gain);
lowestRestore=collectCut(cur);
sort(lowestRestore.begin(),lowestRestore.end(),
[&](int a,int b){return restore[a]<restore[b];});
highestGain.clear();
for(int v:trees)if(!cur[v])highestGain.push_back(v);
sort(highestGain.begin(),highestGain.end(),
[&](int a,int b){return gain[a]>gain[b];});
if(highestGain.size()>150)highestGain.resize(150);
lastRefresh=acceptCount;
};
while(elapsed()<annealEnd) {
++iter;
if((iter&63)==0)used=collectCut(cur);
if((lastRefresh<0 || acceptCount-lastRefresh>=16) &&
elapsed()+0.15<annealEnd) refreshRanking();
double p=min(1.0,max(0.0,(elapsed()-annealStart)/duration));
// Keep structural moves separately: these move coherent long
// separators, rather than an arbitrary fixed number of trees.
bool done=false;
int r=rng.integer(100);
if(r<8) done=mutate(candidate,cur, r<4?2:3,used);
if(!done) {
int k=sampleSwapSize(p);
if(k==0)continue;
if(!mutateSwapRange(candidate,cur,used,k,lowestRestore,highestGain))
continue;
}
ll val=eval(candidate,0,SA_S);
ll diff=val-curVal;
double temperature=T0*pow(T1/T0,p)*SA_S;
if(diff<=0 || ((double)diff<100.0*temperature &&
rng.real()<exp(-(double)diff/temperature))) {
cur.swap(candidate);
curVal=val;
++acceptCount;
used=collectCut(cur);
}
if(curVal<saBest) {
saBest=curVal;
// Avoid selecting a solution solely because of sampling noise.
ll v=eval(cur,VALID_START,VALID_SIZE);
if(v<bestVal) {bestVal=v;best=cur;}
}
}
if(timeLeft(0.03)) consider(cur);
}
#ifdef PROFILE
cerr << "finish " << elapsed() << " best " << bestVal/192.0 << "\n";
#endif
return best;
}
};
int main() {
ios::sync_with_stdio(false);
cin.tie(nullptr);
int N,K;
if(!(cin>>N>>K))return 0;
vector<int> d(N*N);
for(int&i:d)cin>>i;
Solver solver(N,K,move(d));
auto cut=solver.run();
vector<int> positions;
for(int v=0;v<N*N;v++)if(cut[v])positions.push_back(v);
cout<<positions.size()<<'\n';
for(int v:positions) cout<<v/N<<' '<<v%N<<'\n';
return 0;
}