pub mod library {
    //! ヒューリスティック問題向けの共通ライブラリ。

    // annealing.rsは未使用で、yukicoderで利用できない外部クレートに依存するため公開しない。
    pub mod diagnostics {
        use itertools::Itertools;
        use std::collections::HashMap;
        use std::{
            borrow::Cow,
            cell::RefCell,
            fmt::Display,
            io::{self, IsTerminal as _},
            rc::Rc,
            time::{Duration, Instant},
        };

        /// プログラムの経過時間を計測する。
        pub struct Timer {
            start: Instant,
        }

        impl Default for Timer {
            fn default() -> Self {
                Self {
                    start: Instant::now(),
                }
            }
        }

        impl Timer {
            pub fn elapsed(&self) -> Duration {
                self.start.elapsed()
            }
        }

        /// 実行時間を計測・集計する構造体
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::diagnostics::Perf;
        /// use cp_lib_rs::perf;
        ///
        /// // 計測グループを作成する
        /// // dropされるときに計測結果を出力する
        /// let mut perf = Perf::new("Group");
        /// let mut _sum = 0u64;
        ///
        /// for i in 0..100000 {
        ///     // start-stop間の処理時間を計測する
        ///     let sw = perf.start("sum");
        ///     _sum += i;
        ///     sw.stop();
        ///
        ///     // 名前は&'static strでもStringでもOK
        ///     let sw = perf.start(format!("no-op"));
        ///     sw.stop();
        /// }
        ///
        /// let mut _sum_sqrt = 0f64;
        ///
        /// for i in 0..100000 {
        ///     // 明示的にstop()を呼ばなくても、スコープを抜ける際に自動でstopする
        ///     let _sw = perf.start("sum sqrt");
        ///     _sum_sqrt += (i as f64).sqrt();
        /// }
        ///
        /// let mut _sum_sq = 0u64;
        ///
        /// for i in 0..100000 {
        ///     // perfのインスタンス化が面倒な場合はシングルトンを使う
        ///     let sw = Perf::start_singleton("sum sq");
        ///     _sum_sq += i * i;
        ///     sw.stop();
        /// }
        ///
        /// let mut _sum_pow = 0f64;
        ///
        /// for i in 0..100000 {
        ///     // マクロでショートハンド化もできる
        ///     _sum_pow += perf!("sum pow", (i as f64).powf(3.1415926));
        /// }
        /// ```
        pub struct Perf {
            name: Option<Cow<'static, str>>,
            measures: HashMap<Cow<'static, str>, Measure>,
        }

        impl Perf {
            thread_local!(static SINGLETON: Rc<RefCell<Perf>> = Rc::new(RefCell::new(Perf::new("Singleton"))));

            /// 新しい Perf インスタンスを作成する
            pub fn new(name: impl Into<Cow<'static, str>>) -> Self {
                Self {
                    name: Some(name.into()),
                    measures: HashMap::default(),
                }
            }

            /// 匿名の Perf インスタンスを作成する
            pub fn new_anonymous() -> Self {
                Self {
                    name: None,
                    measures: HashMap::default(),
                }
            }

            /// 計測を開始する
            pub fn start(&mut self, name: impl Into<Cow<'static, str>>) -> StopWatch<&mut Perf> {
                let name = name.into();

                StopWatch {
                    start: Instant::now(),
                    perf: self,
                    name,
                }
            }

            /// 計測を開始する（スレッドごとのシングルトン）
            pub fn start_singleton(
                name: impl Into<Cow<'static, str>>,
            ) -> StopWatch<Rc<RefCell<Perf>>> {
                Self::SINGLETON.with(|perf| StopWatch {
                    start: Instant::now(),
                    perf: perf.clone(),
                    name: name.into(),
                })
            }
        }

        impl Drop for Perf {
            fn drop(&mut self) {
                if self.measures.is_empty() {
                    return;
                }

                // ターミナルかどうかで色を付けるかどうか分岐
                let is_tty = io::stderr().is_terminal();
                let name = self.name.as_deref().unwrap_or("Anonymous");

                if is_tty {
                    // コンソール → 色付き
                    eprintln!("\x1b[35m[{name}] Performance measures\x1b[0m");
                } else {
                    // リダイレクト → 色なし
                    eprintln!("[{name}] Performance measures");
                }

                for (name, measure) in self
                    .measures
                    .iter()
                    .sorted_unstable_by(|(a, _), (b, _)| a.as_ref().cmp(b.as_ref()))
                {
                    eprintln!("{name}: {measure}");
                }
            }
        }

        #[macro_export]
        macro_rules! perf {
            ($name: expr, $body: expr) => {{
                let sw = Perf::start_singleton($name);
                let result = $body;
                sw.stop();
                result
            }};
        }

        pub struct StopWatch<M: WithMut<Perf>> {
            start: Instant,
            perf: M,
            name: Cow<'static, str>,
        }

        impl<M: WithMut<Perf>> StopWatch<M> {
            pub fn stop(self) {
                std::mem::drop(self);
            }
        }

        impl<M: WithMut<Perf>> Drop for StopWatch<M> {
            fn drop(&mut self) {
                let duration = self.start.elapsed();
                let name = std::mem::take(&mut self.name);

                self.perf.with_mut(|perf| {
                    perf.measures
                        .entry(name)
                        .or_default()
                        .add_measure(&duration)
                })
            }
        }

        #[derive(Debug, Clone, Copy, Default)]
        struct Measure {
            sum: f64,
            sum_sq: f64,
            cnt: usize,
        }

        impl Measure {
            fn add_measure(&mut self, duration: &Duration) {
                let sec = duration.as_secs_f64();
                self.sum += sec;
                self.sum_sq += sec * sec;
                self.cnt += 1;
            }

            fn mean(&self) -> Duration {
                assert_ne!(self.cnt, 0);
                Duration::from_secs_f64(self.sum / self.cnt as f64)
            }

            fn std_dev(&self) -> Duration {
                assert_ne!(self.cnt, 0);

                // Durationへの丸めは最後に行い、浮動小数点の丸め誤差も吸収する。
                let mean = self.sum / self.cnt as f64;
                let variance = ((self.sum_sq / self.cnt as f64) - (mean * mean)).max(0.0);
                Duration::from_secs_f64(variance.sqrt())
            }
        }

        impl Display for Measure {
            fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
                write!(
                    f,
                    "{:?} ± {:?} ({} samples)",
                    self.mean(),
                    self.std_dev(),
                    self.cnt
                )
            }
        }

        /// `&mut T` と `Rc<RefCell<T>>` を統一的に扱って処理を行うためのトレイト
        pub trait WithMut<T> {
            fn with_mut<R>(&mut self, f: impl FnOnce(&mut T) -> R) -> R;
        }

        impl<T> WithMut<T> for &mut T {
            fn with_mut<R>(&mut self, f: impl FnOnce(&mut T) -> R) -> R {
                f(*self)
            }
        }

        impl<T> WithMut<T> for Rc<RefCell<T>> {
            fn with_mut<R>(&mut self, f: impl FnOnce(&mut T) -> R) -> R {
                let mut guard = self.borrow_mut();
                f(&mut *guard)
            }
        }

        #[cfg(test)]
        mod tests {
            use super::*;

            #[test]
            fn standard_deviation_handles_subnanosecond_mean() {
                let mut measure = Measure::default();
                for nanos in [1000, 1001, 1001] {
                    measure.add_measure(&Duration::from_nanos(nanos));
                }
                assert_eq!(measure.mean(), Duration::from_nanos(1001));
                assert_eq!(measure.std_dev(), Duration::ZERO);
                assert!(!measure.to_string().is_empty());
            }

            #[test]
            fn standard_deviation_matches_known_values() {
                for (samples, expected) in
                    [(vec![1000, 2000], 500), (vec![1001; 1000], 0), (vec![0], 0)]
                {
                    let mut measure = Measure::default();
                    for nanos in samples {
                        measure.add_measure(&Duration::from_nanos(nanos));
                    }
                    assert_eq!(measure.std_dev(), Duration::from_nanos(expected));
                }
            }
        }
    }
    pub mod graph {
        use std::ops::Index;

        /// 隣接リスト形式のグラフを行圧縮したグラフ
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::graph::{RowCompressedGraph, UnweightedEdge};
        /// use itertools::Itertools;
        ///
        /// let nodes = vec!["hoge", "fuga", "piyo"];
        /// let edges = vec![
        ///     (1, UnweightedEdge::new(2)),
        ///     (2, UnweightedEdge::new(0)),
        ///     (0, UnweightedEdge::new(1)),
        /// ];
        /// let graph = RowCompressedGraph::new(3, nodes, edges);
        ///
        /// assert_eq!(graph.len(), 3);
        /// assert_eq!(graph.nodes()[0], "hoge");
        /// assert_eq!(graph[0].iter().sorted().copied().collect_vec(), [UnweightedEdge::new(1)]);
        /// ```
        #[derive(Debug, Clone)]
        pub struct RowCompressedGraph<V, E> {
            nodes: Vec<V>,
            edges: Vec<E>,
            pos: Vec<usize>,
            len: usize,
        }

        impl<V, E> RowCompressedGraph<V, E> {
            pub fn new(n: usize, nodes: Vec<V>, edges: Vec<(usize, E)>) -> Self {
                assert!(nodes.len() == n);
                let mut raw_edges = edges;
                raw_edges.sort_unstable_by_key(|(u, _)| *u);
                let mut edges = Vec::with_capacity(raw_edges.len());
                let mut pos = vec![0; n + 1];
                let mut index = 0;

                for (i, (u, e)) in raw_edges.into_iter().enumerate() {
                    assert!(u < n);

                    while index < u {
                        index += 1;
                        pos[index] = i;
                    }

                    edges.push(e);
                }

                pos[index + 1..].fill(edges.len());

                Self {
                    nodes,
                    edges,
                    pos,
                    len: n,
                }
            }

            #[allow(clippy::len_without_is_empty)]
            pub fn len(&self) -> usize {
                self.len
            }

            pub fn nodes(&self) -> &[V] {
                &self.nodes
            }

            pub fn nodes_mut(&mut self) -> &mut [V] {
                &mut self.nodes
            }
        }

        impl<E> RowCompressedGraph<(), E> {
            pub fn new_edges(n: usize, edges: Vec<(usize, E)>) -> Self {
                Self::new(n, vec![(); n], edges)
            }
        }

        impl<V, E> Index<usize> for RowCompressedGraph<V, E> {
            type Output = [E];

            fn index(&self, index: usize) -> &Self::Output {
                let start = self.pos[index];
                let end = self.pos[index + 1];
                &self.edges[start..end]
            }
        }

        pub trait EdgeTo {
            fn to(&self) -> usize;
        }

        pub trait EdgeWeight<T>: EdgeTo {
            fn weight(&self) -> &T;
        }

        #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash)]
        pub struct UnweightedEdge {
            to: usize,
        }

        impl UnweightedEdge {
            pub fn new(to: usize) -> Self {
                Self { to }
            }
        }

        impl EdgeTo for UnweightedEdge {
            fn to(&self) -> usize {
                self.to
            }
        }

        #[derive(Debug, Clone, Copy)]
        pub struct WeightedEdge<T> {
            to: usize,
            weight: T,
        }

        impl<T> WeightedEdge<T> {
            pub fn new(to: usize, weight: T) -> Self {
                Self { to, weight }
            }
        }

        impl<T> EdgeTo for WeightedEdge<T> {
            fn to(&self) -> usize {
                self.to
            }
        }

        impl<T> EdgeWeight<T> for WeightedEdge<T> {
            fn weight(&self) -> &T {
                &self.weight
            }
        }
    }
    pub mod grid {
        use std::{
            fmt::Display,
            ops::{
                Add, AddAssign, Div, DivAssign, Index, IndexMut, Mul, MulAssign, Neg, Sub,
                SubAssign,
            },
        };

        #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash, Default)]
        pub struct Coord {
            row: u8,
            col: u8,
        }

        impl Coord {
            pub const fn new(row: usize, col: usize) -> Self {
                Self {
                    row: row as u8,
                    col: col as u8,
                }
            }

            pub const fn row(&self) -> usize {
                self.row as usize
            }

            pub const fn col(&self) -> usize {
                self.col as usize
            }

            pub fn in_map(&self, size: usize) -> bool {
                self.row < size as u8 && self.col < size as u8
            }

            pub const fn to_index(&self, size: usize) -> CoordIndex {
                CoordIndex(self.row as usize * size + self.col as usize)
            }

            pub const fn dist(&self, other: &Self) -> usize {
                Self::dist_1d(self.row, other.row) + Self::dist_1d(self.col, other.col)
            }

            const fn dist_1d(x0: u8, x1: u8) -> usize {
                x0.abs_diff(x1) as usize
            }
        }

        impl Display for Coord {
            fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
                write!(f, "({}, {})", self.row, self.col)
            }
        }

        #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash, Default)]
        pub struct CoordIndex(pub usize);

        impl CoordIndex {
            pub const fn to_coord(&self, size: usize) -> Coord {
                Coord::new(self.0 / size, self.0 % size)
            }
        }

        #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash, Default)]
        pub struct CoordDiff {
            dr: i8,
            dc: i8,
        }

        impl CoordDiff {
            pub const fn new(dr: i32, dc: i32) -> Self {
                Self {
                    dr: dr as i8,
                    dc: dc as i8,
                }
            }

            pub const fn invert(&self) -> Self {
                Self {
                    dr: -self.dr,
                    dc: -self.dc,
                }
            }

            pub const fn dr(&self) -> i32 {
                self.dr as i32
            }

            pub const fn dc(&self) -> i32 {
                self.dc as i32
            }
        }

        impl Display for CoordDiff {
            fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
                write!(f, "({}, {})", self.dr, self.dc)
            }
        }

        impl Add<CoordDiff> for Coord {
            type Output = Coord;

            fn add(self, rhs: CoordDiff) -> Self::Output {
                Coord {
                    row: self.row.wrapping_add(rhs.dr as u8),
                    col: self.col.wrapping_add(rhs.dc as u8),
                }
            }
        }

        impl AddAssign<CoordDiff> for Coord {
            fn add_assign(&mut self, rhs: CoordDiff) {
                self.row = self.row.wrapping_add(rhs.dr as u8);
                self.col = self.col.wrapping_add(rhs.dc as u8);
            }
        }

        impl Add<CoordDiff> for CoordDiff {
            type Output = CoordDiff;

            fn add(self, rhs: CoordDiff) -> Self::Output {
                CoordDiff {
                    dr: self.dr + rhs.dr,
                    dc: self.dc + rhs.dc,
                }
            }
        }

        impl Add<Coord> for CoordDiff {
            type Output = Coord;

            fn add(self, rhs: Coord) -> Self::Output {
                rhs + self
            }
        }

        impl Sub<CoordDiff> for Coord {
            type Output = Coord;

            fn sub(self, rhs: CoordDiff) -> Self::Output {
                Coord {
                    row: self.row.wrapping_sub(rhs.dr as u8),
                    col: self.col.wrapping_sub(rhs.dc as u8),
                }
            }
        }

        impl SubAssign<CoordDiff> for Coord {
            fn sub_assign(&mut self, rhs: CoordDiff) {
                *self = *self - rhs;
            }
        }

        // Differences outside the i8 range wrap, as in CoordDiff::new.
        impl Sub<Coord> for Coord {
            type Output = CoordDiff;

            fn sub(self, rhs: Coord) -> Self::Output {
                CoordDiff::new(
                    self.row as i32 - rhs.row as i32,
                    self.col as i32 - rhs.col as i32,
                )
            }
        }

        impl Sub<CoordDiff> for CoordDiff {
            type Output = CoordDiff;

            fn sub(self, rhs: CoordDiff) -> Self::Output {
                CoordDiff {
                    dr: self.dr - rhs.dr,
                    dc: self.dc - rhs.dc,
                }
            }
        }

        impl SubAssign<CoordDiff> for CoordDiff {
            fn sub_assign(&mut self, rhs: CoordDiff) {
                *self = *self - rhs;
            }
        }

        // Scalar operations use i32 arithmetic, then narrow as in CoordDiff::new.
        impl Mul<i32> for CoordDiff {
            type Output = CoordDiff;

            fn mul(self, rhs: i32) -> Self::Output {
                CoordDiff::new(self.dr() * rhs, self.dc() * rhs)
            }
        }

        impl Mul<CoordDiff> for i32 {
            type Output = CoordDiff;

            fn mul(self, rhs: CoordDiff) -> Self::Output {
                rhs * self
            }
        }

        impl MulAssign<i32> for CoordDiff {
            fn mul_assign(&mut self, rhs: i32) {
                *self = *self * rhs;
            }
        }

        impl Div<i32> for CoordDiff {
            type Output = CoordDiff;

            fn div(self, rhs: i32) -> Self::Output {
                CoordDiff::new(self.dr() / rhs, self.dc() / rhs)
            }
        }

        impl DivAssign<i32> for CoordDiff {
            fn div_assign(&mut self, rhs: i32) {
                *self = *self / rhs;
            }
        }

        impl AddAssign<CoordDiff> for CoordDiff {
            fn add_assign(&mut self, rhs: CoordDiff) {
                self.dr += rhs.dr;
                self.dc += rhs.dc;
            }
        }

        impl Neg for CoordDiff {
            type Output = Self;

            fn neg(self) -> Self::Output {
                Self {
                    dr: -self.dr,
                    dc: -self.dc,
                }
            }
        }

        pub const ADJACENTS: [CoordDiff; 4] = [
            CoordDiff::new(-1, 0),
            CoordDiff::new(0, 1),
            CoordDiff::new(1, 0),
            CoordDiff::new(0, -1),
        ];

        pub const U: usize = 0;
        pub const R: usize = 1;
        pub const D: usize = 2;
        pub const L: usize = 3;
        pub const DIRECTIONS: [char; 4] = ['U', 'R', 'D', 'L'];

        #[derive(Debug, Clone)]
        pub struct Map2d<T> {
            size: usize,
            map: Vec<T>,
        }

        impl<T> Map2d<T> {
            pub fn new(map: Vec<T>, size: usize) -> Self {
                debug_assert!(size * size == map.len());
                Self { size, map }
            }

            pub fn from_fn(mut f: impl FnMut(Coord) -> T, size: usize) -> Self {
                let mut map = Vec::with_capacity(size * size);

                for row in 0..size {
                    for col in 0..size {
                        map.push(f(Coord::new(row, col)));
                    }
                }

                Self { size, map }
            }

            pub fn iter(&self) -> impl Iterator<Item = &T> {
                self.map.iter()
            }
        }

        impl<T: Default + Clone> Map2d<T> {
            pub fn with_default(size: usize) -> Self {
                let map = vec![T::default(); size * size];
                Self { size, map }
            }
        }

        impl<T> Index<Coord> for Map2d<T> {
            type Output = T;

            #[inline]
            fn index(&self, coordinate: Coord) -> &Self::Output {
                &self.map[coordinate.to_index(self.size).0]
            }
        }

        impl<T> IndexMut<Coord> for Map2d<T> {
            #[inline]
            fn index_mut(&mut self, coordinate: Coord) -> &mut Self::Output {
                &mut self.map[coordinate.to_index(self.size).0]
            }
        }

        impl<T> Index<&Coord> for Map2d<T> {
            type Output = T;

            #[inline]
            fn index(&self, coordinate: &Coord) -> &Self::Output {
                &self.map[coordinate.to_index(self.size).0]
            }
        }

        impl<T> IndexMut<&Coord> for Map2d<T> {
            #[inline]
            fn index_mut(&mut self, coordinate: &Coord) -> &mut Self::Output {
                &mut self.map[coordinate.to_index(self.size).0]
            }
        }

        impl<T> Index<usize> for Map2d<T> {
            type Output = [T];

            #[inline]
            fn index(&self, row: usize) -> &Self::Output {
                let begin = row * self.size;
                let end = begin + self.size;
                &self.map[begin..end]
            }
        }

        impl<T> IndexMut<usize> for Map2d<T> {
            #[inline]
            fn index_mut(&mut self, row: usize) -> &mut Self::Output {
                let begin = row * self.size;
                let end = begin + self.size;
                &mut self.map[begin..end]
            }
        }

        impl<T> Index<CoordIndex> for Map2d<T> {
            type Output = T;

            #[inline]
            fn index(&self, index: CoordIndex) -> &Self::Output {
                &self.map[index.0]
            }
        }

        impl<T> IndexMut<CoordIndex> for Map2d<T> {
            #[inline]
            fn index_mut(&mut self, index: CoordIndex) -> &mut Self::Output {
                &mut self.map[index.0]
            }
        }

        #[derive(Debug, Clone)]
        pub struct ConstMap2d<T, const N: usize, const LEN: usize> {
            map: [T; LEN],
        }

        impl<T, const N: usize, const LEN: usize> ConstMap2d<T, N, LEN> {
            pub fn new(map: Vec<T>) -> Self {
                assert_eq!(LEN, N * N);
                assert_eq!(map.len(), LEN);
                let Ok(map) = map.try_into() else {
                    unreachable!("長さ検証済みのVecから配列への変換に失敗しました。");
                };

                Self { map }
            }

            pub fn from_fn(mut f: impl FnMut(Coord) -> T) -> Self {
                assert_eq!(LEN, N * N);
                let map = std::array::from_fn(|index| f(Coord::new(index / N, index % N)));

                Self { map }
            }

            pub fn iter(&self) -> impl Iterator<Item = &T> {
                self.map.iter()
            }
        }

        impl<T: Default, const N: usize, const LEN: usize> ConstMap2d<T, N, LEN> {
            pub fn with_default() -> Self {
                assert_eq!(LEN, N * N);
                let map = std::array::from_fn(|_| T::default());

                Self { map }
            }
        }

        impl<T, const N: usize, const LEN: usize> Index<Coord> for ConstMap2d<T, N, LEN> {
            type Output = T;

            #[inline]
            fn index(&self, coordinate: Coord) -> &Self::Output {
                &self.map[coordinate.to_index(N).0]
            }
        }

        impl<T, const N: usize, const LEN: usize> IndexMut<Coord> for ConstMap2d<T, N, LEN> {
            #[inline]
            fn index_mut(&mut self, coordinate: Coord) -> &mut Self::Output {
                &mut self.map[coordinate.to_index(N).0]
            }
        }

        impl<T, const N: usize, const LEN: usize> Index<&Coord> for ConstMap2d<T, N, LEN> {
            type Output = T;

            #[inline]
            fn index(&self, coordinate: &Coord) -> &Self::Output {
                &self.map[coordinate.to_index(N).0]
            }
        }

        impl<T, const N: usize, const LEN: usize> IndexMut<&Coord> for ConstMap2d<T, N, LEN> {
            #[inline]
            fn index_mut(&mut self, coordinate: &Coord) -> &mut Self::Output {
                &mut self.map[coordinate.to_index(N).0]
            }
        }

        impl<T, const N: usize, const LEN: usize> Index<usize> for ConstMap2d<T, N, LEN> {
            type Output = [T];

            #[inline]
            fn index(&self, row: usize) -> &Self::Output {
                let begin = row * N;
                let end = begin + N;
                &self.map[begin..end]
            }
        }

        impl<T, const N: usize, const LEN: usize> IndexMut<usize> for ConstMap2d<T, N, LEN> {
            #[inline]
            fn index_mut(&mut self, row: usize) -> &mut Self::Output {
                let begin = row * N;
                let end = begin + N;
                &mut self.map[begin..end]
            }
        }

        impl<T, const N: usize, const LEN: usize> Index<CoordIndex> for ConstMap2d<T, N, LEN> {
            type Output = T;

            #[inline]
            fn index(&self, index: CoordIndex) -> &Self::Output {
                &self.map[index.0]
            }
        }

        impl<T, const N: usize, const LEN: usize> IndexMut<CoordIndex> for ConstMap2d<T, N, LEN> {
            #[inline]
            fn index_mut(&mut self, index: CoordIndex) -> &mut Self::Output {
                &mut self.map[index.0]
            }
        }

        #[cfg(test)]
        mod test {
            use super::{ConstMap2d, Coord, CoordDiff, Map2d};

            #[test]
            fn coord_add() {
                let c = Coord::new(2, 4);
                let d = CoordDiff::new(-3, 5);
                let actual = c + d;

                let expected = Coord::new(!0, 9);
                assert_eq!(expected, actual);
            }

            #[test]
            fn coord_add_assign() {
                let mut c = Coord::new(2, 4);
                let d = CoordDiff::new(-3, 5);
                c += d;

                let expected = Coord::new(!0, 9);
                assert_eq!(expected, c);
            }

            #[test]
            fn coord_translation_round_trip() {
                for row in 0..=255 {
                    let c = Coord::new(row, 255 - row);
                    for dr in -128..=127 {
                        let d = CoordDiff::new(dr, -1 - dr);
                        assert_eq!(c + d, d + c);
                        assert_eq!((c + d) - d, c);
                        assert_eq!((c + d) - c, d);

                        let mut translated = c;
                        translated -= d;
                        assert_eq!(translated + d, c);
                    }
                }
            }

            #[test]
            fn coord_difference() {
                assert_eq!(Coord::new(2, 9) - Coord::new(5, 4), CoordDiff::new(-3, 5));
                assert_eq!(
                    Coord::new(255, 0) - Coord::new(0, 255),
                    CoordDiff::new(-1, 1)
                );
            }

            #[test]
            fn coord_diff_arithmetic() {
                let a = CoordDiff::new(-3, 5);
                let b = CoordDiff::new(7, -2);
                assert_eq!(a + b, CoordDiff::new(4, 3));
                assert_eq!(a - b, CoordDiff::new(-10, 7));
                assert_eq!(-a, a.invert());
                assert_eq!(a + (-a), CoordDiff::default());

                let mut d = a;
                d += b;
                assert_eq!(d, a + b);
                d -= b;
                assert_eq!(d, a);
            }

            #[test]
            fn coord_diff_scalar_arithmetic() {
                let d = CoordDiff::new(-3, 5);
                assert_eq!(d * -2, CoordDiff::new(6, -10));
                assert_eq!(-2 * d, d * -2);
                assert_eq!(d * 0, CoordDiff::default());
                assert_eq!(d / 2, CoordDiff::new(-1, 2));
                assert_eq!(d / -2, CoordDiff::new(1, -2));
                assert_eq!(CoordDiff::new(64, -128) * 2, CoordDiff::new(-128, 0));
                assert_eq!(CoordDiff::new(-128, 0) / -1, CoordDiff::new(-128, 0));
                assert_eq!(d / 256, CoordDiff::default());

                let mut scaled = d;
                scaled *= -2;
                assert_eq!(scaled, CoordDiff::new(6, -10));
                scaled /= -2;
                assert_eq!(scaled, d);
            }

            #[test]
            #[should_panic]
            fn coord_diff_divide_by_zero() {
                let _ = CoordDiff::new(1, 2) / 0;
            }

            #[test]
            fn map_new() {
                let map = Map2d::new(vec![0, 1, 2, 3], 2);
                let actual = map[Coord::new(1, 0)];
                let expected = 2;
                assert_eq!(expected, actual);
            }

            #[test]
            fn map_default() {
                let map = Map2d::with_default(2);
                let actual = map[Coord::new(1, 0)];
                let expected = 0;
                assert_eq!(expected, actual);
            }

            #[test]
            fn const_map_new() {
                let map = ConstMap2d::<_, 2, 4>::new(vec![0, 1, 2, 3]);
                let actual = map[Coord::new(1, 0)];
                let expected = 2;
                assert_eq!(expected, actual);
            }

            #[test]
            fn const_map_default() {
                let map = ConstMap2d::<_, 2, 4>::with_default();
                let actual = map[Coord::new(1, 0)];
                let expected = 0;
                assert_eq!(expected, actual);
            }
        }
    }
    pub mod random {
        use num::{FromPrimitive, PrimInt};
        use std::ops::{Range, RangeInclusive};

        /// 外部クレートに依存しないSplitMix64乱数生成器。
        #[derive(Debug, Clone)]
        pub struct Random {
            state: u64,
        }

        impl Random {
            pub fn new(seed: u64) -> Self {
                Self { state: seed }
            }

            pub fn next_u64(&mut self) -> u64 {
                self.state = self.state.wrapping_add(0x9e3779b97f4a7c15);
                let mut z = self.state;
                z = (z ^ (z >> 30)).wrapping_mul(0xbf58476d1ce4e5b9);
                z = (z ^ (z >> 27)).wrapping_mul(0x94d049bb133111eb);
                z ^ (z >> 31)
            }

            /// 棄却法により剰余の偏りを除き、0..upperから一様に選ぶ。
            pub fn index(&mut self, upper: usize) -> usize {
                assert!(upper > 0);
                let bound = upper as u64;
                let threshold = bound.wrapping_neg() % bound;
                loop {
                    let value = self.next_u64();
                    if value >= threshold {
                        return (value % bound) as usize;
                    }
                }
            }
        }

        /// 高速な乱数生成の拡張トレイト
        pub trait RandExtension {
            fn fast_gen_range_u64x1<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
                &mut self,
                range: R,
            ) -> T;

            fn fast_gen_range_u32x1<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
                &mut self,
                range: R,
            ) -> T;

            fn fast_gen_range_u32x2<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
            ) -> (T0, T1);

            fn fast_gen_range_u16x1<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
                &mut self,
                range: R,
            ) -> T;

            fn fast_gen_range_u16x2<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
            ) -> (T0, T1);

            fn fast_gen_range_u16x3<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                T2: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
                R2: BoundedRange<T2>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
                range2: R2,
            ) -> (T0, T1, T2);

            fn fast_gen_range_u16x4<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                T2: PrimInt + FromPrimitive,
                T3: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
                R2: BoundedRange<T2>,
                R3: BoundedRange<T3>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
                range2: R2,
                range3: R3,
            ) -> (T0, T1, T2, T3);
        }

        impl RandExtension for Random {
            fn fast_gen_range_u64x1<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
                &mut self,
                range: R,
            ) -> T {
                let rand_value = self.next_u64();
                gen_range_u64(range, rand_value)
            }

            fn fast_gen_range_u32x1<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
                &mut self,
                range: R,
            ) -> T {
                let rand_value = self.next_u64();
                gen_range_u32(range, rand_value as u32)
            }

            fn fast_gen_range_u32x2<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
            ) -> (T0, T1) {
                let rand_value = self.next_u64();
                let v0 = gen_range_u32(range0, rand_value as u32);
                let v1 = gen_range_u32(range1, (rand_value >> 32) as u32);
                (v0, v1)
            }

            fn fast_gen_range_u16x1<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
                &mut self,
                range: R,
            ) -> T {
                let rand_value = self.next_u64();
                gen_range_u16(range, rand_value as u16)
            }

            fn fast_gen_range_u16x2<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
            ) -> (T0, T1) {
                let rand_value = self.next_u64();
                let v0 = gen_range_u16(range0, rand_value as u16);
                let v1 = gen_range_u16(range1, (rand_value >> 16) as u16);
                (v0, v1)
            }

            fn fast_gen_range_u16x3<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                T2: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
                R2: BoundedRange<T2>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
                range2: R2,
            ) -> (T0, T1, T2) {
                let rand_value = self.next_u64();
                let v0 = gen_range_u16(range0, rand_value as u16);
                let v1 = gen_range_u16(range1, (rand_value >> 16) as u16);
                let v2 = gen_range_u16(range2, (rand_value >> 32) as u16);
                (v0, v1, v2)
            }

            fn fast_gen_range_u16x4<
                T0: PrimInt + FromPrimitive,
                T1: PrimInt + FromPrimitive,
                T2: PrimInt + FromPrimitive,
                T3: PrimInt + FromPrimitive,
                R0: BoundedRange<T0>,
                R1: BoundedRange<T1>,
                R2: BoundedRange<T2>,
                R3: BoundedRange<T3>,
            >(
                &mut self,
                range0: R0,
                range1: R1,
                range2: R2,
                range3: R3,
            ) -> (T0, T1, T2, T3) {
                let rand_value = self.next_u64();
                let v0 = gen_range_u16(range0, rand_value as u16);
                let v1 = gen_range_u16(range1, (rand_value >> 16) as u16);
                let v2 = gen_range_u16(range2, (rand_value >> 32) as u16);
                let v3 = gen_range_u16(range3, (rand_value >> 48) as u16);
                (v0, v1, v2, v3)
            }
        }

        fn gen_range_u64<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
            range: R,
            rand_value: u64,
        ) -> T {
            assert!(!range.is_empty(), "cannot sample empty range");
            let width = range.width().to_u64().expect("width must fit in u64");
            let start = range.start();
            let value = (((width as u128) * (rand_value as u128)) >> 64) as u64;
            T::from_u64(value).expect("value must fit in T") + start
        }

        fn gen_range_u32<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
            range: R,
            rand_value: u32,
        ) -> T {
            assert!(!range.is_empty(), "cannot sample empty range");
            let width = range.width().to_u32().expect("width must fit in u32");
            let start = range.start();
            let value = (((width as u64) * (rand_value as u64)) >> 32) as u32;
            T::from_u32(value).expect("value must fit in T") + start
        }

        fn gen_range_u16<T: PrimInt + FromPrimitive, R: BoundedRange<T>>(
            range: R,
            rand_value: u16,
        ) -> T {
            assert!(!range.is_empty(), "cannot sample empty range");
            let width = range.width().to_u16().expect("width must fit in u16");
            let start = range.start();
            let value = (((width as u32) * (rand_value as u32)) >> 16) as u16;
            T::from_u16(value).expect("value must fit in T") + start
        }

        pub trait BoundedRange<T> {
            fn start(&self) -> T;
            fn width(&self) -> T;
            fn is_empty(&self) -> bool;
        }

        impl<T: PrimInt> BoundedRange<T> for Range<T> {
            fn start(&self) -> T {
                self.start
            }

            fn width(&self) -> T {
                self.end - self.start
            }

            fn is_empty(&self) -> bool {
                self.start >= self.end
            }
        }

        impl<T: PrimInt> BoundedRange<T> for RangeInclusive<T> {
            fn start(&self) -> T {
                *self.start()
            }

            fn width(&self) -> T {
                *self.end() - *self.start() + T::one()
            }

            fn is_empty(&self) -> bool {
                self.start() > self.end()
            }
        }

        #[cfg(test)]
        mod tests {
            use super::*;
            use std::collections::HashSet;

            #[test]
            fn splitmix64_known_sequence() {
                let mut rng = Random::new(0);
                assert_eq!(rng.next_u64(), 0xe220a8397b1dcdaf);
                assert_eq!(rng.next_u64(), 0x6e789e6aa1b965f4);
            }

            #[test]
            fn index_stays_in_range() {
                let mut rng = Random::new(42);
                for upper in [1, 2, 10, 2500] {
                    for _ in 0..10000 {
                        assert!(rng.index(upper) < upper);
                    }
                }
            }

            fn check_range<R: BoundedRange<usize> + Clone, F: Fn(&mut Random, R) -> usize>(
                range: R,
                generator: F,
            ) {
                let mut rng = Random::new(42);
                let mut values = HashSet::new();
                let start = range.start();
                let width = range.width();
                let end = start + width;

                for _ in 0..10000 {
                    let v = generator(&mut rng, range.clone());
                    assert!(start <= v && v < end, "value {v:?} out of range");
                    values.insert(v);
                }

                for v in start..end {
                    assert!(values.contains(&v), "value {v:?} not generated");
                }
            }

            #[test]
            fn test_fast_gen_range_u64x1() {
                check_range(0usize..10usize, |rng, r| rng.fast_gen_range_u64x1(r));
            }

            #[test]
            #[should_panic(expected = "cannot sample empty range")]
            fn test_fast_gen_range_u64x1_empty_range() {
                let mut rng = Random::new(42);
                let _ = rng.fast_gen_range_u64x1(5usize..5usize);
            }

            #[test]
            fn test_fast_gen_range_u32x1() {
                check_range(0usize..10usize, |rng, r| rng.fast_gen_range_u32x1(r));
            }

            #[test]
            fn test_fast_gen_range_u32x2() {
                let mut rng = Random::new(123);
                let mut set0 = HashSet::new();
                let mut set1 = HashSet::new();

                for _ in 0..10000 {
                    let (v0, v1) = rng.fast_gen_range_u32x2(0u32..5u32, 10u32..=15u32);
                    assert!((0..5).contains(&v0));
                    assert!((10..=15).contains(&v1));
                    set0.insert(v0);
                    set1.insert(v1);
                }

                for v in 0u32..5u32 {
                    assert!(set0.contains(&v));
                }

                for v in 10u32..=15u32 {
                    assert!(set1.contains(&v));
                }
            }

            #[test]
            fn test_fast_gen_range_u16x1() {
                check_range(0usize..10usize, |rng, r| rng.fast_gen_range_u16x1(r));
            }

            #[test]
            fn test_fast_gen_range_u16x2() {
                let mut rng = Random::new(456);
                let mut set0 = HashSet::new();
                let mut set1 = HashSet::new();

                for _ in 0..10000 {
                    let (v0, v1) = rng.fast_gen_range_u16x2(0u16..5u16, 10u16..15u16);
                    assert!((0..5).contains(&v0));
                    assert!((10..15).contains(&v1));
                    set0.insert(v0);
                    set1.insert(v1);
                }

                for v in 0u16..5u16 {
                    assert!(set0.contains(&v));
                }

                for v in 10u16..15u16 {
                    assert!(set1.contains(&v));
                }
            }

            #[test]
            fn test_fast_gen_range_u16x3() {
                let mut rng = Random::new(789);
                let mut set0 = HashSet::new();
                let mut set1 = HashSet::new();
                let mut set2 = HashSet::new();

                for _ in 0..10000 {
                    let (v0, v1, v2) =
                        rng.fast_gen_range_u16x3(0u16..3u16, 10u16..13u16, 20u16..=23u16);
                    assert!((0..3).contains(&v0));
                    assert!((10..13).contains(&v1));
                    assert!((20..=23).contains(&v2));
                    set0.insert(v0);
                    set1.insert(v1);
                    set2.insert(v2);
                }

                for v in 0u16..3u16 {
                    assert!(set0.contains(&v));
                }

                for v in 10u16..13u16 {
                    assert!(set1.contains(&v));
                }

                for v in 20u16..=23u16 {
                    assert!(set2.contains(&v));
                }
            }

            #[test]
            fn test_fast_gen_range_u16x4() {
                let mut rng = Random::new(321);
                let mut set0 = HashSet::new();
                let mut set1 = HashSet::new();
                let mut set2 = HashSet::new();
                let mut set3 = HashSet::new();

                for _ in 0..10000 {
                    let (v0, v1, v2, v3) = rng.fast_gen_range_u16x4(
                        0u16..2u16,
                        10u16..12u16,
                        20u16..=22u16,
                        30u16..32u16,
                    );
                    assert!((0..2).contains(&v0));
                    assert!((10..12).contains(&v1));
                    assert!((20..=22).contains(&v2));
                    assert!((30..32).contains(&v3));
                    set0.insert(v0);
                    set1.insert(v1);
                    set2.insert(v2);
                    set3.insert(v3);
                }

                for v in 0u16..2u16 {
                    assert!(set0.contains(&v));
                }

                for v in 10u16..12u16 {
                    assert!(set1.contains(&v));
                }

                for v in 20u16..=22u16 {
                    assert!(set2.contains(&v));
                }

                for v in 30u16..32u16 {
                    assert!(set3.contains(&v));
                }
            }
        }
    }
    pub mod util {
        //! よく使われるユーティリティ関数をまとめたモジュール

        use std::{
            fmt::Display,
            io::{self, BufWriter, Write as _},
        };

        use num::PrimInt;

        /// 最小値と最大値を更新するトレイト
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::util::ChangeMinMax;
        ///
        /// let mut x = 10;
        /// assert!(x.change_min(3));
        /// assert_eq!(x, 3);
        /// ```
        pub trait ChangeMinMax {
            fn change_min(&mut self, v: Self) -> bool;
            fn change_max(&mut self, v: Self) -> bool;
        }

        impl<T: PartialOrd> ChangeMinMax for T {
            fn change_min(&mut self, v: T) -> bool {
                *self > v && {
                    *self = v;
                    true
                }
            }

            fn change_max(&mut self, v: T) -> bool {
                *self < v && {
                    *self = v;
                    true
                }
            }
        }

        /// 条件に従ってYes/Noを出力する
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::yesno;
        ///
        /// let n = 3;
        /// yesno!(3 % 2 == 0)
        /// ```
        #[macro_export]
        macro_rules! yesno {
            ($p:expr) => {
                if $p {
                    println!("Yes");
                } else {
                    println!("No");
                }
            };
        }

        /// 標準出力・標準エラー出力に出力するトレイト
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::util::PrintLine as _;
        ///
        /// let x = 5;
        /// x.println();
        /// x.eprintln();
        ///
        /// let y = [1, 2, 3];
        /// y.println();
        /// y.eprintln();
        /// ```
        pub trait PrintLine {
            fn println(&self);
            fn eprintln(&self);
        }

        /// 単体値版
        impl<T: Display> PrintLine for T {
            fn println(&self) {
                println!("{self}");
            }

            fn eprintln(&self) {
                eprintln!("{self}");
            }
        }

        /// スライス版
        impl<T: Display> PrintLine for [T] {
            fn println(&self) {
                let stdout = io::stdout();
                let mut out = BufWriter::new(stdout.lock());

                let mut first = true;
                for x in self {
                    if !first {
                        out.write_all(b" ").unwrap();
                    }
                    write!(out, "{x}").unwrap();
                    first = false;
                }

                out.write_all(b"\n").unwrap();
                // drop(out)でflush
            }

            fn eprintln(&self) {
                let stderr = io::stderr();
                let mut out = BufWriter::new(stderr.lock());

                let mut first = true;
                for x in self {
                    if !first {
                        out.write_all(b" ").unwrap();
                    }
                    write!(out, "{x}").unwrap();
                    first = false;
                }

                out.write_all(b"\n").unwrap();
            }
        }

        /// 多次元配列を作成する
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::mat;
        ///
        /// let a = mat![0; 4; 3];
        /// assert_eq!(a, vec![vec![0; 3]; 4]);
        /// ```
        #[macro_export]
        macro_rules! mat {
	($($e:expr),*) => { vec![$($e),*] };
	($($e:expr,)*) => { vec![$($e),*] };
	($e:expr; $d:expr) => { vec![$e; $d] };
	($e:expr; $d:expr $(; $ds:expr)+) => { vec![mat![$e $(; $ds)*]; $d] };
}

        /// 整数の二分探索を行う
        ///
        /// # Examples
        ///
        /// ```
        /// use cp_lib_rs::util::binary_search;
        ///
        /// let result = binary_search(0, 10, |x| x * x <= 5);
        /// assert_eq!(result, 2);
        /// ```
        pub fn binary_search<T: PrimInt>(ok: T, ng: T, f: impl Fn(T) -> bool) -> T {
            let mut ok = ok;
            let mut ng = ng;

            while ok.max(ng) - ok.min(ng) > T::one() {
                let mid = (ok + ng) / (T::one() << 1);

                if f(mid) {
                    ok = mid;
                } else {
                    ng = mid;
                }
            }

            ok
        }

        #[cfg(test)]
        mod test {
            use super::*;

            #[test]
            fn binary_search_test() {
                assert_eq!(binary_search(0, 10, |x| x * x <= 5), 2);
                assert_eq!(binary_search(10, 0, |x| x * x >= 5), 3);
            }
        }
    }
}
pub mod problem {
    use std::fmt::{self, Display, Formatter};
    use std::io::{self, BufRead};

    use proconio::{input, marker::Bytes, source::once::OnceSource};

    use crate::library::grid::{Coord, Map2d};

    #[derive(Debug, Clone)]
    pub struct Input {
        /// 各マスの乾燥度。0は木のないマス、1〜9は木の乾燥度を表す。
        pub d: Map2d<u8>,
    }

    impl Input {
        /// 森の一辺のマス数。
        pub const N: usize = 50;
        /// 伐採できる木の本数の上限。
        pub const K: usize = 100;

        pub fn read() -> Self {
            Self::read_from(io::stdin().lock())
        }

        fn read_from(reader: impl BufRead) -> Self {
            let mut source = OnceSource::new(reader);
            input! {
                from &mut source,
                n: usize,
                k: usize,
            }
            assert_eq!(n, Self::N, "入力のNが定数と一致しません");
            assert_eq!(k, Self::K, "入力のKが定数と一致しません");
            let mut cells = Vec::with_capacity(Self::N * Self::N);
            for _ in 0..Self::N {
                input! { from &mut source, row: Bytes }
                assert_eq!(row.len(), Self::N, "乾燥度の行の長さがNと一致しません");
                for digit in row {
                    assert!(digit.is_ascii_digit(), "乾燥度は0〜9の数字で指定します");
                    cells.push(digit - b'0');
                }
            }
            Self {
                d: Map2d::new(cells, Self::N),
            }
        }
    }

    #[derive(Debug, Clone, Default)]
    pub struct Output {
        /// 伐採する木のマス（0始まりの行・列）。相異なる木のあるマスをK個以下指定する。
        pub cuts: Vec<Coord>,
    }

    impl Display for Output {
        fn fmt(&self, f: &mut Formatter<'_>) -> fmt::Result {
            // 本数は座標列から求め、出力行数との不一致を防ぐ。
            writeln!(f, "{}", self.cuts.len())?;
            for coord in &self.cuts {
                writeln!(f, "{} {}", coord.row(), coord.col())?;
            }
            Ok(())
        }
    }

    #[cfg(test)]
    mod tests {
        use super::*;

        #[test]
        fn reads_digits_in_row_major_order_including_leading_zero() {
            let mut text = String::from("50 100\n");
            for row in 0..50 {
                for col in 0..50 {
                    text.push(char::from(b'0' + ((row + col) % 10) as u8));
                }
                text.push('\n');
            }
            let input = Input::read_from(text.as_bytes());
            assert_eq!(input.d.iter().count(), 2500);
            for row in 0..50 {
                for col in 0..50 {
                    assert_eq!(input.d[Coord::new(row, col)], ((row + col) % 10) as u8);
                }
            }
        }

        #[test]
        fn formats_sample_output() {
            let output = Output {
                cuts: vec![Coord::new(0, 32), Coord::new(1, 33), Coord::new(2, 34)],
            };
            assert_eq!(output.to_string(), "3\n0 32\n1 33\n2 34\n");
        }

        #[test]
        fn formats_empty_output() {
            assert_eq!(Output::default().to_string(), "0\n");
        }
    }
}
pub mod solver {
    mod vertex_cut {
        mod evaluation {
            use crate::library::{
                grid::{ADJACENTS, Coord, Map2d},
                random::Random,
            };
            use crate::problem::Input;

            const SAMPLES: usize = 100;

            pub(super) struct Evaluator {
                samples: Vec<Map2d<bool>>,
                baseline: u64,
            }

            fn sample_seed(input: &Input) -> u64 {
                // 配布ツールと同じ入力由来のseedで、先頭100サンプルを固定する。
                let mut hash = 1469598103934665603_u64;
                for &d in input.d.iter() {
                    hash ^= u64::from(b'0' + d);
                    hash = hash.wrapping_mul(1099511628211);
                }
                hash ^ (Input::N as u64)
                    .wrapping_mul(1000003)
                    .wrapping_add((Input::K as u64).wrapping_mul(998244353))
            }

            fn component_square_sum(dry: &Map2d<bool>, cut: &Map2d<bool>) -> u64 {
                let mut seen = Map2d::with_default(Input::N);
                let mut queue = Vec::with_capacity(Input::N * Input::N);
                let mut sum = 0;
                for row in 0..Input::N {
                    for col in 0..Input::N {
                        let start = Coord::new(row, col);
                        if seen[start] || !dry[start] || cut[start] {
                            continue;
                        }
                        seen[start] = true;
                        queue.clear();
                        queue.push(start);
                        let mut head = 0;
                        while head < queue.len() {
                            let p = queue[head];
                            head += 1;
                            for diff in ADJACENTS {
                                let next = p + diff;
                                if next.in_map(Input::N) && !seen[next] && dry[next] && !cut[next] {
                                    seen[next] = true;
                                    queue.push(next);
                                }
                            }
                        }
                        sum += (queue.len() as u64).pow(2);
                    }
                }
                sum
            }

            impl Evaluator {
                pub(super) fn new(input: &Input) -> Self {
                    let mut rng = Random::new(sample_seed(input));
                    let thresholds: Vec<u64> = input
                        .d
                        .iter()
                        .map(|&d| u64::from(d) * (1_u64 << 32) / 10)
                        .collect();
                    let samples: Vec<_> = (0..SAMPLES)
                        .map(|_| {
                            let mut cells = Vec::with_capacity(Input::N * Input::N);
                            for pair in thresholds.chunks(2) {
                                let bits = rng.next_u64();
                                cells.push((bits & 0xffff_ffff) < pair[0]);
                                if pair.len() == 2 {
                                    cells.push((bits >> 32) < pair[1]);
                                }
                            }
                            Map2d::new(cells, Input::N)
                        })
                        .collect();
                    let no_cuts = Map2d::with_default(Input::N);
                    let baseline = samples
                        .iter()
                        .map(|sample| component_square_sum(sample, &no_cuts))
                        .sum();
                    Self { samples, baseline }
                }

                pub(super) fn score(&self, cut: &Map2d<bool>) -> i64 {
                    let total: u64 = self
                        .samples
                        .iter()
                        .map(|sample| component_square_sum(sample, cut))
                        .sum();
                    let denominator = (SAMPLES * Input::N * Input::N) as f64;
                    let e0 = self.baseline as f64 / denominator;
                    let e = total as f64 / denominator;
                    // ゼロ近傍の処理・丸めを配布ツールと揃える。
                    (1e6 * (1.0 + (e0 / e.max(1e-12)).log2())).round() as i64
                }
            }

            #[cfg(test)]
            mod tests {
                use super::*;

                #[test]
                fn square_sum_and_cut_splitting() {
                    let mut dry = Map2d::with_default(Input::N);
                    for col in 0..3 {
                        dry[Coord::new(0, col)] = true;
                    }
                    dry[Coord::new(2, 2)] = true;
                    let mut cut = Map2d::with_default(Input::N);
                    assert_eq!(component_square_sum(&dry, &cut), 10);
                    cut[Coord::new(0, 1)] = true;
                    assert_eq!(component_square_sum(&dry, &cut), 3);
                }

                #[test]
                fn no_cut_score_is_one_million_and_samples_are_fixed() {
                    let input = Input {
                        d: Map2d::from_fn(|p| ((p.row() + p.col()) % 10) as u8, Input::N),
                    };
                    let first = Evaluator::new(&input);
                    let second = Evaluator::new(&input);
                    assert_eq!(first.score(&Map2d::with_default(Input::N)), 1_000_000);
                    for (a, b) in first.samples.iter().zip(&second.samples) {
                        assert!(a.iter().eq(b.iter()));
                        for row in 0..Input::N {
                            for col in 0..Input::N {
                                let p = Coord::new(row, col);
                                if input.d[p] == 0 {
                                    assert!(!a[p]);
                                }
                            }
                        }
                    }
                }
            }
        }

        use std::time::Duration;

        use ac_library::MfGraph;

        use crate::library::{
            diagnostics::Timer,
            grid::{ADJACENTS, Coord, Map2d},
            random::Random,
        };
        use crate::problem::{Input, Output};

        use evaluation::Evaluator;

        const D_THRESHOLD: u8 = 5;
        const PAIR_TRIALS: usize = 1;
        const MIN_DISTANCE: usize = 15;
        const TERMINAL_RADIUS: usize = 7;
        const TARGET_CUTS: usize = 95;
        const TIME_LIMIT: Duration = Duration::from_millis(1800);

        fn components(input: &Input, cut: &Map2d<bool>) -> Vec<Vec<Coord>> {
            let mut seen = Map2d::with_default(Input::N);
            let mut result = Vec::new();
            for row in 0..Input::N {
                for col in 0..Input::N {
                    let start = Coord::new(row, col);
                    if seen[start] || cut[start] || input.d[start] < D_THRESHOLD {
                        continue;
                    }
                    seen[start] = true;
                    let mut vertices = vec![start];
                    let mut head = 0;
                    while head < vertices.len() {
                        let p = vertices[head];
                        head += 1;
                        for diff in ADJACENTS {
                            let next = p + diff;
                            if next.in_map(Input::N)
                                && !seen[next]
                                && !cut[next]
                                && input.d[next] >= D_THRESHOLD
                            {
                                seen[next] = true;
                                vertices.push(next);
                            }
                        }
                    }
                    result.push(vertices);
                }
            }
            result
        }

        fn choose_component<'a>(
            components: &'a [Vec<Coord>],
            rng: &mut Random,
        ) -> Option<&'a [Coord]> {
            let total: usize = components.iter().map(Vec::len).sum();
            if total == 0 {
                return None;
            }
            // 全頂点から1個を選ぶことで、成分サイズに比例した確率になる。
            let mut index = rng.index(total);
            for component in components {
                if index < component.len() {
                    return Some(component);
                }
                index -= component.len();
            }
            unreachable!()
        }

        fn choose_pair(vertices: &[Coord], rng: &mut Random) -> (Coord, Coord) {
            let mut best = (vertices[0], vertices[0]);
            for _ in 0..PAIR_TRIALS {
                let pair = (
                    vertices[rng.index(vertices.len())],
                    vertices[rng.index(vertices.len())],
                );
                if pair.0.dist(&pair.1) > best.0.dist(&best.1) {
                    best = pair;
                }
            }
            best
        }

        fn minimum_vertex_cut(
            input: &Input,
            vertices: &[Coord],
            p1: Coord,
            p2: Coord,
        ) -> Vec<Coord> {
            let len = vertices.len();
            let mut indices = Map2d::from_fn(|_| usize::MAX, Input::N);
            for (i, &p) in vertices.iter().enumerate() {
                indices[p] = i;
            }
            let source = 2 * len;
            let sink = source + 1;
            let mut graph = MfGraph::<i64>::new(sink + 1);
            // 全頂点の有限コストの総和より大きく、隣接辺はカットされない。
            let infinity = (len as i64) * 10009 + 1;
            for (i, &p) in vertices.iter().enumerate() {
                let vin = 2 * i;
                let vout = vin + 1;
                graph.add_edge(vin, vout, 10009 - i64::from(input.d[p]));
                // source/sink内の頂点も有限コストのvin→voutを通るので伐採可能。
                if p.dist(&p1) <= TERMINAL_RADIUS {
                    graph.add_edge(source, vin, infinity);
                }
                if p.dist(&p2) <= TERMINAL_RADIUS {
                    graph.add_edge(vout, sink, infinity);
                }
                for diff in ADJACENTS {
                    let next = p + diff;
                    if next.in_map(Input::N) && indices[next] != usize::MAX {
                        graph.add_edge(vout, 2 * indices[next], infinity);
                    }
                }
            }
            graph.flow(source, sink);
            let reachable = graph.min_cut(source);
            vertices
                .iter()
                .enumerate()
                .filter_map(|(i, &p)| (reachable[2 * i] && !reachable[2 * i + 1]).then_some(p))
                .collect()
        }

        fn trial(input: &Input, rng: &mut Random) -> (Output, Map2d<bool>) {
            let mut output = Output::default();
            let mut cut = Map2d::with_default(Input::N);
            let mut current = components(input, &cut);
            // 指定仕様により、95本に達するまで再抽選する。時間では中断しない。
            while output.cuts.len() < TARGET_CUTS {
                let Some(vertices) = choose_component(&current, rng) else {
                    continue;
                };
                let (p1, p2) = choose_pair(vertices, rng);
                if p1.dist(&p2) <= MIN_DISTANCE {
                    continue;
                }
                let added = minimum_vertex_cut(input, vertices, p1, p2);
                if output.cuts.len() + added.len() > Input::K {
                    continue;
                }
                for p in added {
                    cut[p] = true;
                    output.cuts.push(p);
                }
                current = components(input, &cut);
            }
            (output, cut)
        }

        pub fn solve(input: &Input, timer: &Timer) -> Output {
            let evaluator = Evaluator::new(input);
            let mut rng = Random::new(42);
            let mut best = Output::default();
            let mut best_score = i64::MIN;
            let mut trials = 0;
            // 最初の試行は必ず完了させ、以後は試行開始時のみ時間を確認する。
            while trials == 0 || timer.elapsed() <= TIME_LIMIT {
                let (output, cut) = trial(input, &mut rng);
                let score = evaluator.score(&cut);
                trials += 1;
                if score > best_score {
                    best_score = score;
                    best = output;
                }
            }
            eprintln!(
                "trials={trials} estimated_score={best_score} cuts={} elapsed_ms={:.3}",
                best.cuts.len(),
                timer.elapsed().as_secs_f64() * 1000.0
            );
            best
        }

        #[cfg(test)]
        mod tests {
            use super::*;

            #[test]
            fn components_exclude_low_dryness_and_cut_cells() {
                let mut input = Input {
                    d: Map2d::with_default(Input::N),
                };
                for (col, d) in [5, 9, 4, 5, 5, 9].into_iter().enumerate() {
                    input.d[Coord::new(0, col)] = d;
                }
                let mut cut = Map2d::with_default(Input::N);
                cut[Coord::new(0, 4)] = true;
                let actual = components(&input, &cut);
                assert_eq!(actual.iter().map(Vec::len).collect::<Vec<_>>(), [2, 1, 1]);
                assert_eq!(actual[0], [Coord::new(0, 0), Coord::new(0, 1)]);
            }

            #[test]
            fn minimum_cut_may_remove_a_source_vertex() {
                let mut input = Input {
                    d: Map2d::from_fn(|_| 5, Input::N),
                };
                let vertices: Vec<_> = (0..12).map(|col| Coord::new(0, col)).collect();
                // 選択成分外にも木があるが、フローには含めない。
                // sourceの右端を最安にし、その頂点自身を切れることも確認する。
                input.d[Coord::new(0, 4)] = 9;
                assert_eq!(
                    minimum_vertex_cut(&input, &vertices, vertices[0], vertices[11]),
                    [Coord::new(0, 4)]
                );
            }

            #[test]
            fn minimum_cut_cost_matches_exhaustive_search() {
                let mut input = Input {
                    d: Map2d::with_default(Input::N),
                };
                let vertices: Vec<_> = (0..12).map(|i| Coord::new(i / 6, i % 6)).collect();
                for (i, &p) in vertices.iter().enumerate() {
                    input.d[p] = 5 + (i % 5) as u8;
                }
                let p1 = vertices[0];
                let p2 = vertices[11];
                let mut optimum = i64::MAX;
                for mask in 0..(1 << vertices.len()) {
                    let cost: i64 = vertices
                        .iter()
                        .enumerate()
                        .filter(|(i, _)| mask & (1 << i) != 0)
                        .map(|(_, &p)| 10009 - i64::from(input.d[p]))
                        .sum();
                    if cost >= optimum {
                        continue;
                    }
                    let mut seen = vec![false; vertices.len()];
                    let mut stack = Vec::new();
                    for (i, &p) in vertices.iter().enumerate() {
                        if mask & (1 << i) == 0 && p.dist(&p1) <= TERMINAL_RADIUS {
                            seen[i] = true;
                            stack.push(i);
                        }
                    }
                    while let Some(i) = stack.pop() {
                        for (j, &p) in vertices.iter().enumerate() {
                            if mask & (1 << j) == 0 && !seen[j] && vertices[i].dist(&p) == 1 {
                                seen[j] = true;
                                stack.push(j);
                            }
                        }
                    }
                    if !vertices
                        .iter()
                        .enumerate()
                        .any(|(i, &p)| seen[i] && p.dist(&p2) <= TERMINAL_RADIUS)
                    {
                        optimum = cost;
                    }
                }
                let actual = minimum_vertex_cut(&input, &vertices, p1, p2);
                let cost: i64 = actual.iter().map(|&p| 10009 - i64::from(input.d[p])).sum();
                assert_eq!(cost, optimum);
            }
        }
    }

    pub use vertex_cut::solve;
}

fn main() {
    let timer = library::diagnostics::Timer::default();
    let input = problem::Input::read();
    let output = solver::solve(&input, &timer);
    print!("{output}");
}