Log encoding for integer decision variables (#260)

- Introduce `Instance::log_encoding` for encoding integer variable into
binary variables
- This method does not substitute yielded polynomial into objective and
constraints. It is done by `Instance::substitute`

rust/ommx/src/convert/instance.rs

@@ -1,11 +1,11 @@
 use crate::{
     sorted_ids::{BinaryIdPair, BinaryIds},
     v1::{
-        decision_variable::Kind, instance::Sense, Equality, Function, Instance, Parameter,
-        ParametricInstance, RemovedConstraint,
+        decision_variable::Kind, instance::Sense, DecisionVariable, Equality, Function, Instance,
+        Linear, Parameter, ParametricInstance, RemovedConstraint,
     },
 };
-use anyhow::{bail, Context, Result};
+use anyhow::{bail, ensure, Context, Result};
 use approx::AbsDiffEq;
 use maplit::hashmap;
 use num::Zero;
@@ -296,6 +296,118 @@ impl Instance {
         }
         Ok((quad, constant))
     }
+
+    /// Encode an integer decision variable into binary decision variables.
+    ///
+    /// Note that this method does not substitute the yielded binary representation into the objective and constraints.
+    /// Call [`Instance::substitute`] with the returned [`Linear`] representation.
+    ///
+    /// Mutability
+    /// ----------
+    /// - This adds new binary decision variables introduced for binary encoding to the instance.
+    ///
+    /// Errors
+    /// ------
+    /// Returns [anyhow::Error] in the following cases:
+    ///
+    /// - The given decision variable ID is not found
+    /// - The specified decision variable is not an integer type.
+    /// - The bound of the decision variable is not set or not finite.
+    ///
+    pub fn log_encoding(&mut self, decision_variable_id: u64) -> Result<Linear> {
+        let v = self
+            .decision_variables
+            .iter()
+            .find(|dv| dv.id == decision_variable_id)
+            .with_context(|| format!("Decision variable ID {} not found", decision_variable_id))?;
+        if v.kind() != Kind::Integer {
+            bail!(
+                "The decision variable is not an integer type: ID={}",
+                decision_variable_id
+            );
+        }
+
+        let bound = v.bound.as_ref().with_context(|| {
+            format!(
+                "Bound must be set and finite for log-encoding: ID={}",
+                decision_variable_id
+            )
+        })?;
+
+        // Bound of integer may be non-integer value
+        let upper = bound.upper.floor();
+        let lower = bound.lower.ceil();
+        let u_l = upper - lower;
+        ensure!(
+            u_l >= 0.0,
+            "No feasible integer found in the bound: ID={}, lower={}, upper={}",
+            decision_variable_id,
+            bound.lower,
+            bound.upper
+        );
+
+        // There is only one feasible integer, and no need to encode
+        if u_l == 0.0 {
+            return Ok(Linear::from(lower));
+        }
+
+        // Log-encoding
+        let n = (u_l + 1.0).log2().ceil() as usize;
+        let id_base = self
+            .defined_ids()
+            .last()
+            .map(|id| id + 1)
+            .expect("At least one decision variable here");
+
+        let mut terms = Vec::new();
+        for i in 0..n {
+            let id = id_base + i as u64;
+            terms.push((
+                id,
+                if i == n - 1 {
+                    u_l - 2.0f64.powi(i as i32) + 1.0
+                } else {
+                    2.0f64.powi(i as i32)
+                },
+            ));
+            self.decision_variables.push(DecisionVariable {
+                id,
+                name: Some("ommx.log_encoding".to_string()),
+                subscripts: vec![decision_variable_id as i64, i as i64],
+                kind: Kind::Binary as i32,
+                bound: Some(crate::v1::Bound {
+                    lower: 0.0,
+                    upper: 1.0,
+                }),
+                ..Default::default()
+            });
+        }
+        Ok(Linear::new(terms.into_iter(), lower))
+    }
+
+    /// Substitute dependent decision variables with given [Function]s.
+    pub fn substitute(&mut self, replacement: HashMap<u64, Function>) -> Result<()> {
+        if let Some(obj) = self.objective.as_mut() {
+            *obj = obj.substitute(&replacement)?;
+        }
+        for c in &mut self.constraints {
+            if let Some(f) = c.function.as_mut() {
+                *f = f.substitute(&replacement)?;
+            }
+        }
+        for c in &mut self.removed_constraints {
+            if let Some(c) = &mut c.constraint {
+                if let Some(f) = c.function.as_mut() {
+                    *f = f.substitute(&replacement)?;
+                }
+            }
+        }
+        for (_id, f) in self.decision_variable_dependency.iter_mut() {
+            *f = f.substitute(&replacement)?;
+        }
+        self.decision_variable_dependency.extend(replacement);
+        Ok(())
+    }
 }
 
 /// Compare two instances as mathematical programming problems. This does not compare the metadata.
@@ -358,7 +470,11 @@ impl AbsDiffEq for Instance {
 #[cfg(test)]
 mod tests {
     use super::*;
-    use crate::{random::InstanceParameters, v1::Parameters};
+    use crate::{
+        random::InstanceParameters,
+        v1::{Parameters, State},
+        Evaluate,
+    };
     use proptest::prelude::*;
 
     proptest! {
@@ -454,5 +570,63 @@ mod tests {
                 prop_assert!(c.abs() > f64::EPSILON);
             }
         }
+
+        #[test]
+        fn log_encoding((lower, upper) in (-10.0_f64..10.0, -10.0_f64..10.0)
+            .prop_filter("At least one integer", |(lower, upper)| lower.ceil() <= upper.floor())
+        ) {
+            let mut instance = Instance::default();
+            instance.decision_variables.push(DecisionVariable {
+                id: 0,
+                name: Some("x".to_string()),
+                kind: Kind::Integer as i32,
+                bound: Some(crate::v1::Bound { lower, upper }),
+                ..Default::default()
+            });
+            let encoded = instance.log_encoding(0).unwrap();
+
+            // Test the ID of yielded decision variables are not duplicated
+            instance.validate().unwrap();
+
+            // Get decision variables introduced for log-encoding
+            let aux_bits = instance
+                .decision_variables
+                .iter()
+                .filter_map(|dv| {
+                    if dv.name == Some("ommx.log_encoding".to_string()) && dv.subscripts[0] == 0 {
+                        Some(dv.id)
+                    } else {
+                        None
+                    }
+                })
+                .collect::<Vec<_>>();
+
+            if lower.ceil() == upper.floor() {
+                // No need to encode
+                prop_assert_eq!(encoded.as_constant().unwrap(), lower.ceil());
+                prop_assert_eq!(aux_bits.len(), 0);
+                return Ok(());
+            }
+
+            let state = State { entries: aux_bits.iter().map(|&id| (id, 0.0)).collect::<HashMap<_, _>>() };
+            let (lower_evaluated, _) = encoded.evaluate(&state).unwrap();
+            prop_assert_eq!(lower_evaluated, lower.ceil());
+
+            let state = State { entries: aux_bits.iter().map(|&id| (id, 1.0)).collect::<HashMap<_, _>>() };
+            let (upper_evaluated, _) = encoded.evaluate(&state).unwrap();
+            prop_assert_eq!(upper_evaluated, upper.floor());
+        }
+
+        /// Compare the result of partial_evaluate and substitute with `Function::Constant`.
+        #[test]
+        fn substitute_fixed_value(instance in Instance::arbitrary(), value in -3.0..3.0) {
+            for id in instance.defined_ids() {
+                let mut partially_evaluated = instance.clone();
+                partially_evaluated.partial_evaluate(&State { entries: hashmap! { id => value } }).unwrap();
+                let mut substituted = instance.clone();
+                substituted.substitute(hashmap! { id => Function::from(value) }).unwrap();
+                prop_assert!(partially_evaluated.abs_diff_eq(&substituted, 1e-10));
+            }
+        }
     }
 }