RsBundle  Diff

// You should have received a copy of the GNU General Public License
// along with this program.  If not, see  <http://www.gnu.org/licenses/>
//

use crate::mcf;
use crate::{Aggregatable, DVector, FirstOrderProblem, Minorant, Real, SimpleEvaluation};

use itertools::izip;
use log::debug;

use std::f64::INFINITY;
use std::fmt;
use std::fs::File;
use std::io::Read;
use std::result;

    ind: usize,
    val: Real,
}

/// A single MMCF subproblem, i.e. one network.
struct Subproblem {
    net: mcf::Solver,
    lhs: Vec<Vec<Elem>>,
    cbase: DVector,
    c: DVector,
}

pub struct MMCFProblem {
    pub multimodel: bool,

    subs: Vec<Subproblem>,
    lhs: Vec<Vec<Vec<Elem>>>,
    rhs: DVector,
    rhsval: Real,
}

impl MMCFProblem {
    pub fn read_mnetgen(basename: &str) -> std::result::Result<MMCFProblem, MMCFReadError> {
        let mut buffer = String::new();

            let mut data = line.split_whitespace();
            let mt = data.next().unwrap().parse::<usize>()? - 1;
            let cap = data.next().unwrap().parse::<Real>()?;
            rhs[mt] = cap;
        }

        // set lhs
        let mut lhs = vec![vec![vec![]; ncom]; ncaps];
        let mut lhs = vec![vec![vec![]; ncaps]; ncom];
        for i in 0..ncaps {
            for fidx in 0..ncom {
                lhs[i][fidx] = lhsidx[i][fidx].iter().map(|&j| Elem { ind: j, val: 1.0 }).collect();
        for fidx in 0..ncom {
            for i in 0..ncaps {
                lhs[fidx][i] = lhsidx[i][fidx].iter().map(|&j| Elem { ind: j, val: 1.0 }).collect();
            }
        }

        let subproblems = nets
            .into_iter()
            .zip(cbase)
            .map(|(net, cbase)| Subproblem { net, cbase, c: dvec![] })
        let subproblems = izip!(nets, cbase, lhs)
            .map(|(net, cbase, lhs)| Subproblem {
                net,
                cbase,
                c: dvec![],
                lhs,
            })
            .collect();

        Ok(MMCFProblem {
            multimodel: false,
            subs: subproblems,
            lhs,
            rhs,
            rhsval: 0.0,
        })
    }

    /// Compute costs for a primal solution.
    pub fn get_primal_costs(&self, fidx: usize, primals: &[DVector]) -> Real {

    type Err = Error;

    type Primal = Vec<DVector>;

    type EvalResult = SimpleEvaluation<Vec<DVector>>;

    fn num_variables(&self) -> usize {
        self.lhs.len()
        self.rhs.len()
    }

    fn lower_bounds(&self) -> Option<Vec<Real>> {
        Some(vec![0.0; self.lhs.len()])
        Some(vec![0.0; self.rhs.len()])
    }

    fn upper_bounds(&self) -> Option<Vec<Real>> {
        None
    }

    fn num_subproblems(&self) -> usize {

            sub.c.clear();
            sub.c.extend(sub.cbase.iter());
        }

        for (i, &y) in y.iter().enumerate().filter(|&(i, &y)| y != 0.0) {
            self.rhsval += self.rhs[i] * y;
            for (fidx, sub) in self.subs.iter_mut().enumerate() {
                for elem in &self.lhs[i][fidx] {
                for elem in &sub.lhs[i] {
                    sub.c[elem.ind] += y * elem.val;
                }
            }
        }

        debug!("y={:?}", y);
        for (i, sub) in self.subs.iter_mut().enumerate() {

                objective = self.rhsval - self.subs[fidx].net.objective()?;
            } else {
                subg = dvec![0.0; self.rhs.len()];
                objective = -self.subs[fidx].net.objective()?;
            }

            let sol = self.subs[fidx].net.get_solution()?;
            for (i, lhs) in self.lhs.iter().enumerate() {
                subg[i] -= lhs[fidx].iter().map(|elem| elem.val * sol[elem.ind]).sum::<Real>();
            for (i, lhs) in self.subs[fidx].lhs.iter().enumerate() {
                subg[i] -= lhs.iter().map(|elem| elem.val * sol[elem.ind]).sum::<Real>();
            }

            Ok(SimpleEvaluation {
                objective,
                minorants: vec![(
                    Minorant {
                        constant: objective,

            let mut sols = Vec::with_capacity(self.subs.len());
            for sub in &self.subs {
                objective -= sub.net.objective()?;
                sols.push(sub.net.get_solution()?);
            }

            let mut subg = self.rhs.clone();
            for (i, lhs) in self.lhs.iter().enumerate() {
                for (fidx, flhs) in lhs.iter().enumerate() {
                    subg[i] -= flhs.iter().map(|elem| elem.val * sols[fidx][elem.ind]).sum::<Real>();
            for (fidx, sub) in self.subs.iter().enumerate() {
                for (i, lhs) in sub.lhs.iter().enumerate() {
                    subg[i] -= lhs.iter().map(|elem| elem.val * sols[fidx][elem.ind]).sum::<Real>();
                }
            }

            Ok(SimpleEvaluation {
                objective,
                minorants: vec![(
                    Minorant {