Skip to main content

problemreductions/models/graph/
minimum_cost_maximum_flow.rs

1//! Minimum-Cost Maximum-Flow problem implementation.
2//!
3//! Given a directed graph `G = (V, A)` with arc capacities `c(a)` and
4//! arc costs `cost(a)`, a source `s`, and a sink `t`, find a flow `f`
5//! that
6//!
7//! 1. first maximizes the flow value `|f| = sum_{a in delta^+(s)} f(a)
8//!    - sum_{a in delta^-(s)} f(a)`, and then
9//! 2. among all maximum-value flows, minimizes the total arc cost
10//!    `sum_a cost(a) * f(a)`.
11//!
12//! The objective is lexicographic: maximum-value first, ties broken by
13//! lower total cost.
14//!
15//! # Integral-flow restriction
16//!
17//! The mathematical CellRouter formulation in issue #1029 uses
18//! continuous flows `f: A -> R_{>= 0}`, but the [`Problem`] trait
19//! requires a discrete configuration space via [`Problem::dims`].
20//! Following the same precedent as
21//! [`MinimumEdgeCostFlow`](super::MinimumEdgeCostFlow) (see
22//! `src/models/graph/minimum_edge_cost_flow.rs`), we therefore restrict
23//! to **integer** flows: each variable `f(a)` ranges over
24//! `{0, 1, ..., c(a)}`. When capacities and costs are integral, the
25//! standard minimum-cost flow theory (see e.g. the MIT 6.854 notes,
26//! Ahuja-Magnanti-Orlin) guarantees that an integral optimum exists, so
27//! this restriction does not change the optimal value on integer
28//! instances.
29//!
30//! # Lexicographic encoding
31//!
32//! The lexicographic objective `(max |f|, min cost(f))` is encoded as a
33//! single scalar score
34//!
35//! `score = M * (max_possible_flow - |f|) + cost(f)`
36//!
37//! where `M = sum_e c(e) * cost(e) + 1` strictly exceeds any feasible
38//! cost. Minimizing this scalar therefore minimizes
39//! `max_possible_flow - |f|` first (i.e. maximizes `|f|`), then breaks
40//! ties by `cost(f)`. The optimum is always non-negative, and a smaller
41//! score is strictly better in the lex order.
42
43use crate::registry::{FieldInfo, ProblemSchemaEntry, ProblemSizeFieldEntry};
44use crate::topology::DirectedGraph;
45use crate::traits::Problem;
46use serde::{Deserialize, Serialize};
47
48inventory::submit! {
49    ProblemSchemaEntry {
50        name: "MinimumCostMaximumFlow",
51        display_name: "Minimum-Cost Maximum-Flow",
52        aliases: &["MCMF"],
53        dimensions: &[],
54        module_path: module_path!(),
55        description: "Integral flow that lexicographically maximizes value then minimizes total arc cost",
56        fields: &[
57            FieldInfo { name: "graph", type_name: "DirectedGraph", description: "Directed graph G = (V, A)" },
58            FieldInfo { name: "source", type_name: "usize", description: "Source vertex s" },
59            FieldInfo { name: "sink", type_name: "usize", description: "Sink vertex t" },
60            FieldInfo { name: "capacities", type_name: "Vec<i64>", description: "Arc capacity c(a) in graph arc order (non-negative)" },
61            FieldInfo { name: "costs", type_name: "Vec<i64>", description: "Arc cost cost(a) in graph arc order (non-negative)" },
62        ],
63    }
64}
65
66inventory::submit! {
67    ProblemSizeFieldEntry {
68        name: "MinimumCostMaximumFlow",
69        fields: &["num_vertices", "num_arcs"],
70    }
71}
72
73/// Minimum-Cost Maximum-Flow problem.
74///
75/// # Variables
76///
77/// `|A|` variables: variable `a` ranges over `{0, ..., c(a)}`
78/// representing the integral flow on arc `a`.
79///
80/// # Example
81///
82/// ```
83/// use problemreductions::models::graph::MinimumCostMaximumFlow;
84/// use problemreductions::topology::DirectedGraph;
85/// use problemreductions::{Problem, Solver, BruteForce};
86///
87/// // Diamond network from the canonical example.
88/// let graph = DirectedGraph::new(4, vec![
89///     (0, 1), (0, 2), (1, 2), (1, 3), (2, 3),
90/// ]);
91/// let problem = MinimumCostMaximumFlow::new(
92///     graph,
93///     0, 3,
94///     vec![2, 1, 1, 1, 2], // capacities
95///     vec![1, 0, 0, 1, 2], // costs
96/// );
97/// let solver = BruteForce::new();
98/// let witness = solver.find_witness(&problem).unwrap();
99/// // Optimal flow has value 3 and cost 7.
100/// assert_eq!(problem.flow_value(&witness), 3);
101/// assert_eq!(problem.total_cost(&witness), 7);
102/// ```
103#[derive(Debug, Clone, Serialize, Deserialize)]
104pub struct MinimumCostMaximumFlow {
105    /// The directed graph G = (V, A).
106    graph: DirectedGraph,
107    /// Source vertex s.
108    source: usize,
109    /// Sink vertex t.
110    sink: usize,
111    /// Capacity c(a) for each arc.
112    capacities: Vec<i64>,
113    /// Cost cost(a) for each arc.
114    costs: Vec<i64>,
115}
116
117impl MinimumCostMaximumFlow {
118    /// Create a new Minimum-Cost Maximum-Flow problem.
119    ///
120    /// # Panics
121    ///
122    /// Panics if any of the following holds:
123    /// - `capacities.len() != graph.num_arcs()`
124    /// - `costs.len() != graph.num_arcs()`
125    /// - `source >= graph.num_vertices()`
126    /// - `sink >= graph.num_vertices()`
127    /// - `source == sink`
128    /// - Any capacity is negative
129    /// - Any cost is negative
130    pub fn new(
131        graph: DirectedGraph,
132        source: usize,
133        sink: usize,
134        capacities: Vec<i64>,
135        costs: Vec<i64>,
136    ) -> Self {
137        let n = graph.num_vertices();
138        let m = graph.num_arcs();
139        assert_eq!(
140            capacities.len(),
141            m,
142            "capacities length ({}) must match num_arcs ({m})",
143            capacities.len()
144        );
145        assert_eq!(
146            costs.len(),
147            m,
148            "costs length ({}) must match num_arcs ({m})",
149            costs.len()
150        );
151        assert!(source < n, "source ({source}) >= num_vertices ({n})");
152        assert!(sink < n, "sink ({sink}) >= num_vertices ({n})");
153        assert_ne!(source, sink, "source and sink must be distinct");
154        for (i, &c) in capacities.iter().enumerate() {
155            assert!(c >= 0, "capacity[{i}] = {c} is negative");
156        }
157        for (i, &c) in costs.iter().enumerate() {
158            assert!(c >= 0, "cost[{i}] = {c} is negative");
159        }
160        Self {
161            graph,
162            source,
163            sink,
164            capacities,
165            costs,
166        }
167    }
168
169    /// Get a reference to the underlying directed graph.
170    pub fn graph(&self) -> &DirectedGraph {
171        &self.graph
172    }
173
174    /// Get the source vertex.
175    pub fn source(&self) -> usize {
176        self.source
177    }
178
179    /// Get the sink vertex.
180    pub fn sink(&self) -> usize {
181        self.sink
182    }
183
184    /// Get a reference to the arc capacities.
185    pub fn capacities(&self) -> &[i64] {
186        &self.capacities
187    }
188
189    /// Get a reference to the arc costs.
190    pub fn costs(&self) -> &[i64] {
191        &self.costs
192    }
193
194    /// Get the number of vertices `|V|`.
195    pub fn num_vertices(&self) -> usize {
196        self.graph.num_vertices()
197    }
198
199    /// Get the number of arcs `|A|`.
200    pub fn num_arcs(&self) -> usize {
201        self.graph.num_arcs()
202    }
203
204    /// Check whether a flow assignment is feasible.
205    ///
206    /// A flow is feasible iff
207    /// 1. `config.len() == num_arcs`,
208    /// 2. each `0 <= f(a) <= c(a)`, and
209    /// 3. flow is conserved at every non-terminal vertex.
210    pub fn is_feasible(&self, config: &[usize]) -> bool {
211        let m = self.graph.num_arcs();
212        if config.len() != m {
213            return false;
214        }
215        // (1) Capacity constraints
216        for (flow, cap) in config.iter().zip(self.capacities.iter()) {
217            if (*flow as i64) > *cap {
218                return false;
219            }
220        }
221        // (2) Flow conservation at non-terminal vertices
222        let n = self.graph.num_vertices();
223        let mut balance = vec![0_i64; n];
224        for (a, &(u, v)) in self.graph.arcs().iter().enumerate() {
225            let flow = config[a] as i64;
226            balance[u] -= flow;
227            balance[v] += flow;
228        }
229        for (v, &bal) in balance.iter().enumerate() {
230            if v != self.source && v != self.sink && bal != 0 {
231                return false;
232            }
233        }
234        true
235    }
236
237    /// Compute the flow value `|f|` = net outflow from the source for a
238    /// feasible configuration. Result is meaningless if `config` is not
239    /// feasible.
240    pub fn flow_value(&self, config: &[usize]) -> i64 {
241        let mut net_out: i64 = 0;
242        for (a, &(u, v)) in self.graph.arcs().iter().enumerate() {
243            let f = config[a] as i64;
244            if u == self.source {
245                net_out += f;
246            }
247            if v == self.source {
248                net_out -= f;
249            }
250        }
251        net_out
252    }
253
254    /// Compute the total cost `sum_a cost(a) * f(a)` of a flow.
255    pub fn total_cost(&self, config: &[usize]) -> i64 {
256        config
257            .iter()
258            .zip(self.costs.iter())
259            .map(|(&f, &c)| (f as i64) * c)
260            .sum()
261    }
262
263    /// Upper bound on the integral flow value: `sum_a c(a)` (a trivial
264    /// but valid bound, since `|f|` is bounded by the total capacity).
265    fn max_possible_flow(&self) -> i64 {
266        self.capacities.iter().sum()
267    }
268
269    /// Strict upper bound on any feasible cost, used as the
270    /// lex-multiplier `M` so that the scalar `score = M * (B - |f|)
271    /// + cost(f)` orders by `(max |f|, min cost(f))`.
272    fn cost_multiplier(&self) -> i64 {
273        self.capacities
274            .iter()
275            .zip(self.costs.iter())
276            .map(|(&c, &k)| c * k)
277            .sum::<i64>()
278            + 1
279    }
280}
281
282impl Problem for MinimumCostMaximumFlow {
283    const NAME: &'static str = "MinimumCostMaximumFlow";
284    type Value = crate::types::Min<i64>;
285
286    fn dims(&self) -> Vec<usize> {
287        self.capacities.iter().map(|&c| (c as usize) + 1).collect()
288    }
289
290    fn evaluate(&self, config: &[usize]) -> crate::types::Min<i64> {
291        if !self.is_feasible(config) {
292            return crate::types::Min(None);
293        }
294        let m = self.cost_multiplier();
295        let value = self.flow_value(config);
296        let cost = self.total_cost(config);
297        let bound = self.max_possible_flow();
298        // score = M * (max_possible_flow - |f|) + cost(f)
299        crate::types::Min(Some(m * (bound - value) + cost))
300    }
301
302    fn variant() -> Vec<(&'static str, &'static str)> {
303        crate::variant_params![]
304    }
305}
306
307crate::declare_variants! {
308    default MinimumCostMaximumFlow => "(num_vertices + num_arcs)^6",
309}
310
311#[cfg(feature = "example-db")]
312pub(crate) fn canonical_model_example_specs() -> Vec<crate::example_db::specs::ModelExampleSpec> {
313    let problem = MinimumCostMaximumFlow::new(
314        crate::topology::DirectedGraph::new(4, vec![(0, 1), (0, 2), (1, 2), (1, 3), (2, 3)]),
315        0,
316        3,
317        vec![2, 1, 1, 1, 2],
318        vec![1, 0, 0, 1, 2],
319    );
320    // Optimal flow has value 3, routed as:
321    //   - 1 unit on 0->1->3        via arcs 0,3 (cost 1 + 1 = 2)
322    //   - 1 unit on 0->1->2->3     via arcs 0,2,4 (cost 1 + 0 + 2 = 3)
323    //   - 1 unit on 0->2->3        via arcs 1,4 (cost 0 + 2 = 2)
324    // Arc flows sum to f = [2, 1, 1, 1, 2]: value = 3,
325    // cost = 2*1 + 1*0 + 1*0 + 1*1 + 2*2 = 7.
326    let optimal_config = vec![2_usize, 1, 1, 1, 2];
327    let optimal_value = problem.evaluate(&optimal_config);
328    let scalar = match optimal_value {
329        crate::types::Min(Some(v)) => v,
330        crate::types::Min(None) => panic!("canonical example must be feasible"),
331    };
332    vec![crate::example_db::specs::ModelExampleSpec {
333        id: "minimum_cost_maximum_flow",
334        instance: Box::new(problem),
335        optimal_config,
336        optimal_value: serde_json::json!(scalar),
337    }]
338}
339
340#[cfg(test)]
341#[path = "../../unit_tests/models/graph/minimum_cost_maximum_flow.rs"]
342mod tests;