Browse Source

Added some test code for affine scaling -- it doesn't work yet.

karmarkar
Thomas Johnson 3 years ago
parent
commit
af53b7145b
  1. 239
      src/affine_scaling.rs
  2. 8
      src/clause.rs
  3. 5
      src/main.rs

239
src/affine_scaling.rs

@ -2,7 +2,7 @@
use core::fmt::{Display, Formatter};
use core::mem;
use core::ops::{Add, AddAssign, Mul, MulAssign, Neg};
use core::ops::{Add, AddAssign, Mul, MulAssign, Neg, Shl};
use std::collections::{HashMap, HashSet};
use nalgebra::{DMatrix, DVector};
@ -65,6 +65,17 @@ impl Add for Term {
}
}
impl Add<&Variable> for Term {
type Output = Expr;
fn add(self, other: &Variable) -> Self::Output {
let mut expr = Expr::new_zero();
expr.add_term(self);
expr.add_term(Self::from(other));
expr
}
}
impl Neg for Term {
type Output = Self;
@ -97,6 +108,15 @@ impl Display for Term {
}
}
impl From<&Variable> for Term {
fn from(var: &Variable) -> Self {
Self {
coeff: 1.0,
var: var.clone(),
}
}
}
pub fn combine_like_terms<I: IntoIterator<Item = Term>>(iter: I) -> impl Iterator<Item = Term> {
let mut map = HashMap::new();
for term in iter {
@ -199,6 +219,29 @@ impl AddAssign<Term> for Expr {
}
}
impl Add<&Variable> for Expr {
type Output = Expr;
fn add(mut self, addend: &Variable) -> Self::Output {
self += addend;
self
}
}
impl AddAssign<&Variable> for Expr {
fn add_assign(&mut self, addend: &Variable) {
self.add_term(Term::from(addend));
}
}
impl Shl<Expr> for Expr {
type Output = LeqInequality;
fn shl(self, rhs: Expr) -> LeqInequality {
LeqInequality { lhs: self, rhs }
}
}
impl Display for Expr {
fn fmt(&self, fmt: &mut Formatter) -> core::fmt::Result {
for term in &self.terms {
@ -208,6 +251,30 @@ impl Display for Expr {
}
}
impl From<Term> for Expr {
fn from(t: Term) -> Self {
Expr {
terms: vec![t],
offset: 0.0,
}
}
}
impl From<Fl> for Expr {
fn from(offset: Fl) -> Self {
Expr {
terms: Vec::new(),
offset,
}
}
}
impl From<&Variable> for Expr {
fn from(var: &Variable) -> Self {
Self::from(Term::from(var))
}
}
#[derive(Debug)]
pub struct NewVariableCounter(usize);
@ -291,6 +358,29 @@ impl LinearCombination {
}
}
impl Display for LinearCombination {
fn fmt(&self, fmt: &mut Formatter) -> core::fmt::Result {
let mut first = true;
for (var, val) in self.dir.iter() {
if first {
if *val >= 0.0 {
write!(fmt, "{}{}", val, var.get_name())?;
} else {
write!(fmt, "-{}{}", -val, var.get_name())?;
}
first = false;
} else {
if *val >= 0.0 {
write!(fmt, " + {}{}", val, var.get_name())?;
} else {
write!(fmt, " - {}{}", -val, var.get_name())?;
}
}
}
Ok(())
}
}
#[derive(Clone, Debug)]
pub struct Equality {
lhs: LinearCombination,
@ -322,10 +412,16 @@ impl Equality {
}
}
impl Display for Equality {
fn fmt(&self, fmt: &mut Formatter) -> core::fmt::Result {
write!(fmt, "{} = {}", self.lhs, self.rhs)
}
}
#[derive(Debug, Clone)]
#[allow(non_snake_case)]
pub struct StandardFormLP {
vars: Vec<Variable>,
pub vars: Vec<Variable>,
c: DVector<Fl>,
A: DMatrix<Fl>,
b: DVector<Fl>,
@ -349,16 +445,16 @@ impl StandardFormLP {
let nvars = vars.len();
let mut c = DVector::zeros(nvars);
#[allow(non_snake_case)]
let mut A = DMatrix::zeros(nvars, ncstr);
let mut A = DMatrix::zeros(ncstr, nvars);
let mut b = DVector::zeros(ncstr);
for (i, var) in vars.iter().enumerate() {
c[i] = objective.get_coeff(var);
for (j, constraint) in constraints.iter().enumerate() {
for (i, constraint) in constraints.iter().enumerate() {
b[i] = constraint.get_rhs();
for (j, var) in vars.iter().enumerate() {
A[(i, j)] = constraint.get_coeff(var);
}
}
for (j, constr) in constraints.into_iter().enumerate() {
b[j] = constr.get_rhs();
for (j, var) in vars.iter().enumerate() {
c[j] = objective.get_coeff(var);
}
StandardFormLP { vars, c, A, b }
}
@ -376,7 +472,7 @@ impl StandardFormLP {
#[allow(non_snake_case)]
let mut A = self.A.insert_column(lastcol, 0.0);
A.set_column(lastcol, &v);
let mut c = DVector::zeros(lastcol);
let mut c = DVector::zeros(lastcol + 1);
c[lastcol] = -1.0; // Try to minimize u
let u = Variable::from_name("u".into());
self.vars.push(u);
@ -384,6 +480,24 @@ impl StandardFormLP {
}
}
impl Display for StandardFormLP {
fn fmt(&self, fmt: &mut Formatter) -> core::fmt::Result {
write!(fmt, "variables: ")?;
for var in &self.vars {
write!(fmt, "{} ", var.get_name())?;
}
writeln!(fmt, "")?;
writeln!(fmt, "A:")?;
write!(fmt, "{}", &self.A)?;
writeln!(fmt, "b:")?;
write!(fmt, "{}", &self.b)?;
writeln!(fmt, "c:")?;
write!(fmt, "{}", &self.c)?;
Ok(())
}
}
#[derive(Debug, Clone, Copy)]
pub enum IterationResult {
Distance(Fl),
Indeterminate(),
@ -437,17 +551,51 @@ impl AffineScaling {
// Returns an upper bound on the error
pub fn iterate(&mut self) -> IterationResult {
use IterationResult::*;
#[cfg(debug_affine_scaling)]
{
print!("current infeasiblity:");
print!("{}", &self.problem.b - &self.problem.A * &self.x);
}
self.ADk2.copy_from(&self.problem.A);
for (i, mut col) in self.ADk2.column_iter_mut().enumerate() {
col *= self.x[i] * self.x[i];
}
#[cfg(debug_affine_scaling)]
{
print!("ADk2:");
print!("{}", &self.ADk2);
}
self.ADk2.mul_to(&self.problem.c, &mut self.w);
self.ADk2.mul_to(&self.AT, &mut self.lm);
#[cfg(debug_affine_scaling)]
let old_rhs = self.w.clone();
self.lm.clone().lu().solve_mut(&mut self.w);
#[cfg(debug_affine_scaling)]
{
print!("error in w:");
print!("{}", &self.lm * &self.w - old_rhs);
print!("w:");
print!("{}", &self.w);
}
self.problem.A.tr_mul_to(&self.w, &mut self.r);
self.r -= &self.problem.c;
#[cfg(debug_affine_scaling)]
{
print!("r:");
print!("{}", &self.r);
}
let min = self.r[self.r.imin()];
self.r.component_mul_assign(&self.x);
self.r.component_mul_assign(&self.x); // Reuse r no less than twice.
let rv = {
if min >= 0.0 {
Distance(self.r.sum())
@ -459,6 +607,13 @@ impl AffineScaling {
};
self.r.component_mul_assign(&self.x);
self.r *= self.β / self.r.norm();
#[cfg(debug_affine_scaling)]
{
print!("step:");
print!("{}", -&self.r);
}
self.x -= &self.r;
rv
}
@ -476,8 +631,70 @@ impl AffineScaling {
}
}
impl Display for AffineScaling {
fn fmt(&self, fmt: &mut Formatter) -> core::fmt::Result {
write!(fmt, "{}", &self.problem)?;
writeln!(fmt, "ε = {}", self.ε)?;
writeln!(fmt, "β = {}", self.β)?;
writeln!(fmt, "x:")?;
write!(fmt, "{}", &self.x)?;
Ok(())
}
}
/// Defines a DSL which allows for the declaration of many variables at once. See the unit tests
/// for examples.
#[macro_export]
macro_rules! linear_vars {
[] => {};
[$varname:ident$(,)?] => {
let $varname = $crate::affine_scaling::Variable::from_name(stringify!($varname).into());
};
[$varname:ident, $($varnames:ident),+] => {
let $varname = $crate::affine_scaling::Variable::from_name(stringify!($varname).into());
linear_vars![$($varnames),+]
};
}
#[cfg(test)]
mod test {
use super::*;
#[test]
fn inequality_test() {}
fn inequality_test() {
linear_vars![x, y, z, w];
let mut constraints = Vec::<LeqInequality>::new();
constraints.push(Expr::from(&x) + &y + &z + &w << 5.0.into());
constraints.push(Expr::from(&x) << 4.0.into());
constraints.push(Expr::from(&y) << 4.0.into());
constraints.push(Expr::from(&z) << 4.0.into());
constraints.push(Expr::from(&w) << 4.0.into());
let mut objective = LinearCombination::new();
objective.insert_variable(x.clone(), 1.0);
objective.insert_variable(y.clone(), 2.0);
objective.insert_variable(z.clone(), 3.0);
objective.insert_variable(w.clone(), 4.0);
println!("Maximize {}", objective);
println!("subject to");
for constraint in constraints.iter() {
println!("{}", constraint);
}
let sf = StandardFormLP::from_ineqs_and_objective(constraints, objective);
println!("");
println!("{}", sf);
let frp = sf.clone().to_feasible_region_problem();
println!("");
println!("{}", frp);
let len = frp.vars.len();
let mut fras =
AffineScaling::from_feasible_point(frp, DVector::from_element(len, 1.0)).unwrap();
println!("");
println!("{}", fras);
for _ in 0..10 {
println!("{:?}", fras.iterate());
println!("{}", fras.x);
}
}
}

8
src/clause.rs

@ -255,17 +255,17 @@ impl core::ops::IndexMut<usize> for ClauseList {
}
}
#[macro_export]
/// Defines a DSL which allows for the declaration of many variables at once. See the unit tests
/// for examples.
#[macro_export]
macro_rules! vars {
macro_rules! clause_vars {
[] => {};
[$varname:ident$(,)?] => {
let $varname = $crate::clause::Variable::from_name(stringify!($varname).into());
};
[$varname:ident, $($varnames:ident),+] => {
let $varname = $crate::clause::Variable::from_name(stringify!($varname).into());
vars![$($varnames),+]
clause_vars![$($varnames),+]
};
}
@ -321,7 +321,7 @@ mod test {
#[test]
fn clause_test() {
vars![x, y, z, w];
clause_vars![x, y, z, w];
println!(
"names: {} {} {} {}",
x.get_name(),

5
src/main.rs

@ -3,6 +3,9 @@
#![feature(specialization)]
#![feature(non_ascii_idents)]
#[macro_use]
/// Variables -- both boolean, and otherwise.
pub mod variable;
#[macro_use]
/// Useful for disjunctive clauses.
pub mod clause;
@ -11,7 +14,5 @@ pub mod affine_scaling;
/// I made this before I knew I didn't need it. :/
pub mod hash_wrapper;
pub mod messy_minsat;
/// Variables -- both boolean, and otherwise.
pub mod variable;
fn main() {}
Loading…
Cancel
Save