RsBundle  Changes On Branch mcf-trait

/*
 * Copyright (c) 2016-2020 Frank Fischer <frank-fischer@shadow-soft.de>
 * Copyright (c) 2016-2021 Frank Fischer <frank-fischer@shadow-soft.de>
 *
 * This program is free software: you can redistribute it and/or
 * modify it under the terms of the GNU General Public License as
 * published by the Free Software Foundation, either version 3 of the
 * License, or (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful, but

use env_logger;
use env_logger::fmt::Color;
use log::{info, Level};
use rustop::opts;
use std::io::Write;

use bundle::master::{Builder, FullMasterBuilder, MasterProblem, MinimalMasterBuilder};
use bundle::mcf::MMCFProblem;
use bundle::mcf::{self, MMCFProblem};
use bundle::solver::asyn::{GuessModelType, Parameters, Solver as AsyncSolver};
use bundle::solver::sync::Solver as SyncSolver;
use bundle::{terminator::StandardTerminator, weighter::HKWeighter};
use bundle::{DVector, Real};

use std::error::Error;
use std::result::Result;

        opt aggregated:bool, desc:"Use aggregated model";
        opt sync:bool, desc:"Use synchronized solver";
        opt bundle_size:usize = 0, name:"NUM", desc:"The bundle size";
        opt nearest_guess:bool, desc:"Use nearest guess model (default: cut-model)";
        opt descent_factor:Real = 0.5, name:"ρ", desc:"Descent parameter in (0,1)";
        opt imprecision_factor:Real = 0.9, name:"α", desc:"Imprecision parameter in (0,1)";
        opt l_guess_factor:Real = 0.1, name:"β", desc:"Lipschitz-guess increase factor in (0,1]";
        opt cpx:bool, desc:"Use Cplex as network flow solver";
        opt netspx:bool, desc:"Use network simplex implementation as network flow solver (default)";
        param file:String, desc:"Input file in MNetGen format";
    }
    .parse_or_exit();

    let filename = args.file;
    info!("Reading instance: {}", filename);

    if !args.minimal {
        let mut mmcf = if args.netspx || !args.cpx {
            info!("Use rs-graph network simplex solver");
        let mut mmcf = MMCFProblem::read_mnetgen(&filename)?;
            MMCFProblem::read_mnetgen_with_solver::<mcf::NetSpxSolver>(&filename)?
        } else {
            info!("Use Cplex network simplex solver");
            MMCFProblem::read_mnetgen_with_solver::<mcf::CplexSolver>(&filename)?
        };
        mmcf.set_separate_constraints(args.separate);
        mmcf.multimodel = true;

        let mut master = FullMasterBuilder::default();
        if args.aggregated {
            master.max_bundle_size(if args.bundle_size <= 1 { 50 } else { args.bundle_size });
            master.use_full_aggregation();

// Copyright (c) 2016 Frank Fischer <frank-fischer@shadow-soft.de>
// Copyright (c) 2016, 2021 Frank Fischer <frank-fischer@shadow-soft.de>
//
// This program is free software: you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation, either version 3 of the
// License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program.  If not, see  <http://www.gnu.org/licenses/>
//

//! Multi-commodity min-cost-flow subproblems.

mod solver;
pub use crate::mcf::solver::CplexSolver;
pub use crate::mcf::solver::NetSpxSolver;
pub use crate::mcf::solver::Solver;

mod problem;
pub use crate::mcf::problem::MMCFProblem;

// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// General Public License for more details.
//
// 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::mcf::{self, Solver};
use crate::problem::{
    FirstOrderProblem as ParallelProblem, ResultSender, UpdateSender, UpdateState as ParallelUpdateState,
};
use crate::{DVector, Minorant, Real};

use log::{debug, warn};
use num_traits::Float;

    ind: usize,
    val: Real,
}

/// A single MMCF subproblem, i.e. one network.
struct Subproblem {
    /// The (net, cost) pair.
    net: mcf::Solver,
    net: mcf::NetSpxSolver,
    c: DVector,
}

/// Constraint data of one subproblem.
struct SubData {
    lhs: Vec<Vec<Elem>>,
    /// The right-hand side ... might be empty.

impl MMCFProblem {
    pub fn num_subproblems(&self) -> usize {
        self.subs.len()
    }

    pub fn read_mnetgen(basename: &str) -> std::result::Result<MMCFProblem, MMCFReadError> {
        MMCFProblem::read_mnetgen_with_solver::<mcf::NetSpxSolver>(basename)
    }

    pub fn read_mnetgen_with_solver<S>(basename: &str) -> std::result::Result<MMCFProblem, MMCFReadError>
    where
        S: Solver + 'static,
    {
        use MMCFReadError::*;

        let mut buffer = String::new();
        {
            let mut f = File::open(&format!("{}.nod", basename))?;
            f.read_to_string(&mut buffer)?;
        }

                return Err(ExtraNodField);
            }
        }

        // read nodes
        let mut nets = Vec::with_capacity(ncom);
        for i in 0..ncom {
            nets.push(mcf::Solver::new(i, nnodes)?)
            nets.push(mcf::NetSpxSolver::new(i, nnodes)?)
        }
        {
            let mut f = File::open(&format!("{}.sup", basename))?;
            buffer.clear();
            f.read_to_string(&mut buffer)?;
        }
        for line in buffer.lines() {

// along with this program.  If not, see  <http://www.gnu.org/licenses/>
//

//! Solver for MCF-problems based on Cplex.

#![allow(unused_unsafe)]

use super::{Error, Result, Solver};
use crate::{DVector, Real};

use c_str_macro::c_str;
use cplex_sys as cpx;
use cplex_sys::trycpx;
use thiserror::Error;
use cplex_sys::CplexError;

use std::ffi::CString;
use std::ptr;
use std::result;

use std::os::raw::{c_char, c_double, c_int};

#[derive(Debug, Error)]
#[non_exhaustive]
pub enum Error {
    #[error("CPLEX error")]
    Solver(#[from] cpx::CplexError),
}

pub type Result<T> = result::Result<T, Error>;

pub struct Solver {
pub struct CplexSolver {
    env: *mut cpx::Env,
    net: *mut cpx::Net,
}

impl From<CplexError> for Error {
    fn from(err: CplexError) -> Error {
        Error::Unknown(Box::new(err))
    }
}

impl Drop for Solver {
impl Drop for CplexSolver {
    fn drop(&mut self) {
        unsafe {
            cpx::NETfreeprob(self.env, &mut self.net);
        }
        unsafe { cpx::closeCPLEX(&mut self.env) };
    }
}

impl Solver {
    pub fn new(nnodes: usize) -> Result<Solver> {
impl Solver for CplexSolver {
    fn new(_id: usize, nnodes: usize) -> Result<CplexSolver> {
        let mut status: c_int;
        let env;
        let mut net = std::ptr::null_mut();

        unsafe {
            trycpx!({
                let mut status = 0;

                    code: status,
                    msg: CString::from_raw(msg).to_string_lossy().into_owned(),
                }
                .into());
            }
        }

        Ok(Solver { env, net })
        Ok(CplexSolver { env, net })
    }

    pub fn num_nodes(&self) -> usize {
    fn num_nodes(&self) -> usize {
        unsafe { cpx::NETgetnumnodes(self.env, self.net) as usize }
    }

    pub fn num_arcs(&self) -> usize {
    fn num_arcs(&self) -> usize {
        unsafe { cpx::NETgetnumarcs(self.env, self.net) as usize }
    }

    pub fn set_balance(&mut self, node: usize, supply: Real) -> Result<()> {
    fn set_balance(&mut self, node: usize, supply: Real) -> Result<()> {
        let n = node as c_int;
        let s = supply as c_double;
        trycpx!(cpx::NETchgsupply(self.env, self.net, 1, &n, &s as *const c_double));
        Ok(())
    }

    pub fn set_objective(&mut self, obj: &DVector) -> Result<()> {
    fn set_objective(&mut self, obj: &DVector) -> Result<()> {
        let inds = (0..obj.len() as c_int).collect::<Vec<_>>();
        trycpx!(cpx::NETchgobj(
            self.env,
            self.net,
            obj.len() as c_int,
            inds.as_ptr(),
            obj.as_ptr()
        ));
        Ok(())
    }

    pub fn add_arc(&mut self, src: usize, snk: usize, cost: Real, cap: Real) -> Result<()> {
    fn add_arc(&mut self, src: usize, snk: usize, cost: Real, cap: Real) -> Result<()> {
        let f = src as c_int;
        let t = snk as c_int;
        let c = cost as c_double;
        let u = cap as c_double;
        let name = CString::new(format!("x{}#{}_{}", self.num_arcs() + 1, f + 1, t + 1)).unwrap();
        let cname = name.as_ptr();
        trycpx!(cpx::NETaddarcs(
            self.env,
            self.net,
            1,
            &f,
            &t,
            ptr::null(),
            &u,
            &c,
            &cname as *const *const c_char
        ));
        Ok(())
    }

    pub fn solve(&mut self) -> Result<()> {
    fn solve(&mut self) -> Result<()> {
        trycpx!(cpx::NETprimopt(self.env, self.net));
        Ok(())
    }

    pub fn objective(&self) -> Result<Real> {
    fn objective(&self) -> Result<Real> {
        let mut objval: c_double = 0.0;
        trycpx!(cpx::NETgetobjval(self.env, self.net, &mut objval as *mut c_double));
        Ok(objval)
    }

    pub fn get_solution(&self) -> Result<DVector> {
    fn get_solution(&self) -> Result<DVector> {
        let mut sol = dvec![0.0; self.num_arcs()];
        let mut stat: c_int = 0;
        let mut objval: c_double = 0.0;
        trycpx!(cpx::NETsolution(
            self.env,
            self.net,
            &mut stat as *mut c_int,
            &mut objval as *mut c_double,
            sol.as_mut_ptr(),
            ptr::null_mut(),
            ptr::null_mut(),
            ptr::null_mut()
        ));
        Ok(sol)
    }

    pub fn writelp(&self, filename: &str) -> Result<()> {
    fn writelp(&self, filename: &str) -> Result<()> {
        let fname = CString::new(filename).unwrap();
        trycpx!(cpx::NETwriteprob(self.env, self.net, fname.as_ptr(), ptr::null_mut()));
        Ok(())
    }
}

 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 *
 * 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 thiserror::Error;

use crate::{DVector, Real};

#[derive(Debug, Error)]
#[non_exhaustive]
pub enum Error {
    #[error("MCF problem is infeasible")]
    Infeasible,
    #[error("MCF problem is unbounded")]
    Unbounded,
    #[error("MCF problem has not been solved")]
    Unsolved,
    #[error("The mcf network cannot be modified after is has been solved")]
    NotInBuildState,
    #[error("An unknown error has been raised by the backend solver")]
    Unknown(Box<dyn std::error::Error + Send + 'static>),
}

pub type Result<T> = std::result::Result<T, Error>;

/// Generic trait for MCF solvers.
pub trait Solver {
    /// Create a new solver.
    fn new(id: usize, nnodes: usize) -> Result<Self>
    where
        Self: Sized;

    /// Return the number of nodes.
    fn num_nodes(&self) -> usize;

    /// Return the number of arcs.
    fn num_arcs(&self) -> usize;

    /// Set the supply (balance) of a node.
    fn set_balance(&mut self, node: usize, supply: Real) -> Result<()>;

    /// Set the objective coefficients for each arg.
    fn set_objective(&mut self, obj: &DVector) -> Result<()>;

    /// Add an arc to the network.
    fn add_arc(&mut self, src: usize, snk: usize, cost: Real, cap: Real) -> Result<()>;

    /// Solve the network-flow problem.
    fn solve(&mut self) -> Result<()>;

    /// Return the optimal objective value.
    fn objective(&self) -> Result<Real>;

    /// Return the optimal flow on each arc.
    fn get_solution(&self) -> Result<DVector>;

    /// Write an LP representation of the problem to the given file.
    fn writelp(&self, _filename: &str) -> Result<()> {
        unimplemented!("writelp is not implemented for this solver")
    }
}

pub mod cpx;
pub mod netspx;
pub use cpx::CplexSolver;

/// The default solver type.
pub type Solver = netspx::Solver;
pub type Error = netspx::Error;
pub mod netspx;
pub use netspx::NetSpxSolver;

use num_traits::{Float, Zero};
use rs_graph::{
    mcf::{MinCostFlow, NetworkSimplex, SolutionState},
    traits::{GraphSize, IndexGraph, Indexable},
    Buildable, Builder, Net,
};
use thiserror::Error;

#[derive(Debug, Error)]
#[non_exhaustive]
pub enum Error {
    #[error("MCF problem is infeasible")]
    Infeasible,
    #[error("MCF problem is unbounded")]
    Unbounded,
    #[error("MCF problem has not been solved")]
    Unsolved,
    #[error("The mcf network cannot be modified after is has been solved")]
    NotInBuildState,
}

use super::{Error, Result, Solver};
pub type Result<T> = std::result::Result<T, Error>;

type Mcf = NetworkSimplex<Rc<Net>, Real>;

enum SolverState {
    Building {
        builder: <Net as Buildable>::Builder,
        costs: Vec<Real>,
        capacities: Vec<Real>,
        balances: Vec<Real>,
    },
    Solving {
        net: Rc<Net>,
        mcf: Box<Mcf>,
    },
}

use SolverState::*;

pub struct Solver {
pub struct NetSpxSolver {
    _id: usize,
    state: SolverState,
}

impl Solver {
    pub fn new(id: usize, nnodes: usize) -> Result<Solver> {
impl Solver for NetSpxSolver {
    fn new(id: usize, nnodes: usize) -> Result<NetSpxSolver> {
        let mut builder = Net::new_builder();
        for _ in 0..nnodes {
            builder.add_node();
        }
        Ok(Solver {
        Ok(NetSpxSolver {
            _id: id,
            state: Building {
                builder,
                costs: vec![],
                capacities: vec![],
                balances: vec![Real::zero(); nnodes],
            },
        })
    }

    pub fn num_nodes(&self) -> usize {
    fn num_nodes(&self) -> usize {
        match &self.state {
            Building { ref builder, .. } => builder.num_nodes(),
            Solving { ref net, .. } => net.num_nodes(),
        }
    }

    pub fn num_arcs(&self) -> usize {
    fn num_arcs(&self) -> usize {
        match &self.state {
            Building { ref builder, .. } => builder.num_edges(),
            Solving { ref net, .. } => net.num_edges(),
        }
    }

    pub fn set_balance(&mut self, node: usize, supply: Real) -> Result<()> {
    fn set_balance(&mut self, node: usize, supply: Real) -> Result<()> {
        match &mut self.state {
            Building { balances, .. } => balances[node] = supply,
            Solving { net, mcf } => mcf.set_balance(net.id2node(node), supply),
        }
        Ok(())
    }

    pub fn set_objective(&mut self, obj: &DVector) -> Result<()> {
    fn set_objective(&mut self, obj: &DVector) -> Result<()> {
        let (net, mcf) = self.state.solving()?;
        for (e, &c) in net.edges().zip(obj.iter()) {
            mcf.set_cost(e, c);
        }
        Ok(())
    }

    pub fn add_arc(&mut self, src: usize, snk: usize, cost: Real, cap: Real) -> Result<()> {
    fn add_arc(&mut self, src: usize, snk: usize, cost: Real, cap: Real) -> Result<()> {
        // The graph must be in build-state, otherwise this won't work.
        match &mut self.state {
            Building {
                builder,
                costs,
                capacities,
                ..
            } => {
                builder.add_edge(builder.id2node(src), builder.id2node(snk));
                costs.push(cost);
                capacities.push(cap);
                Ok(())
            }
            _ => Err(Error::NotInBuildState),
        }
    }

    pub fn solve(&mut self) -> Result<()> {
    fn solve(&mut self) -> Result<()> {
        let (_, mcf) = self.state.solving()?;
        let state = mcf.solve();
        match state {
            SolutionState::Optimal => Ok(()),
            SolutionState::Infeasible => Err(Error::Infeasible),
            SolutionState::Unbounded => Err(Error::Unbounded),
            _ => unreachable!(),
        }
    }

    pub fn objective(&self) -> Result<Real> {
    fn objective(&self) -> Result<Real> {
        if let Solving { mcf, .. } = &self.state {
            Ok(mcf.value())
        } else {
            Err(Error::Unsolved)
        }
    }

    pub fn get_solution(&self) -> Result<DVector> {
    fn get_solution(&self) -> Result<DVector> {
        use SolutionState::*;
        if let Solving { net, mcf } = &self.state {
            match mcf.solution_state() {
                Optimal => Ok(net.edges().map(|e| mcf.flow(e)).collect()),
                Infeasible => Err(Error::Infeasible),
                Unbounded => Err(Error::Unbounded),
                Unknown => Err(Error::Unsolved),
            }
        } else {
            Err(Error::Unsolved)
        }
    }

    pub fn writelp(&self, _filename: &str) -> Result<()> {
        // let fname = CString::new(filename).unwrap();
        // trycpx!(cpx::NETwriteprob(self.env, self.net, fname.as_ptr(), ptr::null_mut()));
        // Ok(())
        todo!()
    }
}

impl SolverState {
    /// Transit into solving state.
    fn solving(&mut self) -> Result<(&Net, &mut Mcf)> {
        replace_with::replace_with_or_abort(self, |state| {
            if let Building {