Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
55 changes: 36 additions & 19 deletions examples/column_range.rs
Original file line number Diff line number Diff line change
Expand Up @@ -18,19 +18,19 @@
//! minimums and column maximums are, as well as the slack variable for each of the constraints that
//! don't already contain a pivot from the column minimum or maximum.
use num_traits::Zero;
use relp_num::RationalBig;
use relp_num::RB;
use relp_num::RationalBig;

use relp::algorithm::OptimizationResult;
use relp::algorithm::two_phase::matrix_provider::matrix_data::MatrixData;
use relp::algorithm::two_phase::matrix_provider::MatrixProvider;
use relp::algorithm::two_phase::matrix_provider::matrix_data::MatrixData;
use relp::algorithm::two_phase::phase_two;
use relp::algorithm::two_phase::strategy::pivot_rule::FirstProfitable;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::basis_inverse_rows::BasisInverseRows;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::Carry;
use relp::algorithm::two_phase::tableau::Tableau;
use relp::algorithm::two_phase::tableau::inverse_maintenance::InverseMaintainer;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::Carry;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::basis_inverse_rows::BasisInverseRows;
use relp::algorithm::two_phase::tableau::kind::non_artificial::NonArtificial;
use relp::algorithm::two_phase::tableau::Tableau;
use relp::data::linear_algebra::matrix::{ColumnMajor, MatrixOrder};
use relp::data::linear_algebra::vector::{DenseVector, Vector};
use relp::data::linear_program::elements::VariableType;
Expand All @@ -42,21 +42,24 @@ fn main() {
// non-negative.
let input_matrix = [
[Some(3), Some(3), Some(3)],
[None, Some(3), Some(3)],
[None, Some(3), Some(3)],
[Some(1), Some(2), Some(3)],
];

let m = input_matrix.len();
let n = input_matrix[0].len();

let column_extreme = |f: fn(_) -> Option<i32>| (0..n)
.map(|j| f((0..m).filter_map(move |i| input_matrix[i][j])).unwrap())
.collect::<Vec<_>>();
let column_extreme = |f: fn(_) -> Option<i32>| {
(0..n)
.map(|j| f((0..m).filter_map(move |i| input_matrix[i][j])).unwrap())
.collect::<Vec<_>>()
};
let column_min = column_extreme(Iterator::min);
let column_max = column_extreme(Iterator::max);

// Variables
let subtraction_amount = (0..m).map(|_| Variable {
let subtraction_amount = (0..m)
.map(|_| Variable {
variable_type: VariableType::Continuous,
cost: RB!(0),
lower_bound: Some(RB!(0)),
Expand All @@ -65,7 +68,8 @@ fn main() {
flipped: false,
})
.collect::<Vec<_>>();
let column_minimum = (0..n).map(|_| Variable {
let column_minimum = (0..n)
.map(|_| Variable {
variable_type: VariableType::Continuous,
// We subtract this variable from a corresponding variable to compute a column range.
cost: RB!(-1),
Expand All @@ -75,7 +79,8 @@ fn main() {
flipped: false,
})
.collect::<Vec<_>>();
let column_maximum = (0..n).map(|_| Variable {
let column_maximum = (0..n)
.map(|_| Variable {
variable_type: VariableType::Continuous,
cost: RB!(1),
lower_bound: Some(RB!(0)),
Expand Down Expand Up @@ -107,9 +112,9 @@ fn main() {
if let Some(value) = input_matrix[i][j] {
// There is a value there, so we add a constraint.
row_major_constraints.push(vec![
(i, RB!(1)), // The "subtraction amount" variable index
(i, RB!(1)), // The "subtraction amount" variable index
(m + j, RB!(1)), // The "column minimum" variable index (there are `m` of the
// "subtraction amount" variables)
// "subtraction amount" variables)
]);
b.push(RB!(value));
basis_columns.push(if value == column_min[j] {
Expand Down Expand Up @@ -163,15 +168,25 @@ fn main() {
constraints[column_index].push((row_index, value));
}
}
let constraints = ColumnMajor::new(constraints, nr_upper_bounded_constraints + nr_lower_bounded_constraints, m + 2 * n);
let b = DenseVector::new(b, nr_upper_bounded_constraints + nr_lower_bounded_constraints);
let constraints = ColumnMajor::new(
constraints,
nr_upper_bounded_constraints + nr_lower_bounded_constraints,
m + 2 * n,
);
let b = DenseVector::new(
b,
nr_upper_bounded_constraints + nr_lower_bounded_constraints,
);

// Create the datastructure that will serve as a `MatrixProvider`
let matrix = MatrixData::new(
&constraints,
&b,
Vec::with_capacity(0),
0, 0, nr_upper_bounded_constraints, nr_lower_bounded_constraints,
0,
0,
nr_upper_bounded_constraints,
nr_lower_bounded_constraints,
&variables,
);

Expand All @@ -182,7 +197,9 @@ fn main() {
// We create a tableau using the constructed matrix. The basis is initialized using the
// specified columns.
let mut tableau = Tableau::<_, NonArtificial<_>>::new_with_inverse_maintainer(
&matrix, inverse_maintainer, basis_columns.into_iter().collect(),
&matrix,
inverse_maintainer,
basis_columns.into_iter().collect(),
);

// We apply primal simplex to improve the solution.
Expand All @@ -197,7 +214,7 @@ fn main() {
})
.collect::<Vec<_>>();
println!("{:?}", shifts);
},
}
_ => panic!("We started with a feasible solution, and has at least value 0."),
}
}
82 changes: 47 additions & 35 deletions examples/max_flow.rs
Original file line number Diff line number Diff line change
Expand Up @@ -8,21 +8,23 @@ use std::ops::{Add, Mul, Range};

use index_utils::remove_sparse_indices;
use num_traits::Zero;
use relp_num::{Rational64, RationalBig};
use relp_num::NonZero;
use relp_num::One;
use relp_num::{Rational64, RationalBig};

use relp::algorithm::{OptimizationResult, SolveRelaxation};
use relp::algorithm::two_phase::matrix_provider::column::{Column as ColumnTrait, SparseSliceIterator};
use relp::algorithm::two_phase::matrix_provider::MatrixProvider;
use relp::algorithm::two_phase::matrix_provider::column::identity::Identity;
use relp::algorithm::two_phase::matrix_provider::column::{
Column as ColumnTrait, SparseSliceIterator,
};
use relp::algorithm::two_phase::matrix_provider::filter::generic_wrapper::IntoFilteredColumn;
use relp::algorithm::two_phase::matrix_provider::MatrixProvider;
use relp::algorithm::two_phase::phase_one::PartialInitialBasis;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::basis_inverse_rows::BasisInverseRows;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::Carry;
use relp::data::linear_algebra::matrix::{ColumnMajor, SparseMatrix};
use relp::data::linear_algebra::matrix::MatrixOrder;
use relp::algorithm::two_phase::tableau::inverse_maintenance::carry::basis_inverse_rows::BasisInverseRows;
use relp::algorithm::{OptimizationResult, SolveRelaxation};
use relp::data::linear_algebra::SparseTuple;
use relp::data::linear_algebra::matrix::MatrixOrder;
use relp::data::linear_algebra::matrix::{ColumnMajor, SparseMatrix};
use relp::data::linear_algebra::traits::{SparseComparator, SparseElement};
use relp::data::linear_algebra::vector::{DenseVector, SparseVector, Vector};
use relp::data::linear_program::elements::BoundDirection;
Expand Down Expand Up @@ -51,11 +53,7 @@ impl<F> Primal<F>
where
F: SparseElement<F> + SparseComparator,
{
pub fn new(
adjacency_matrix: SparseMatrix<F, F, ColumnMajor>,
s: usize,
t: usize,
) -> Self {
pub fn new(adjacency_matrix: SparseMatrix<F, F, ColumnMajor>, s: usize, t: usize) -> Self {
let nr_vertices = adjacency_matrix.nr_columns();
debug_assert!(s < nr_vertices && t < nr_vertices);

Expand Down Expand Up @@ -91,7 +89,7 @@ struct Column {

impl ColumnTrait for Column {
type F = ArcDirection;
type Iter<'a> = impl Iterator<Item=SparseTuple<&'a Self::F>> + Clone;
type Iter<'a> = impl Iterator<Item = SparseTuple<&'a Self::F>> + Clone;

fn iter(&self) -> Self::Iter<'_> {
SparseSliceIterator::new(&self.constraint_values)
Expand All @@ -114,7 +112,7 @@ impl Identity for Column {
fn identity(i: usize, _len: usize) -> Self {
Self {
constraint_values: Vec::with_capacity(0),
slack: (i, ArcDirection::Incoming)
slack: (i, ArcDirection::Incoming),
}
}
}
Expand All @@ -131,7 +129,7 @@ impl IntoFilteredColumn for Column {

impl IntoIterator for Column {
type Item = SparseTuple<ArcDirection>;
type IntoIter = impl Iterator<Item=Self::Item>;
type IntoIter = impl Iterator<Item = Self::Item>;

fn into_iter(self) -> Self::IntoIter {
self.constraint_values.into_iter().chain(once(self.slack))
Expand All @@ -143,7 +141,10 @@ where
F: SparseElement<F> + Zero + Eq + NonZero,
{
type Column = Column;
type Cost<'a> = Cost where Self: 'a;
type Cost<'a>
= Cost
where
Self: 'a;
type Rhs = F;

fn column(&self, j: usize) -> Self::Column {
Expand All @@ -157,7 +158,10 @@ where
} else {
Column {
constraint_values: Vec::with_capacity(0),
slack: (self.nr_constraints() + j - self.nr_edges(), ArcDirection::Incoming)
slack: (
self.nr_constraints() + j - self.nr_edges(),
ArcDirection::Incoming,
),
}
}
}
Expand All @@ -183,11 +187,13 @@ where

match bound_type {
BoundDirection::Lower => None,
BoundDirection::Upper => if j < self.nr_edges() {
Some(self.nr_constraints() + j)
} else {
None
},
BoundDirection::Upper => {
if j < self.nr_edges() {
Some(self.nr_constraints() + j)
} else {
None
}
}
}
}

Expand All @@ -214,7 +220,9 @@ where
F: SparseElement<F> + Zero + Eq + NonZero,
{
fn pivot_element_indices(&self) -> Vec<(usize, usize)> {
(0..self.nr_edges()).map(|j| (j + self.nr_constraints(), self.nr_edges() + j)).collect()
(0..self.nr_edges())
.map(|j| (j + self.nr_constraints(), self.nr_edges() + j))
.collect()
}

fn nr_initial_elements(&self) -> usize {
Expand All @@ -228,9 +236,7 @@ impl Add<Cost> for RationalBig {
fn add(self, rhs: Cost) -> Self::Output {
match rhs {
Cost::Zero => self,
Cost::MinusOne => {
self - One
}
Cost::MinusOne => self - One,
}
}
}
Expand Down Expand Up @@ -262,22 +268,28 @@ fn main() {
type S = RationalBig;

// Example from Papadimitriou's Combinatorial Optimization.
let data = ColumnMajor::from_test_data::<T, T, _>(&vec![
// Directed; from is top, to is on the right
// s a b t
vec![0, 0, 0, 0], // s
vec![2, 0, 0, 0], // a
vec![1, 1, 0, 0], // b
vec![0, 1, 2, 0], // t
], 4);
let data = ColumnMajor::from_test_data::<T, T, _>(
&vec![
// Directed; from is top, to is on the right
// s a b t
vec![0, 0, 0, 0], // s
vec![2, 0, 0, 0], // a
vec![1, 1, 0, 0], // b
vec![0, 1, 2, 0], // t
],
4,
);
let problem = Primal::new(data, 0, 3);

SolveRelaxation::solve_relaxation::<Carry<S, BasisInverseRows<S>>>(&problem);

assert_eq!(
problem.solve_relaxation::<Carry<S, BasisInverseRows<S>>>(),
OptimizationResult::FiniteOptimum(
[2, 1, 1, 1, 2, 0, 0, 0, 0, 0].iter().map(|&v| RationalBig::from(v)).collect()
[2, 1, 1, 1, 2, 0, 0, 0, 0, 0]
.iter()
.map(|&v| RationalBig::from(v))
.collect()
),
);
}
Loading