#include <bits/stdc++.h>
using namespace std;

// C++17, no external libraries. N = 50, K = 100.
// This solver minimizes the expected squared connected component sizes.
// Additional search: weighted 8-neighbor minimum-cost separators,
// exact re-ranked swaps, annealed 1-15-point exchange,
// route-preserving shortest-path rerouting, reverse DSU restoration greedy,
// and 3 cooling cycles with 2 reheats on common samples.
// SA neighbor: randomly restore k cut trees, rank ALL previously uncut trees
// by Monte-Carlo Lowlink marginal savings, and cut the best k (k=1..15).
// 2 reheats yield 3 SA cooling cycles. Geometric separators are
// retained for initial solutions; no LNS is used.
// Judges should compile with: g++ -std=c++17 -O2
#ifndef GUIDED_MODE
#define GUIDED_MODE 2
#endif
#ifndef GUIDE_WIDTH_MODE
#define GUIDE_WIDTH_MODE 4
#endif
// BESTK_MODE=0: original search
// BESTK_MODE=1: use exact top-k target selection for the broad (35%) neighbor
// BESTK_MODE=2: use top-k target selection as the only SA neighbor
#ifndef BESTK_MODE
#define BESTK_MODE 2
#endif
#ifndef BESTK_SAMPLES
#define BESTK_SAMPLES 8
#endif
// Minimax-focused DSU initial solution:
// Alongside expected sum of squared component sizes, penalize growth of
// the largest dry connected component in each Monte-Carlo sample.
// The first DSU candidate uses the original score, while the second
// generates a structurally different max-fragmentation-oriented candidate.
// 0 reproduces the original greedy evaluation.
#ifndef RESTORE_MAX_WEIGHT
#define RESTORE_MAX_WEIGHT 12
#endif
#ifndef REHEATS
#define REHEATS 2
#endif
using ll = long long;
using Clock = chrono::steady_clock;

struct Random {
    uint64_t state;
    explicit Random(uint64_t s) : state(s) {}
    uint64_t next() {
        uint64_t z = (state += 0x9e3779b97f4a7c15ULL);
        z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL;
        z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL;
        return z ^ (z >> 31);
    }
    int integer(int n) { return (int)(next() % (uint64_t)n); }
    double real() { return (next() >> 11) * 0x1.0p-53; }
};

class Solver {
    static constexpr int NS = 192;
    int N, K, V;
    bool sparseEnough = false;
    vector<int> d, trees;
    vector<array<int,4>> adj;
    vector<uint8_t> degree;
    vector<vector<uint8_t>> dry;
    Random rnd;
    Clock::time_point start;
    const double TIME_LIMIT = 1.92; // Reserve ~0.11s before the 2s judge limit.
    static constexpr int COOLING_CYCLES = REHEATS + 1; // Number of independent cooling schedules.

    vector<int> tin, low, sub, separateSum, separateSq, dfsOrder;
    int timer = 0;
    const vector<uint8_t>* currentCut = nullptr;
    const vector<uint8_t>* currentDry = nullptr;

    vector<int> mark, queueBuffer, comp, sizes;
    int stamp = 0;

    double elapsed() const {
        return chrono::duration<double>(Clock::now() - start).count();
    }
    bool timeOkay(double deadline) const { return elapsed() < deadline; }

    ll evalOne(const vector<uint8_t>& cut, const vector<uint8_t>& sample) {
        if (++stamp == INT_MAX) { fill(mark.begin(), mark.end(), 0); stamp = 1; }
        ll result = 0;
        for (int root = 0; root < V; ++root) {
            if (!sample[root] || cut[root] || mark[root] == stamp) continue;
            int h = 0, t = 1;
            queueBuffer[0] = root;
            mark[root] = stamp;
            while (h < t) {
                int u = queueBuffer[h++];
                for (int z = 0; z < degree[u]; ++z) {
                    int w = adj[u][z];
                    if (sample[w] && !cut[w] && mark[w] != stamp) {
                        mark[w] = stamp;
                        queueBuffer[t++] = w;
                    }
                }
            }
            result += (ll)t*t;
        }
        return result;
    }

    ll evaluate(const vector<uint8_t>& cut, int from, int count) {
        ll val = 0;
        for (int s = from; s < from + count; ++s) val += evalOne(cut, dry[s]);
        return val;
    }

    // Exact local objective updates for the fixed SA sample set.
    // Label the currently retained dry components once after every accepted move.
    // For a trial only components touching changed cells can change their sizes.
    static constexpr int SA_SAMPLES = 32;
    vector<vector<int16_t>> cachedId;
    vector<vector<int16_t>> cachedSize;
    ll cachedScore[SA_SAMPLES] = {};
    vector<int> touchedComp;
    int touchedGeneration=0;
    bool saCacheReady=false;

    void invalidateSACache() { saCacheReady=false; }

    void buildSACache(const vector<uint8_t>& cut) {
        if((int)cachedId.size()!=SA_SAMPLES) {
            cachedId.assign(SA_SAMPLES,vector<int16_t>(V));
            cachedSize.resize(SA_SAMPLES);
            touchedComp.assign(V,0);
        }
        for(int s=0;s<SA_SAMPLES;s++) {
            auto &ids=cachedId[s];
            auto &sizesHere=cachedSize[s];
            fill(ids.begin(),ids.end(),int16_t(-1));
            sizesHere.clear();
            ll total=0;
            const auto &sample=dry[s];
            for(int r=0;r<V;r++) {
                if(!sample[r] || cut[r] || ids[r]>=0)continue;
                int id=(int)sizesHere.size();
                int h=0,t=1;queueBuffer[0]=r;ids[r]=id;
                while(h<t) {
                    int u=queueBuffer[h++];
                    for(int j=0;j<degree[u];j++) {
                        int v=adj[u][j];
                        if(ids[v]<0 && sample[v] && !cut[v]) {
                            ids[v]=id;
                            queueBuffer[t++]=v;
                        }
                    }
                }
                sizesHere.push_back((int16_t)t);
                total+=(ll)t*t;
            }
            cachedScore[s]=total;
        }
        saCacheReady=true;
    }

    ll evaluateSATrial(const vector<uint8_t>& current,
                       const vector<uint8_t>& trial) {
        if(!saCacheReady)buildSACache(current);
        int changes[2500], nc=0;
        for(int v:trees)if(current[v]!=trial[v])changes[nc++]=v;
        if(nc==0){ll sum=0;for(ll val:cachedScore)sum+=val;return sum;}
        // For very large modifications, a full scan is faster than
        // building an affected-component frontier.
        if(nc>45)return evaluate(trial,0,SA_SAMPLES);
        ll total=0;
        for(int s=0;s<SA_SAMPLES;s++) {
            const auto &sample=dry[s];
            const auto &ids=cachedId[s];
            const auto &sizesHere=cachedSize[s];
            ++touchedGeneration;
            ll oldSq=0;
            auto include=[&](int u) {
                int id=ids[u];
                if(id<0 || touchedComp[id]==touchedGeneration)return;
                touchedComp[id]=touchedGeneration;
                ll sz=sizesHere[id];
                oldSq+=sz*sz;
            };
            for(int i=0;i<nc;i++) {
                int u=changes[i];
                if(!sample[u])continue;
                if(trial[u])include(u); // a previously retained dry vertex is cut
                else for(int j=0;j<degree[u];j++)include(adj[u][j]);
            }
            if(++stamp==INT_MAX){fill(mark.begin(),mark.end(),0);stamp=1;}
            ll newSq=0;
            auto bfs=[&](int root) {
                if(!sample[root] || trial[root] || mark[root]==stamp)return;
                int h=0,t=1;
                queueBuffer[0]=root;mark[root]=stamp;
                while(h<t) {
                    int u=queueBuffer[h++];
                    for(int j=0;j<degree[u];j++) {
                        int v=adj[u][j];
                        if(sample[v] && !trial[v] && mark[v]!=stamp) {
                            mark[v]=stamp;queueBuffer[t++]=v;
                        }
                    }
                }
                newSq+=(ll)t*t;
            };
            for(int i=0;i<nc;i++) {
                int u=changes[i];
                if(!sample[u])continue;
                if(!trial[u])bfs(u);
                else for(int j=0;j<degree[u];j++)bfs(adj[u][j]);
            }
            total+=cachedScore[s]-oldSq+newSq;
        }
#ifdef CHECK_INCREMENTAL
        ll full=evaluate(trial,0,SA_SAMPLES);
        if(full!=total){cerr<<"INVALID incremental="<<total<<" full="<<full<<"\n";abort();}
#endif
        return total;
    }

    void dfs(int u, int parent) {
        tin[u] = low[u] = timer++;
        sub[u] = 1;
        separateSum[u] = separateSq[u] = 0;
        dfsOrder.push_back(u);
        for (int z = 0; z < degree[u]; ++z) {
            int v = adj[u][z];
            if ((*currentCut)[v] || !(*currentDry)[v]) continue;
            if (tin[v] == -1) {
                dfs(v, u);
                sub[u] += sub[v];
                low[u] = min(low[u], low[v]);
                if (low[v] >= tin[u]) {
                    separateSum[u] += sub[v];
                    separateSq[u] += sub[v]*sub[v];
                }
            } else if (v != parent) {
                low[u] = min(low[u], tin[v]);
            }
        }
    }

    void removalGains(const vector<uint8_t>& cut, int S, vector<ll>& gain, int startSample=0) {
        fill(gain.begin(), gain.end(), 0);
        currentCut = &cut;
        for (int s = 0; s < S; ++s) {
            currentDry = &dry[(s+startSample)%NS];
            fill(tin.begin(), tin.end(), -1);
            dfsOrder.clear();
            timer = 0;
            for (int r = 0; r < V; ++r) {
                if (cut[r] || !(*currentDry)[r] || tin[r] != -1) continue;
                int first = (int)dfsOrder.size();
                dfs(r, -1);
                int componentSize = (int)dfsOrder.size() - first;
                int totalSq = componentSize*componentSize;
                for (int i = first; i < (int)dfsOrder.size(); ++i) {
                    int u = dfsOrder[i];
                    int leftover = componentSize - 1 - separateSum[u];
                    gain[u] += (ll)totalSq - separateSq[u] - leftover*leftover;
                }
            }
        }
    }

    void restoreCosts(const vector<uint8_t>& cut, int S, vector<ll>& cost) {
        fill(cost.begin(), cost.end(), 0LL);
        for (int s = 0; s < S; ++s) {
            fill(comp.begin(), comp.end(), -1);
            sizes.clear();
            for (int r = 0; r < V; ++r) {
                if (cut[r] || !dry[s][r] || comp[r] != -1) continue;
                int id = (int)sizes.size();
                int h = 0, t = 1;
                queueBuffer[0] = r;
                comp[r] = id;
                while (h < t) {
                    int u = queueBuffer[h++];
                    for (int z = 0; z < degree[u]; ++z) {
                        int w = adj[u][z];
                        if (!cut[w] && dry[s][w] && comp[w] == -1) {
                            comp[w] = id;
                            queueBuffer[t++] = w;
                        }
                    }
                }
                sizes.push_back(t);
            }
            for (int v : trees) {
                if (!cut[v] || !dry[s][v]) continue;
                int uniqueIds[4], count = 0;
                ll joined = 1, previous = 0;
                for (int z = 0; z < degree[v]; ++z) {
                    int id = comp[adj[v][z]];
                    if (id == -1) continue;
                    bool duplicate = false;
                    for (int j = 0; j < count; ++j) duplicate |= (uniqueIds[j] == id);
                    if (duplicate) continue;
                    uniqueIds[count++] = id;
                    ll q = sizes[id];
                    joined += q;
                    previous += q*q;
                }
                cost[v] += joined*joined - previous;
            }
        }
    }

    // Reverse construction: start with every tree removed, then restore
    // V-K vertices with minimum marginal increase in sampled fire damage.
    // A lazy heap avoids rescanning every candidate after each insertion.
    // In cyclic graphs marginal costs are NOT necessarily monotone (an
    // alternative bridge can lower a vertex's cost). Periodic full heap
    // refreshes are therefore essential; lazy steps are approximations.
    bool reverseRestoration(vector<uint8_t>& answer,
                            int samplesCount, double deadline,
                            const vector<uint8_t>* forwardGuide=nullptr,
                            int largestPenaltyWeight=0) {
        const int nt=(int)trees.size();
        const int need=nt-K;
        if(need<=0) return false;
        const int nslots=samplesCount*V;
        vector<int> parent(nslots,-1);
        vector<int16_t> componentSize(nslots,0);
        vector<uint8_t> retained(V,0);
        vector<int> maxDryComponent(samplesCount,0);
        auto findRoot=[&](int p)->int {
            while(parent[p]!=p){
                int q=parent[p];
                parent[p]=parent[q];
                p=parent[p];
            }
            return p;
        };
        auto roots=[&](int sample,int vertex,int *out)->int {
            int n=0,base=sample*V;
            for(int j=0;j<degree[vertex];j++){
                int v=adj[vertex][j];
                if(parent[base+v]<0)continue;
                int root=findRoot(base+v);
                bool dup=false;
                for(int k=0;k<n;k++)if(out[k]==root){dup=true;break;}
                if(!dup)out[n++]=root;
            }
            return n;
        };
        auto calculate=[&](int vertex)->ll {
            ll total=0;
            int nearby[4];
            for(int s=0;s<samplesCount;s++){
                if(!dry[s][vertex])continue;
                int count=roots(s,vertex,nearby);
                ll sum=1,prior=0;
                for(int j=0;j<count;j++){
                    ll z=componentSize[nearby[j]];
                    sum+=z;prior+=z*z;
                }
                total+=sum*sum-prior;
                // Minimax fragmentation penalty: discourage forming a new
                // largest dry component. Sums are over independent dry samples.
                if(largestPenaltyWeight && sum>maxDryComponent[s]) {
                    ll m=maxDryComponent[s];
                    total+=(ll)largestPenaltyWeight*(sum*sum-m*m);
                }
            }
            return total;
        };
        // Penalize restoring a cut identified as important by the
        // independent forward lowlink construction.
        auto guideCost=[&](int v,ll raw)->ll {
            return (forwardGuide && (*forwardGuide)[v])
                ? raw + raw/5 + 1 : raw;
        };
        auto add=[&](int vertex) {
            retained[vertex]=1;
            int nearby[4];
            for(int s=0;s<samplesCount;s++){
                if(!dry[s][vertex])continue;
                int base=s*V;
                int count=roots(s,vertex,nearby);
                int r=base+vertex;
                parent[r]=r;
                int total=1;
                for(int j=0;j<count;j++){
                    int q=nearby[j];
                    total+=componentSize[q];
                    parent[q]=r;
                }
                componentSize[r]=total;
                if(total>maxDryComponent[s])maxDryComponent[s]=total;
            }
        };
        struct Entry { ll cost; uint64_t tie; int vertex; };
        struct Greater {bool operator()(const Entry&a,const Entry&b)const {
            if(a.cost!=b.cost)return a.cost>b.cost;
            return a.tie>b.tie;
        }};
        priority_queue<Entry,vector<Entry>,Greater> pq;
        // Keep the existing SA's RNG trajectory unchanged whenever this
        // additional initial candidate is rejected.
        Random reverseRnd(0x9f4a7c159e3779b9ULL);
        // Precompute singleton cost exactly, avoiding 2,500 DSU scans.
        for(int v:trees){
            int singleton=0;
            for(int s=0;s<samplesCount;s++)singleton+=dry[s][v];
            pq.push({guideCost(v,(ll)singleton*(1+largestPenaltyWeight)),reverseRnd.next(),v});
        }
        int restored=0,checks=0,lastRefresh=-1;
        while(restored<need && !pq.empty()){
            if((++checks&511)==0 && !timeOkay(deadline))return false;
            if(restored>0 && restored%96==0 && restored!=lastRefresh) {
                // Fresh priorities also account for marginal costs that
                // DECREASE when other bridges connect adjacent components.
                priority_queue<Entry,vector<Entry>,Greater> refreshed;
                for(int v:trees)if(!retained[v])
                    refreshed.push({guideCost(v,calculate(v)),reverseRnd.next(),v});
                pq.swap(refreshed);
                lastRefresh=restored;
            }
            Entry e=pq.top();pq.pop();
            if(retained[e.vertex])continue;
            ll actual=guideCost(e.vertex,calculate(e.vertex));
            // Compare the refreshed cost to the best remaining *cached*
            // cost. Periodic global refresh keeps stale priorities in check.
            if(!pq.empty() &&
               (actual>pq.top().cost ||
                (actual==pq.top().cost && e.tie>pq.top().tie))){
                pq.push({actual,e.tie,e.vertex});
                continue;
            }
            add(e.vertex);
            ++restored;
        }
        answer.assign(V,0);
        for(int v:trees)answer[v]=!retained[v];
        return restored==need;
    }

    int numberCut(const vector<uint8_t>& cut) const {
        int c = 0;
        for (int v : trees) c += cut[v];
        return c;
    }

    void greedyComplete(vector<uint8_t>& cut, int samplesCount, double deadline) {
        vector<ll> gain(V);
        int count = numberCut(cut);
        while (count < K && timeOkay(deadline)) {
            removalGains(cut, samplesCount, gain);
            ll bestVal = -1;
            int bestVertex = -1;
            for (int v : trees) {
                if (!cut[v] && gain[v] > bestVal) {
                    bestVal = gain[v];
                    bestVertex = v;
                }
            }
            if (bestVertex < 0) break;
            cut[bestVertex] = 1;
            ++count;
        }
    }

    // Cheap single pass to fill a geometric candidate to exactly K cuts.
    void fillOnce(vector<uint8_t>& cut) {
        if (numberCut(cut) >= K) return;
        vector<ll> gains(V);
        removalGains(cut, 12, gains);
        vector<int> remaining;
        for (int v : trees) if (!cut[v]) remaining.push_back(v);
        sort(remaining.begin(), remaining.end(), [&](int a, int b) {
            return gains[a] > gains[b];
        });
        int current = numberCut(cut);
        for (int v : remaining) {
            if (current >= K) break;
            cut[v] = 1;
            ++current;
        }
    }

    // Water-aware thin blocking barriers. The resulting path can bend.
    vector<int> getBarrier(bool horizontal, int target, int seed) {
        const int INF = 1000000000;
        vector<int> dp(V, INF), parent(V, -1);
        auto at = [&](int progress, int pos) {
            return horizontal ? pos*N + progress : progress*N + pos;
        };
        for (int pos = 0; pos < N; ++pos) {
            int v = at(0, pos);
            dp[pos] = 500*!!d[v] + 6*d[v] + 9*abs(pos-target);
        }
        for (int x = 1; x < N; ++x) {
            for (int p = 0; p < N; ++p) {
                int best = INF, idx = -1;
                for (int delta = -1; delta <= 1; ++delta) {
                    int prev = p+delta;
                    if (prev < 0 || prev >= N) continue;
                    int trial = dp[(x-1)*N+prev];
                    if (trial < best) { best = trial; idx = prev; }
                }
                int v = at(x,p);
                int waviness = (seed&1) ? 5 : 12;
                dp[x*N+p] = best + 500*!!d[v] + 6*d[v] + waviness*abs(p-target);
                parent[x*N+p] = idx;
            }
        }
        int p = 0, best = INF;
        for (int j = 0; j < N; ++j) {
            int score = dp[(N-1)*N+j] + 9*abs(j-target);
            if (score < best) best = score, p = j;
        }
        vector<int> result;
        result.reserve(N);
        for (int x = N-1; x >= 0; --x) {
            result.push_back(at(x,p));
            if (x) p = parent[x*N+p];
        }
        return result;
    }


    // Node-weighted shortest 8-connected separator path. Unlike the old
    // monotone DP path, this can detour backwards through large lakes.
    vector<int> minCostBarrier(bool horizontal, int target, int kind) {
        const int INF = 1000000000;
        vector<int> dist(V,INF),previous(V,-1);
        // We split the grid into sides by removing an 8-connected path.
        // On a 4-neighbor primal grid, 8-connectivity is the correct dual.
        priority_queue<pair<int,int>, vector<pair<int,int>>, greater<pair<int,int>>> q;
        auto cost=[&](int v)->int {
            int y=v/N,x=v%N;
            int perp=horizontal?y:x;
            int off=abs(perp-target);
            int width=kind==0?2:(kind==1?4:7);
            int bias=off*off/width;
            if(d[v]==0)return 2+bias;
            int w=110+(kind==2?0:d[v]*2);
            return w+bias;
        };
        for(int v=0;v<V;v++) {
            int y=v/N,x=v%N;
            if((horizontal?x:y)!=0)continue;
            int edgePerp=horizontal?y:x;
            int startCost=cost(v)+40*abs(edgePerp-target);
            dist[v]=startCost;
            q.push({startCost,v});
        }
        int dest=-1, destCost=INF;
        while(!q.empty()){
            auto [ds,u]=q.top();q.pop();
            if(ds!=dist[u])continue;
            int y=u/N,x=u%N;
            int progress=horizontal?x:y;
            if(progress==N-1){
                int finalCost=ds+40*abs((horizontal?y:x)-target);
                if(finalCost<destCost) destCost=finalCost,dest=u;
                // still explore; water may provide a better target end
                if(ds>destCost)break;
            }
            for(int dy=-1;dy<=1;dy++)for(int dx=-1;dx<=1;dx++){
                if(!dy&&!dx)continue;
                int ny=y+dy,nx=x+dx;
                if(ny<0||ny>=N||nx<0||nx>=N)continue;
                int v=ny*N+nx;
                int nd=ds+cost(v)+((dy&&dx)?5:0);
                if(nd<dist[v]) {dist[v]=nd;previous[v]=u;q.push({nd,v});}
            }
        }
        vector<int> path;
        for(int u=dest;u>=0;u=previous[u]) path.push_back(u);
        return path;
    }

    // Steepest-exchange search: approximate every single swap with exact
    // marginal restoration/removal costs, check top-ranked swaps with the
    // full stochastic objective. This is not a random SA neighbor.
    bool steepestExchange(vector<uint8_t>& current, ll& currentScore,
                          vector<uint8_t>& best, ll& bestValidation,
                          double deadline, int sampleOffset) {
        vector<ll> gain(V), cost(V);
        removalGains(current,12,gain,sampleOffset);
        restoreCosts(current,12,cost);
        vector<int> a,b;
        for(int v:trees){if(current[v])a.push_back(v);else b.push_back(v);}
        if(a.empty()||b.empty())return false;
        int ac=min(14,(int)a.size()),bc=min(24,(int)b.size());
        partial_sort(a.begin(),a.begin()+ac,a.end(),[&](int x,int y){return cost[x]<cost[y];});
        partial_sort(b.begin(),b.begin()+bc,b.end(),[&](int x,int y){return gain[x]>gain[y];});
        struct Move{ll approx;int a,b;};
        vector<Move> cand;
        for(int i=0;i<ac;i++)for(int j=0;j<bc;j++)cand.push_back({cost[a[i]]-gain[b[j]],a[i],b[j]});
        int sz=min(42,(int)cand.size());
        partial_sort(cand.begin(),cand.begin()+sz,cand.end(),[&](auto&x,auto&y){return x.approx<y.approx;});
        vector<uint8_t> trial=current,chosen;
        ll minScore=currentScore;
        for(int i=0;i<sz&&timeOkay(deadline);i++){
            int x=cand[i].a,y=cand[i].b;
            trial[x]=0;trial[y]=1;
            ll score=evaluate(trial,0,32);
            if(score<minScore){minScore=score;chosen=trial;}
            trial[x]=1;trial[y]=0;
        }
        if(chosen.empty())return false;
        current.swap(chosen);currentScore=minScore;
        ll val=evaluate(current,96,96);
        if(val<bestValidation){bestValidation=val;best=current;}
        return true;
    }

    vector<uint8_t> makeCut(const vector<int>& a, const vector<int>& b) {
        vector<uint8_t> cut(V, 0);
        for (int v : a) if (d[v]) cut[v] = 1;
        for (int v : b) if (d[v]) cut[v] = 1;
        return cut;
    }

    // At the beginning of EACH cooling cycle, sample all swap sizes 1..15
    // uniformly. Gradually concentrate the distribution towards 1-point swaps.
    // P(k|p) = exp(-2*p*(k-1)) / sum_j exp(-2*p*(j-1)), where p in [0,1].
    // Reheating resets p to 0, restoring the uniform distribution.
    int swapSize(double progress) {
        const int maxK = min({15, K, (int)trees.size() - K});
        if (maxK == 0) return 0;
        const double lambda = 2.0 * max(0.0, min(1.0, progress));
        double weightSum = 0.0;
        double weights[16];
        for (int k=1; k<=maxK; ++k) {
            weights[k] = exp(-lambda * (k - 1));
            weightSum += weights[k];
        }
        double choice = rnd.real() * weightSum;
        for (int k=1; k<=maxK; ++k) {
            choice -= weights[k];
            if (choice < 0) return k;
        }
        return maxK;
    }

    // Path-preserving nonlocal move. Reroute an 8-connected segment of an
    // existing firebreak through an alternative narrow winding corridor.
    // Restrict to a local box around both endpoints for speed and locality.
    bool reroute(const vector<uint8_t>& cut, vector<uint8_t>& result,
                 const vector<int>& selected, const vector<int>& lowRest) {
        if(selected.size()<12)return false;
        int src=selected[rnd.integer((int)selected.size())];
        // Find nearby 8-connected cut vertices.
        vector<int> seen(V,-1),q(V),pre(V,-1);
        int h=0,t=1;seen[src]=0;q[0]=src;
        while(h<t&&t<200){
            int u=q[h++],y=u/N,x=u%N;
            for(int dy=-1;dy<=1;dy++)for(int dx=-1;dx<=1;dx++){
                if(!dx&&!dy)continue;
                int ny=y+dy,nx=x+dx;
                if(ny<0||ny>=N||nx<0||nx>=N)continue;
                int v=ny*N+nx;
                if(cut[v]&&seen[v]<0){seen[v]=seen[u]+1;pre[v]=u;q[t++]=v;}
            }
        }
        vector<int> ends;
        for(int i=0;i<t;i++)if(seen[q[i]]>=7&&seen[q[i]]<=20)ends.push_back(q[i]);
        if(ends.empty())return false;
        int dst=ends[rnd.integer((int)ends.size())];
        vector<int> oldPath;
        for(int u=dst;u!=src&&u>=0;u=pre[u])oldPath.push_back(u);
        oldPath.push_back(src);
        if(oldPath.size()<7)return false;
        int minY=N,maxY=0,minX=N,maxX=0;
        for(int v:oldPath){minY=min(minY,v/N);maxY=max(maxY,v/N);minX=min(minX,v%N);maxX=max(maxX,v%N);}
        minY=max(0,minY-3);maxY=min(N-1,maxY+3);
        minX=max(0,minX-3);maxX=min(N-1,maxX+3);
        // Early anchors are harder to remove; unused trees in the narrow
        // corridor have a small penalty to produce a true alternative path.
        result=cut;
        for(int v:oldPath)if(v!=src&&v!=dst) result[v]=0;
        const int INF=1000000000;
        vector<int> dist(V,INF),par(V,-1);
        priority_queue<pair<int,int>,vector<pair<int,int>>,greater<pair<int,int>>> pq;
        dist[src]=0;pq.push({0,src});
        int wander=rnd.integer(13);
        while(!pq.empty()){
            auto [cd,u]=pq.top();pq.pop();
            if(cd!=dist[u])continue;
            if(u==dst)break;
            int y=u/N,x=u%N;
            for(int dy=-1;dy<=1;dy++)for(int dx=-1;dx<=1;dx++){
                if(!dx&&!dy)continue;
                int ny=y+dy,nx=x+dx;
                if(ny<minY||ny>maxY||nx<minX||nx>maxX)continue;
                int v=ny*N+nx;
                int cost=result[v]?0:(d[v]?90+3*d[v]:1);
                // perturb path tie breaking and discourage the original path
                if(!result[v]&&d[v])cost+=int((uint64_t(v)*1836311903ULL+uint64_t(wander)*5039ULL)%17);
                int nd=cd+cost;
                if(nd<dist[v]){dist[v]=nd;par[v]=u;pq.push({nd,v});}
            }
        }
        if(dist[dst]>=INF)return false;
        for(int u=dst;u>=0;u=par[u])if(d[u])result[u]=1;
        int delta=numberCut(result)-K;
        if(delta>0){
            // Remove noncritical cuts, preferring those with the lowest
            // restore marginal gain. Avoid destroying the new path.
            vector<int> pathFlag(V,0);
            for(int u=dst;u>=0;u=par[u])pathFlag[u]=1;
            for(int v:lowRest)if(delta>0&&result[v]&&!pathFlag[v])result[v]=0,--delta;
            if(delta>0)return false;
        } else if(delta<0){
            vector<ll> gain(V);
            removalGains(result,8,gain,rnd.integer(NS));
            vector<int> remaining;
            for(int v:trees)if(!result[v])remaining.push_back(v);
            int k=min((int)remaining.size(),-delta);
            if(k>0){
                partial_sort(remaining.begin(),remaining.begin()+k,remaining.end(),
                    [&](int a,int b){return gain[a]>gain[b];});
                for(int i=0;i<k;i++)result[remaining[i]]=1;
            }
        }
        return numberCut(result)==K && result!=cut;
    }

    // Guided 1..7-point swaps. At the start of every SA cooling cycle,
    // each k is equiprobable; close to the end, prefer k=1.
    // Guided destinations are sampled WITHOUT replacement from a
    // temperature-dependent top fraction of latest Lowlink removal gains.
    // Original cut locations are uniform as requested.
    int guidedSwapSize(double progress) {
        const int m = min({7, K, (int)trees.size()-K});
        if (m <= 0) return 0;
        const double lambda = 2.0 * clamp(progress, 0.0, 1.0);
        // One exp() instead of m calls, using a geometric series.
        double decay = exp(-lambda);
        double sum = 0.0, w = 1.0;
        for (int j=0; j<m; ++j) {sum+=w; w*=decay;}
        double r = rnd.real()*sum;
        w=1.0;
        for (int k=1; k<=m; ++k) {
            r-=w;
            if (r < 0.0) return k;
            w*=decay;
        }
        return m;
    }

    // Sample destinations from a temperature-dependent top ranked fraction.
    // Mode 0: baseline fixed top 20%.
    // Mode 1: 100% -> 6%, gently curved schedule.
    // Mode 2: 60% -> 10%.
    // Mode 3: 100% -> 4%, sharp early narrowing.
    // Mode 4: top 20% -> top 1 (or top k for k-point swaps).
    double guidedWidthFraction(double progress) const {
        double q = clamp(1.0-progress, 0.0, 1.0);
#if GUIDE_WIDTH_MODE == 0
        return 0.20;
#elif GUIDE_WIDTH_MODE == 1
        return 0.06 + 0.94 * pow(q, 1.4);
#elif GUIDE_WIDTH_MODE == 2
        return 0.10 + 0.50 * pow(q, 1.3);
#elif GUIDE_WIDTH_MODE == 3
        return 0.04 + 0.96 * q*q;
#elif GUIDE_WIDTH_MODE == 4
        // Ends with a fraction equivalent to exactly one ranked position.
        // guidedMultiSwap enforces pool >= k for k distinct destinations.
        return (1.0 + (0.20 * max(0, (int)trees.size()-K) - 1.0) * pow(q, 1.3)) / max(1, (int)trees.size()-K);
#else
#error "Invalid GUIDE_WIDTH_MODE"
#endif
    }

    bool guidedMultiSwap(vector<uint8_t>& out,
                         const vector<uint8_t>& cut,
                         const vector<int>& cutList,
                         const vector<int>& topGain,
                         double progress) {
        if (cutList.empty() || topGain.empty()) return false;
        int k=guidedSwapSize(progress);
        if (k<=0)return false;
        // Keep the FULL ranking and restrict its sampling width each step.
        // At the cold endpoint the pool is 1 for k=1, or k for k>1;
        // if a ranked cell is already cut, extend the prefix to k live cells.
        // Stale entries are ignored, and the prefix is enlarged if needed.
        int pool=min((int)topGain.size(),max(k,(int)ceil(guidedWidthFraction(progress)*topGain.size())));
        int live=0;
        for(int i=0;i<pool;i++) live+=(!cut[topGain[i]] && d[topGain[i]]>0);
        while(live<k && pool<(int)topGain.size()) {
            int v=topGain[pool++];
            live+=(!cut[v] && d[v]>0);
        }
        k=min(k,live);
        if (k<=0)return false;
        out=cut;
        // First remove k distinct old cuts. out[a] makes duplicates impossible.
        for(int j=0;j<k;j++) {
            int a=-1;
            for(int trial=0;trial<24;trial++) {
                int v=cutList[rnd.integer((int)cutList.size())];
                if(out[v]) {a=v;break;}
            }
            if(a<0) for(int v:cutList) if(out[v]) {a=v;break;}
            if(a<0) return false;
            out[a]=0;
        }
        // Add k distinct new cuts. `cut[v]` (not out[v]) ensures a restored
        // old cut is never accidentally selected again as a destination.
        for(int j=0;j<k;j++) {
            int b=-1;
            for(int trial=0;trial<24;trial++) {
                int v=topGain[rnd.integer(pool)];
                if(!cut[v] && !out[v] && d[v]>0) {b=v;break;}
            }
            if(b<0) for(int i=0;i<pool;i++) {
                int v=topGain[i];
                if(!cut[v] && !out[v] && d[v]>0) {b=v;break;}
            }
            if(b<0)return false;
            out[b]=1;
        }
        return true;
    }

    // Replace exactly k RANDOM existing cuts with the k highest-marginal-
    // gain NEW sites, scored AFTER the k cuts have been restored.
    // The newly cut sites must be uncut in the ORIGINAL solution: this is a
    // genuine k-exchange, not a no-op that picks back the restored bridges.
    bool bestKSwap(vector<uint8_t>& out,
                   const vector<uint8_t>& cut,
                   const vector<int>& cutList,
                   double progress) {
        int available = (int)trees.size() - (int)cutList.size();
        int maxK = min({15, (int)cutList.size(), available});
        if (maxK <= 0) return false;
        int k = swapSize(progress);
        k = min(maxK, k);
        if(k <= 0) return false;

        out = cut;
        // A uniform random k-subset, sampled without replacement. O(k)
        // choices from at most 100 cuts; never prefer low restoration costs.
        int used[16], nused=0;
        for(int j=0;j<k;++j) {
            int idx, tries=0;
            do {
                idx=rnd.integer((int)cutList.size());
                bool seen=false;
                for(int x=0;x<nused;++x)seen |= (used[x]==idx);
                if(!seen)break;
            } while(++tries<100);
            if(tries==100) {
                idx=0;
                for(;idx<(int)cutList.size();++idx) {
                    bool seen=false;
                    for(int x=0;x<nused;++x)seen|=(used[x]==idx);
                    if(!seen)break;
                }
            }
            used[nused++]=idx;
            out[cutList[idx]]=0;
        }
        // Tarjan lowlink computes every uncut site's *exact* one-step
        // marginal gain for a dry-state sample in one DFS. Summing across
        // randomized dry states gives an approximate expected gain.
        vector<ll> gains(V,0);
        removalGains(out, BESTK_SAMPLES, gains,
                     rnd.integer(NS-BESTK_SAMPLES+1));
        vector<int> choices;
        choices.reserve(trees.size()-cutList.size());
        for(int v:trees)if(!cut[v])choices.push_back(v);
        if((int)choices.size()<k) return false;
        // Deterministic tie break keeps each trial reproducible; top-k
        // candidates selected without replacement and cut simultaneously.
        partial_sort(choices.begin(),choices.begin()+k,choices.end(),
            [&](int a,int b){
                if(gains[a]!=gains[b]) return gains[a]>gains[b];
                return a<b;
            });
        for(int j=0;j<k;j++)out[choices[j]]=1;
        return true;
    }

    bool unbiasedSwap(vector<uint8_t>& out,
                      const vector<uint8_t>& cut,
                      const vector<int>& cutList) {
        if(cutList.empty() || (int)trees.size() == (int)cutList.size()) return false;
        out=cut;
        int a=cutList[rnd.integer((int)cutList.size())];
        int b=-1;
        for(int k=0;k<35;k++){
            int v=trees[rnd.integer((int)trees.size())];
            if(!cut[v]) { b=v; break; }
        }
        if(b<0)for(int v:trees)if(!cut[v]){b=v;break;}
        if(b<0)return false;
        out[a]=0;
        out[b]=1;
        return true;
    }

    // Every successful neighbor retains exactly the same number of cuts.
    bool propose(vector<uint8_t>& out, const vector<uint8_t>& cut,
                 const vector<int>& cutList, const vector<int>& lowRest,
                 const vector<int>& highGain, double progress) {
        if (cutList.empty()) return false;
#if BESTK_MODE == 2
        return bestKSwap(out,cut,cutList,progress);
#endif
#if GUIDED_MODE == 1
        // Pure guided 1..7-point neighbor, no LNS.
        return guidedMultiSwap(out,cut,cutList,highGain,progress);
#elif GUIDED_MODE == 2
        // 60% temperature-width guided 1..7 swap; 35% 1..15 swap,
        // 5% completely random single swap.
        int modeRoll=rnd.integer(100);
        if(modeRoll<60) return guidedMultiSwap(out,cut,cutList,highGain,progress);
        if(modeRoll<65) return unbiasedSwap(out,cut,cutList);
#if BESTK_MODE == 1
        return bestKSwap(out,cut,cutList,progress);
#endif
#elif GUIDED_MODE == 3
        // Temperature-dependent mixture: broad early, guided late.
        // p=0: 20% guided, 75% broad, 5% uniform.
        // p=1: 85% guided, 10% broad, 5% uniform.
        int modeRoll=rnd.integer(10000);
        int guidedThreshold=(int)(2000+6500*progress);
        if(modeRoll<guidedThreshold) return guidedMultiSwap(out,cut,cutList,highGain,progress);
        if(modeRoll<guidedThreshold+500) return unbiasedSwap(out,cut,cutList);
#if BESTK_MODE == 1
        return bestKSwap(out,cut,cutList,progress);
#endif
#endif
        out = cut;
        int roll = rnd.integer(100);
        if(roll<8 && reroute(cut,out,cutList,lowRest))return true;
        if (roll>=8 && roll<14) {
            // Move all cut sites on one (popular) row or column by one cell.
            bool row = (rnd.integer(2) == 0);
            int counts[50] = {};
            for (int u : cutList) ++counts[row ? u/N : u%N];
            int maxCount = *max_element(counts, counts+N);
            if (maxCount >= 5) {
                vector<int> lines;
                for (int i=0;i<N;i++) if (counts[i] >= max(5,maxCount-3)) lines.push_back(i);
                int line = lines[rnd.integer((int)lines.size())];
                int off = rnd.integer(2) ? 1 : -1;
                if (line+off >= 0 && line+off < N) {
                    vector<pair<int,int>> changes;
                    for (int j = 0; j < N; ++j) {
                        int old = row ? line*N+j : j*N+line;
                        int nxt = row ? (line+off)*N+j : j*N+line+off;
                        if (cut[old] && d[nxt] && !cut[nxt]) changes.push_back({old,nxt});
                    }
                    if (changes.size() >= 3) {
                        for (auto [a,b]:changes) { out[a]=0; out[b]=1; }
                        return true;
                    }
                }
            }
        }
        int k = swapSize(progress);
        if (k == 0) return false;
        vector<int> removed;
        removed.reserve(k);
        bool cluster = k>=2 && rnd.integer(100)<45;
        int anchor = cutList[rnd.integer((int)cutList.size())];
        for (int j=0;j<k;j++) {
            int a = -1;
            if (cluster) {
                for (int trial=0;trial<30;trial++) {
                    int q=cutList[rnd.integer((int)cutList.size())];
                    int distance=abs(q/N-anchor/N)+abs(q%N-anchor%N);
                    if (out[q] && distance <= k+2) {a=q;break;}
                }
            }
            if (a < 0 && !lowRest.empty() && rnd.integer(100)<60) {
                int top = min((int)lowRest.size(), 30+7*k);
                for (int trial=0;trial<20;trial++) {
                    int q=lowRest[rnd.integer(top)];
                    if (out[q]) { a=q;break; }
                }
            }
            if (a < 0) {
                for (int trial=0;trial<20;trial++) {
                    int q=cutList[rnd.integer((int)cutList.size())];
                    if (out[q]) {a=q;break;}
                }
            }
            if (a < 0) for (int q:cutList) if (out[q]) {a=q;break;}
            if (a < 0) return false;
            removed.push_back(a);
            out[a]=0;
        }
        for (int j=0;j<k;j++) {
            int a=removed[j], b=-1;
            if (rnd.integer(100)<75) {
                int radius=rnd.integer(100)<80?2:5;
                for (int tries=0;tries<35;tries++) {
                    int y=a/N + rnd.integer(2*radius+1)-radius;
                    int x=a%N + rnd.integer(2*radius+1)-radius;
                    if (x<0||x>=N||y<0||y>=N) continue;
                    int q=y*N+x;
                    if (d[q] && !cut[q] && !out[q]) {b=q;break;}
                }
            }
            if (b < 0 && !highGain.empty()) {
                int top=min((int)highGain.size(), 100+15*k);
                for (int tries=0;tries<35;tries++) {
                    int q=highGain[rnd.integer(top)];
                    if (!cut[q] && !out[q]) {b=q;break;}
                }
            }
            if (b<0) {
                for (int tries=0;tries<35;tries++) {
                    int q=trees[rnd.integer((int)trees.size())];
                    if (!cut[q] && !out[q]) {b=q;break;}
                }
            }
            if (b<0) for (int q:trees) if (!cut[q] && !out[q]) {b=q;break;}
            if (b<0) return false;
            out[b]=1;
        }
        return true;
    }

public:
    Solver(int n, int k, vector<int> dd)
        : N(n), K(k), V(n*n), d(move(dd)), adj(V), degree(V,0),
          rnd(0x31aa707f682e1cb9ULL), start(Clock::now()),
          tin(V), low(V), sub(V), separateSum(V), separateSq(V),
          mark(V), queueBuffer(V), comp(V) {
        for (int v=0;v<V;v++) {
            if (d[v]>0) trees.push_back(v);
            int y=v/N,x=v%N;
            if (y) adj[v][degree[v]++]=v-N;
            if (y+1<N) adj[v][degree[v]++]=v+N;
            if (x) adj[v][degree[v]++]=v-1;
            if (x+1<N) adj[v][degree[v]++]=v+1;
        }
        // Guard against log(0) in scorers that don't handle E==0.
        // If <=K trees exist, leaving one is always legal and ensures E>0.
        if (!trees.empty() && (int)trees.size()<=K) {
            sparseEnough = true;
            K=(int)trees.size()-1;
        }
        else K=min(K,(int)trees.size());
        dry.assign(NS,vector<uint8_t>(V));
        for (int s=0;s<NS;s++) {
            for (int v=0;v<V;v++) dry[s][v]=(d[v] && rnd.integer(10)<d[v]);
        }
    }

    vector<uint8_t> solve() {
        vector<uint8_t> empty(V,0);
        if (sparseEnough) {
            // Keep exactly one lowest-d tree: E remains positive, so log(E)
            // is defined even if the checker has no special E==0 handling.
            int survivor=trees[0];
            for (int v:trees) if (d[v]<d[survivor]) survivor=v;
            for (int v:trees) if (v!=survivor) empty[v]=1;
            return empty;
        }
        if (K==0) return empty;
        if (K==(int)trees.size()) {
            for(int v:trees) empty[v]=1;
            return empty;
        }
        vector<uint8_t> best=empty;
        ll bestValidation = LLONG_MAX;
        auto registerCandidate = [&](const vector<uint8_t>& cand) {
            ll value = evaluate(cand,96,96);
            if (value < bestValidation) {bestValidation=value;best=cand;}
        };

        // Reverse-only initial construction. Forward lowlink greedy has
        // been removed entirely; its ~0.5s budget is now available for SA.
        // Keep a small deterministic geometric fallback, and compare it
        // against reverse-DSU candidates on the same validation samples.
        vector<uint8_t> reverseReference;
        if (N == 50 && K >= 2*N-1) {
            for (int kind=0;kind<3;kind++) {
                vector<uint8_t> simple(V,0);
                if (kind==0) {
                    for(int j=0;j<N;j++) if(d[(N/2)*N+j])simple[(N/2)*N+j]=1;
                    for(int i=0;i<N;i++) if(d[i*N+(N/2)])simple[i*N+(N/2)]=1;
                } else if(kind==1) {
                    for(int r:{N/3,2*N/3})for(int j=0;j<N;j++)if(d[r*N+j])simple[r*N+j]=1;
                } else {
                    for(int c:{N/3,2*N/3})for(int i=0;i<N;i++)if(d[i*N+c])simple[i*N+c]=1;
                }
                fillOnce(simple);
                registerCandidate(simple);
            }
        }
        // Always attempt the independent DSU restoration before spending
        // time on optional shape generation.
        if(timeOkay(0.45) && (int)trees.size()>K) {
            vector<uint8_t> reverse;
            if(reverseRestoration(reverse,24,0.69)) {
                reverseReference=reverse;
                registerCandidate(reverse);
#ifdef PROFILE
                cerr<<"reverse "<<elapsed()<<" best="<<bestValidation<<"\n";
#endif
            }
        }
        // Guided second DSU restoration: the first reverse solution is
        // a structural prior; min-max penalization diversifies the second
        // solution by minimizing the largest sampled dry component. There is no forward solution involved.
        if(!reverseReference.empty() && timeOkay(0.69) && (int)trees.size()>K+20) {
            vector<uint8_t> guidedReverse;
            if(reverseRestoration(guidedReverse,32,0.84,&reverseReference,RESTORE_MAX_WEIGHT)) {
                registerCandidate(guidedReverse);
#ifdef PROFILE
                cerr<<"guided-reverse "<<elapsed()<<" best="<<bestValidation<<"\n";
#endif
            }
        }
        // The optional Dijkstra/DP separators are a separate, non-greedy
        // family of initial solutions, retained only when useful.
        double fraction = double(bestValidation) / max(1LL,evaluate(empty,96,96));
        if (timeOkay(0.80) && K>=50 && fraction>0.40) {
            vector<vector<int>> horizontal, vertical;
            for (int t : {10,15,20,25,30,35,40}) {
                horizontal.push_back(getBarrier(true,t,0));
                vertical.push_back(getBarrier(false,t,0));
            }
            for(int t:{12,20,25,30,38})for(int kind=0;kind<2;kind++){
                horizontal.push_back(minCostBarrier(true,t,kind));
                vertical.push_back(minCostBarrier(false,t,kind));
            }
            vector<pair<ll,vector<uint8_t>>> candidates;
            auto addCandidate=[&](const vector<int>&a,const vector<int>&b) {
                vector<uint8_t> c=makeCut(a,b);
                if (numberCut(c)>K) return;
                candidates.emplace_back(evaluate(c,0,6),move(c));
            };
            for (const auto &h:horizontal) for (const auto &v:vertical) {
                if (!timeOkay(0.91)) break;
                addCandidate(h,v);
            }
            for(int i=0;i<(int)horizontal.size();i++) for(int j=i+1;j<(int)horizontal.size();j++) {
                if (!timeOkay(0.91)) break;
                addCandidate(horizontal[i],horizontal[j]);
                addCandidate(vertical[i],vertical[j]);
            }
            sort(candidates.begin(),candidates.end(),
                 [](const auto&a,const auto&b){return a.first<b.first;});
            if(candidates.size()>5) candidates.resize(5);
            for(auto& item:candidates) {
                if(!timeOkay(0.99)) break;
                fillOnce(item.second);
                registerCandidate(item.second);
            }
        }
        // For sparse cases where no geometric/reverse candidate could be
        // completed by the deadline, a one-pass lowlink fill stays valid.
        if(bestValidation==LLONG_MAX) {
            vector<uint8_t> fallback(V,0);
            fillOnce(fallback);
            registerCandidate(fallback);
        }

#ifdef PROFILE
        cerr<<"init "<<elapsed()<<" score="<<bestValidation<<"\n";
#endif
#ifdef INITIAL_ONLY
        return best;
#endif
        // Simulated annealing on common random numbers, with independent validation.
        if(!timeOkay(1.02)) return best;
        vector<uint8_t> current=best, trial;
        const int SA_S=SA_SAMPLES;
        ll currentScore=evaluate(current,0,SA_S);
        ll bestSA=currentScore;
        double saStart=elapsed(), saFinish=TIME_LIMIT-0.030;
        double T0=max(10.0,(double)currentScore/SA_S*0.010);
        double T1=max(0.3,(double)currentScore/SA_S*0.000010);
        int accepted=0, iter=0;
#if BESTK_MODE != 2
        int lastRankUpdate=-10000;
#endif
        int coolingCycle=0, reheats=0;
        vector<int> cycleIterations(COOLING_CYCLES,0);
        vector<int> cutList, lowRest, highGain;
        vector<ll> gains(V), costs(V);
        auto refreshCutList=[&]() {
            cutList.clear();
            for(int v:trees) if(current[v]) cutList.push_back(v);
        };
        // BESTK_MODE=2 computes marginal gains inside every proposed best-K
        // exchange. Precomputing stale cut/uncut rankings is pure overhead.
#if BESTK_MODE != 2
        auto refreshRanks=[&]() {
            removalGains(current,12,gains);
            restoreCosts(current,12,costs);
            lowRest=cutList;
            sort(lowRest.begin(),lowRest.end(),[&](int a,int b){return costs[a]<costs[b];});
            highGain.clear();
            for(int v:trees) if(!current[v]) highGain.push_back(v);
            sort(highGain.begin(),highGain.end(),[&](int a,int b){return gains[a]>gains[b];});
#if GUIDED_MODE == 0
            if(highGain.size()>160)highGain.resize(160);
#endif
            lastRankUpdate=accepted;
        };
#endif
        refreshCutList();
        while(timeOkay(saFinish)) {
            ++iter;
            // Keep the earlier cooling schedule through 1.58s, but let
            // the FINAL (third) SA continue until saFinish (~1.825s).
            // After the third round has cooled completely, hold T1 instead
            // of restarting temperature, so the final extra time is spent
            // on refinement without disrupting the first two rounds.
            const double now = elapsed();
            const double ordinaryFinish = min(saFinish, 1.58);
            const double ordinarySpan = max(0.001, ordinaryFinish-saStart);
            const double ordinaryProgress = max(0.0, (now-saStart)/ordinarySpan);
            int targetCycle;
            double localProgress;
            if constexpr (COOLING_CYCLES == 3) {
                if(ordinaryProgress < 1.0/3.0) {
                    targetCycle = 0;
                    localProgress = ordinaryProgress*3.0;
                } else if (ordinaryProgress < 2.0/3.0) {
                    targetCycle = 1;
                    localProgress = ordinaryProgress*3.0-1.0;
                } else {
                    targetCycle = 2;
                    localProgress = ordinaryProgress*3.0-2.0;
                }
            } else {
                targetCycle = min(COOLING_CYCLES-1, int(ordinaryProgress*COOLING_CYCLES));
                localProgress = ordinaryProgress*COOLING_CYCLES-targetCycle;
            }
            if (targetCycle != coolingCycle) {
                reheats += targetCycle-coolingCycle;
                coolingCycle = targetCycle;
#if BESTK_MODE != 2
                lastRankUpdate = -10000;
#endif
            }
            ++cycleIterations[coolingCycle];
            localProgress = clamp(localProgress,0.0,1.0);
            double temp = max(1e-9, T0 * pow(T1/T0,localProgress) * SA_S);
            // Periodically run exact-guided one-exchange alongside guided 1..7 and 1..15-point swaps.
            if(iter%97==0 && timeOkay(saFinish-0.12)){
                if(steepestExchange(current,currentScore,best,bestValidation,
                                    min(saFinish-0.10,elapsed()+0.045),rnd.integer(NS))){
                    accepted++;
                    refreshCutList();
                    invalidateSACache();
                }
                continue;
            }
#if BESTK_MODE != 2
            if (lastRankUpdate < 0 || (accepted-lastRankUpdate>=24 && timeOkay(saFinish-0.08))) {
                refreshRanks();
            }
#endif
            if(!propose(trial,current,cutList,lowRest,highGain,localProgress))continue;
            ll newVal=evaluateSATrial(current,trial);
            ll delta=newVal-currentScore;
            if(delta<=0 || (rnd.real()<exp(-double(delta)/max(1e-9,temp)))) {
                current.swap(trial);
                currentScore=newVal;
                invalidateSACache();
                ++accepted;
                refreshCutList();
            }
            if(currentScore<bestSA) {
                bestSA=currentScore;
                if(timeOkay(saFinish-0.005)) registerCandidate(current);
            }
        }
        // Validate the final state if possible.
        if(timeOkay(TIME_LIMIT-0.006)) registerCandidate(current);
#ifdef PROFILE
        cerr<<"done "<<elapsed()<<" iter="<<iter<<" accepted="<<accepted
            <<" reheats="<<reheats<<" cycles=";
        for (int i=0;i<COOLING_CYCLES;i++) cerr<<(i?",":"")<<cycleIterations[i];
        cerr<<" score="<<bestValidation<<"\n";
#endif
        return best;
    }
};

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    int N, K;
    if(!(cin >> N >> K)) return 0;
    if(N <= 0 || N > 50 || K < 0) { cout << "0\n"; return 0; }
    // Each dry level is a SINGLE DIGIT (0..9).
    // Accept both contiguous digit rows and whitespace-separated cells.
    // Read characters rather than extracting multi-digit ints.
    vector<int> d(N*N);
    int filled = 0;
    char ch;
    while (filled < N*N && cin.get(ch)) {
        if (ch >= '0' && ch <= '9') d[filled++] = ch - '0';
    }
    if (filled != N*N) {
        cout << "0\n";
        return 0;
    }
    Solver solver(N,K,d);
    auto cut=solver.solve();
    // Output sanitizer: never emit an invalid cell, repeated coordinate or >K cuts.
    vector<int> output;
    for (int v=0;v<N*N;v++) if (cut[v] && d[v]>0 && (int)output.size()<K)
        output.push_back(v);
    cout << output.size() << '\n';
    for (int v : output) cout << v/N << ' ' << v%N << '\n';
    return 0;
}