Skip to content

Commit

Permalink
Merge pull request #25 from kehlert/dev
Browse files Browse the repository at this point in the history 
Dev
  • Loading branch information
@kehlert
kehlert authored Jun 11, 2021
2 parents 3427b3e + a691ea4 commit 348dd51
Show file tree
Hide file tree
Showing 8 changed files with 126 additions and 264 deletions.
80 changes: 40 additions & 40 deletions README.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,60 +12,60 @@ Here is example code that sets up a linear program, and then solves it with both


```rust
use ellp::*;
use ellp::*;

let mut prob = Problem::new();
let mut prob = Problem::new();

let x1 = prob
.add_var(2., Bound::TwoSided(-1., 1.), Some("x1".to_string()))
.unwrap();
let x1 = prob
.add_var(2., Bound::TwoSided(-1., 1.), Some("x1".to_string()))
.unwrap();

let x2 = prob
.add_var(10., Bound::Upper(6.), Some("x2".to_string()))
.unwrap();
let x2 = prob
.add_var(10., Bound::Upper(6.), Some("x2".to_string()))
.unwrap();

let x3 = prob
.add_var(0., Bound::Lower(0.), Some("x3".to_string()))
.unwrap();
let x3 = prob
.add_var(0., Bound::Lower(0.), Some("x3".to_string()))
.unwrap();

let x4 = prob
.add_var(1., Bound::Fixed(0.), Some("x4".to_string()))
.unwrap();
let x4 = prob
.add_var(1., Bound::Fixed(0.), Some("x4".to_string()))
.unwrap();

let x5 = prob
.add_var(0., Bound::Free, Some("x5".to_string()))
.unwrap();
let x5 = prob
.add_var(0., Bound::Free, Some("x5".to_string()))
.unwrap();

prob.add_constraint(vec![(x1, 2.5), (x2, 3.5)], ConstraintOp::Gte, 5.)
.unwrap();
prob.add_constraint(vec![(x1, 2.5), (x2, 3.5)], ConstraintOp::Gte, 5.)
.unwrap();

prob.add_constraint(vec![(x2, 2.5), (x1, 4.5)], ConstraintOp::Lte, 1.)
.unwrap();
prob.add_constraint(vec![(x2, 2.5), (x1, 4.5)], ConstraintOp::Lte, 1.)
.unwrap();

prob.add_constraint(vec![(x3, -1.), (x4, -3.), (x5, -4.)], ConstraintOp::Eq, 2.)
.unwrap();
prob.add_constraint(vec![(x3, -1.), (x4, -3.), (x5, -4.)], ConstraintOp::Eq, 2.)
.unwrap();

println!("{}", prob);
println!("{}", prob);

let primal_solver = PrimalSimplexSolver::default();
let dual_solver = DualSimplexSolver::default();
let primal_solver = PrimalSimplexSolver::default();
let dual_solver = DualSimplexSolver::default();

let primal_result = primal_solver.solve(&prob).unwrap();
let dual_result = dual_solver.solve(&prob).unwrap();
let primal_result = primal_solver.solve(prob.clone()).unwrap();
let dual_result = dual_solver.solve(prob).unwrap();

if let SolverResult::Optimal(sol) = primal_result {
println!("primal obj: {}", sol.obj());
println!("primal opt point: {}", sol.x());
} else {
panic!("should have an optimal point");
}
if let SolverResult::Optimal(sol) = primal_result {
println!("primal obj: {}", sol.obj());
println!("primal opt point: {}", sol.x());
} else {
panic!("should have an optimal point");
}

if let SolverResult::Optimal(sol) = dual_result {
println!("dual obj: {}", sol.obj());
println!("dual opt point: {}", sol.x());
} else {
panic!("should have an optimal point");
}
if let SolverResult::Optimal(sol) = dual_result {
println!("dual obj: {}", sol.obj());
println!("dual opt point: {}", sol.x());
} else {
panic!("should have an optimal point");
}
```

The output is
Expand Down
93 changes: 21 additions & 72 deletions src/dual/dual_problem.rs
View file
Original file line number Diff line number Diff line change
@@ -1,22 +1,17 @@
#![allow(non_snake_case)]

use crate::problem::{Bound, ConstraintOp, Problem, VariableId};
use crate::standard_form::{Basic, BasicPoint, Nonbasic, NonbasicBound, Point, StandardForm};
use crate::problem::{Bound, ConstraintOp, Problem};
use crate::standard_form::{
Basic, BasicPoint, Nonbasic, NonbasicBound, Point, StandardForm, StandardizedProblem,
};
use crate::util::EPS;

use std::collections::HashSet;

pub trait DualProblem {
fn obj(&self) -> f64;
fn std_form(&self) -> &StandardForm;
fn pt(&self) -> &DualFeasiblePoint;
fn pt_mut(&mut self) -> &mut DualFeasiblePoint;
fn unpack(&mut self) -> (&StandardForm, &mut DualFeasiblePoint);
}

#[derive(Debug, Clone)]
pub struct DualFeasiblePoint {
pub y: nalgebra::DVector<f64>,
pub d: nalgebra::DVector<f64>,
pub point: Point,
}

Expand Down Expand Up @@ -57,28 +52,16 @@ impl DualPhase1 {
}
}

impl DualProblem for DualPhase1 {
fn obj(&self) -> f64 {
self.std_form.obj(&self.point.x)
}

#[inline]
fn std_form(&self) -> &StandardForm {
&self.std_form
}

#[inline]
fn pt(&self) -> &DualFeasiblePoint {
&self.point
}
impl StandardizedProblem for DualPhase1 {
type FeasiblePoint = DualFeasiblePoint;

#[inline]
fn pt_mut(&mut self) -> &mut DualFeasiblePoint {
&mut self.point
fn obj(&self) -> f64 {
self.std_form.dual_obj(&self.point.y, &self.point.d)
}

#[inline]
fn unpack(&mut self) -> (&StandardForm, &mut DualFeasiblePoint) {
fn unpack(&mut self) -> (&StandardForm, &mut Self::FeasiblePoint) {
(&self.std_form, &mut self.point)
}
}
Expand All @@ -89,36 +72,22 @@ pub struct DualPhase2 {
pub point: DualFeasiblePoint,
}

impl DualProblem for DualPhase2 {
fn obj(&self) -> f64 {
self.std_form.obj(&self.point.x)
}

#[inline]
fn std_form(&self) -> &StandardForm {
&self.std_form
}

#[inline]
fn pt(&self) -> &DualFeasiblePoint {
&self.point
}
impl StandardizedProblem for DualPhase2 {
type FeasiblePoint = DualFeasiblePoint;

#[inline]
fn pt_mut(&mut self) -> &mut DualFeasiblePoint {
&mut self.point
fn obj(&self) -> f64 {
self.std_form.dual_obj(&self.point.y, &self.point.d)
}

#[inline]
fn unpack(&mut self) -> (&StandardForm, &mut DualFeasiblePoint) {
fn unpack(&mut self) -> (&StandardForm, &mut Self::FeasiblePoint) {
(&self.std_form, &mut self.point)
}
}

impl std::convert::From<Problem> for Option<DualPhase1> {
fn from(prob: Problem) -> Self {
println!("orig prob:\n{}", prob);

let orig_std_form: StandardForm = match prob.into() {
Some(std_form) => std_form,
None => return None,
Expand Down Expand Up @@ -162,8 +131,6 @@ impl std::convert::From<Problem> for Option<DualPhase1> {
.unwrap();
}

println!("{}", phase_1_prob);

let std_form: StandardForm = match phase_1_prob.into() {
Some(std_form) => std_form,
None => return None,
Expand Down Expand Up @@ -192,12 +159,6 @@ impl std::convert::From<Problem> for Option<DualPhase1> {

let c_B = nalgebra::DVector::from_iterator(B.len(), B.iter().map(|i| std_form.c[i.index]));
let A_B = std_form.A.select_columns(B.iter().map(|b| &b.index));

println!("c:{}", std_form.c);
println!("A:{}", std_form.A);
println!("A_B:{}", A_B);
println!("B:{:?}", B);

let A_B_lu = A_B.lu();

if !B.is_empty() {
Expand All @@ -207,8 +168,6 @@ impl std::convert::From<Problem> for Option<DualPhase1> {
A_B_lu.p().inv_permute_rows(&mut y);

let d = &std_form.c - std_form.A.tr_mul(&y);
println!("y: {}", y);
println!("d: {}", d);

let mut x = nalgebra::DVector::zeros(std_form.bounds.len());

Expand Down Expand Up @@ -252,17 +211,9 @@ impl std::convert::From<Problem> for Option<DualPhase1> {
x[i.index] = *val;
}

println!("N:{:?}", N);

println!("{:?}", std_form.bounds);
println!("y: {}", y);
println!("x: {}", x);
println!("d: {}", d);
println!("Ax:{}", &std_form.A * &x);
println!("b:{}", std_form.b);

let point = DualFeasiblePoint {
y,
d,
point: Point { x, N, B },
};

Expand All @@ -287,11 +238,9 @@ impl std::convert::From<Problem> for Option<DualPhase1> {
}
}

println!("HEREERERE: {}", x);
println!("N: {:?}", N);

let point = DualFeasiblePoint {
y: empty_vec,
d: std_form.c.clone(),
point: Point { x, N, B },
};

Expand Down Expand Up @@ -365,8 +314,6 @@ impl std::convert::From<DualPhase1> for DualPhase2 {
Bound::Fixed(val) => (*val, NonbasicBound::Lower),
};

println!("i: {}, x_i: {}, d_i: {}", i, x_i, d_i);

Some((x_i, Nonbasic::new(i, bound)))
} else {
None
Expand Down Expand Up @@ -394,10 +341,11 @@ impl std::convert::From<DualPhase1> for DualPhase2 {

let point = DualFeasiblePoint {
y,
d,
point: Point { x, N, B },
};

DualPhase2 { point, std_form }
DualPhase2 { std_form, point }
} else {
let empty_vec = nalgebra::DVector::zeros(0);
let number_nonbasic = std_form.A.ncols();
Expand Down Expand Up @@ -444,10 +392,11 @@ impl std::convert::From<DualPhase1> for DualPhase2 {

let point = DualFeasiblePoint {
y: empty_vec,
d: std_form.c.clone(),
point: Point { x: x_N, N, B },
};

DualPhase2 { point, std_form }
DualPhase2 { std_form, point }
}
}
}
Loading

0 comments on commit 348dd51

Please sign in to comment.