routik-solver 0.1.0

//! Local-search representation and neighbourhood moves.
//!
//! The public [`Solution`] is convenient to hand back but awkward to mutate
//! millions of times, so the search runs on an internal [`State`]: routes as
//! sequences of stop **location indices** (`1..=n`, depot implicit at both
//! ends), with cached per-route load/distance/duration so the rolled-up totals
//! stay correct without re-summing the whole solution.
//!
//! Four neighbourhoods are offered, each intra- **and** inter-route:
//! *relocate* (move one stop), *swap* (exchange two stops), *2-opt* (reverse a
//! segment, intra only), and *or-opt* (move a 2–3 stop segment, optionally
//! reversed). Every move reports its **distance delta in O(1)** from only the
//! edges it touches (the solver rule); feasibility (capacity + time windows) is
//! re-checked with [`walk_route`] over the one or two affected routes, and is
//! only paid when a move is actually applied.

use crate::feasibility::{build_solution, walk_route};
use crate::matrix::CostMatrix;
use crate::model::{LocationId, Problem, Solution};
use crate::objective::Objective;
use rand::Rng;
use std::collections::HashMap;

/// A candidate neighbourhood move. Positions index the *current* routes.
#[derive(Debug, Clone, Copy)]
pub(crate) enum Move {
    /// Move the stop at `(from_r, from_pos)` to slot `to_slot` of `to_r`.
    Relocate {
        from_r: usize,
        from_pos: usize,
        to_r: usize,
        to_slot: usize,
    },
    /// Exchange the stops at `(ra, pa)` and `(rb, pb)`.
    Swap {
        ra: usize,
        pa: usize,
        rb: usize,
        pb: usize,
    },
    /// Reverse `route[i..=j]` within route `r` (intra-route 2-opt).
    TwoOpt { r: usize, i: usize, j: usize },
    /// Move the `len`-stop segment at `(from_r, from_pos)` to slot `to_slot` of
    /// `to_r`, reversed iff `reversed`.
    OrOpt {
        from_r: usize,
        from_pos: usize,
        len: usize,
        to_r: usize,
        to_slot: usize,
        reversed: bool,
    },
}

/// Location at signed slot `k` of a route, with the depot (`0`) at both ends.
fn loc_at(route: &[usize], k: isize) -> usize {
    if k < 0 {
        return 0;
    }
    let k = k as usize;
    if k >= route.len() { 0 } else { route[k] }
}

/// Mutable search state over a single [`Problem`]/[`CostMatrix`] pairing.
pub(crate) struct State<'a, M: CostMatrix> {
    problem: &'a Problem,
    matrix: &'a M,
    capacity: u32,
    /// One slot per vehicle; empty slots are spare (unused) vehicles.
    routes: Vec<Vec<usize>>,
    loads: Vec<u32>,
    route_dist: Vec<f64>,
    route_dur: Vec<f64>,
    /// Route index currently holding each stop (`stop_route[loc - 1]`).
    stop_route: Vec<usize>,
    total_distance: f64,
    total_duration: f64,
    n: usize,
}

impl<'a, M: CostMatrix> State<'a, M> {
    /// Build a state from an existing (feasible) solution, padding the route
    /// list up to the fleet size so moves can open spare vehicles.
    pub(crate) fn from_solution(
        problem: &'a Problem,
        matrix: &'a M,
        capacity: u32,
        solution: &Solution,
    ) -> Self {
        let n = problem.stops.len();
        let loc_of: HashMap<_, _> = problem
            .stops
            .iter()
            .enumerate()
            .map(|(i, s)| (s.id, i + 1))
            .collect();

        // A few spare (empty) route slots let the search open new vehicles, but
        // padding to the *whole* fleet would make most random insertion targets
        // empty and waste iterations — so cap the spares to a small margin.
        const SPARE_ROUTES: usize = 3;
        let used = solution.routes.len();
        let fleet = (used + SPARE_ROUTES).min(problem.vehicles.len().max(used));

        let mut routes: Vec<Vec<usize>> = Vec::with_capacity(fleet);
        for route in &solution.routes {
            let seq: Vec<usize> = route
                .stop_ids
                .iter()
                .filter_map(|sid| loc_of.get(sid).copied())
                .collect();
            routes.push(seq);
        }
        routes.resize_with(fleet, Vec::new);

        let mut state = Self {
            problem,
            matrix,
            capacity,
            routes,
            loads: vec![0; fleet],
            route_dist: vec![0.0; fleet],
            route_dur: vec![0.0; fleet],
            stop_route: vec![0; n],
            total_distance: 0.0,
            total_duration: 0.0,
            n,
        };
        for r in 0..fleet {
            state.recompute_route(r);
        }
        state.total_distance = state.route_dist.iter().sum();
        state.total_duration = state.route_dur.iter().sum();
        state
    }

    /// Refresh the cached metrics and stop→route map for one route.
    fn recompute_route(&mut self, r: usize) {
        let (metrics, _) = walk_route(self.problem, self.matrix, self.capacity, &self.routes[r]);
        self.loads[r] = metrics.load;
        self.route_dist[r] = metrics.distance;
        self.route_dur[r] = metrics.duration;
        for &loc in &self.routes[r] {
            self.stop_route[loc - 1] = r;
        }
    }

    fn d(&self, u: usize, v: usize) -> f64 {
        self.matrix.distance(LocationId(u), LocationId(v))
    }

    #[cfg(test)]
    pub(crate) fn total_distance(&self) -> f64 {
        self.total_distance
    }

    fn active_routes(&self) -> usize {
        self.routes.iter().filter(|r| !r.is_empty()).count()
    }

    /// Current objective value from cached totals.
    pub(crate) fn cost(&self, objective: &Objective) -> f64 {
        objective.score(
            self.total_distance,
            self.active_routes(),
            self.total_duration,
        )
    }

    /// Convert back to a public solution (drops empty routes).
    pub(crate) fn to_solution(&self) -> Solution {
        build_solution(self.problem, self.matrix, self.capacity, &self.routes)
    }

    /// Snapshot just the route sequences (for keeping the best solution).
    pub(crate) fn snapshot(&self) -> Vec<Vec<usize>> {
        self.routes.clone()
    }

    /// Restore route sequences from a snapshot and refresh all caches.
    pub(crate) fn restore(&mut self, routes: Vec<Vec<usize>>) {
        self.routes = routes;
        let fleet = self.routes.len();
        self.loads.resize(fleet, 0);
        self.route_dist.resize(fleet, 0.0);
        self.route_dur.resize(fleet, 0.0);
        for r in 0..fleet {
            self.recompute_route(r);
        }
        self.total_distance = self.route_dist.iter().sum();
        self.total_duration = self.route_dur.iter().sum();
    }

    // --- O(1) distance delta per move -------------------------------------

    /// Change in total distance if `mv` were applied, from touched edges only.
    pub(crate) fn distance_delta(&self, mv: &Move) -> f64 {
        match *mv {
            Move::Relocate {
                from_r,
                from_pos,
                to_r,
                to_slot,
            } => {
                let from = &self.routes[from_r];
                let s = from[from_pos];
                let a = loc_at(from, from_pos as isize - 1);
                let b = loc_at(from, from_pos as isize + 1);
                let remove = self.d(a, b) - self.d(a, s) - self.d(s, b);

                let to = &self.routes[to_r];
                let c = loc_at(to, to_slot as isize - 1);
                let e = loc_at(to, to_slot as isize);
                let insert = self.d(c, s) + self.d(s, e) - self.d(c, e);
                remove + insert
            }
            Move::Swap { ra, pa, rb, pb } => {
                let s1 = self.routes[ra][pa];
                let s2 = self.routes[rb][pb];
                let ap = loc_at(&self.routes[ra], pa as isize - 1);
                let an = loc_at(&self.routes[ra], pa as isize + 1);
                let bp = loc_at(&self.routes[rb], pb as isize - 1);
                let bn = loc_at(&self.routes[rb], pb as isize + 1);
                let removed = self.d(ap, s1) + self.d(s1, an) + self.d(bp, s2) + self.d(s2, bn);
                let added = self.d(ap, s2) + self.d(s2, an) + self.d(bp, s1) + self.d(s1, bn);
                added - removed
            }
            Move::TwoOpt { r, i, j } => {
                let route = &self.routes[r];
                let pre = loc_at(route, i as isize - 1);
                let si = route[i];
                let sj = route[j];
                let post = loc_at(route, j as isize + 1);
                self.d(pre, sj) + self.d(si, post) - self.d(pre, si) - self.d(sj, post)
            }
            Move::OrOpt {
                from_r,
                from_pos,
                len,
                to_r,
                to_slot,
                reversed,
            } => {
                let from = &self.routes[from_r];
                let first = from[from_pos];
                let last = from[from_pos + len - 1];
                let a = loc_at(from, from_pos as isize - 1);
                let b = loc_at(from, (from_pos + len) as isize);
                let remove = self.d(a, b) - self.d(a, first) - self.d(last, b);

                let to = &self.routes[to_r];
                let c = loc_at(to, to_slot as isize - 1);
                let e = loc_at(to, to_slot as isize);
                let insert = if reversed {
                    self.d(c, last) + self.d(first, e) - self.d(c, e)
                } else {
                    self.d(c, first) + self.d(last, e) - self.d(c, e)
                };
                remove + insert
            }
        }
    }

    /// Change in the number of used vehicles (empties a source / fills a
    /// spare). O(1).
    fn route_count_delta(&self, mv: &Move) -> i32 {
        let (from_r, removed, to_r) = match *mv {
            Move::Relocate { from_r, to_r, .. } => (from_r, 1usize, to_r),
            Move::OrOpt {
                from_r, len, to_r, ..
            } => (from_r, len, to_r),
            Move::Swap { .. } | Move::TwoOpt { .. } => return 0,
        };
        if from_r == to_r {
            return 0;
        }
        let mut delta = 0;
        if self.routes[from_r].len() == removed {
            delta -= 1; // source emptied
        }
        if self.routes[to_r].is_empty() {
            delta += 1; // spare opened
        }
        delta
    }

    /// Objective-weighted delta. Distance and vehicle terms are O(1); the time
    /// term (only when weighted) costs an O(len) re-walk of the touched routes.
    pub(crate) fn cost_delta(&self, mv: &Move, objective: &Objective) -> f64 {
        let dist = self.distance_delta(mv);
        let veh = f64::from(self.route_count_delta(mv));
        let time = if objective.time == 0.0 {
            0.0
        } else {
            self.duration_delta(mv)
        };
        objective.distance * dist + objective.vehicles * veh + objective.time * time
    }

    /// Duration delta via materialising the touched routes (used only when the
    /// objective weights time).
    fn duration_delta(&self, mv: &Move) -> f64 {
        let touched = self.materialize(mv);
        let mut delta = 0.0;
        for (r, seq) in &touched {
            let (m, _) = walk_route(self.problem, self.matrix, self.capacity, seq);
            delta += m.duration - self.route_dur[*r];
        }
        delta
    }

    // --- materialise + apply ----------------------------------------------

    /// Produce the new sequences for the routes a move touches, without
    /// mutating the state. One entry for intra-route moves, two for inter.
    fn materialize(&self, mv: &Move) -> Vec<(usize, Vec<usize>)> {
        match *mv {
            Move::Relocate {
                from_r,
                from_pos,
                to_r,
                to_slot,
            } => {
                if from_r == to_r {
                    let mut seq = self.routes[from_r].clone();
                    let s = seq.remove(from_pos);
                    let at = if to_slot > from_pos {
                        to_slot - 1
                    } else {
                        to_slot
                    };
                    seq.insert(at, s);
                    vec![(from_r, seq)]
                } else {
                    let mut from = self.routes[from_r].clone();
                    let s = from.remove(from_pos);
                    let mut to = self.routes[to_r].clone();
                    to.insert(to_slot, s);
                    vec![(from_r, from), (to_r, to)]
                }
            }
            Move::Swap { ra, pa, rb, pb } => {
                if ra == rb {
                    let mut seq = self.routes[ra].clone();
                    seq.swap(pa, pb);
                    vec![(ra, seq)]
                } else {
                    let mut a = self.routes[ra].clone();
                    let mut b = self.routes[rb].clone();
                    std::mem::swap(&mut a[pa], &mut b[pb]);
                    vec![(ra, a), (rb, b)]
                }
            }
            Move::TwoOpt { r, i, j } => {
                let mut seq = self.routes[r].clone();
                seq[i..=j].reverse();
                vec![(r, seq)]
            }
            Move::OrOpt {
                from_r,
                from_pos,
                len,
                to_r,
                to_slot,
                reversed,
            } => {
                let mut segment: Vec<usize> =
                    self.routes[from_r][from_pos..from_pos + len].to_vec();
                if reversed {
                    segment.reverse();
                }
                if from_r == to_r {
                    let mut seq = self.routes[from_r].clone();
                    seq.drain(from_pos..from_pos + len);
                    let at = if to_slot > from_pos {
                        to_slot - len
                    } else {
                        to_slot
                    };
                    seq.splice(at..at, segment);
                    vec![(from_r, seq)]
                } else {
                    let mut from = self.routes[from_r].clone();
                    from.drain(from_pos..from_pos + len);
                    let mut to = self.routes[to_r].clone();
                    to.splice(to_slot..to_slot, segment);
                    vec![(from_r, from), (to_r, to)]
                }
            }
        }
    }

    /// Apply `mv` iff the routes it touches stay feasible. Returns whether it
    /// was applied; on `false` the state is unchanged.
    pub(crate) fn try_apply(&mut self, mv: &Move) -> bool {
        let touched = self.materialize(mv);
        // All touched routes must be feasible before committing anything.
        for (_, seq) in &touched {
            if !walk_route(self.problem, self.matrix, self.capacity, seq).1 {
                return false;
            }
        }
        for (r, seq) in touched {
            self.routes[r] = seq;
            self.recompute_route(r);
        }
        self.total_distance = self.route_dist.iter().sum();
        self.total_duration = self.route_dur.iter().sum();
        true
    }

    // --- random move generation -------------------------------------------

    /// Sample a random candidate move, or `None` if the draw was degenerate
    /// (e.g. a no-op intra move); the caller treats `None` as a skipped step.
    pub(crate) fn random_move(&self, rng: &mut impl Rng) -> Option<Move> {
        if self.n == 0 {
            return None;
        }
        match rng.random_range(0..4) {
            0 => self.gen_relocate(rng),
            1 => self.gen_swap(rng),
            2 => self.gen_two_opt(rng),
            _ => self.gen_or_opt(rng),
        }
    }

    /// A uniformly random occupied `(route, pos)` slot.
    fn random_stop_slot(&self, rng: &mut impl Rng) -> Option<(usize, usize)> {
        let loc = rng.random_range(1..=self.n);
        let r = self.stop_route[loc - 1];
        let pos = self.routes[r].iter().position(|&x| x == loc)?;
        Some((r, pos))
    }

    fn gen_relocate(&self, rng: &mut impl Rng) -> Option<Move> {
        let (from_r, from_pos) = self.random_stop_slot(rng)?;
        let to_r = rng.random_range(0..self.routes.len());
        let to_slot = rng.random_range(0..=self.routes[to_r].len());
        if from_r == to_r && (to_slot == from_pos || to_slot == from_pos + 1) {
            return None; // no-op
        }
        Some(Move::Relocate {
            from_r,
            from_pos,
            to_r,
            to_slot,
        })
    }

    fn gen_swap(&self, rng: &mut impl Rng) -> Option<Move> {
        let (ra, pa) = self.random_stop_slot(rng)?;
        let (rb, pb) = self.random_stop_slot(rng)?;
        if ra == rb {
            // Need two distinct, non-adjacent positions; order them.
            let (lo, hi) = (pa.min(pb), pa.max(pb));
            if hi < lo + 2 {
                return None;
            }
            Some(Move::Swap {
                ra,
                pa: lo,
                rb,
                pb: hi,
            })
        } else {
            Some(Move::Swap { ra, pa, rb, pb })
        }
    }

    fn gen_two_opt(&self, rng: &mut impl Rng) -> Option<Move> {
        let (r, _) = self.random_stop_slot(rng)?;
        let len = self.routes[r].len();
        if len < 2 {
            return None;
        }
        let i = rng.random_range(0..len - 1);
        let j = rng.random_range(i + 1..len);
        Some(Move::TwoOpt { r, i, j })
    }

    fn gen_or_opt(&self, rng: &mut impl Rng) -> Option<Move> {
        let (from_r, _) = self.random_stop_slot(rng)?;
        let route_len = self.routes[from_r].len();
        if route_len < 2 {
            return None;
        }
        let len = rng.random_range(2..=3usize.min(route_len));
        let from_pos = rng.random_range(0..=route_len - len);
        let to_r = rng.random_range(0..self.routes.len());
        let to_len = self.routes[to_r].len();
        let to_slot = rng.random_range(0..=to_len);
        if from_r == to_r && to_slot >= from_pos && to_slot <= from_pos + len {
            return None; // insertion inside the removed segment
        }
        let reversed = rng.random_bool(0.5);
        Some(Move::OrOpt {
            from_r,
            from_pos,
            len,
            to_r,
            to_slot,
            reversed,
        })
    }
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::clarke_wright::clarke_wright;
    use crate::matrix::EuclideanMatrix;
    use crate::model::{Coord, Stop, StopId, TimeWindow, Vehicle, VehicleId};
    use rand::SeedableRng;
    use rand_chacha::ChaCha8Rng;

    /// A small VRPTW instance with loose windows so most moves stay feasible.
    fn problem() -> Problem {
        let stops = (1..=10u32)
            .map(|i| {
                let angle = f64::from(i) * 0.6;
                Stop {
                    id: StopId(i),
                    coord: Coord::new(10.0 * angle.sin(), 10.0 * angle.cos()),
                    demand: 5,
                    time_window: Some(TimeWindow {
                        start: 0.0,
                        end: 10_000.0,
                    }),
                    service_time: 1.0,
                }
            })
            .collect();
        let vehicles = (1..=5u32)
            .map(|i| Vehicle {
                id: VehicleId(i),
                capacity: 30,
            })
            .collect();
        Problem {
            depot: Coord::new(0.0, 0.0),
            stops,
            vehicles,
            depot_window: Some(TimeWindow {
                start: 0.0,
                end: 100_000.0,
            }),
        }
    }

    /// Full from-scratch distance of the current routes (independent of caches).
    fn recompute_distance(state: &State<EuclideanMatrix>) -> f64 {
        state
            .routes
            .iter()
            .map(|seq| {
                walk_route(state.problem, state.matrix, state.capacity, seq)
                    .0
                    .distance
            })
            .sum()
    }

    #[test]
    fn delta_matches_full_recompute_over_many_random_moves() {
        let p = problem();
        let m = EuclideanMatrix::from_problem(&p);
        let cw = clarke_wright(&p, &m).expect("feasible");
        let mut state = State::from_solution(&p, &m, 30, &cw);
        let mut rng = ChaCha8Rng::seed_from_u64(7);

        let mut applied = 0;
        for _ in 0..5_000 {
            let Some(mv) = state.random_move(&mut rng) else {
                continue;
            };
            let predicted = state.distance_delta(&mv);
            let before = state.total_distance();
            if state.try_apply(&mv) {
                applied += 1;
                let after = state.total_distance();
                assert!(
                    (after - before - predicted).abs() < 1e-6,
                    "delta mismatch: predicted {predicted}, actual {}",
                    after - before
                );
                // Cache must agree with an independent recompute.
                assert!((recompute_distance(&state) - after).abs() < 1e-6);
            }
        }
        assert!(applied > 50, "suspiciously few moves applied: {applied}");
    }

    #[test]
    fn search_state_roundtrips_and_preserves_stops() {
        let p = problem();
        let m = EuclideanMatrix::from_problem(&p);
        let cw = clarke_wright(&p, &m).expect("feasible");
        let state = State::from_solution(&p, &m, 30, &cw);
        let sol = state.to_solution();
        let served: usize = sol.routes.iter().map(|r| r.stop_ids.len()).sum();
        assert_eq!(served, p.stops.len());
        assert!(sol.feasible);
    }
}