problemreductions/solvers/ilp/
solver.rs

1//! ILP solver implementation using HiGHS.
2
3use crate::models::algebraic::{Comparison, ObjectiveSense, VariableDomain, ILP};
4use crate::solvers::registry::solver_capability_registry;
5use crate::solvers::ExactProblemKey;
6use crate::traits::Problem;
7use crate::types::{i64_to_exact_f64, MAX_EXACT_F64_INTEGER};
8use good_lp::highs;
9use good_lp::solvers::highs::HighsParallelType;
10use good_lp::{
11    variable, ProblemVariables, ResolutionError, Solution, SolutionStatus, SolverModel, Variable,
12};
13
14/// A failure to produce an ILP solution optimal within backend numerical tolerances.
15#[derive(Clone, Debug, PartialEq, Eq, thiserror::Error)]
16pub enum ILPSolveError {
17    /// The constraints have no feasible assignment.
18    #[error("the ILP is infeasible")]
19    Infeasible,
20    /// A target witness did not establish the source decision threshold.
21    #[error(
22        "the ILP witness does not meet the decision threshold for {0}; the decision is unresolved"
23    )]
24    UnresolvedDecision(String),
25    /// The objective is unbounded.
26    #[error("the ILP objective is unbounded")]
27    Unbounded,
28    /// The configured time limit was reached before optimality was proven.
29    #[error("the ILP solver reached its time limit before proving optimality")]
30    Timeout,
31    /// The selected backend failed for another reason.
32    #[error("the ILP backend failed: {0}")]
33    BackendFailure(String),
34    /// Type-erased dispatch received a value other than a supported ILP variant.
35    #[error("the ILP backend requires bool/i64 variables and f64 coefficients")]
36    UnsupportedProblemType,
37    /// No ILP pipeline is registered for the exact problem variant.
38    #[error("no ILP pipeline is registered for {0}")]
39    MissingPipeline(String),
40    /// The solver capability registry is invalid.
41    #[error("solver capability registry is invalid: {0}")]
42    InvalidRegistry(String),
43    /// A registered pipeline returned a solution for a different source type.
44    #[error("registered ILP pipeline returned the wrong solution type for {0}")]
45    PipelineTypeMismatch(String),
46    /// HiGHS reported an optimal solution that is invalid after integer rounding.
47    #[error("the ILP backend returned an invalid rounded solution: {0}")]
48    InvalidSolution(String),
49    /// An exact integer in the model cannot be transported through the f64 backend API.
50    #[error("the ILP backend cannot represent an exact model integer: {0}")]
51    InexactTransport(#[from] crate::types::ExactI64ToF64Error),
52    /// A target witness could not be mapped back to the source problem.
53    #[error(transparent)]
54    Extraction(#[from] crate::rules::ExtractionError),
55    /// A registered reduction could not construct its target instance.
56    #[error(transparent)]
57    Reduction(#[from] crate::rules::ReductionError),
58}
59
60fn classify_backend_error(error: ResolutionError, time_limit: Option<f64>) -> ILPSolveError {
61    match error {
62        ResolutionError::Infeasible => ILPSolveError::Infeasible,
63        ResolutionError::Unbounded => ILPSolveError::Unbounded,
64        ResolutionError::Other("NoSolutionFound") if time_limit.is_some() => ILPSolveError::Timeout,
65        other => ILPSolveError::BackendFailure(other.to_string()),
66    }
67}
68
69/// An ILP solver using the HiGHS backend.
70///
71/// Registered reductions map a source problem to an `ILP<V, f64>` terminal,
72/// which this solver sends to HiGHS before extracting the source solution.
73/// Optimality and infeasibility are assessed within HiGHS numerical tolerances.
74/// Zero MIP gaps do not make floating-point solving mathematically exact.
75///
76/// # Example
77///
78/// ```rust
79/// use problemreductions::models::algebraic::{ILP, LinearConstraint, ObjectiveSense};
80/// use problemreductions::solvers::ILPSolver;
81///
82/// // Create a simple binary ILP: maximize x0 + 2*x1 subject to x0 + x1 <= 1
83/// let ilp = ILP::<bool, f64>::new(
84///     2,
85///     vec![LinearConstraint::le(vec![(0, 1.0), (1, 1.0)], 1.0)],
86///     vec![(0, 1.0), (1, 2.0)],
87///     ObjectiveSense::Maximize,
88/// )?;
89///
90/// let solver = ILPSolver::new();
91/// let solution = solver.solve(&ilp)?;
92/// println!("Solution: {:?}", solution);
93/// # Ok::<(), Box<dyn std::error::Error>>(())
94/// ```
95#[derive(Debug, Clone, Default)]
96pub struct ILPSolver {
97    /// Time limit in seconds (None = no limit).
98    pub time_limit: Option<f64>,
99}
100
101impl ILPSolver {
102    /// Create a new ILP solver with default settings.
103    pub fn new() -> Self {
104        Self::default()
105    }
106
107    /// Create an ILP solver with a time limit.
108    pub fn with_time_limit(seconds: f64) -> Self {
109        Self {
110            time_limit: Some(seconds),
111        }
112    }
113
114    /// Solve a problem through its registered ILP pipeline.
115    ///
116    /// Returns a classified error when the problem is infeasible, the time
117    /// limit is reached, the pipeline is missing, or the backend fails.
118    pub fn solve<P>(&self, problem: &P) -> Result<P::Solution, ILPSolveError>
119    where
120        P: Problem + 'static,
121        P::Solution: 'static,
122    {
123        let key = ExactProblemKey::new(P::NAME, crate::export::variant_to_map(P::variant()));
124        let registry = solver_capability_registry()
125            .map_err(|error| ILPSolveError::InvalidRegistry(error.to_string()))?;
126        let pipeline = registry
127            .lookup(&key)
128            .ilp
129            .ok_or_else(|| ILPSolveError::MissingPipeline(key.label()))?;
130        pipeline.solve_typed(problem, self)
131    }
132
133    fn solve_backend<V>(&self, problem: &ILP<V, f64>) -> Result<Vec<i64>, ILPSolveError>
134    where
135        V: VariableDomain,
136    {
137        self.solve_with_objective(problem, problem.objective())
138    }
139
140    fn solve_with_objective<V>(
141        &self,
142        problem: &ILP<V, f64>,
143        objective_terms: &[(usize, f64)],
144    ) -> Result<Vec<i64>, ILPSolveError>
145    where
146        V: VariableDomain,
147    {
148        let n = problem.num_vars();
149        if n == 0 {
150            return if problem
151                .is_feasible(&[])
152                .map_err(|error| ILPSolveError::InvalidSolution(error.to_string()))?
153            {
154                Ok(vec![])
155            } else {
156                Err(ILPSolveError::Infeasible)
157            };
158        }
159
160        let mut vars_builder = ProblemVariables::new();
161        let vars: Vec<Variable> = problem
162            .variables()
163            .iter()
164            .map(|variable_bounds| {
165                let mut definition = variable().integer();
166                if let Some(lower) = variable_bounds.lower_bound() {
167                    definition = definition.min(i64_to_exact_f64(lower)?);
168                }
169                if let Some(upper) = variable_bounds.upper_bound() {
170                    definition = definition.max(i64_to_exact_f64(upper)?);
171                }
172                Ok(vars_builder.add(definition))
173            })
174            .collect::<Result<_, ILPSolveError>>()?;
175
176        // Build objective expression
177        let objective: good_lp::Expression = objective_terms
178            .iter()
179            .map(|&(var_idx, coefficient)| coefficient * vars[var_idx])
180            .sum();
181
182        // Build the model with objective
183        let unsolved = match problem.sense() {
184            ObjectiveSense::Maximize => vars_builder.maximise(&objective),
185            ObjectiveSense::Minimize => vars_builder.minimise(&objective),
186        };
187
188        // Create the solver model
189        let mut model = {
190            let mut model = unsolved
191                .using(highs)
192                .set_option("random_seed", 0i32)
193                .set_option("mip_rel_gap", 0.0)
194                .set_option("mip_abs_gap", 0.0)
195                .set_parallel(HighsParallelType::Off)
196                .set_threads(1);
197            if let Some(seconds) = self.time_limit {
198                model = model.set_time_limit(seconds);
199            }
200            model
201        };
202
203        // Add constraints
204        for constraint in problem.constraints() {
205            // Build left-hand side expression
206            let lhs: good_lp::Expression = constraint
207                .terms()
208                .iter()
209                .map(|&(var_idx, coefficient)| coefficient * vars[var_idx])
210                .sum();
211
212            let rhs = constraint.rhs();
213
214            // Create the constraint based on comparison type
215            let good_lp_constraint = match constraint.comparison() {
216                Comparison::Le => lhs.leq(rhs),
217                Comparison::Ge => lhs.geq(rhs),
218                Comparison::Eq => lhs.eq(rhs),
219            };
220
221            model = model.with(good_lp_constraint);
222        }
223
224        // Solve
225        let solution = match model.solve() {
226            Ok(solution) => solution,
227            Err(ResolutionError::Infeasible)
228                if !objective_terms.is_empty()
229                    && problem.variables().iter().any(|variable| {
230                        variable.lower_bound().is_none() || variable.upper_bound().is_none()
231                    }) =>
232            {
233                // A zero objective cannot be unbounded, so feasibility distinguishes the two states.
234                self.solve_with_objective(problem, &[])?;
235                return Err(ILPSolveError::Unbounded);
236            }
237            Err(error) => return Err(classify_backend_error(error, self.time_limit)),
238        };
239
240        match solution.status() {
241            SolutionStatus::Optimal => {}
242            SolutionStatus::TimeLimit => return Err(ILPSolveError::Timeout),
243            SolutionStatus::GapLimit => {
244                return Err(ILPSolveError::BackendFailure(
245                    "the backend stopped at its gap limit before proving optimality".to_string(),
246                ));
247            }
248        }
249
250        let result: Vec<i64> = vars
251            .iter()
252            .enumerate()
253            .map(|(index, v)| {
254                let value = solution.value(*v);
255                if !value.is_finite() {
256                    return Err(ILPSolveError::InvalidSolution(format!(
257                        "variable {index} is non-finite"
258                    )));
259                }
260                let rounded = value.round();
261                if (value - rounded).abs() > 1e-6 {
262                    return Err(ILPSolveError::InvalidSolution(format!(
263                        "variable {index} has non-integral value {value}"
264                    )));
265                }
266                if rounded.abs() > MAX_EXACT_F64_INTEGER as f64 {
267                    return Err(ILPSolveError::InvalidSolution(format!(
268                        "variable {index} value {rounded} exceeds exact f64 integer transport"
269                    )));
270                }
271                Ok(rounded as i64)
272            })
273            .collect::<Result<_, _>>()?;
274
275        if !problem
276            .is_feasible(&result)
277            .map_err(|error| ILPSolveError::InvalidSolution(error.to_string()))?
278        {
279            return Err(ILPSolveError::InvalidSolution(
280                "the rounded assignment violates the ILP".into(),
281            ));
282        }
283
284        Ok(result)
285    }
286
287    /// Solve a type-erased supported ILP variant directly.
288    pub(crate) fn solve_dyn(&self, any: &dyn std::any::Any) -> Result<Vec<i64>, ILPSolveError> {
289        if let Some(ilp) = any.downcast_ref::<ILP<bool, f64>>() {
290            return self.solve_backend(ilp);
291        }
292        if let Some(ilp) = any.downcast_ref::<ILP<i64, f64>>() {
293            return self.solve_backend(ilp);
294        }
295        Err(ILPSolveError::UnsupportedProblemType)
296    }
297}
298
299#[cfg(test)]
300#[path = "../../unit_tests/solvers/ilp/solver.rs"]
301mod tests;
problemreductions/solvers/ilp/solver.rs

problemreductions/solvers/ilp/
solver.rs
