RsBundle  Changes On Branch minorant-trait

Many hyperlinks are disabled.
Use anonymous login to enable hyperlinks.

Changes In Branch minorant-trait Excluding Merge-Ins

This is equivalent to a diff from 98774dd6e6 to 76f1c7f25f

2020-07-20
18:11
Merge minorant-trait check-in: fe82c5787c user: fifr tags: trunk
18:03
master: allow only adding of new variables (no moving) Closed-Leaf check-in: 76f1c7f25f user: fifr tags: minorant-trait
17:43
minorant: add a debug assertion to `move_center` check-in: 3a337c0eec user: fifr tags: minorant-trait
2020-07-19
20:47
Make `Minorant` a trait check-in: 095240b3e2 user: fifr tags: minorant-trait
2020-07-18
18:49
Merge trunk check-in: e12bd4bd5d user: fifr tags: async
14:59
Merge minorant-primal check-in: 98774dd6e6 user: fifr tags: trunk
14:56
master: remove callback from `compress` method Closed-Leaf check-in: baa49a9473 user: fifr tags: minorant-primal
2020-07-15
19:52
Update ordered-float to version 2.0 check-in: f53f2d0005 user: fifr tags: trunk

Changes to examples/cflp.rs.
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47


48
49
50
51
52
53
54

use env_logger::{self, fmt::Color};
use ordered_float::NotNan;
use threadpool::ThreadPool;

use bundle::problem::{FirstOrderProblem as ParallelProblem, ResultSender};
use bundle::solver::sync::{DefaultSolver, NoBundleSolver};
use bundle::{dvec, DVector, Minorant, Real};

const Nfac: usize = 3;
const Ncus: usize = 5;
const F: [Real; Nfac] = [1000.0, 1000.0, 1000.0];
const CAP: [Real; Nfac] = [500.0, 500.0, 500.0];
const C: [[Real; Ncus]; Nfac] = [
    [4.0, 5.0, 6.0, 8.0, 10.0], //
    [6.0, 4.0, 3.0, 5.0, 8.0],  //
    [9.0, 7.0, 4.0, 3.0, 4.0],  //
];
const DEMAND: [Real; 5] = [80.0, 270.0, 250.0, 160.0, 180.0];



#[derive(Debug)]
enum EvalError {
    Customer(usize),
}

impl std::fmt::Display for EvalError {







|











>
>







29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56

use env_logger::{self, fmt::Color};
use ordered_float::NotNan;
use threadpool::ThreadPool;

use bundle::problem::{FirstOrderProblem as ParallelProblem, ResultSender};
use bundle::solver::sync::{DefaultSolver, NoBundleSolver};
use bundle::{dvec, DVector, Real};

const Nfac: usize = 3;
const Ncus: usize = 5;
const F: [Real; Nfac] = [1000.0, 1000.0, 1000.0];
const CAP: [Real; Nfac] = [500.0, 500.0, 500.0];
const C: [[Real; Ncus]; Nfac] = [
    [4.0, 5.0, 6.0, 8.0, 10.0], //
    [6.0, 4.0, 3.0, 5.0, 8.0],  //
    [9.0, 7.0, 4.0, 3.0, 4.0],  //
];
const DEMAND: [Real; 5] = [80.0, 270.0, 250.0, 160.0, 180.0];

type Minorant = (Real, DVector, DVector);

#[derive(Debug)]
enum EvalError {
    Customer(usize),
}

impl std::fmt::Display for EvalError {
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132

impl CFLProblem {
    fn new() -> CFLProblem {
        CFLProblem { pool: None }
    }
}

fn solve_facility_problem(j: usize, lambda: &[Real]) -> Result<(Real, Minorant<DVector>), EvalError> {
    // facility subproblem max\{-F[j] + CAP[j] * lambda[j] * y[j] | y_j \in \{0,1\}\}
    let objective;
    let mut subg = dvec![0.0; lambda.len()];
    let primal;
    if -F[j] + CAP[j] * lambda[j] > 0.0 {
        objective = -F[j] + CAP[j] * lambda[j];
        subg[j] = CAP[j];
        primal = dvec![1.0; 1];
    } else {
        objective = 0.0;
        primal = dvec![0.0; 1];
    }

    Ok((
        objective,
        Minorant {
            constant: objective,
            linear: subg,
            primal,
        },
    ))
}

fn solve_customer_problem(i: usize, lambda: &[Real]) -> Result<(Real, Minorant<DVector>), EvalError> {
    // customer subproblem
    let opt = (0..Nfac)
        .map(|j| (j, -DEMAND[i] * C[j][i] - DEMAND[i] * lambda[j]))
        .flat_map(|x| NotNan::new(x.1).map(|y| (x.0, y)))
        .max_by_key(|x| x.1)
        .ok_or(EvalError::Customer(i))?;

    let objective = opt.1.into_inner();
    let mut subg = dvec![0.0; lambda.len()];
    subg[opt.0] = -DEMAND[i];
    let mut primal = dvec![0.0; Nfac];
    primal[opt.0] = 1.0;

    Ok((
        objective,
        Minorant {
            constant: objective,
            linear: subg,
            primal,
        },
    ))
}

impl ParallelProblem for CFLProblem {
    type Err = EvalError;

    type Primal = DVector;

    fn num_variables(&self) -> usize {
        Nfac
    }

    fn lower_bounds(&self) -> Option<Vec<Real>> {
        Some(vec![0.0; ParallelProblem::num_variables(self)])







|













<
|
<
<
<
<
<
<


|













<
|
<
<
<
<
<
<





|







69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89

90






91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106

107






108
109
110
111
112
113
114
115
116
117
118
119
120

impl CFLProblem {
    fn new() -> CFLProblem {
        CFLProblem { pool: None }
    }
}

fn solve_facility_problem(j: usize, lambda: &[Real]) -> Result<(Real, Minorant), EvalError> {
    // facility subproblem max\{-F[j] + CAP[j] * lambda[j] * y[j] | y_j \in \{0,1\}\}
    let objective;
    let mut subg = dvec![0.0; lambda.len()];
    let primal;
    if -F[j] + CAP[j] * lambda[j] > 0.0 {
        objective = -F[j] + CAP[j] * lambda[j];
        subg[j] = CAP[j];
        primal = dvec![1.0; 1];
    } else {
        objective = 0.0;
        primal = dvec![0.0; 1];
    }


    Ok((objective, (objective, subg, primal)))






}

fn solve_customer_problem(i: usize, lambda: &[Real]) -> Result<(Real, Minorant), EvalError> {
    // customer subproblem
    let opt = (0..Nfac)
        .map(|j| (j, -DEMAND[i] * C[j][i] - DEMAND[i] * lambda[j]))
        .flat_map(|x| NotNan::new(x.1).map(|y| (x.0, y)))
        .max_by_key(|x| x.1)
        .ok_or(EvalError::Customer(i))?;

    let objective = opt.1.into_inner();
    let mut subg = dvec![0.0; lambda.len()];
    subg[opt.0] = -DEMAND[i];
    let mut primal = dvec![0.0; Nfac];
    primal[opt.0] = 1.0;


    Ok((objective, (objective, subg, primal)))






}

impl ParallelProblem for CFLProblem {
    type Err = EvalError;

    type Minorant = Minorant;

    fn num_variables(&self) -> usize {
        Nfac
    }

    fn lower_bounds(&self) -> Option<Vec<Real>> {
        Some(vec![0.0; ParallelProblem::num_variables(self)])
Changes to examples/mmcf.rs.
21
22
23
24
25
26
27
28
29

30
31
32


33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
use rustop::opts;
use std::io::Write;

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


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



fn solve_parallel<M>(master: M, mmcf: MMCFProblem) -> Result<(), Box<dyn Error>>
where
    M: Builder<DVector>,
    M::MasterProblem: MasterProblem<DVector, MinorantIndex = usize>,
{
    let mut slv = AsyncSolver::<_, StandardTerminator, HKWeighter, M>::with_master(mmcf, master);
    slv.weighter.set_weight_bounds(1e-1, 100.0);
    slv.terminator.termination_precision = 1e-6;
    slv.solve()?;

    let costs: f64 = (0..slv.problem().num_subproblems())
        .map(|i| {
            let aggr_primals = slv.aggregated_primal(i).unwrap();
            slv.problem().get_primal_costs(i, &[aggr_primals])
        })
        .sum();
    info!("Primal costs: {}", costs);
    Ok(())
}

fn solve_sync<M>(master: M, mmcf: MMCFProblem) -> Result<(), Box<dyn Error>>
where
    M: Builder<DVector>,
    M::MasterProblem: MasterProblem<DVector, MinorantIndex = usize>,
{
    let mut slv = SyncSolver::<_, StandardTerminator, HKWeighter, M>::with_master(mmcf, master);
    slv.weighter.set_weight_bounds(1e-1, 100.0);
    slv.terminator.termination_precision = 1e-6;
    slv.solve()?;

    let costs: f64 = (0..slv.problem().num_subproblems())







<

>



>
>



|
|


















|
|







21
22
23
24
25
26
27

28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
use rustop::opts;
use std::io::Write;

use bundle::master::{Builder, FullMasterBuilder, MasterProblem, MinimalMasterBuilder};
use bundle::mcf::MMCFProblem;
use bundle::solver::asyn::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;

type Minorant = (Real, DVector, DVector);

fn solve_parallel<M>(master: M, mmcf: MMCFProblem) -> Result<(), Box<dyn Error>>
where
    M: Builder<Minorant>,
    M::MasterProblem: MasterProblem<Minorant, MinorantIndex = usize>,
{
    let mut slv = AsyncSolver::<_, StandardTerminator, HKWeighter, M>::with_master(mmcf, master);
    slv.weighter.set_weight_bounds(1e-1, 100.0);
    slv.terminator.termination_precision = 1e-6;
    slv.solve()?;

    let costs: f64 = (0..slv.problem().num_subproblems())
        .map(|i| {
            let aggr_primals = slv.aggregated_primal(i).unwrap();
            slv.problem().get_primal_costs(i, &[aggr_primals])
        })
        .sum();
    info!("Primal costs: {}", costs);
    Ok(())
}

fn solve_sync<M>(master: M, mmcf: MMCFProblem) -> Result<(), Box<dyn Error>>
where
    M: Builder<Minorant>,
    M::MasterProblem: MasterProblem<Minorant, MinorantIndex = usize>,
{
    let mut slv = SyncSolver::<_, StandardTerminator, HKWeighter, M>::with_master(mmcf, master);
    slv.weighter.set_weight_bounds(1e-1, 100.0);
    slv.terminator.termination_precision = 1e-6;
    slv.solve()?;

    let costs: f64 = (0..slv.problem().num_subproblems())
Changes to examples/quadratic.rs.
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52

53
54
55
56
57
58
59
use std::error::Error;
use std::io::Write;
use std::sync::Arc;
use std::thread;

use bundle::problem::{FirstOrderProblem as ParallelProblem, ResultSender};
use bundle::solver::sync::{DefaultSolver, NoBundleSolver};
use bundle::{DVector, Minorant, Real};

#[derive(Clone)]
struct QuadraticProblem {
    a: [[Real; 2]; 2],
    b: [Real; 2],
    c: Real,
}

impl QuadraticProblem {
    fn new() -> QuadraticProblem {
        QuadraticProblem {
            a: [[5.0, 1.0], [1.0, 4.0]],
            b: [-12.0, -10.0],
            c: 3.0,
        }
    }
}

impl ParallelProblem for QuadraticProblem {
    type Err = Box<dyn Error + Send>;
    type Primal = ();


    fn num_variables(&self) -> usize {
        2
    }

    fn num_subproblems(&self) -> usize {
        1







|




















|
>







24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
use std::error::Error;
use std::io::Write;
use std::sync::Arc;
use std::thread;

use bundle::problem::{FirstOrderProblem as ParallelProblem, ResultSender};
use bundle::solver::sync::{DefaultSolver, NoBundleSolver};
use bundle::{DVector, Real};

#[derive(Clone)]
struct QuadraticProblem {
    a: [[Real; 2]; 2],
    b: [Real; 2],
    c: Real,
}

impl QuadraticProblem {
    fn new() -> QuadraticProblem {
        QuadraticProblem {
            a: [[5.0, 1.0], [1.0, 4.0]],
            b: [-12.0, -10.0],
            c: 3.0,
        }
    }
}

impl ParallelProblem for QuadraticProblem {
    type Err = Box<dyn Error + Send>;

    type Minorant = (Real, DVector);

    fn num_variables(&self) -> usize {
        2
    }

    fn num_subproblems(&self) -> usize {
        1
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
            }

            debug!("Evaluation at {:?}", x);
            debug!("  objective={}", objective);
            debug!("  subgradient={}", g);

            tx.objective(objective).unwrap();
            tx.minorant(Minorant {
                constant: objective,
                linear: g,
                primal: (),
            })
            .unwrap();
        });
        Ok(())
    }
}

fn main() -> Result<(), Box<dyn std::error::Error>> {
    better_panic::install();







|
<
<
<
<
<







82
83
84
85
86
87
88
89





90
91
92
93
94
95
96
            }

            debug!("Evaluation at {:?}", x);
            debug!("  objective={}", objective);
            debug!("  subgradient={}", g);

            tx.objective(objective).unwrap();
            tx.minorant((objective, g)).unwrap();





        });
        Ok(())
    }
}

fn main() -> Result<(), Box<dyn std::error::Error>> {
    better_panic::install();
Changes to src/data/aggregatable.rs.
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50

//! Objects that can be combined linearly.

use super::Real;
use std::borrow::Borrow;

/// An aggregatable object.
pub trait Aggregatable: Default {
    /// Return a scaled version of `other`, i.e. `alpha * other`.
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>;

    /// Add a scaled version of `other` to `self`.
    ///
    /// This sets `self = self + alpha * other`.
    fn add_scaled<A>(&mut self, alpha: Real, other: A)
    where
        A: Borrow<Self>;

    /// Return $\sum\_{i=1}\^n alpha_i m_i$.
    ///
    /// If `aggregates` is empty return the default value.
    fn combine<I, A>(aggregates: I) -> Self
    where
        I: IntoIterator<Item = (Real, A)>,
        A: Borrow<Self>,
    {
        let mut it = aggregates.into_iter();
        let mut x;
        if let Some((alpha, y)) = it.next() {
            x = Self::new_scaled(alpha, y);
        } else {
            return Self::default();







|



|






|







|







17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50

//! Objects that can be combined linearly.

use super::Real;
use std::borrow::Borrow;

/// An aggregatable object.
pub trait Aggregatable<T = Self>: Default {
    /// Return a scaled version of `other`, i.e. `alpha * other`.
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<T>;

    /// Add a scaled version of `other` to `self`.
    ///
    /// This sets `self = self + alpha * other`.
    fn add_scaled<A>(&mut self, alpha: Real, other: A)
    where
        A: Borrow<T>;

    /// Return $\sum\_{i=1}\^n alpha_i m_i$.
    ///
    /// If `aggregates` is empty return the default value.
    fn combine<I, A>(aggregates: I) -> Self
    where
        I: IntoIterator<Item = (Real, A)>,
        A: Borrow<T>,
    {
        let mut it = aggregates.into_iter();
        let mut x;
        if let Some((alpha, y)) = it.next() {
            x = Self::new_scaled(alpha, y);
        } else {
            return Self::default();
68
69
70
71
72
73
74






















































75
76
77
78
79
80
81

    fn add_scaled<A>(&mut self, _alpha: Real, _other: A)
    where
        A: Borrow<Self>,
    {
    }
}























































/// Implement for scalar values.
impl Aggregatable for Real {
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>,
    {







>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>







68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135

    fn add_scaled<A>(&mut self, _alpha: Real, _other: A)
    where
        A: Borrow<Self>,
    {
    }
}

/// Implement for pairs.
impl<T1, T2> Aggregatable for (T1, T2)
where
    T1: Aggregatable,
    T2: Aggregatable,
{
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>,
    {
        let x = other.borrow();
        (T1::new_scaled(alpha, &x.0), T2::new_scaled(alpha, &x.1))
    }

    fn add_scaled<A>(&mut self, alpha: Real, other: A)
    where
        A: Borrow<Self>,
    {
        let x = other.borrow();
        self.0.add_scaled(alpha, &x.0);
        self.1.add_scaled(alpha, &x.1);
    }
}

/// Implement for triples.
impl<T1, T2, T3> Aggregatable for (T1, T2, T3)
where
    T1: Aggregatable,
    T2: Aggregatable,
    T3: Aggregatable,
{
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>,
    {
        let x = other.borrow();
        (
            T1::new_scaled(alpha, &x.0),
            T2::new_scaled(alpha, &x.1),
            T3::new_scaled(alpha, &x.2),
        )
    }

    fn add_scaled<A>(&mut self, alpha: Real, other: A)
    where
        A: Borrow<Self>,
    {
        let x = other.borrow();
        self.0.add_scaled(alpha, &x.0);
        self.1.add_scaled(alpha, &x.1);
        self.2.add_scaled(alpha, &x.2);
    }
}

/// Implement for scalar values.
impl Aggregatable for Real {
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>,
    {
Changes to src/data/minorant.rs.
13
14
15
16
17
18
19
20
21
22

23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38

39
40
41
42
43
44
45
46
47
48
49
50
51
52
53


54
55
56
57
58
59
60
61
62
63
64
65
66

67

68
69
70
71






72
73
74
75
76




77

78

79
80


81
82
83
84
85

86
87
88
89
90
91
92
93
94
95
96
97

98
99
100
101
102
103
104
105
106
107
108
109


110


111


112
113



114
115


116

117

118



119






120

121

122




123



124
125


126



127



128



129



130



131








132

133




134


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

//! A linear minorant.

use super::{Aggregatable, DVector, Real};

use std::borrow::Borrow;
use std::fmt;


/// A linear minorant of a convex function.
///
/// A linear minorant of a convex function $f \colon \mathbb{R}\^n \to
/// \mathbb{R}$ is a linear function of the form
///
///   \\[ l \colon \mathbb{R}\^n \to \mathbb{R}, x \mapsto \langle g, x
///   \rangle + c \\]
///
/// such that $l(x) \le f(x)$ for all $x \in \mathbb{R}\^n$.
///
/// Each minorant may have an additional "primal" data, that is
/// linearly combined along with the minorant itself. This data can
/// be used for separation algorithms.
#[derive(Clone)]
pub struct Minorant<P: Aggregatable> {

    /// The constant term.
    pub constant: Real,

    /// The linear term.
    pub linear: DVector,

    /// The associated primal data.
    pub primal: P,
}

impl<P: Aggregatable> fmt::Debug for Minorant<P> {
    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
        write!(f, "{:?}", (self.constant, &self.linear))?;
        Ok(())
    }


}

impl<P: Aggregatable> fmt::Display for Minorant<P> {
    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
        write!(f, "{} + y * {}", self.constant, self.linear)?;
        Ok(())
    }
}

impl<P: Aggregatable> Default for Minorant<P> {
    fn default() -> Minorant<P> {
        Minorant {
            constant: 0.0,

            linear: dvec![],

            primal: P::default(),
        }
    }
}







impl<P: Aggregatable> Minorant<P> {
    /// Return a new 0 minorant.
    pub fn new(constant: Real, linear: Vec<Real>, primal: P) -> Minorant<P> {
        Minorant {




            constant,

            linear: DVector(linear),

            primal,
        }


    }

    /**
     * Evaluate minorant at some point.
     *

     * This function computes $c + \langle g, x \rangle$ for this minorant
     *   \\[\ell \colon \mathbb{R}\^n \to \mathbb{R}, x \mapsto c + \langle g, x \rangle\\]
     * and the given point $x \in \mathbb{R}\^n$.
     */
    pub fn eval(&self, x: &DVector) -> Real {
        self.constant + self.linear.dot(x)
    }

    /**
     * Move the center of the minorant.
     */
    pub fn move_center(&mut self, alpha: Real, d: &DVector) {

        self.constant += alpha * self.linear.dot(d);
    }

    /// Move the center of the minorant.
    ///
    /// In contrast to [`move_center`] the vector `d` might be shorter than the
    /// linear term of the minorant. The missing elements are assumed to be
    /// zero.
    ///
    /// [`move_center`]: Minorant::move_center
    pub fn move_center_begin(&mut self, alpha: Real, d: &DVector) {
        debug_assert!(self.linear.len() >= d.len());


        self.constant += alpha * self.linear.dot_begin(d);


    }


}




impl<P: Aggregatable> Aggregatable for Minorant<P> {
    fn new_scaled<A>(alpha: Real, other: A) -> Self


    where

        A: Borrow<Self>,

    {



        let m = other.borrow();






        Minorant {

            constant: alpha * m.constant,

            linear: DVector::scaled(&m.linear, alpha),




            primal: P::new_scaled(alpha, &m.primal),



        }
    }






    fn add_scaled<A>(&mut self, alpha: Real, other: A)



    where



        A: Borrow<Self>,



    {



        let m = other.borrow();








        self.constant += alpha * m.constant;

        self.linear.add_scaled(alpha, &m.linear);




        self.primal.add_scaled(alpha, &m.primal);


    }
}

/// The subgradient extender is a callback used to update existing minorants
/// given their associated primal data.
pub type SubgradientExtender<P, E> = dyn FnMut(usize, &P, &[usize]) -> Result<DVector, E> + Send;







|

|
>














|
|
>
|
|

|
<
|
<
<
|
|
<
|
<
<
|
>
>
|
|
<
|
<
<
|
<
|
<
<
<
<
>
|
>
<
|
<
<
>
>
>
>
>
>
|
<
|
<
|
>
>
>
>
|
>
|
>
|
|
>
>
|
|
|
|
<
>
|
|
|
|
|
|


<
|
<
|
>
|









|
|
>
>
|
>
>
|
>
>
|

>
>
>
|
|
>
>
|
>
|
>
|
>
>
>
|
>
>
>
>
>
>
|
>
|
>
|
>
>
>
>
|
>
>
>
|
|
>
>

>
>
>
|
>
>
>
|
>
>
>
|
>
>
>
|
>
>
>
|
>
>
>
>
>
>
>
>
|
>
|
>
>
>
>
|
>
>



|
<
|
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44

45


46
47

48


49
50
51
52
53

54


55

56




57
58
59

60


61
62
63
64
65
66
67

68

69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85

86
87
88
89
90
91
92
93
94

95

96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201

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

//! A linear minorant.

use super::{Aggregatable, DVector, Real};
use num_traits::Zero;
use std::borrow::Borrow;

use crate::data::raw::BLAS1;

/// A linear minorant of a convex function.
///
/// A linear minorant of a convex function $f \colon \mathbb{R}\^n \to
/// \mathbb{R}$ is a linear function of the form
///
///   \\[ l \colon \mathbb{R}\^n \to \mathbb{R}, x \mapsto \langle g, x
///   \rangle + c \\]
///
/// such that $l(x) \le f(x)$ for all $x \in \mathbb{R}\^n$.
///
/// Each minorant may have an additional "primal" data, that is
/// linearly combined along with the minorant itself. This data can
/// be used for separation algorithms.
pub trait Minorant: Aggregatable + Default + Send + 'static {
    type Primal: Aggregatable + Send;

    /// Return the constant value.
    fn constant(&self) -> Real;

    /// Return the dimension of the linear term.

    fn dim(&self) -> usize;



    /// Evaluate the linear term at some point.

    fn linear(&self, x: &DVector) -> Real;



    /// Return the associated primal data.
    fn primal(&self) -> &Self::Primal;

    /// Compute the inner product with another minorant.

    fn product(&self, other: &Self) -> Real;




    /// Add a scaled `self.linear` to `target`.




    ///
    /// This sets `target = target + alpha * self.linear`.
    fn add_scaled_to(&self, alpha: Real, target: &mut [Real]);




    /// Compute linear combination of (linear terms of) the given minorants and
    /// store in a dense vector.
    fn combine_to_vec<I, A>(minorants: I, target: &mut (Real, DVector))
    where
        I: IntoIterator<Item = (Real, A)>,
        A: Borrow<Self>,
    {

        target.0 = Real::zero();

        let mut it = minorants.into_iter();
        if let Some((alpha, m)) = it.next() {
            let m = m.borrow();
            target.1.init0(m.dim());
            m.add_scaled_to(alpha, &mut target.1);
            target.0 += alpha * m.constant();
            for (alpha, m) in it {
                let m = m.borrow();
                m.add_scaled_to(alpha, &mut target.1);
                target.0 += alpha * m.constant();
            }
        } else {
            target.1.clear();
        }
    }

    /// Evaluate minorant at some point.

    ///
    /// This function computes $c + \langle g, x \rangle$ for this minorant
    ///   \\[\ell \colon \mathbb{R}\^n \to \mathbb{R}, x \mapsto c + \langle g, x \rangle\\]
    /// and the given point $x \in \mathbb{R}\^n$.
    ///
    fn eval(&self, x: &DVector) -> Real {
        self.constant() + self.linear(x)
    }


    /// Move the center of the minorant along direction `d`.

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        debug_assert_eq!(d.len(), self.dim(), "Minorant and direction must have same dimension");
        self.move_center_begin(alpha, d)
    }

    /// Move the center of the minorant.
    ///
    /// In contrast to [`move_center`] the vector `d` might be shorter than the
    /// linear term of the minorant. The missing elements are assumed to be
    /// zero.
    ///
    /// [`move_center`]: Minorant::move_center
    fn move_center_begin(&mut self, alpha: Real, d: &DVector);

    fn debug(&self) -> &dyn std::fmt::Debug;
}

impl<P: Aggregatable + Send + 'static> Minorant for (Real, DVector, P) {
    type Primal = P;

    fn dim(&self) -> usize {
        self.1.len()
    }

    fn constant(&self) -> Real {
        self.0
    }

    fn linear(&self, x: &DVector) -> Real {
        self.1.dot(x)
    }

    fn primal(&self) -> &P {
        &self.2
    }

    fn product(&self, other: &Self) -> Real {
        self.1.dot(&other.1)
    }

    fn add_scaled_to(&self, alpha: Real, target: &mut [Real]) {
        debug_assert_eq!(self.1.len(), target.len());
        if !alpha.is_zero() {
            BLAS1::add_scaled(target, alpha, &self.1);
        }
    }

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        self.0 += alpha * self.1.dot(d);
    }

    fn move_center_begin(&mut self, alpha: Real, d: &DVector) {
        debug_assert!(self.1.len() >= d.len());
        self.0 += alpha * self.1.dot_begin(d);
    }

    fn debug(&self) -> &dyn std::fmt::Debug {
        &self.1
    }
}

impl Minorant for (Real, DVector) {
    type Primal = ();

    fn dim(&self) -> usize {
        self.1.len()
    }

    fn constant(&self) -> Real {
        self.0
    }

    fn linear(&self, x: &DVector) -> Real {
        self.1.dot(x)
    }

    fn primal(&self) -> &Self::Primal {
        &()
    }

    fn product(&self, other: &Self) -> Real {
        self.1.dot(&other.1)
    }

    fn add_scaled_to(&self, alpha: Real, target: &mut [Real]) {
        debug_assert_eq!(self.1.len(), target.len());
        if !alpha.is_zero() {
            BLAS1::add_scaled(target, alpha, &self.1);
        }
    }

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        self.0 += alpha * self.1.dot(d);
    }

    fn move_center_begin(&mut self, alpha: Real, d: &DVector) {
        debug_assert!(self.1.len() >= d.len());
        self.0 += alpha * self.1.dot_begin(d);
    }

    fn debug(&self) -> &dyn std::fmt::Debug {
        self
    }
}

/// The subgradient extender is a callback used to update an existing minorant.

pub type SubgradientExtender<M, E> = dyn FnMut(usize, &mut M) -> Result<(), E> + Send;
Changes to src/data/mod.rs.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18


19
20
21
22
23
24
25
26
27
28
29
30
/*
 * Copyright (c) 2019 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/>
 */

//! General types and data structures.



pub mod vector;
pub use vector::{DVector, Vector};

pub mod aggregatable;
pub use aggregatable::Aggregatable;

pub mod minorant;
pub use minorant::Minorant;

/// Type used for real numbers throughout the library.
pub type Real = f64;

|
















>
>












1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
/*
 * Copyright (c) 2019, 2020 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/>
 */

//! General types and data structures.

pub(crate) mod raw;

pub mod vector;
pub use vector::{DVector, Vector};

pub mod aggregatable;
pub use aggregatable::Aggregatable;

pub mod minorant;
pub use minorant::Minorant;

/// Type used for real numbers throughout the library.
pub type Real = f64;
Added src/data/raw.rs.


































































































































































>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
/*
 * Copyright (c) 2020 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/>
 */

///! BLAS-like low-level vector routines.

#[cfg(feature = "blas")]
use {openblas_src as _, rs_blas as blas, std::os::raw::c_int};

pub trait BLAS1<T> {
    /// Compute the inner product.
    fn dot(&self, y: &[T]) -> T;

    /// Compute the inner product.
    ///
    /// The inner product is computed on the smaller of the two
    /// dimensions. All other elements are assumed to be zero.
    fn dot_begin(&self, y: &[T]) -> T;

    /// Compute `self = self + alpha * y`.
    fn add_scaled(&mut self, alpha: T, y: &[T]);

    /// Return the 2-norm of this vector.
    fn norm2(&self) -> T;
}

impl BLAS1<f64> for [f64] {
    fn dot(&self, other: &[f64]) -> f64 {
        debug_assert_eq!(self.len(), other.len(), "Vectors must have the same size");
        Self::dot_begin(self, other)
    }

    fn dot_begin(&self, other: &[f64]) -> f64 {
        #[cfg(feature = "blas")]
        unsafe {
            blas::ddot(self.len().min(other.len()) as c_int, &self, 1, &other, 1)
        }
        #[cfg(not(feature = "blas"))]
        {
            self.iter().zip(other.iter()).map(|(x, y)| x * y).sum::<f64>()
        }
    }

    fn add_scaled(&mut self, alpha: f64, y: &[f64]) {
        assert_eq!(self.len(), y.len());
        #[cfg(feature = "blas")]
        unsafe {
            blas::daxpy(self.len() as c_int, alpha, &y, 1, &mut self[..], 1)
        }
        #[cfg(not(feature = "blas"))]
        {
            for (x, y) in self.iter_mut().zip(y.iter()) {
                *x += alpha * y;
            }
        }
    }

    fn norm2(&self) -> f64 {
        #[cfg(feature = "blas")]
        unsafe {
            blas::dnrm2(self.len() as c_int, &self, 1)
        }
        #[cfg(not(feature = "blas"))]
        {
            self.iter().map(|x| x * x).sum::<f64>().sqrt()
        }
    }
}
Changes to src/data/vector.rs.
1
2
3
4
5
6
7
8
// Copyright (c) 2016, 2017, 2018, 2019 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
|







1
2
3
4
5
6
7
8
// Copyright (c) 2016-2020 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
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
use num_traits::Zero;
use std::borrow::Borrow;
use std::fmt;
use std::iter::FromIterator;
use std::ops::{Deref, DerefMut};
use std::vec::IntoIter;

#[cfg(feature = "blas")]
use {openblas_src as _, rs_blas as blas, std::os::raw::c_int};

/// Type of dense vectors.
#[derive(Debug, Clone, PartialEq, Default)]
pub struct DVector(pub Vec<Real>);

impl Deref for DVector {
    type Target = Vec<Real>;







<
|







20
21
22
23
24
25
26

27
28
29
30
31
32
33
34
use num_traits::Zero;
use std::borrow::Borrow;
use std::fmt;
use std::iter::FromIterator;
use std::ops::{Deref, DerefMut};
use std::vec::IntoIter;


use super::raw::BLAS1;

/// Type of dense vectors.
#[derive(Debug, Clone, PartialEq, Default)]
pub struct DVector(pub Vec<Real>);

impl Deref for DVector {
    type Target = Vec<Real>;
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
    }

    /// Return the inner product with another vector.
    ///
    /// The inner product is computed on the smaller of the two
    /// dimensions. All other elements are assumed to be zero.
    pub fn dot_begin(&self, other: &DVector) -> Real {
        #[cfg(feature = "blas")]
        unsafe {
            blas::ddot(self.len().min(other.len()) as c_int, &self, 1, &other, 1)
        }
        #[cfg(not(feature = "blas"))]
        {
            self.iter().zip(other.iter()).map(|(x, y)| x * y).sum::<Real>()
        }
    }

    /// Add two vectors and store result in this vector.
    pub fn add(&mut self, x: &DVector, y: &DVector) {
        assert_eq!(x.len(), y.len());
        self.clear();
        self.extend(x.iter().zip(y.iter()).map(|(a, b)| a + b));
    }

    /// Add two vectors and store result in this vector.
    pub fn add_scaled(&mut self, alpha: Real, y: &DVector) {
        assert_eq!(self.len(), y.len());
        #[cfg(feature = "blas")]
        unsafe {
            blas::daxpy(self.len() as c_int, alpha, &y, 1, &mut self[..], 1)
        }
        #[cfg(not(feature = "blas"))]
        {
            for (x, y) in self.iter_mut().zip(y.iter()) {
                *x += alpha * y;
            }
        }
    }

    /// Add two vectors and store result in this vector.
    ///
    /// In contrast to `add_scaled`, the two vectors might have
    /// different sizes. The size of the resulting vector is the
    /// larger of the two vector sizes and the remaining entries of
    /// the smaller vector are assumed to be 0.0.
    pub fn add_scaled_begin(&mut self, alpha: Real, y: &DVector) {
        #[cfg(feature = "blas")]
        unsafe {
            let n = self.len();
            blas::daxpy(n.min(y.len()) as c_int, alpha, &y, 1, &mut self[..], 1);
        }
        #[cfg(not(feature = "blas"))]
        {
            for (x, y) in self.iter_mut().zip(y.iter()) {
                *x += alpha * y;
            }
        }
        let n = self.len();
        if n < y.len() {
            self.extend(y[n..].iter().map(|y| alpha * y));
        }
    }

    /// Return the 2-norm of this vector.
    pub fn norm2(&self) -> Real {
        #[cfg(feature = "blas")]
        unsafe {
            blas::dnrm2(self.len() as c_int, &self, 1)
        }
        #[cfg(not(feature = "blas"))]
        {
            self.iter().map(|x| x * x).sum::<Real>().sqrt()
        }
    }
}

impl Aggregatable for DVector {
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>,







<
<
|
<
<
<
<
<











<
<
<
|
<
<
<
<
<
<
<









<
<
|
<
|
<
<
<
|
|
<
<
|






<
<
|
<
<
<
<
<







134
135
136
137
138
139
140


141





142
143
144
145
146
147
148
149
150
151
152



153







154
155
156
157
158
159
160
161
162


163

164



165
166


167
168
169
170
171
172
173


174





175
176
177
178
179
180
181
    }

    /// Return the inner product with another vector.
    ///
    /// The inner product is computed on the smaller of the two
    /// dimensions. All other elements are assumed to be zero.
    pub fn dot_begin(&self, other: &DVector) -> Real {


        BLAS1::dot_begin(&self[..], other)





    }

    /// Add two vectors and store result in this vector.
    pub fn add(&mut self, x: &DVector, y: &DVector) {
        assert_eq!(x.len(), y.len());
        self.clear();
        self.extend(x.iter().zip(y.iter()).map(|(a, b)| a + b));
    }

    /// Add two vectors and store result in this vector.
    pub fn add_scaled(&mut self, alpha: Real, y: &DVector) {



        BLAS1::add_scaled(&mut self[..], alpha, y);







    }

    /// Add two vectors and store result in this vector.
    ///
    /// In contrast to `add_scaled`, the two vectors might have
    /// different sizes. The size of the resulting vector is the
    /// larger of the two vector sizes and the remaining entries of
    /// the smaller vector are assumed to be 0.0.
    pub fn add_scaled_begin(&mut self, alpha: Real, y: &DVector) {


        let n = self.len().min(y.len());





        BLAS1::add_scaled(&mut self[..n], alpha, &y[..n]);



        if self.len() < y.len() {
            self.extend(y[n..].iter().map(|y| alpha * y));
        }
    }

    /// Return the 2-norm of this vector.
    pub fn norm2(&self) -> Real {


        BLAS1::norm2(&self[..])





    }
}

impl Aggregatable for DVector {
    fn new_scaled<A>(alpha: Real, other: A) -> Self
    where
        A: Borrow<Self>,
Changes to src/master/boxed.rs.
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
//

pub mod unconstrained;
use self::unconstrained::UnconstrainedMasterProblem;

use super::MasterProblem;
use crate::problem::SubgradientExtender;
use crate::{Aggregatable, DVector, Minorant, Real};

use either::Either;
use log::debug;
use std::f64::{EPSILON, INFINITY, NEG_INFINITY};
use std::marker::PhantomData;

/**
 * Turn unconstrained master problem into box-constrained one.
 *
 * This master problem adds box constraints to an unconstrained
 * master problem implementation. The box constraints are enforced by
 * an additional outer optimization loop.
 */
pub struct BoxedMasterProblem<P, M>
where
    P: Aggregatable,
{
    /// The unconstrained master problem solver.
    pub master: M,

    /// Maximal number of updates of box multipliers.
    pub max_updates: usize,








|













|

|







15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
//

pub mod unconstrained;
use self::unconstrained::UnconstrainedMasterProblem;

use super::MasterProblem;
use crate::problem::SubgradientExtender;
use crate::{DVector, Minorant, Real};

use either::Either;
use log::debug;
use std::f64::{EPSILON, INFINITY, NEG_INFINITY};
use std::marker::PhantomData;

/**
 * Turn unconstrained master problem into box-constrained one.
 *
 * This master problem adds box constraints to an unconstrained
 * master problem implementation. The box constraints are enforced by
 * an additional outer optimization loop.
 */
pub struct BoxedMasterProblem<Mn, M>
where
    Mn: Minorant,
{
    /// The unconstrained master problem solver.
    pub master: M,

    /// Maximal number of updates of box multipliers.
    pub max_updates: usize,

65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
    model_eps: Real,

    need_new_candidate: bool,

    /// Current number of updates.
    cnt_updates: usize,

    phantom: PhantomData<P>,
}

impl<P, M> BoxedMasterProblem<P, M>
where
    P: Aggregatable,
    M: UnconstrainedMasterProblem<P>,
{
    pub fn with_master(master: M) -> BoxedMasterProblem<P, M> {
        BoxedMasterProblem {
            master,
            max_updates: 50,
            lb: dvec![],
            ub: dvec![],
            eta: dvec![],
            primopt: dvec![],







|


|

|
|

|







65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
    model_eps: Real,

    need_new_candidate: bool,

    /// Current number of updates.
    cnt_updates: usize,

    phantom: PhantomData<Mn>,
}

impl<Mn, M> BoxedMasterProblem<Mn, M>
where
    Mn: Minorant,
    M: UnconstrainedMasterProblem<Mn>,
{
    pub fn with_master(master: M) -> BoxedMasterProblem<Mn, M> {
        BoxedMasterProblem {
            master,
            max_updates: 50,
            lb: dvec![],
            ub: dvec![],
            eta: dvec![],
            primopt: dvec![],
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
     * the current box-multipliers $\eta$.
     */
    fn get_norm_subg2(&self) -> Real {
        self.eta.dot(self.master.dualopt())
    }
}

impl<P, M> MasterProblem<P> for BoxedMasterProblem<P, M>
where
    P: Aggregatable + Send + 'static,
    M: UnconstrainedMasterProblem<P>,
{
    type MinorantIndex = M::MinorantIndex;

    type Err = M::Err;

    fn new() -> Result<Self, Self::Err> {
        M::new().map(BoxedMasterProblem::with_master)







|

|
|







187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
     * the current box-multipliers $\eta$.
     */
    fn get_norm_subg2(&self) -> Real {
        self.eta.dot(self.master.dualopt())
    }
}

impl<Mn, M> MasterProblem<Mn> for BoxedMasterProblem<Mn, M>
where
    Mn: Minorant,
    M: UnconstrainedMasterProblem<Mn>,
{
    type MinorantIndex = M::MinorantIndex;

    type Err = M::Err;

    fn new() -> Result<Self, Self::Err> {
        M::new().map(BoxedMasterProblem::with_master)
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251

252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
        Ok(())
    }

    fn num_minorants(&self, fidx: usize) -> usize {
        self.master.num_minorants(fidx)
    }

    fn add_minorant(&mut self, fidx: usize, minorant: Minorant<P>) -> Result<Self::MinorantIndex, Self::Err> {
        self.master.add_minorant(fidx, minorant)
    }

    fn weight(&self) -> Real {
        self.master.weight()
    }

    fn set_weight(&mut self, weight: Real) -> Result<(), Self::Err> {
        self.master.set_weight(weight)
    }

    fn add_vars<SErr>(
        &mut self,
        bounds: &[(Option<usize>, Real, Real)],
        extend_subgradient: &mut SubgradientExtender<P, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>> {

        if !bounds.is_empty() {
            for (index, l, u) in bounds.iter().filter_map(|v| v.0.map(|i| (i, v.1, v.2))) {
                self.lb[index] = l;
                self.ub[index] = u;
            }
            self.lb.extend(bounds.iter().filter(|v| v.0.is_none()).map(|x| x.1));
            self.ub.extend(bounds.iter().filter(|v| v.0.is_none()).map(|x| x.2));
            self.eta.resize(self.lb.len(), 0.0);
            self.need_new_candidate = true;
            let nnew = bounds.iter().filter(|v| v.0.is_none()).count();
            let changed = bounds.iter().filter_map(|v| v.0).collect::<Vec<_>>();
            self.master.add_vars(nnew, &changed, extend_subgradient)
        } else {
            Ok(())
        }
    }

    #[allow(clippy::cognitive_complexity)]
    fn solve(&mut self, center_value: Real) -> Result<(), Self::Err> {







|













|
|

>

<
<
<
<
|
|


|
<
|







228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253




254
255
256
257
258

259
260
261
262
263
264
265
266
        Ok(())
    }

    fn num_minorants(&self, fidx: usize) -> usize {
        self.master.num_minorants(fidx)
    }

    fn add_minorant(&mut self, fidx: usize, minorant: Mn) -> Result<Self::MinorantIndex, Self::Err> {
        self.master.add_minorant(fidx, minorant)
    }

    fn weight(&self) -> Real {
        self.master.weight()
    }

    fn set_weight(&mut self, weight: Real) -> Result<(), Self::Err> {
        self.master.set_weight(weight)
    }

    fn add_vars<SErr>(
        &mut self,
        bounds: &[(Real, Real)],
        extend_subgradient: &mut SubgradientExtender<Mn, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>> {
        debug_assert!(bounds.iter().all(|&(l, u)| l <= u));
        if !bounds.is_empty() {




            self.lb.extend(bounds.iter().map(|x| x.0));
            self.ub.extend(bounds.iter().map(|x| x.1));
            self.eta.resize(self.lb.len(), 0.0);
            self.need_new_candidate = true;
            let nnew = bounds.len();

            self.master.add_vars(nnew, extend_subgradient)
        } else {
            Ok(())
        }
    }

    #[allow(clippy::cognitive_complexity)]
    fn solve(&mut self, center_value: Real) -> Result<(), Self::Err> {
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
        self.dualoptnorm2
    }

    fn multiplier(&self, min: Self::MinorantIndex) -> Real {
        self.master.multiplier(min)
    }

    fn aggregated_primal(&self, fidx: usize) -> Result<P, Self::Err> {
        self.master.aggregated_primal(fidx)
    }

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        self.need_new_candidate = true;
        self.master.move_center(alpha, d);
        self.lb.add_scaled_begin(-alpha, d);







|







352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
        self.dualoptnorm2
    }

    fn multiplier(&self, min: Self::MinorantIndex) -> Real {
        self.master.multiplier(min)
    }

    fn aggregated_primal(&self, fidx: usize) -> Result<Mn::Primal, Self::Err> {
        self.master.aggregated_primal(fidx)
    }

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        self.need_new_candidate = true;
        self.master.move_center(alpha, d);
        self.lb.add_scaled_begin(-alpha, d);
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404

/// Builder for `BoxedMasterProblem`.
///
/// `B` is a builder of the underlying `UnconstrainedMasterProblem`.
#[derive(Default)]
pub struct Builder<B>(B);

impl<P, B> super::Builder<P> for Builder<B>
where
    P: Aggregatable + Send + 'static,
    B: unconstrained::Builder<P>,
    B::MasterProblem: UnconstrainedMasterProblem<P>,
{
    type MasterProblem = BoxedMasterProblem<P, B::MasterProblem>;

    fn build(&mut self) -> Result<Self::MasterProblem, <Self::MasterProblem as MasterProblem<P>>::Err> {
        self.0.build().map(BoxedMasterProblem::with_master)
    }
}

impl<B> std::ops::Deref for Builder<B> {
    type Target = B;








|

|
|
|

|

|







378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400

/// Builder for `BoxedMasterProblem`.
///
/// `B` is a builder of the underlying `UnconstrainedMasterProblem`.
#[derive(Default)]
pub struct Builder<B>(B);

impl<Mn, B> super::Builder<Mn> for Builder<B>
where
    Mn: Minorant,
    B: unconstrained::Builder<Mn>,
    B::MasterProblem: UnconstrainedMasterProblem<Mn>,
{
    type MasterProblem = BoxedMasterProblem<Mn, B::MasterProblem>;

    fn build(&mut self) -> Result<Self::MasterProblem, <Self::MasterProblem as MasterProblem<Mn>>::Err> {
        self.0.build().map(BoxedMasterProblem::with_master)
    }
}

impl<B> std::ops::Deref for Builder<B> {
    type Target = B;

Changes to src/master/boxed/unconstrained.rs.
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
// along with this program.  If not, see  <http://www.gnu.org/licenses/>
//

pub mod cpx;
pub mod minimal;

use crate::problem::SubgradientExtender;
use crate::{Aggregatable, DVector, Minorant, Real};

use either::Either;
use std::error::Error;

/**
 * Trait for master problems without box constraints.
 *







|







14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
// along with this program.  If not, see  <http://www.gnu.org/licenses/>
//

pub mod cpx;
pub mod minimal;

use crate::problem::SubgradientExtender;
use crate::{DVector, Minorant, Real};

use either::Either;
use std::error::Error;

/**
 * Trait for master problems without box constraints.
 *
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
 * \alpha_i = 1 \right\\}$. Note, the unconstrained solver is expected
 * to compute *dual* optimal solutions, i.e. the solver must compute
 * optimal coefficients $\bar{\alpha}$ for the dual problem
 *
 * \\[ \max_{\alpha \in \Delta} \min_{d \in \mathbb{R}\^n}
 *     \sum_{i=1}\^k \alpha_i \ell_i(d) + \frac{u}{2} \\| d \\|\^2. \\]
 */
pub trait UnconstrainedMasterProblem<P: Aggregatable>: Send + 'static {
    /// Unique index for a minorant.
    type MinorantIndex: Copy + Eq;

    /// Error type for this master problem.
    type Err: Error + Send + Sync;

    /// Return a new instance of the unconstrained master problem.
    fn new() -> Result<Self, Self::Err>







|
<







43
44
45
46
47
48
49
50

51
52
53
54
55
56
57
 * \alpha_i = 1 \right\\}$. Note, the unconstrained solver is expected
 * to compute *dual* optimal solutions, i.e. the solver must compute
 * optimal coefficients $\bar{\alpha}$ for the dual problem
 *
 * \\[ \max_{\alpha \in \Delta} \min_{d \in \mathbb{R}\^n}
 *     \sum_{i=1}\^k \alpha_i \ell_i(d) + \frac{u}{2} \\| d \\|\^2. \\]
 */
pub trait UnconstrainedMasterProblem<M: Minorant>: Send + 'static {

    type MinorantIndex: Copy + Eq;

    /// Error type for this master problem.
    type Err: Error + Send + Sync;

    /// Return a new instance of the unconstrained master problem.
    fn new() -> Result<Self, Self::Err>
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
    ///
    /// When some minorants are compressed, the callback is called with the
    /// coefficients and indices of the compressed minorants and the index of
    /// the new minorant. The callback may be called several times.
    fn compress(&mut self) -> Result<(), Self::Err>;

    /// Add a new minorant to the model.
    fn add_minorant(&mut self, fidx: usize, minorant: Minorant<P>) -> Result<Self::MinorantIndex, Self::Err>;

    /// Add or move some variables.
    ///
    /// The variables in `changed` have been changed, so the subgradient
    /// information must be updated. Furthermore, `nnew` new variables
    /// are added.
    fn add_vars<SErr>(
        &mut self,
        nnew: usize,
        changed: &[usize],
        extend_subgradient: &mut SubgradientExtender<P, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>>;

    /// Solve the master problem.
    fn solve(&mut self, eta: &DVector, fbound: Real, augbound: Real, relprec: Real) -> Result<(), Self::Err>;

    /// Return the current dual optimal solution.
    fn dualopt(&self) -> &DVector;

    /// Return the current dual optimal solution value.
    fn dualopt_cutval(&self) -> Real;

    /// Return the multiplier associated with a minorant.
    fn multiplier(&self, min: Self::MinorantIndex) -> Real;

    /// Return the value of the current model at the given point.
    fn eval_model(&self, y: &DVector) -> Real;

    /// Return the aggregated primal.
    fn aggregated_primal(&self, fidx: usize) -> Result<P, Self::Err>;

    /// Move the center of the master problem along $\alpha \cdot d$.
    fn move_center(&mut self, alpha: Real, d: &DVector);
}

/// A builder for creating unconstrained master problem solvers.
pub trait Builder<P: Aggregatable> {
    /// The master problem to be build.
    type MasterProblem: UnconstrainedMasterProblem<P>;

    /// Create a new master problem.
    fn build(&mut self) -> Result<Self::MasterProblem, <Self::MasterProblem as UnconstrainedMasterProblem<P>>::Err>;
}







|

<
<
<
|
<



<
|


















|






|

|


|

77
78
79
80
81
82
83
84
85



86

87
88
89

90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
    ///
    /// When some minorants are compressed, the callback is called with the
    /// coefficients and indices of the compressed minorants and the index of
    /// the new minorant. The callback may be called several times.
    fn compress(&mut self) -> Result<(), Self::Err>;

    /// Add a new minorant to the model.
    fn add_minorant(&mut self, fidx: usize, minorant: M) -> Result<Self::MinorantIndex, Self::Err>;




    /// Add `nnew` new variables.

    fn add_vars<SErr>(
        &mut self,
        nnew: usize,

        extend_subgradient: &mut SubgradientExtender<M, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>>;

    /// Solve the master problem.
    fn solve(&mut self, eta: &DVector, fbound: Real, augbound: Real, relprec: Real) -> Result<(), Self::Err>;

    /// Return the current dual optimal solution.
    fn dualopt(&self) -> &DVector;

    /// Return the current dual optimal solution value.
    fn dualopt_cutval(&self) -> Real;

    /// Return the multiplier associated with a minorant.
    fn multiplier(&self, min: Self::MinorantIndex) -> Real;

    /// Return the value of the current model at the given point.
    fn eval_model(&self, y: &DVector) -> Real;

    /// Return the aggregated primal.
    fn aggregated_primal(&self, fidx: usize) -> Result<M::Primal, Self::Err>;

    /// Move the center of the master problem along $\alpha \cdot d$.
    fn move_center(&mut self, alpha: Real, d: &DVector);
}

/// A builder for creating unconstrained master problem solvers.
pub trait Builder<Mn: Minorant> {
    /// The master problem to be build.
    type MasterProblem: UnconstrainedMasterProblem<Mn>;

    /// Create a new master problem.
    fn build(&mut self) -> Result<Self::MasterProblem, <Self::MasterProblem as UnconstrainedMasterProblem<Mn>>::Err>;
}
Changes to src/master/boxed/unconstrained/cpx.rs.
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
        CplexMasterError::Cplex(err)
    }
}

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

/// A minorant and its unique index.
struct MinorantInfo<P: Aggregatable> {
    minorant: Minorant<P>,
    index: usize,
}

impl<P: Aggregatable> Deref for MinorantInfo<P> {
    type Target = Minorant<P>;
    fn deref(&self) -> &Minorant<P> {
        &self.minorant
    }
}

impl<P: Aggregatable> DerefMut for MinorantInfo<P> {
    fn deref_mut(&mut self) -> &mut Minorant<P> {
        &mut self.minorant
    }
}

/// Maps a global index to a minorant.
#[derive(Clone, Copy)]
struct MinorantIdx {
    /// The function (subproblem) index.
    fidx: usize,
    /// The minorant index within the subproblem.
    idx: usize,
}

/// A submodel.
///
/// A submodel is a list of subproblems being aggregated in that model.
#[derive(Debug, Clone, Default)]
struct SubModel<P: Aggregatable> {
    /// The list of subproblems.
    subproblems: Vec<usize>,

    /// The number of minorants in this subproblem.
    ///
    /// The is the minimal number of minorants in each contained subproblem.
    num_mins: usize,
    /// The aggregated minorants of this submodel.
    ///
    /// This is just the sum of the corresponding minorants of the single
    /// functions contained in this submodel.
    minorants: Vec<Minorant<P>>,
}

impl<P: Aggregatable> Deref for SubModel<P> {
    type Target = Vec<usize>;
    fn deref(&self) -> &Vec<usize> {
        &self.subproblems
    }
}

impl<P: Aggregatable> DerefMut for SubModel<P> {
    fn deref_mut(&mut self) -> &mut Vec<usize> {
        &mut self.subproblems
    }
}

pub struct CplexMaster<P: Aggregatable> {
    env: *mut cpx::Env,

    lp: *mut cpx::Lp,

    /// True if the QP must be updated.
    force_update: bool,








|
|



|
<
<
<
<
<
|
<
<
|
















|











|


|






|





|







67
68
69
70
71
72
73
74
75
76
77
78
79





80


81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
        CplexMasterError::Cplex(err)
    }
}

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

/// A minorant and its unique index.
struct MinorantInfo<Mn: Minorant> {
    minorant: Mn,
    index: usize,
}

impl<Mn: Minorant> MinorantInfo<Mn> {





    fn new(minorant: Mn, index: usize) -> Self {


        MinorantInfo { minorant, index }
    }
}

/// Maps a global index to a minorant.
#[derive(Clone, Copy)]
struct MinorantIdx {
    /// The function (subproblem) index.
    fidx: usize,
    /// The minorant index within the subproblem.
    idx: usize,
}

/// A submodel.
///
/// A submodel is a list of subproblems being aggregated in that model.
#[derive(Debug, Clone, Default)]
struct SubModel {
    /// The list of subproblems.
    subproblems: Vec<usize>,

    /// The number of minorants in this subproblem.
    ///
    /// The is the minimal number of minorants in each contained subproblem.
    num_mins: usize,
    /// The aggregated minorants of this submodel.
    ///
    /// This is just the sum of the corresponding minorants of the single
    /// functions contained in this submodel.
    minorants: Vec<(Real, DVector)>,
}

impl Deref for SubModel {
    type Target = Vec<usize>;
    fn deref(&self) -> &Vec<usize> {
        &self.subproblems
    }
}

impl DerefMut for SubModel {
    fn deref_mut(&mut self) -> &mut Vec<usize> {
        &mut self.subproblems
    }
}

pub struct CplexMaster<Mn: Minorant> {
    env: *mut cpx::Env,

    lp: *mut cpx::Lp,

    /// True if the QP must be updated.
    force_update: bool,

153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
    /// The additional diagonal term to ensure positive definiteness.
    qdiag: Real,

    /// The weight of the quadratic term.
    weight: Real,

    /// The minorants for each subproblem in the model.
    minorants: Vec<Vec<MinorantInfo<P>>>,

    /// The callback for the submodel for each subproblem.
    select_model: Arc<dyn Fn(usize) -> usize>,

    /// For each submodel the list of subproblems contained in that model.
    submodels: Vec<SubModel<P>>,
    /// For each subproblem the submodel it is contained in.
    in_submodel: Vec<usize>,

    /// Optimal multipliers for each subproblem in the model.
    opt_mults: Vec<DVector>,
    /// Optimal aggregated minorant.
    opt_minorant: Minorant<P>,

    /// Maximal bundle size.
    pub max_bundle_size: usize,
}

unsafe impl<P: Aggregatable> Send for CplexMaster<P> {}

impl<P: Aggregatable> Drop for CplexMaster<P> {
    fn drop(&mut self) {
        unsafe { cpx::freeprob(self.env, &mut self.lp) };
        unsafe { cpx::closeCPLEX(&mut self.env) };
    }
}

impl<P: Aggregatable + 'static> UnconstrainedMasterProblem<P> for CplexMaster<P> {
    type MinorantIndex = usize;

    type Err = CplexMasterError;

    fn new() -> Result<CplexMaster<P>> {
        let env;

        trycpx!({
            let mut status = 0;
            env = cpx::openCPLEX(&mut status);
            status
        });







|





|






|





|

|






|




|







146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
    /// The additional diagonal term to ensure positive definiteness.
    qdiag: Real,

    /// The weight of the quadratic term.
    weight: Real,

    /// The minorants for each subproblem in the model.
    minorants: Vec<Vec<MinorantInfo<Mn>>>,

    /// The callback for the submodel for each subproblem.
    select_model: Arc<dyn Fn(usize) -> usize>,

    /// For each submodel the list of subproblems contained in that model.
    submodels: Vec<SubModel>,
    /// For each subproblem the submodel it is contained in.
    in_submodel: Vec<usize>,

    /// Optimal multipliers for each subproblem in the model.
    opt_mults: Vec<DVector>,
    /// Optimal aggregated minorant.
    opt_minorant: (Real, DVector),

    /// Maximal bundle size.
    pub max_bundle_size: usize,
}

unsafe impl<Mn: Minorant> Send for CplexMaster<Mn> {}

impl<Mn: Minorant> Drop for CplexMaster<Mn> {
    fn drop(&mut self) {
        unsafe { cpx::freeprob(self.env, &mut self.lp) };
        unsafe { cpx::closeCPLEX(&mut self.env) };
    }
}

impl<Mn: Minorant> UnconstrainedMasterProblem<Mn> for CplexMaster<Mn> {
    type MinorantIndex = usize;

    type Err = CplexMasterError;

    fn new() -> Result<CplexMaster<Mn>> {
        let env;

        trycpx!({
            let mut status = 0;
            env = cpx::openCPLEX(&mut status);
            status
        });
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
            qdiag: 0.0,
            weight: 1.0,
            minorants: vec![],
            select_model: Arc::new(|i| i),
            submodels: vec![],
            in_submodel: vec![],
            opt_mults: vec![],
            opt_minorant: Minorant::default(),
            max_bundle_size: 50,
        })
    }

    fn num_subproblems(&self) -> usize {
        self.minorants.len()
    }







|







203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
            qdiag: 0.0,
            weight: 1.0,
            minorants: vec![],
            select_model: Arc::new(|i| i),
            submodels: vec![],
            in_submodel: vec![],
            opt_mults: vec![],
            opt_minorant: Default::default(),
            max_bundle_size: 50,
        })
    }

    fn num_subproblems(&self) -> usize {
        self.minorants.len()
    }
268
269
270
271
272
273
274
275
276



277


278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
                    .collect::<Vec<_>>();
                self.aggregate(i, &inds)?;
            }
        }
        Ok(())
    }

    fn add_minorant(&mut self, fidx: usize, minorant: Minorant<P>) -> Result<usize> {
        debug!("Add minorant");



        debug!("  fidx={} index={}: {}", fidx, self.minorants[fidx].len(), minorant);



        debug_assert!(
            self.minorants[fidx]
                .last()
                .map_or(true, |m| m.linear.len() == minorant.linear.len()),
            "Newly added minorant has wrong size (got: {}, expected: {})",
            minorant.linear.len(),
            self.minorants[fidx].last().map_or(0, |m| m.linear.len())
        );

        // find new unique minorant index
        let min_idx = self.minorants[fidx].len();
        let index = if let Some(index) = self.freeinds.pop() {
            self.index2min[index] = MinorantIdx { fidx, idx: min_idx };
            index
        } else {
            self.index2min.push(MinorantIdx { fidx, idx: min_idx });
            self.index2min.len() - 1
        };
        self.updateinds.push(index);

        // store minorant
        self.minorants[fidx].push(MinorantInfo { minorant, index });
        self.opt_mults[fidx].push(0.0);

        Ok(index)
    }

    fn add_vars<SErr>(
        &mut self,
        nnew: usize,
        changed: &[usize],
        extend_subgradient: &mut SubgradientExtender<P, SErr>,
    ) -> std::result::Result<(), Either<CplexMasterError, SErr>> {
        debug_assert!(!self.minorants[0].is_empty());
        if changed.is_empty() && nnew == 0 {
            return Ok(());
        }
        let noldvars = self.minorants[0][0].linear.len();
        let nnewvars = noldvars + nnew;

        let mut changedvars = vec![];
        changedvars.extend_from_slice(changed);
        changedvars.extend(noldvars..nnewvars);
        for (fidx, mins) in self.minorants.iter_mut().enumerate() {
            for m in &mut mins[..] {
                let new_subg = (extend_subgradient)(fidx, &m.primal, &changedvars).map_err(Either::Right)?;
                for (&j, &g) in changed.iter().zip(new_subg.iter()) {
                    m.linear[j] = g;
                }
                m.linear.extend_from_slice(&new_subg[changed.len()..]);
            }
        }

        // update qterm because all minorants have changed
        self.force_update = true;

        Ok(())







|

>
>
>
|
>
>




|

|
|














|








<
|


|


<
<

<
<
<


|
<
<
<
<







261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306

307
308
309
310
311
312


313



314
315
316




317
318
319
320
321
322
323
                    .collect::<Vec<_>>();
                self.aggregate(i, &inds)?;
            }
        }
        Ok(())
    }

    fn add_minorant(&mut self, fidx: usize, minorant: Mn) -> Result<usize> {
        debug!("Add minorant");
        debug!(
            "  fidx={} index={}: {:?}",
            fidx,
            self.minorants[fidx].len(),
            minorant.debug()
        );

        debug_assert!(
            self.minorants[fidx]
                .last()
                .map_or(true, |m| m.minorant.dim() == minorant.dim()),
            "Newly added minorant has wrong size (got: {}, expected: {})",
            minorant.dim(),
            self.minorants[fidx].last().map_or(0, |m| m.minorant.dim())
        );

        // find new unique minorant index
        let min_idx = self.minorants[fidx].len();
        let index = if let Some(index) = self.freeinds.pop() {
            self.index2min[index] = MinorantIdx { fidx, idx: min_idx };
            index
        } else {
            self.index2min.push(MinorantIdx { fidx, idx: min_idx });
            self.index2min.len() - 1
        };
        self.updateinds.push(index);

        // store minorant
        self.minorants[fidx].push(MinorantInfo::new(minorant, index));
        self.opt_mults[fidx].push(0.0);

        Ok(index)
    }

    fn add_vars<SErr>(
        &mut self,
        nnew: usize,

        extend_subgradient: &mut SubgradientExtender<Mn, SErr>,
    ) -> std::result::Result<(), Either<CplexMasterError, SErr>> {
        debug_assert!(!self.minorants[0].is_empty());
        if nnew == 0 {
            return Ok(());
        }






        for (fidx, mins) in self.minorants.iter_mut().enumerate() {
            for m in &mut mins[..] {
                (extend_subgradient)(fidx, &mut m.minorant).map_err(Either::Right)?;




            }
        }

        // update qterm because all minorants have changed
        self.force_update = true;

        Ok(())
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
        // update linear costs
        {
            let mut c = Vec::with_capacity(nvars);
            let mut inds = Vec::with_capacity(nvars);
            for submodel in &self.submodels {
                for i in 0..submodel.num_mins {
                    let m = &submodel.minorants[i];
                    let cost = -m.constant * self.weight - m.linear.dot(eta);
                    inds.push(c.len() as c_int);
                    c.push(cost);
                }
            }
            debug_assert_eq!(inds.len(), nvars);
            trycpx!(cpx::chgobj(
                self.env,







|







342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
        // update linear costs
        {
            let mut c = Vec::with_capacity(nvars);
            let mut inds = Vec::with_capacity(nvars);
            for submodel in &self.submodels {
                for i in 0..submodel.num_mins {
                    let m = &submodel.minorants[i];
                    let cost = -m.constant() * self.weight - m.linear(eta);
                    inds.push(c.len() as c_int);
                    c.push(cost);
                }
            }
            debug_assert_eq!(inds.len(), nvars);
            trycpx!(cpx::chgobj(
                self.env,
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
            for &fidx in submodel.iter() {
                for mult in &mut self.opt_mults[fidx][submodel.num_mins..] {
                    *mult = 0.0;
                }
            }
        }

        self.opt_minorant = Aggregatable::combine(mults.into_iter().zip(mins));

        Ok(())
    }

    fn dualopt(&self) -> &DVector {
        &self.opt_minorant.linear
    }

    fn dualopt_cutval(&self) -> Real {
        self.opt_minorant.constant
    }

    fn multiplier(&self, min: usize) -> Real {
        let MinorantIdx { fidx, idx } = self.index2min[min];
        self.opt_mults[fidx][idx]
    }

    fn eval_model(&self, y: &DVector) -> Real {
        let mut result = 0.0;
        for submodel in &self.submodels {
            let mut this_val = NEG_INFINITY;
            for m in &submodel.minorants[0..submodel.num_mins] {
                this_val = this_val.max(m.eval(y));
            }
            result += this_val;
        }
        result
    }

    fn aggregated_primal(&self, fidx: usize) -> Result<P> {
        Ok(P::combine(
            self.opt_mults[fidx]
                .iter()
                .enumerate()
                .map(|(i, alpha)| (*alpha, &self.minorants[fidx][i].primal)),
        ))
    }

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        for mins in &mut self.minorants {
            for m in mins.iter_mut() {
                m.move_center_begin(alpha, d);
            }
        }
        for submod in &mut self.submodels {
            for m in &mut submod.minorants {
                m.move_center_begin(alpha, d);
            }
        }
    }
}

impl<P: Aggregatable + 'static> CplexMaster<P> {
    /// Set a custom submodel selector.
    ///
    /// For each subproblem index the selector should return a submodel index.
    /// All subproblems with the same submodel index are aggregated in a single
    /// cutting plane model.
    fn set_submodel_selection<F>(&mut self, selector: F)
    where







|





|



|



















|
|



|






|










|







386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
            for &fidx in submodel.iter() {
                for mult in &mut self.opt_mults[fidx][submodel.num_mins..] {
                    *mult = 0.0;
                }
            }
        }

        Mn::combine_to_vec(mults.into_iter().zip(mins), &mut self.opt_minorant);

        Ok(())
    }

    fn dualopt(&self) -> &DVector {
        &self.opt_minorant.1
    }

    fn dualopt_cutval(&self) -> Real {
        self.opt_minorant.constant()
    }

    fn multiplier(&self, min: usize) -> Real {
        let MinorantIdx { fidx, idx } = self.index2min[min];
        self.opt_mults[fidx][idx]
    }

    fn eval_model(&self, y: &DVector) -> Real {
        let mut result = 0.0;
        for submodel in &self.submodels {
            let mut this_val = NEG_INFINITY;
            for m in &submodel.minorants[0..submodel.num_mins] {
                this_val = this_val.max(m.eval(y));
            }
            result += this_val;
        }
        result
    }

    fn aggregated_primal(&self, fidx: usize) -> Result<Mn::Primal> {
        Ok(Aggregatable::combine(
            self.opt_mults[fidx]
                .iter()
                .enumerate()
                .map(|(i, alpha)| (*alpha, self.minorants[fidx][i].minorant.primal())),
        ))
    }

    fn move_center(&mut self, alpha: Real, d: &DVector) {
        for mins in &mut self.minorants {
            for m in mins.iter_mut() {
                m.minorant.move_center_begin(alpha, d);
            }
        }
        for submod in &mut self.submodels {
            for m in &mut submod.minorants {
                m.move_center_begin(alpha, d);
            }
        }
    }
}

impl<Mn: Minorant> CplexMaster<Mn> {
    /// Set a custom submodel selector.
    ///
    /// For each subproblem index the selector should return a submodel index.
    /// All subproblems with the same submodel index are aggregated in a single
    /// cutting plane model.
    fn set_submodel_selection<F>(&mut self, selector: F)
    where
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
        }

        let minorants = &self.minorants;

        // Compute the number of minorants in each submodel.
        for submodel in self.submodels.iter_mut() {
            submodel.num_mins = submodel.iter().map(|&fidx| minorants[fidx].len()).min().unwrap_or(0);
            submodel.minorants.resize_with(submodel.num_mins, Minorant::default);
        }

        // Only minorants belonging to the first subproblem of each submodel
        // must be updated.
        //
        // We filter all indices that
        // 1. belong to a subproblem being the first in its model







|







498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
        }

        let minorants = &self.minorants;

        // Compute the number of minorants in each submodel.
        for submodel in self.submodels.iter_mut() {
            submodel.num_mins = submodel.iter().map(|&fidx| minorants[fidx].len()).min().unwrap_or(0);
            submodel.minorants.resize_with(submodel.num_mins, Default::default);
        }

        // Only minorants belonging to the first subproblem of each submodel
        // must be updated.
        //
        // We filter all indices that
        // 1. belong to a subproblem being the first in its model
541
542
543
544
545
546
547

548


549


550
551
552
553
554
555
556
                }
            })
            .collect::<Vec<_>>();

        // Compute the aggregated minorants.
        for &(_, mod_i, i) in &updateinds {
            let submodel = &mut self.submodels[mod_i];

            submodel.minorants[i] =


                Aggregatable::combine(submodel.iter().map(|&fidx| (1.0, &minorants[fidx][i].minorant)));


        }

        let submodels = &self.submodels;

        // Build quadratic term, this is <g_i, g_j> for all pairs of minorants where each g_i
        // is an aggregated minorant of a submodel.
        //







>
|
>
>
|
>
>







529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
                }
            })
            .collect::<Vec<_>>();

        // Compute the aggregated minorants.
        for &(_, mod_i, i) in &updateinds {
            let submodel = &mut self.submodels[mod_i];
            Mn::combine_to_vec(
                submodel
                    .subproblems
                    .iter()
                    .map(|&fidx| (1.0, &minorants[fidx][i].minorant)),
                &mut submodel.minorants[i],
            );
        }

        let submodels = &self.submodels;

        // Build quadratic term, this is <g_i, g_j> for all pairs of minorants where each g_i
        // is an aggregated minorant of a submodel.
        //
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
        // - j is the number of the minorant within the submodel submodel_j
        // - idx_j is the unique index of that minorant (of the first subproblem)
        for submodel_i in submodels.iter() {
            // Compute the minorant g_i for each i
            for i in 0..submodel_i.num_mins {
                // Store the computed values at the index of the first subproblem in this model.
                let idx_i = minorants[submodel_i[0]][i].index;
                let g_i = &submodel_i.minorants[i].linear;
                // Now compute the inner product with each other minorant
                // that has to be updated.
                for &(idx_j, mod_j, j) in updateinds.iter() {
                    let x = submodels[mod_j].minorants[j].linear.dot(g_i);
                    self.qterm[idx_i][idx_j] = x;
                    self.qterm[idx_j][idx_i] = x;
                }
            }
        }

        // We verify that the qterm is correct
        if cfg!(debug_assertions) {
            for submod_i in submodels.iter() {
                for i in 0..submod_i.num_mins {
                    let idx_i = self.minorants[submod_i[0]][i].index;
                    for submod_j in submodels.iter() {
                        for j in 0..submod_j.num_mins {
                            let idx_j = self.minorants[submod_j[0]][j].index;
                            let x = submod_i.minorants[i].linear.dot(&submod_j.minorants[j].linear);
                            debug_assert!((x - self.qterm[idx_i][idx_j]) < 1e-6);
                        }
                    }
                }
            }
        }








|



|














|







562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
        // - j is the number of the minorant within the submodel submodel_j
        // - idx_j is the unique index of that minorant (of the first subproblem)
        for submodel_i in submodels.iter() {
            // Compute the minorant g_i for each i
            for i in 0..submodel_i.num_mins {
                // Store the computed values at the index of the first subproblem in this model.
                let idx_i = minorants[submodel_i[0]][i].index;
                let g_i = &submodel_i.minorants[i].1;
                // Now compute the inner product with each other minorant
                // that has to be updated.
                for &(idx_j, mod_j, j) in updateinds.iter() {
                    let x = submodels[mod_j].minorants[j].linear(g_i);
                    self.qterm[idx_i][idx_j] = x;
                    self.qterm[idx_j][idx_i] = x;
                }
            }
        }

        // We verify that the qterm is correct
        if cfg!(debug_assertions) {
            for submod_i in submodels.iter() {
                for i in 0..submod_i.num_mins {
                    let idx_i = self.minorants[submod_i[0]][i].index;
                    for submod_j in submodels.iter() {
                        for j in 0..submod_j.num_mins {
                            let idx_j = self.minorants[submod_j[0]][j].index;
                            let x = submod_i.minorants[i].linear(&submod_j.minorants[j].1);
                            debug_assert!((x - self.qterm[idx_i][idx_j]) < 1e-6);
                        }
                    }
                }
            }
        }

776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
        Builder {
            max_bundle_size: 50,
            select_model: Arc::new(|i| i),
        }
    }
}

impl<P: Aggregatable + 'static> super::Builder<P> for Builder {
    type MasterProblem = CplexMaster<P>;

    fn build(&mut self) -> Result<CplexMaster<P>> {
        let mut cpx = CplexMaster::new()?;
        cpx.max_bundle_size = self.max_bundle_size;
        cpx.select_model = self.select_model.clone();
        cpx.update_submodels();
        Ok(cpx)
    }
}







|
|

|







769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
        Builder {
            max_bundle_size: 50,
            select_model: Arc::new(|i| i),
        }
    }
}

impl<Mn: Minorant> super::Builder<Mn> for Builder {
    type MasterProblem = CplexMaster<Mn>;

    fn build(&mut self) -> Result<CplexMaster<Mn>> {
        let mut cpx = CplexMaster::new()?;
        cpx.max_bundle_size = self.max_bundle_size;
        cpx.select_model = self.select_model.clone();
        cpx.update_submodels();
        Ok(cpx)
    }
}
Changes to src/master/boxed/unconstrained/minimal.rs.
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130


131




132


133
134
135
136
137
138
139
140
141



142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
 * is that this model can be solved explicitely and very quickly, but
 * it is only a very loose approximation of the objective function.
 *
 * Because of its properties, it can only be used if the problem to be
 * solved has a maximal number of minorants of two and only one
 * subproblem.
 */
pub struct MinimalMaster<P: Aggregatable> {
    /// The weight of the quadratic term.
    weight: Real,

    /// The minorants in the model.
    ///
    /// There are up to two minorants, each minorant consists of one part for
    /// each subproblem.
    ///
    /// The minorants for the i-th subproblem have the indices `2*i` and
    /// `2*i+1`.
    minorants: [Vec<Minorant<P>>; 2],
    /// The number of minorants. Only the minorants with index less than this
    /// number are valid.
    num_minorants: usize,
    /// The number of minorants for each subproblem.
    num_minorants_of: Vec<usize>,
    /// The number of subproblems.
    num_subproblems: usize,
    /// The number of subproblems with at least 1 minorant.
    num_subproblems_with_1: usize,
    /// The number of subproblems with at least 2 minorants.
    num_subproblems_with_2: usize,
    /// Optimal multipliers.
    opt_mult: [Real; 2],
    /// Optimal aggregated minorant.
    opt_minorant: Minorant<P>,
}

impl<P: Aggregatable + Send + 'static> UnconstrainedMasterProblem<P> for MinimalMaster<P> {
    type MinorantIndex = usize;

    type Err = MinimalMasterError;

    fn new() -> Result<MinimalMaster<P>, Self::Err> {
        Ok(MinimalMaster {
            weight: 1.0,
            num_minorants: 0,
            num_minorants_of: vec![],
            num_subproblems: 0,
            num_subproblems_with_1: 0,
            num_subproblems_with_2: 0,
            minorants: [vec![], vec![]],
            opt_mult: [0.0, 0.0],
            opt_minorant: Minorant::default(),
        })
    }

    fn num_subproblems(&self) -> usize {
        self.num_subproblems
    }

    fn set_num_subproblems(&mut self, n: usize) -> Result<(), Self::Err> {
        self.num_subproblems = n;
        self.num_minorants = 0;
        self.num_minorants_of = vec![0; n];
        self.num_subproblems_with_1 = 0;
        self.num_subproblems_with_2 = 0;
        self.minorants = [
            (0..n).map(|_| Minorant::default()).collect(),
            (0..n).map(|_| Minorant::default()).collect(),
        ];
        Ok(())
    }

    fn compress(&mut self) -> Result<(), Self::Err> {
        if self.num_minorants == 2 {
            debug!("Aggregate");


            debug!("  {} * {:?}", self.opt_mult[0], self.minorants[0]);




            debug!("  {} * {:?}", self.opt_mult[1], self.minorants[1]);



            self.minorants[0] = Aggregatable::combine(self.opt_mult.iter().cloned().zip(&self.minorants));
            self.opt_mult[0] = 1.0;
            self.num_minorants = 1;
            self.num_minorants_of.clear();
            self.num_minorants_of.resize(self.num_subproblems, 1);
            self.num_subproblems_with_2 = 0;

            debug!("  {:?}", self.minorants[0]);



        }
        Ok(())
    }

    fn weight(&self) -> Real {
        self.weight
    }

    fn set_weight(&mut self, weight: Real) -> Result<(), Self::Err> {
        assert!(weight > 0.0);
        self.weight = weight;
        Ok(())
    }

    fn num_minorants(&self, fidx: usize) -> usize {
        self.num_minorants_of[fidx]
    }

    fn add_minorant(&mut self, fidx: usize, minorant: Minorant<P>) -> Result<usize, Self::Err> {
        if self.num_minorants_of[fidx] >= 2 {
            return Err(MinimalMasterError::MaxMinorants { subproblem: fidx });
        }

        let minidx = self.num_minorants_of[fidx];
        self.num_minorants_of[fidx] += 1;
        self.minorants[minidx][fidx] = minorant;







|










|














|


|




|









|














|
|







>
>
|
>
>
>
>
|
>
>








|
>
>
>


















|







56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
 * is that this model can be solved explicitely and very quickly, but
 * it is only a very loose approximation of the objective function.
 *
 * Because of its properties, it can only be used if the problem to be
 * solved has a maximal number of minorants of two and only one
 * subproblem.
 */
pub struct MinimalMaster<Mn: Minorant> {
    /// The weight of the quadratic term.
    weight: Real,

    /// The minorants in the model.
    ///
    /// There are up to two minorants, each minorant consists of one part for
    /// each subproblem.
    ///
    /// The minorants for the i-th subproblem have the indices `2*i` and
    /// `2*i+1`.
    minorants: [Vec<Mn>; 2],
    /// The number of minorants. Only the minorants with index less than this
    /// number are valid.
    num_minorants: usize,
    /// The number of minorants for each subproblem.
    num_minorants_of: Vec<usize>,
    /// The number of subproblems.
    num_subproblems: usize,
    /// The number of subproblems with at least 1 minorant.
    num_subproblems_with_1: usize,
    /// The number of subproblems with at least 2 minorants.
    num_subproblems_with_2: usize,
    /// Optimal multipliers.
    opt_mult: [Real; 2],
    /// Optimal aggregated minorant.
    opt_minorant: (Real, DVector),
}

impl<Mn: Minorant> UnconstrainedMasterProblem<Mn> for MinimalMaster<Mn> {
    type MinorantIndex = usize;

    type Err = MinimalMasterError;

    fn new() -> Result<MinimalMaster<Mn>, Self::Err> {
        Ok(MinimalMaster {
            weight: 1.0,
            num_minorants: 0,
            num_minorants_of: vec![],
            num_subproblems: 0,
            num_subproblems_with_1: 0,
            num_subproblems_with_2: 0,
            minorants: [vec![], vec![]],
            opt_mult: [0.0, 0.0],
            opt_minorant: Default::default(),
        })
    }

    fn num_subproblems(&self) -> usize {
        self.num_subproblems
    }

    fn set_num_subproblems(&mut self, n: usize) -> Result<(), Self::Err> {
        self.num_subproblems = n;
        self.num_minorants = 0;
        self.num_minorants_of = vec![0; n];
        self.num_subproblems_with_1 = 0;
        self.num_subproblems_with_2 = 0;
        self.minorants = [
            (0..n).map(|_| Default::default()).collect(),
            (0..n).map(|_| Default::default()).collect(),
        ];
        Ok(())
    }

    fn compress(&mut self) -> Result<(), Self::Err> {
        if self.num_minorants == 2 {
            debug!("Aggregate");
            debug!(
                "  {} * {:?}",
                self.opt_mult[0],
                self.minorants[0].iter().map(|m| m.debug()).collect::<Vec<_>>()
            );
            debug!(
                "  {} * {:?}",
                self.opt_mult[1],
                self.minorants[1].iter().map(|m| m.debug()).collect::<Vec<_>>()
            );

            self.minorants[0] = Aggregatable::combine(self.opt_mult.iter().cloned().zip(&self.minorants));
            self.opt_mult[0] = 1.0;
            self.num_minorants = 1;
            self.num_minorants_of.clear();
            self.num_minorants_of.resize(self.num_subproblems, 1);
            self.num_subproblems_with_2 = 0;

            debug!(
                "  {:?}",
                self.minorants[0].iter().map(|m| m.debug()).collect::<Vec<_>>()
            );
        }
        Ok(())
    }

    fn weight(&self) -> Real {
        self.weight
    }

    fn set_weight(&mut self, weight: Real) -> Result<(), Self::Err> {
        assert!(weight > 0.0);
        self.weight = weight;
        Ok(())
    }

    fn num_minorants(&self, fidx: usize) -> usize {
        self.num_minorants_of[fidx]
    }

    fn add_minorant(&mut self, fidx: usize, minorant: Mn) -> Result<usize, Self::Err> {
        if self.num_minorants_of[fidx] >= 2 {
            return Err(MinimalMasterError::MaxMinorants { subproblem: fidx });
        }

        let minidx = self.num_minorants_of[fidx];
        self.num_minorants_of[fidx] += 1;
        self.minorants[minidx][fidx] = minorant;
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252

253


254
255


256

257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
            _ => unreachable!("Invalid number of minorants in subproblem {}", fidx),
        }
    }

    fn add_vars<SErr>(
        &mut self,
        nnew: usize,
        changed: &[usize],
        extend_subgradient: &mut SubgradientExtender<P, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>> {
        if self.num_subproblems_with_1 == 0 {
            return Ok(());
        }

        let noldvars = self.minorants[0][self
            .num_minorants_of
            .iter()
            .position(|&n| n > 0)
            .ok_or(MinimalMasterError::NoMinorants)
            .map_err(Either::Left)?]
        .linear
        .len();
        let mut changedvars = vec![];
        changedvars.extend_from_slice(changed);
        changedvars.extend(noldvars..noldvars + nnew);

        for fidx in 0..self.num_subproblems {
            for i in 0..self.num_minorants_of[fidx] {
                let new_subg =
                    (extend_subgradient)(fidx, &self.minorants[i][fidx].primal, &changedvars).map_err(Either::Right)?;
                let m = &mut self.minorants[i][fidx];
                for (&j, &g) in changed.iter().zip(new_subg.iter()) {
                    m.linear[j] = g;
                }
                m.linear.extend_from_slice(&new_subg[changed.len()..]);
            }
        }

        Ok(())
    }

    #[allow(unused_variables)]
    fn solve(&mut self, eta: &DVector, fbound: Real, augbound: Real, relprec: Real) -> Result<(), Self::Err> {
        for fidx in 0..self.num_subproblems {
            for i in 0..self.num_minorants_of[fidx] {
                debug!("  min(fidx:{}, i:{}) = {}", fidx, i, self.minorants[i][fidx]);
            }
        }

        if self.num_minorants == 2 {
            let min0 = Minorant::combine((0..self.num_subproblems).map(|fidx| (1.0, &self.minorants[0][fidx])));
            let min1 = Minorant::combine((0..self.num_subproblems).map(|fidx| (1.0, &self.minorants[1][fidx])));
            let xx = min0.linear.dot(&min0.linear);
            let yy = min1.linear.dot(&min1.linear);
            let xy = min0.linear.dot(&min1.linear);
            let xeta = min0.linear.dot(eta);
            let yeta = min1.linear.dot(eta);
            let a = yy - 2.0 * xy + xx;
            let b = xy - xx - yeta + xeta;

            let mut alpha2 = 0.0;
            if a > 0.0 {
                alpha2 = ((min1.constant - min0.constant) * self.weight - b) / a;
                alpha2 = alpha2.max(0.0).min(1.0);
            }
            self.opt_mult[0] = 1.0 - alpha2;
            self.opt_mult[1] = alpha2;

            self.opt_minorant = Aggregatable::combine(self.opt_mult.iter().cloned().zip([min0, min1].iter()));


        } else if self.num_minorants == 1 {
            let min0 = Aggregatable::combine((0..self.num_subproblems).map(|fidx| (1.0, &self.minorants[0][fidx])));


            self.opt_minorant = min0;

            self.opt_mult[0] = 1.0;
        } else {
            return Err(MinimalMasterError::NoMinorants);
        }

        debug!("Unrestricted");
        debug!("  opt_minorant={}", self.opt_minorant);
        debug!("  opt_mult={:?}", &self.opt_mult[0..self.num_minorants]);

        Ok(())
    }

    fn dualopt(&self) -> &DVector {
        &self.opt_minorant.linear
    }

    fn dualopt_cutval(&self) -> Real {
        self.opt_minorant.constant
    }

    fn multiplier(&self, min: usize) -> Real {
        self.opt_mult[min % 2]
    }

    fn aggregated_primal(&self, fidx: usize) -> Result<P, Self::Err> {
        Ok(P::combine(
            self.opt_mult
                .iter()
                .take(self.num_minorants_of[fidx])
                .enumerate()
                .map(|(i, alpha)| (*alpha, &self.minorants[i][fidx].primal)),
        ))
    }

    fn eval_model(&self, y: &DVector) -> Real {
        let mut result = NEG_INFINITY;
        for mins in &self.minorants[0..self.num_minorants] {
            result = result.max(mins.iter().map(|m| m.eval(y)).sum());







<
|

|



<
<
<
<
<
<
<
<
<
<
<
<


<
<
|
<
<
<
<










|




|
|
|
|
|
|
|





|




>
|
>
>

|
>
>
|
>






|






|



|






|
|




|







197
198
199
200
201
202
203

204
205
206
207
208
209












210
211


212




213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
            _ => unreachable!("Invalid number of minorants in subproblem {}", fidx),
        }
    }

    fn add_vars<SErr>(
        &mut self,
        nnew: usize,

        extend_subgradient: &mut SubgradientExtender<Mn, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>> {
        if nnew == 0 || self.num_subproblems_with_1 == 0 {
            return Ok(());
        }













        for fidx in 0..self.num_subproblems {
            for i in 0..self.num_minorants_of[fidx] {


                (extend_subgradient)(fidx, &mut self.minorants[i][fidx]).map_err(Either::Right)?;




            }
        }

        Ok(())
    }

    #[allow(unused_variables)]
    fn solve(&mut self, eta: &DVector, fbound: Real, augbound: Real, relprec: Real) -> Result<(), Self::Err> {
        for fidx in 0..self.num_subproblems {
            for i in 0..self.num_minorants_of[fidx] {
                debug!("  min(fidx:{}, i:{}) = {:?}", fidx, i, self.minorants[i][fidx].debug());
            }
        }

        if self.num_minorants == 2 {
            let min0: Mn = Aggregatable::combine((0..self.num_subproblems).map(|fidx| (1.0, &self.minorants[0][fidx])));
            let min1: Mn = Aggregatable::combine((0..self.num_subproblems).map(|fidx| (1.0, &self.minorants[1][fidx])));
            let xx = min0.product(&min0);
            let yy = min1.product(&min1);
            let xy = min0.product(&min1);
            let xeta = min0.linear(eta);
            let yeta = min1.linear(eta);
            let a = yy - 2.0 * xy + xx;
            let b = xy - xx - yeta + xeta;

            let mut alpha2 = 0.0;
            if a > 0.0 {
                alpha2 = ((min1.constant() - min0.constant()) * self.weight - b) / a;
                alpha2 = alpha2.max(0.0).min(1.0);
            }
            self.opt_mult[0] = 1.0 - alpha2;
            self.opt_mult[1] = alpha2;
            Mn::combine_to_vec(
                self.opt_mult.iter().cloned().zip([min0, min1].iter()),
                &mut self.opt_minorant,
            );
        } else if self.num_minorants == 1 {
            let mins0 = &self.minorants[0];
            Mn::combine_to_vec(
                std::iter::repeat(1.0).zip(mins0).take(self.num_subproblems),
                &mut self.opt_minorant,
            );
            self.opt_mult[0] = 1.0;
        } else {
            return Err(MinimalMasterError::NoMinorants);
        }

        debug!("Unrestricted");
        debug!("  opt_minorant={:?}", self.opt_minorant.debug());
        debug!("  opt_mult={:?}", &self.opt_mult[0..self.num_minorants]);

        Ok(())
    }

    fn dualopt(&self) -> &DVector {
        &self.opt_minorant.1
    }

    fn dualopt_cutval(&self) -> Real {
        self.opt_minorant.0
    }

    fn multiplier(&self, min: usize) -> Real {
        self.opt_mult[min % 2]
    }

    fn aggregated_primal(&self, fidx: usize) -> Result<Mn::Primal, Self::Err> {
        Ok(Aggregatable::combine(
            self.opt_mult
                .iter()
                .take(self.num_minorants_of[fidx])
                .enumerate()
                .map(|(i, alpha)| (*alpha, self.minorants[i][fidx].primal())),
        ))
    }

    fn eval_model(&self, y: &DVector) -> Real {
        let mut result = NEG_INFINITY;
        for mins in &self.minorants[0..self.num_minorants] {
            result = result.max(mins.iter().map(|m| m.eval(y)).sum());
309
310
311
312
313
314
315
316
317
318
319
320
321
322

impl Default for Builder {
    fn default() -> Self {
        Builder
    }
}

impl<P: Aggregatable + Send + 'static> super::Builder<P> for Builder {
    type MasterProblem = MinimalMaster<P>;

    fn build(&mut self) -> Result<MinimalMaster<P>, MinimalMasterError> {
        MinimalMaster::new()
    }
}







|
|

|



307
308
309
310
311
312
313
314
315
316
317
318
319
320

impl Default for Builder {
    fn default() -> Self {
        Builder
    }
}

impl<Mn: Minorant> super::Builder<Mn> for Builder {
    type MasterProblem = MinimalMaster<Mn>;

    fn build(&mut self) -> Result<MinimalMaster<Mn>, MinimalMasterError> {
        MinimalMaster::new()
    }
}
Changes to src/master/mod.rs.
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
//!   bounds), i.e. replacing $\hat{f}$ by $d \mapsto \hat{f}(d -
//!   \hat{d})$ for some given $\hat{d} \in \mathbb{R}\^n$.

pub mod boxed;
pub use self::boxed::BoxedMasterProblem;

use crate::problem::SubgradientExtender;
use crate::{Aggregatable, DVector, Minorant, Real};
use either::Either;
use std::error::Error;
use std::result::Result;

/// Trait for master problems of a proximal bundle methods.
///
/// Note that solvers are allowed to create multiple master problems potentially
/// running them in different threads. Because of this they must implement `Send
/// + 'static`, i.e. they must be quite self-contained and independent from
/// everything else.
pub trait MasterProblem<P: Aggregatable>: Send + 'static {
    /// Unique index for a minorant.
    type MinorantIndex: Copy + Eq;

    /// The error returned by the master problem.
    type Err: Error + Send + Sync;

    /// Create a new master problem.







|










|







36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
//!   bounds), i.e. replacing $\hat{f}$ by $d \mapsto \hat{f}(d -
//!   \hat{d})$ for some given $\hat{d} \in \mathbb{R}\^n$.

pub mod boxed;
pub use self::boxed::BoxedMasterProblem;

use crate::problem::SubgradientExtender;
use crate::{DVector, Minorant, Real};
use either::Either;
use std::error::Error;
use std::result::Result;

/// Trait for master problems of a proximal bundle methods.
///
/// Note that solvers are allowed to create multiple master problems potentially
/// running them in different threads. Because of this they must implement `Send
/// + 'static`, i.e. they must be quite self-contained and independent from
/// everything else.
pub trait MasterProblem<M: Minorant>: Send + 'static {
    /// Unique index for a minorant.
    type MinorantIndex: Copy + Eq;

    /// The error returned by the master problem.
    type Err: Error + Send + Sync;

    /// Create a new master problem.
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115

    /// Set the maximal number of inner iterations.
    fn set_max_updates(&mut self, max_updates: usize) -> Result<(), Self::Err>;

    /// Return the current number of inner iterations.
    fn cnt_updates(&self) -> usize;

    /// Add or move some variables with bounds.
    ///
    /// If an index is specified, existing variables are moved,
    /// otherwise new variables are generated.
    fn add_vars<SErr>(
        &mut self,
        bounds: &[(Option<usize>, Real, Real)],
        extend_subgradient: &mut SubgradientExtender<P, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>>;

    /// Add a new minorant to the model.
    ///
    /// The function returns a unique (among all minorants of all
    /// subproblems) index of the minorant. This index must remain
    /// valid until the minorant is aggregated.
    fn add_minorant(&mut self, fidx: usize, minorant: Minorant<P>) -> Result<Self::MinorantIndex, Self::Err>;

    /// Solve the master problem.
    fn solve(&mut self, cur_value: Real) -> Result<(), Self::Err>;

    /// Compress the bundle.
    ///
    /// When some minorants are compressed, the callback is called with the







|
<
<
<


|
|







|







86
87
88
89
90
91
92
93



94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112

    /// Set the maximal number of inner iterations.
    fn set_max_updates(&mut self, max_updates: usize) -> Result<(), Self::Err>;

    /// Return the current number of inner iterations.
    fn cnt_updates(&self) -> usize;

    /// Add some variables with bounds.



    fn add_vars<SErr>(
        &mut self,
        bounds: &[(Real, Real)],
        extend_subgradient: &mut SubgradientExtender<M, SErr>,
    ) -> Result<(), Either<Self::Err, SErr>>;

    /// Add a new minorant to the model.
    ///
    /// The function returns a unique (among all minorants of all
    /// subproblems) index of the minorant. This index must remain
    /// valid until the minorant is aggregated.
    fn add_minorant(&mut self, fidx: usize, minorant: M) -> Result<Self::MinorantIndex, Self::Err>;

    /// Solve the master problem.
    fn solve(&mut self, cur_value: Real) -> Result<(), Self::Err>;

    /// Compress the bundle.
    ///
    /// When some minorants are compressed, the callback is called with the
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
    /// $g\^*$ is the optimal aggregated subgradient.
    fn get_dualoptnorm2(&self) -> Real;

    /// Return the multiplier associated with a minorant.
    fn multiplier(&self, min: Self::MinorantIndex) -> Real;

    /// Return the aggregated primal.
    fn aggregated_primal(&self, fidx: usize) -> Result<P, Self::Err>;

    /// Move the center of the master problem to $\alpha \cdot d$.
    ///
    /// The vector `d` is allowed to be shorter than the current number of
    /// variables in the master. The missing elements are assumed to be zero.
    fn move_center(&mut self, alpha: Real, d: &DVector);
}

/// A builder for creating a master problem solvers.
pub trait Builder<P: Aggregatable> {
    /// The master problem to be build.
    type MasterProblem: MasterProblem<P>;

    /// Create a new master problem.
    fn build(&mut self) -> Result<Self::MasterProblem, <Self::MasterProblem as MasterProblem<P>>::Err>;
}

/// The full bundle builder.
pub type FullMasterBuilder = boxed::Builder<boxed::unconstrained::cpx::Builder>;

/// The minimal bundle builder.
pub type MinimalMasterBuilder = boxed::Builder<boxed::unconstrained::minimal::Builder>;







|









|

|


|







129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
    /// $g\^*$ is the optimal aggregated subgradient.
    fn get_dualoptnorm2(&self) -> Real;

    /// Return the multiplier associated with a minorant.
    fn multiplier(&self, min: Self::MinorantIndex) -> Real;

    /// Return the aggregated primal.
    fn aggregated_primal(&self, fidx: usize) -> Result<M::Primal, Self::Err>;

    /// Move the center of the master problem to $\alpha \cdot d$.
    ///
    /// The vector `d` is allowed to be shorter than the current number of
    /// variables in the master. The missing elements are assumed to be zero.
    fn move_center(&mut self, alpha: Real, d: &DVector);
}

/// A builder for creating a master problem solvers.
pub trait Builder<M: Minorant> {
    /// The master problem to be build.
    type MasterProblem: MasterProblem<M>;

    /// Create a new master problem.
    fn build(&mut self) -> Result<Self::MasterProblem, <Self::MasterProblem as MasterProblem<M>>::Err>;
}

/// The full bundle builder.
pub type FullMasterBuilder = boxed::Builder<boxed::unconstrained::cpx::Builder>;

/// The minimal bundle builder.
pub type MinimalMasterBuilder = boxed::Builder<boxed::unconstrained::minimal::Builder>;
Changes to src/mcf/problem.rs.
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
        aggr
    }
}

impl ParallelProblem for MMCFProblem {
    type Err = Error;

    type Primal = DVector;

    fn num_variables(&self) -> usize {
        self.active_constraints.read().unwrap().len()
    }

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







|







417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
        aggr
    }
}

impl ParallelProblem for MMCFProblem {
    type Err = Error;

    type Minorant = (Real, DVector, DVector);

    fn num_variables(&self) -> usize {
        self.active_constraints.read().unwrap().len()
    }

    fn lower_bounds(&self) -> Option<Vec<Real>> {
        Some(vec![0.0; self.num_variables()])
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
        let active_constraints = self.active_constraints.read().unwrap()[0..y.len()].to_vec();
        self.pool
            .as_ref()
            .unwrap()
            .execute(move || match sub.write().unwrap().evaluate(&y, active_constraints) {
                Ok((objective, subg, primal)) => {
                    tx.objective(objective).unwrap();
                    tx.minorant(Minorant {
                        constant: objective,
                        linear: subg,
                        primal,
                    })
                    .unwrap();
                }
                Err(err) => tx.error(err).unwrap(),
            });
        Ok(())
    }

    fn update<U, S>(&mut self, state: U, tx: S) -> Result<()>
    where
        U: ParallelUpdateState<Self::Primal>,
        S: UpdateSender<Self> + 'static,
    {
        if self.inactive_constraints.read().unwrap().is_empty() {
            return Ok(());
        }

        let inactive_constraints = self.inactive_constraints.clone();







|
<
<
<
<
<








|







458
459
460
461
462
463
464
465





466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
        let active_constraints = self.active_constraints.read().unwrap()[0..y.len()].to_vec();
        self.pool
            .as_ref()
            .unwrap()
            .execute(move || match sub.write().unwrap().evaluate(&y, active_constraints) {
                Ok((objective, subg, primal)) => {
                    tx.objective(objective).unwrap();
                    tx.minorant((objective, subg, primal)).unwrap();





                }
                Err(err) => tx.error(err).unwrap(),
            });
        Ok(())
    }

    fn update<U, S>(&mut self, state: U, tx: S) -> Result<()>
    where
        U: ParallelUpdateState<<Self::Minorant as Minorant>::Primal>,
        S: UpdateSender<Self> + 'static,
    {
        if self.inactive_constraints.read().unwrap().is_empty() {
            return Ok(());
        }

        let inactive_constraints = self.inactive_constraints.clone();
525
526
527
528
529
530
531
532
533
534
535

536
537
538
539
540
541
542
543
544
545
546
547
            //     }));

            let nnew = active.len() - nold;
            if nnew > 0 {
                let actives = active.clone();
                if let Err(err) = tx.add_variables(
                    vec![(0.0, Real::infinity()); nnew],
                    Box::new(move |fidx: usize, primal: &DVector, inds: &[usize]| {
                        let sub = subs[fidx].read().unwrap();
                        Ok(inds
                            .iter()

                            .map(|&i| sub.evaluate_constraint(primal, actives[i]))
                            .collect())
                    }),
                ) {
                    warn!("Error sending problem update: {}", err);
                }
            }
        });

        Ok(())
    }
}







|

<
|
>
|
|










520
521
522
523
524
525
526
527
528

529
530
531
532
533
534
535
536
537
538
539
540
541
542
            //     }));

            let nnew = active.len() - nold;
            if nnew > 0 {
                let actives = active.clone();
                if let Err(err) = tx.add_variables(
                    vec![(0.0, Real::infinity()); nnew],
                    Box::new(move |fidx: usize, m: &mut Self::Minorant| {
                        let sub = subs[fidx].read().unwrap();

                        let n = m.dim();
                        let primal = &m.2;
                        m.1.extend((n..n + nnew).map(|i| sub.evaluate_constraint(primal, actives[i])));
                        Ok(())
                    }),
                ) {
                    warn!("Error sending problem update: {}", err);
                }
            }
        });

        Ok(())
    }
}
Changes to src/problem.rs.
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see  <http://www.gnu.org/licenses/>
 */

//! An asynchronous first-order oracle.

pub use crate::data::minorant::SubgradientExtender;
use crate::{Aggregatable, DVector, Minorant, Real};
use std::sync::Arc;

/// Channel to send evaluation results to.
pub trait ResultSender<P>: Send
where
    P: FirstOrderProblem,
{
    type Err: std::error::Error + Send;

    /// Send a new objective `value`.
    fn objective(&self, value: Real) -> Result<(), Self::Err>;

    /// Send a new `minorant` with associated `primal`.
    fn minorant(&self, minorant: Minorant<P::Primal>) -> Result<(), Self::Err>;

    /// Send an error message.
    fn error(&self, err: P::Err) -> Result<(), Self::Err>;
}

/// Channel to send problem updates to.
pub trait UpdateSender<P>: Send
where
    P: FirstOrderProblem,
{
    type Err: std::error::Error + Send;

    /// Add new variables with given bounds.
    ///
    /// The oracle must provide a "subgradient extension callback", that
    /// computes subgradients for the given variables based on the associated
    /// primals.
    fn add_variables(
        &self,
        bounds: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<P::Primal, P::Err>>,
    ) -> Result<(), Self::Err>;

    /// Send an error message.
    fn error(&self, err: P::Err) -> Result<(), Self::Err>;
}

/// This trait provides information available in the







|













|




















|







14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see  <http://www.gnu.org/licenses/>
 */

//! An asynchronous first-order oracle.

pub use crate::data::minorant::SubgradientExtender;
use crate::{DVector, Minorant, Real};
use std::sync::Arc;

/// Channel to send evaluation results to.
pub trait ResultSender<P>: Send
where
    P: FirstOrderProblem,
{
    type Err: std::error::Error + Send;

    /// Send a new objective `value`.
    fn objective(&self, value: Real) -> Result<(), Self::Err>;

    /// Send a new `minorant` with associated `primal`.
    fn minorant(&self, minorant: P::Minorant) -> Result<(), Self::Err>;

    /// Send an error message.
    fn error(&self, err: P::Err) -> Result<(), Self::Err>;
}

/// Channel to send problem updates to.
pub trait UpdateSender<P>: Send
where
    P: FirstOrderProblem,
{
    type Err: std::error::Error + Send;

    /// Add new variables with given bounds.
    ///
    /// The oracle must provide a "subgradient extension callback", that
    /// computes subgradients for the given variables based on the associated
    /// primals.
    fn add_variables(
        &self,
        bounds: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<P::Minorant, P::Err>>,
    ) -> Result<(), Self::Err>;

    /// Send an error message.
    fn error(&self, err: P::Err) -> Result<(), Self::Err>;
}

/// This trait provides information available in the
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
///
/// Note that all computations made by an implementation are supposed to be
/// asynchronous.
pub trait FirstOrderProblem {
    /// Error raised by this oracle.
    type Err;

    /// The primal information associated with a minorant.
    type Primal: Aggregatable + Send + 'static;

    /// Return the number of variables.
    fn num_variables(&self) -> usize;

    /// Return the lower bounds on the variables.
    ///
    /// If no lower bounds are specified, $-\infty$ is assumed.







|
|







90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
///
/// Note that all computations made by an implementation are supposed to be
/// asynchronous.
pub trait FirstOrderProblem {
    /// Error raised by this oracle.
    type Err;

    /// The minorant information.
    type Minorant: Minorant;

    /// Return the number of variables.
    fn num_variables(&self) -> usize;

    /// Return the lower bounds on the variables.
    ///
    /// If no lower bounds are specified, $-\infty$ is assumed.
173
174
175
176
177
178
179
180
181
182
183
184
185
186
    /// information (e.g. adding new variables) to the provided channel.
    ///
    /// The updates might be generated asynchronously.
    ///
    /// The default implementation does nothing.
    fn update<U, S>(&mut self, _state: U, _tx: S) -> Result<(), Self::Err>
    where
        U: UpdateState<Self::Primal>,
        S: UpdateSender<Self> + 'static,
        Self: Sized,
    {
        Ok(())
    }
}







|






173
174
175
176
177
178
179
180
181
182
183
184
185
186
    /// information (e.g. adding new variables) to the provided channel.
    ///
    /// The updates might be generated asynchronously.
    ///
    /// The default implementation does nothing.
    fn update<U, S>(&mut self, _state: U, _tx: S) -> Result<(), Self::Err>
    where
        U: UpdateState<<Self::Minorant as Minorant>::Primal>,
        S: UpdateSender<Self> + 'static,
        Self: Sized,
    {
        Ok(())
    }
}
Changes to src/solver/asyn.rs.
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
use num_cpus;
use num_traits::Float;
use std::iter::repeat;
use std::sync::Arc;
use std::time::Instant;
use threadpool::ThreadPool;

use crate::{DVector, Real};

use super::channels::{
    ChannelResultSender, ChannelUpdateSender, ClientReceiver, ClientSender, EvalResult, Message, Update,
};
use super::masterprocess::{self, MasterConfig, MasterProcess, MasterResponse, Response};
use crate::master::{Builder as MasterBuilder, MasterProblem};
use crate::problem::{FirstOrderProblem, UpdateState};







|







26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
use num_cpus;
use num_traits::Float;
use std::iter::repeat;
use std::sync::Arc;
use std::time::Instant;
use threadpool::ThreadPool;

use crate::{DVector, Minorant, Real};

use super::channels::{
    ChannelResultSender, ChannelUpdateSender, ClientReceiver, ClientSender, EvalResult, Message, Update,
};
use super::masterprocess::{self, MasterConfig, MasterProcess, MasterResponse, Response};
use crate::master::{Builder as MasterBuilder, MasterProblem};
use crate::problem::{FirstOrderProblem, UpdateState};
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
///
/// - `T` is the value type,
/// - `P` is the `FirstOrderProblem` associated with the solver,
/// - `M` is the `MasterBuilder` associated with the solver.
pub type Result<T, P, M> = std::result::Result<
    T,
    Error<
        <<M as MasterBuilder<<P as FirstOrderProblem>::Primal>>::MasterProblem as MasterProblem<
            <P as FirstOrderProblem>::Primal,
        >>::Err,
        <P as FirstOrderProblem>::Err,
    >,
>;

impl<MErr, PErr> std::fmt::Display for Error<MErr, PErr>
where







|
|







85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
///
/// - `T` is the value type,
/// - `P` is the `FirstOrderProblem` associated with the solver,
/// - `M` is the `MasterBuilder` associated with the solver.
pub type Result<T, P, M> = std::result::Result<
    T,
    Error<
        <<M as MasterBuilder<<P as FirstOrderProblem>::Minorant>>::MasterProblem as MasterProblem<
            <P as FirstOrderProblem>::Minorant,
        >>::Err,
        <P as FirstOrderProblem>::Err,
    >,
>;

impl<MErr, PErr> std::fmt::Display for Error<MErr, PErr>
where
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
}

/// Implementation of a parallel bundle method.
pub struct Solver<P, T = StandardTerminator, W = HKWeighter, M = crate::master::FullMasterBuilder>
where
    P: FirstOrderProblem,
    P::Err: 'static,
    M: MasterBuilder<P::Primal>,
{
    /// Parameters for the solver.
    pub params: Parameters,

    /// Termination predicate.
    pub terminator: T,








|







395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
}

/// Implementation of a parallel bundle method.
pub struct Solver<P, T = StandardTerminator, W = HKWeighter, M = crate::master::FullMasterBuilder>
where
    P: FirstOrderProblem,
    P::Err: 'static,
    M: MasterBuilder<P::Minorant>,
{
    /// Parameters for the solver.
    pub params: Parameters,

    /// Termination predicate.
    pub terminator: T,

448
449
450
451
452
453
454
455
456
457
458
459
460
461
462

impl<P, T, W, M> Solver<P, T, W, M>
where
    P: FirstOrderProblem,
    P::Err: Send + 'static,
    T: Terminator<SolverData> + Default,
    W: Weighter<SolverData> + Default,
    M: MasterBuilder<P::Primal>,
{
    /// Create a new parallel bundle solver.
    pub fn new(problem: P) -> Self
    where
        M: Default,
    {
        Self::with_master(problem, M::default())







|







448
449
450
451
452
453
454
455
456
457
458
459
460
461
462

impl<P, T, W, M> Solver<P, T, W, M>
where
    P: FirstOrderProblem,
    P::Err: Send + 'static,
    T: Terminator<SolverData> + Default,
    W: Weighter<SolverData> + Default,
    M: MasterBuilder<P::Minorant>,
{
    /// Create a new parallel bundle solver.
    pub fn new(problem: P) -> Self
    where
        M: Default,
    {
        Self::with_master(problem, M::default())
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
    /// variables), and the master problem is started to compute a new
    /// candidate.
    ///
    /// The function returns
    ///   - `Ok(true)` if the final iteration count has been reached,
    ///   - `Ok(false)` if the final iteration count has not been reached,
    ///   - `Err(_)` on error.
    fn handle_client_response(&mut self, msg: EvalResult<usize, P::Primal, P::Err>) -> Result<bool, P, M> {
        let master = self.master_proc.as_mut().ok_or(Error::NotInitialized)?;
        match msg {
            EvalResult::ObjectiveValue { index, value } => {
                debug!("Receive objective from subproblem {}: {}", index, value);
                if self.data.nxt_ubs[index].is_infinite() {
                    self.data.cnt_remaining_ubs -= 1;
                }
                self.data.nxt_ubs[index] = self.data.nxt_ubs[index].min(value);
            }
            EvalResult::Minorant { index, mut minorant } => {
                debug!("Receive minorant from subproblem {}", index);
                if self.data.nxt_cutvals[index].is_infinite() {
                    self.data.cnt_remaining_mins -= 1;
                }
                // move center of minorant to cur_y
                minorant.move_center(-1.0, &self.data.nxt_d);
                self.data.nxt_cutvals[index] = self.data.nxt_cutvals[index].max(minorant.constant);
                // add minorant to master problem
                master.add_minorant(index, minorant)?;
            }
            EvalResult::Done { .. } => return Ok(false), // currently we do nothing ...
            EvalResult::Error { err, .. } => return Err(Error::Evaluation(err)),
        }








|
















|







678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
    /// variables), and the master problem is started to compute a new
    /// candidate.
    ///
    /// The function returns
    ///   - `Ok(true)` if the final iteration count has been reached,
    ///   - `Ok(false)` if the final iteration count has not been reached,
    ///   - `Err(_)` on error.
    fn handle_client_response(&mut self, msg: EvalResult<usize, P::Minorant, P::Err>) -> Result<bool, P, M> {
        let master = self.master_proc.as_mut().ok_or(Error::NotInitialized)?;
        match msg {
            EvalResult::ObjectiveValue { index, value } => {
                debug!("Receive objective from subproblem {}: {}", index, value);
                if self.data.nxt_ubs[index].is_infinite() {
                    self.data.cnt_remaining_ubs -= 1;
                }
                self.data.nxt_ubs[index] = self.data.nxt_ubs[index].min(value);
            }
            EvalResult::Minorant { index, mut minorant } => {
                debug!("Receive minorant from subproblem {}", index);
                if self.data.nxt_cutvals[index].is_infinite() {
                    self.data.cnt_remaining_mins -= 1;
                }
                // move center of minorant to cur_y
                minorant.move_center(-1.0, &self.data.nxt_d);
                self.data.nxt_cutvals[index] = self.data.nxt_cutvals[index].max(minorant.constant());
                // add minorant to master problem
                master.add_minorant(index, minorant)?;
            }
            EvalResult::Done { .. } => return Ok(false), // currently we do nothing ...
            EvalResult::Error { err, .. } => return Err(Error::Evaluation(err)),
        }

926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
    }

    /// Handles an update response `update` of the problem.
    ///
    /// The method is called if the problem informs the solver about a change of
    /// the problem, e.g. adding a new variable. This method updates the other
    /// parts of the solver, e.g. the master problem, about the modification.
    fn handle_update_response(&mut self, update: Update<usize, P::Primal, P::Err>) -> Result<(), P, M> {
        match update {
            Update::AddVariables { bounds, sgext, .. } => {
                if !bounds.is_empty() {
                    // add new variables
                    self.master_proc
                        .as_mut()
                        .unwrap()
                        .add_vars(bounds.into_iter().map(|(l, u)| (None, l, u)).collect(), sgext)?;
                }
            }
            Update::Done { .. } => self.data.update_in_progress = false, // currently we do nothing ...
            Update::Error { err, .. } => {
                self.data.update_in_progress = false;
                return Err(Error::Update(err));
            }







|







|







926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
    }

    /// Handles an update response `update` of the problem.
    ///
    /// The method is called if the problem informs the solver about a change of
    /// the problem, e.g. adding a new variable. This method updates the other
    /// parts of the solver, e.g. the master problem, about the modification.
    fn handle_update_response(&mut self, update: Update<usize, P::Minorant, P::Err>) -> Result<(), P, M> {
        match update {
            Update::AddVariables { bounds, sgext, .. } => {
                if !bounds.is_empty() {
                    // add new variables
                    self.master_proc
                        .as_mut()
                        .unwrap()
                        .add_vars(bounds.into_iter().map(|(l, u)| (l, u)).collect(), sgext)?;
                }
            }
            Update::Done { .. } => self.data.update_in_progress = false, // currently we do nothing ...
            Update::Error { err, .. } => {
                self.data.update_in_progress = false;
                return Err(Error::Update(err));
            }
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
            self.data.nxt_mod,
            self.data.nxt_val,
            self.data.cur_val
        );
    }

    /// Return the aggregated primal of the given subproblem.
    pub fn aggregated_primal(&self, subproblem: usize) -> Result<P::Primal, P, M> {
        Ok(self
            .master_proc
            .as_ref()
            .ok_or(Error::NotSolved)?
            .get_aggregated_primal(subproblem)?)
    }
}







|







983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
            self.data.nxt_mod,
            self.data.nxt_val,
            self.data.cur_val
        );
    }

    /// Return the aggregated primal of the given subproblem.
    pub fn aggregated_primal(&self, subproblem: usize) -> Result<<P::Minorant as Minorant>::Primal, P, M> {
        Ok(self
            .master_proc
            .as_ref()
            .ok_or(Error::NotSolved)?
            .get_aggregated_primal(subproblem)?)
    }
}
Changes to src/solver/channels.rs.
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98



99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
 */

//! Implementation for `ResultSender` using channels.

use super::masterprocess::Response as MasterResponse;
use crate::master::MasterProblem;
use crate::problem::{FirstOrderProblem, ResultSender, SubgradientExtender, UpdateSender};
use crate::{Aggregatable, Minorant, Real};

#[cfg(feature = "crossbeam")]
use rs_crossbeam::channel::{Receiver, SendError, Sender};
#[cfg(not(feature = "crossbeam"))]
use std::sync::mpsc::{Receiver, SendError, Sender};

/// A message received from a client.
///
/// The client might be either a subproblem or a master problem process.
///
/// A client may send one of two possible messages: either an evaluation result
/// or a problem update.
pub enum Message<I, P, E, MErr>
where
    P: Aggregatable,
{
    /// An evaluation result.
    Eval(EvalResult<I, P, E>),
    /// A problem update.
    Update(Update<I, P, E>),
    /// A master problem result.
    Master(MasterResponse<MErr, E>),
}

/// The sender for client messages.
pub type ClientSender<P, M> = Sender<
    Message<
        usize,
        <P as FirstOrderProblem>::Primal,
        <P as FirstOrderProblem>::Err,
        <M as MasterProblem<<P as FirstOrderProblem>::Primal>>::Err,
    >,
>;

/// The receiver for client messages.
pub type ClientReceiver<P, M> = Receiver<
    Message<
        usize,
        <P as FirstOrderProblem>::Primal,
        <P as FirstOrderProblem>::Err,
        <M as MasterProblem<<P as FirstOrderProblem>::Primal>>::Err,
    >,
>;

/// Parameters for tuning the solver.

/// Evaluation result.
///
/// The result of an evaluation is new information to be made
/// available to the solver and the master problem. There are
/// essentially two types of information:
///
///    1. The (exact) function value of a sub-function at some point.
///    2. A minorant of some sub-function.
#[derive(Debug)]
pub enum EvalResult<I, P, E>
where
    P: Aggregatable,
{
    /// The objective value at some point.
    ObjectiveValue { index: I, value: Real },
    /// A minorant with an associated primal.
    Minorant { index: I, minorant: Minorant<P> },
    /// An error occurred.
    Error { index: I, err: E },
    /// Done, no further message will be received from this evaluation.
    Done { index: I },
}

/// Problem update information.
///
/// The solver calls the `update` method of the problem regularly.
/// This method can modify the problem by adding (or moving)
/// variables. The possible updates are encoded in this type.
pub enum Update<I, P, E> {



    /// Add new variables with bounds.
    ///
    /// The initial value of each variable will be the feasible value
    /// closest to 0.
    AddVariables {
        index: I,
        bounds: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<P, E>>,
    },
    /// An error occurred.
    Error { index: I, err: E },
    /// Done, no further message will be received from this evaluation.
    Done { index: I },
}

pub struct ChannelResultSender<I, P, E, M>
where
    I: Clone,
    P: Aggregatable,
{
    index: I,
    tx: Sender<Message<I, P, E, M>>,
}

impl<I, P, E, M> ChannelResultSender<I, P, E, M>
where
    I: Clone,
    P: Aggregatable,
{
    pub fn new(index: I, tx: Sender<Message<I, P, E, M>>) -> Self {
        ChannelResultSender { index, tx }
    }
}

impl<I, P, M> ResultSender<P> for ChannelResultSender<I, P::Primal, P::Err, M>
where
    I: Clone + Send,
    P: FirstOrderProblem,
    P::Err: Send,
    M: Send,
{
    type Err = SendError<Message<I, P::Primal, P::Err, M>>;

    fn objective(&self, value: Real) -> Result<(), Self::Err> {
        self.tx.send(Message::Eval(EvalResult::ObjectiveValue {
            index: self.index.clone(),
            value,
        }))
    }

    fn minorant(&self, minorant: Minorant<P::Primal>) -> Result<(), Self::Err> {
        self.tx.send(Message::Eval(EvalResult::Minorant {
            index: self.index.clone(),
            minorant,
        }))
    }

    fn error(&self, err: P::Err) -> Result<(), Self::Err> {
        self.tx.send(Message::Eval(EvalResult::Error {
            index: self.index.clone(),
            err,
        }))
    }
}

impl<I, P, E, M> Drop for ChannelResultSender<I, P, E, M>
where
    I: Clone,
    P: Aggregatable,
{
    fn drop(&mut self) {
        self.tx
            .send(Message::Eval(EvalResult::Done {
                index: self.index.clone(),
            }))
            .unwrap()
    }
}

pub struct ChannelUpdateSender<I, P, E, M>
where
    I: Clone,
    P: Aggregatable,
{
    index: I,
    tx: Sender<Message<I, P, E, M>>,
}

impl<I, P, E, M> ChannelUpdateSender<I, P, E, M>
where
    I: Clone,
    P: Aggregatable,
{
    pub fn new(index: I, tx: Sender<Message<I, P, E, M>>) -> Self {
        ChannelUpdateSender { index, tx }
    }
}

impl<I, P, M> UpdateSender<P> for ChannelUpdateSender<I, P::Primal, P::Err, M>
where
    I: Clone + Send,
    P: FirstOrderProblem,
    P::Err: Send,
    M: Send,
{
    type Err = SendError<Message<I, P::Primal, P::Err, M>>;

    fn add_variables(
        &self,
        bounds: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<P::Primal, P::Err>>,
    ) -> Result<(), Self::Err> {
        self.tx.send(Message::Update(Update::AddVariables {
            index: self.index.clone(),
            bounds,
            sgext,
        }))
    }

    fn error(&self, err: P::Err) -> Result<(), Self::Err> {
        self.tx.send(Message::Update(Update::Error {
            index: self.index.clone(),
            err,
        }))
    }
}

impl<I, P, E, M> Drop for ChannelUpdateSender<I, P, E, M>
where
    I: Clone,
    P: Aggregatable,
{
    fn drop(&mut self) {
        self.tx
            .send(Message::Update(Update::Done {
                index: self.index.clone(),
            }))
            .unwrap()
    }
}







|












|

|


|

|








|

|







|

|














|

|




|











|
>
>
>







|







|


|


|


|


|

|




|






|








|














|


|










|


|


|


|


|

|




|






|




|
















|


|









16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
 */

//! Implementation for `ResultSender` using channels.

use super::masterprocess::Response as MasterResponse;
use crate::master::MasterProblem;
use crate::problem::{FirstOrderProblem, ResultSender, SubgradientExtender, UpdateSender};
use crate::{Minorant, Real};

#[cfg(feature = "crossbeam")]
use rs_crossbeam::channel::{Receiver, SendError, Sender};
#[cfg(not(feature = "crossbeam"))]
use std::sync::mpsc::{Receiver, SendError, Sender};

/// A message received from a client.
///
/// The client might be either a subproblem or a master problem process.
///
/// A client may send one of two possible messages: either an evaluation result
/// or a problem update.
pub enum Message<I, Mn, E, MErr>
where
    Mn: Minorant,
{
    /// An evaluation result.
    Eval(EvalResult<I, Mn, E>),
    /// A problem update.
    Update(Update<I, Mn, E>),
    /// A master problem result.
    Master(MasterResponse<MErr, E>),
}

/// The sender for client messages.
pub type ClientSender<P, M> = Sender<
    Message<
        usize,
        <P as FirstOrderProblem>::Minorant,
        <P as FirstOrderProblem>::Err,
        <M as MasterProblem<<P as FirstOrderProblem>::Minorant>>::Err,
    >,
>;

/// The receiver for client messages.
pub type ClientReceiver<P, M> = Receiver<
    Message<
        usize,
        <P as FirstOrderProblem>::Minorant,
        <P as FirstOrderProblem>::Err,
        <M as MasterProblem<<P as FirstOrderProblem>::Minorant>>::Err,
    >,
>;

/// Parameters for tuning the solver.

/// Evaluation result.
///
/// The result of an evaluation is new information to be made
/// available to the solver and the master problem. There are
/// essentially two types of information:
///
///    1. The (exact) function value of a sub-function at some point.
///    2. A minorant of some sub-function.
#[derive(Debug)]
pub enum EvalResult<I, Mn, E>
where
    Mn: Minorant,
{
    /// The objective value at some point.
    ObjectiveValue { index: I, value: Real },
    /// A minorant with an associated primal.
    Minorant { index: I, minorant: Mn },
    /// An error occurred.
    Error { index: I, err: E },
    /// Done, no further message will be received from this evaluation.
    Done { index: I },
}

/// Problem update information.
///
/// The solver calls the `update` method of the problem regularly.
/// This method can modify the problem by adding (or moving)
/// variables. The possible updates are encoded in this type.
pub enum Update<I, Mn, E>
where
    Mn: Minorant,
{
    /// Add new variables with bounds.
    ///
    /// The initial value of each variable will be the feasible value
    /// closest to 0.
    AddVariables {
        index: I,
        bounds: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<Mn, E>>,
    },
    /// An error occurred.
    Error { index: I, err: E },
    /// Done, no further message will be received from this evaluation.
    Done { index: I },
}

pub struct ChannelResultSender<I, Mn, E, M>
where
    I: Clone,
    Mn: Minorant,
{
    index: I,
    tx: Sender<Message<I, Mn, E, M>>,
}

impl<I, Mn, E, M> ChannelResultSender<I, Mn, E, M>
where
    I: Clone,
    Mn: Minorant,
{
    pub fn new(index: I, tx: Sender<Message<I, Mn, E, M>>) -> Self {
        ChannelResultSender { index, tx }
    }
}

impl<I, P, M> ResultSender<P> for ChannelResultSender<I, P::Minorant, P::Err, M>
where
    I: Clone + Send,
    P: FirstOrderProblem,
    P::Err: Send,
    M: Send,
{
    type Err = SendError<Message<I, P::Minorant, P::Err, M>>;

    fn objective(&self, value: Real) -> Result<(), Self::Err> {
        self.tx.send(Message::Eval(EvalResult::ObjectiveValue {
            index: self.index.clone(),
            value,
        }))
    }

    fn minorant(&self, minorant: P::Minorant) -> Result<(), Self::Err> {
        self.tx.send(Message::Eval(EvalResult::Minorant {
            index: self.index.clone(),
            minorant,
        }))
    }

    fn error(&self, err: P::Err) -> Result<(), Self::Err> {
        self.tx.send(Message::Eval(EvalResult::Error {
            index: self.index.clone(),
            err,
        }))
    }
}

impl<I, Mn, E, M> Drop for ChannelResultSender<I, Mn, E, M>
where
    I: Clone,
    Mn: Minorant,
{
    fn drop(&mut self) {
        self.tx
            .send(Message::Eval(EvalResult::Done {
                index: self.index.clone(),
            }))
            .unwrap()
    }
}

pub struct ChannelUpdateSender<I, Mn, E, M>
where
    I: Clone,
    Mn: Minorant,
{
    index: I,
    tx: Sender<Message<I, Mn, E, M>>,
}

impl<I, Mn, E, M> ChannelUpdateSender<I, Mn, E, M>
where
    I: Clone,
    Mn: Minorant,
{
    pub fn new(index: I, tx: Sender<Message<I, Mn, E, M>>) -> Self {
        ChannelUpdateSender { index, tx }
    }
}

impl<I, P, M> UpdateSender<P> for ChannelUpdateSender<I, P::Minorant, P::Err, M>
where
    I: Clone + Send,
    P: FirstOrderProblem,
    P::Err: Send,
    M: Send,
{
    type Err = SendError<Message<I, P::Minorant, P::Err, M>>;

    fn add_variables(
        &self,
        bounds: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<P::Minorant, P::Err>>,
    ) -> Result<(), Self::Err> {
        self.tx.send(Message::Update(Update::AddVariables {
            index: self.index.clone(),
            bounds,
            sgext,
        }))
    }

    fn error(&self, err: P::Err) -> Result<(), Self::Err> {
        self.tx.send(Message::Update(Update::Error {
            index: self.index.clone(),
            err,
        }))
    }
}

impl<I, Mn, E, M> Drop for ChannelUpdateSender<I, Mn, E, M>
where
    I: Clone,
    Mn: Minorant,
{
    fn drop(&mut self) {
        self.tx
            .send(Message::Update(Update::Done {
                index: self.index.clone(),
            }))
            .unwrap()
    }
}
Changes to src/solver/masterprocess.rs.
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
use log::{debug, warn};
use std::sync::Arc;
use threadpool::ThreadPool;

use super::channels::{ClientSender, Message};
use crate::master::MasterProblem;
use crate::problem::{FirstOrderProblem, SubgradientExtender};
use crate::{Aggregatable, DVector, Minorant, Real};

#[derive(Debug)]
pub enum Error<MErr, PErr> {
    /// The communication channel for sending requests has been disconnected.
    DisconnectedSender,
    /// The communication channel for sending responds has been disconnected.
    DisconnectedReceiver,







|







26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
use log::{debug, warn};
use std::sync::Arc;
use threadpool::ThreadPool;

use super::channels::{ClientSender, Message};
use crate::master::MasterProblem;
use crate::problem::{FirstOrderProblem, SubgradientExtender};
use crate::{DVector, Minorant, Real};

#[derive(Debug)]
pub enum Error<MErr, PErr> {
    /// The communication channel for sending requests has been disconnected.
    DisconnectedSender,
    /// The communication channel for sending responds has been disconnected.
    DisconnectedReceiver,
104
105
106
107
108
109
110
111
112

113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
    /// The lower bounds on the variables.
    pub lower_bounds: Option<DVector>,
    /// The lower bounds on the variables.
    pub upper_bounds: Option<DVector>,
}

/// A task for the master problem.
enum MasterTask<Pr, PErr, M>
where

    M: MasterProblem<Pr>,
    Pr: Aggregatable,
{
    /// Add new variables to the master problem.
    AddVariables(Vec<(Option<usize>, Real, Real)>, Box<SubgradientExtender<Pr, PErr>>),

    /// Add a new minorant for a subfunction to the master problem.
    AddMinorant(usize, Minorant<Pr>),

    /// Move the center of the master problem in the given direction.
    MoveCenter(Real, Arc<DVector>),

    /// Start a new computation of the master problem.
    Solve { center_value: Real },

    /// Compress the bundle.
    Compress,

    /// Set the weight parameter of the master problem.
    SetWeight { weight: Real },

    /// Return the current aggregated primal.
    GetAggregatedPrimal {
        subproblem: usize,
        tx: Sender<Result<Pr, M::Err>>,
    },
}

/// The response send from a master process.
///
/// The response contains the evaluation results of the latest call to `solve`.
pub enum Response<MErr, PErr> {







|

>
|
<

|
|


|
















|







104
105
106
107
108
109
110
111
112
113
114

115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
    /// The lower bounds on the variables.
    pub lower_bounds: Option<DVector>,
    /// The lower bounds on the variables.
    pub upper_bounds: Option<DVector>,
}

/// A task for the master problem.
enum MasterTask<Mn, PErr, M>
where
    Mn: Minorant,
    M: MasterProblem<Mn>,

{
    /// Add new variables with bounds to the master problem.
    AddVariables(Vec<(Real, Real)>, Box<SubgradientExtender<Mn, PErr>>),

    /// Add a new minorant for a subfunction to the master problem.
    AddMinorant(usize, Mn),

    /// Move the center of the master problem in the given direction.
    MoveCenter(Real, Arc<DVector>),

    /// Start a new computation of the master problem.
    Solve { center_value: Real },

    /// Compress the bundle.
    Compress,

    /// Set the weight parameter of the master problem.
    SetWeight { weight: Real },

    /// Return the current aggregated primal.
    GetAggregatedPrimal {
        subproblem: usize,
        tx: Sender<Result<Mn::Primal, M::Err>>,
    },
}

/// The response send from a master process.
///
/// The response contains the evaluation results of the latest call to `solve`.
pub enum Response<MErr, PErr> {
152
153
154
155
156
157
158
159
160
161
162
163

164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
    },
    /// An error occurred.
    Error(Error<MErr, PErr>),
}

/// Convenient type for `Response`.
pub type MasterResponse<P, M> =
    Response<<M as MasterProblem<<P as FirstOrderProblem>::Primal>>::Err, <P as FirstOrderProblem>::Err>;

type ToMasterSender<P, M> = Sender<MasterTask<<P as FirstOrderProblem>::Primal, <P as FirstOrderProblem>::Err, M>>;

type ToMasterReceiver<P, M> = Receiver<MasterTask<<P as FirstOrderProblem>::Primal, <P as FirstOrderProblem>::Err, M>>;


pub struct MasterProcess<P, M>
where
    P: FirstOrderProblem,
    M: MasterProblem<P::Primal>,
{
    /// The channel to transmit new tasks to the master problem.
    tx: ToMasterSender<P, M>,

    phantom: std::marker::PhantomData<M>,
}

impl<P, M> MasterProcess<P, M>
where
    P: FirstOrderProblem,
    P::Primal: Send + 'static,
    P::Err: Send + 'static,
    M: MasterProblem<P::Primal> + Send + 'static,
    M::Err: Send + 'static,
{
    pub fn start(master: M, master_config: MasterConfig, tx: ClientSender<P, M>, threadpool: &mut ThreadPool) -> Self {
        // Create a pair of communication channels.
        let (to_master_tx, to_master_rx) = channel();

        // The the master process thread.







|

|

|
>




|










|

|







152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
    },
    /// An error occurred.
    Error(Error<MErr, PErr>),
}

/// Convenient type for `Response`.
pub type MasterResponse<P, M> =
    Response<<M as MasterProblem<<P as FirstOrderProblem>::Minorant>>::Err, <P as FirstOrderProblem>::Err>;

type ToMasterSender<P, M> = Sender<MasterTask<<P as FirstOrderProblem>::Minorant, <P as FirstOrderProblem>::Err, M>>;

type ToMasterReceiver<P, M> =
    Receiver<MasterTask<<P as FirstOrderProblem>::Minorant, <P as FirstOrderProblem>::Err, M>>;

pub struct MasterProcess<P, M>
where
    P: FirstOrderProblem,
    M: MasterProblem<P::Minorant>,
{
    /// The channel to transmit new tasks to the master problem.
    tx: ToMasterSender<P, M>,

    phantom: std::marker::PhantomData<M>,
}

impl<P, M> MasterProcess<P, M>
where
    P: FirstOrderProblem,
    P::Minorant: Send + 'static,
    P::Err: Send + 'static,
    M: MasterProblem<P::Minorant> + Send + 'static,
    M::Err: Send + 'static,
{
    pub fn start(master: M, master_config: MasterConfig, tx: ClientSender<P, M>, threadpool: &mut ThreadPool) -> Self {
        // Create a pair of communication channels.
        let (to_master_tx, to_master_rx) = channel();

        // The the master process thread.
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
            phantom: std::marker::PhantomData,
        }
    }

    /// Add new variables to the master problem.
    pub fn add_vars(
        &mut self,
        vars: Vec<(Option<usize>, Real, Real)>,
        sgext: Box<SubgradientExtender<P::Primal, P::Err>>,
    ) -> Result<(), Error<M::Err, P::Err>>
    where
        P::Err: 'static,
    {
        Ok(self.tx.send(MasterTask::AddVariables(vars, sgext))?)
    }

    /// Add a new minorant to the master problem model.
    ///
    /// This adds the specified `minorant` with associated `primal` data to the
    /// model of subproblem `i`.
    pub fn add_minorant(&mut self, i: usize, minorant: Minorant<P::Primal>) -> Result<(), Error<M::Err, P::Err>> {
        Ok(self.tx.send(MasterTask::AddMinorant(i, minorant))?)
    }

    /// Move the center of the master problem.
    ///
    /// This moves the master problem's center in direction $\\alpha \\cdot d$.
    pub fn move_center(&mut self, alpha: Real, d: Arc<DVector>) -> Result<(), Error<M::Err, P::Err>> {







|
|











|







205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
            phantom: std::marker::PhantomData,
        }
    }

    /// Add new variables to the master problem.
    pub fn add_vars(
        &mut self,
        vars: Vec<(Real, Real)>,
        sgext: Box<SubgradientExtender<P::Minorant, P::Err>>,
    ) -> Result<(), Error<M::Err, P::Err>>
    where
        P::Err: 'static,
    {
        Ok(self.tx.send(MasterTask::AddVariables(vars, sgext))?)
    }

    /// Add a new minorant to the master problem model.
    ///
    /// This adds the specified `minorant` with associated `primal` data to the
    /// model of subproblem `i`.
    pub fn add_minorant(&mut self, i: usize, minorant: P::Minorant) -> Result<(), Error<M::Err, P::Err>> {
        Ok(self.tx.send(MasterTask::AddMinorant(i, minorant))?)
    }

    /// Move the center of the master problem.
    ///
    /// This moves the master problem's center in direction $\\alpha \\cdot d$.
    pub fn move_center(&mut self, alpha: Real, d: Arc<DVector>) -> Result<(), Error<M::Err, P::Err>> {
248
249
250
251
252
253
254
255



256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274

    /// Sets the new weight of the proximal term in the master problem.
    pub fn set_weight(&mut self, weight: Real) -> Result<(), Error<M::Err, P::Err>> {
        Ok(self.tx.send(MasterTask::SetWeight { weight })?)
    }

    /// Get the current aggregated primal for a certain subproblem.
    pub fn get_aggregated_primal(&self, subproblem: usize) -> Result<P::Primal, Error<M::Err, P::Err>> {



        let (tx, rx) = channel();
        self.tx.send(MasterTask::GetAggregatedPrimal { subproblem, tx })?;
        rx.recv()?.map_err(Error::Aggregation)
    }

    /// The main loop of the master process.
    fn master_main(
        master: M,
        master_config: MasterConfig,
        tx: &mut ClientSender<P, M>,
        rx: ToMasterReceiver<P, M>,
        subgradient_extender: &mut Option<Box<SubgradientExtender<P::Primal, P::Err>>>,
    ) -> Result<(), Error<M::Err, P::Err>> {
        let mut master = master;

        // Initialize the master problem.
        master
            .set_num_subproblems(master_config.num_subproblems)
            .map_err(Error::Master)?;







|
>
>
>











|







249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278

    /// Sets the new weight of the proximal term in the master problem.
    pub fn set_weight(&mut self, weight: Real) -> Result<(), Error<M::Err, P::Err>> {
        Ok(self.tx.send(MasterTask::SetWeight { weight })?)
    }

    /// Get the current aggregated primal for a certain subproblem.
    pub fn get_aggregated_primal(
        &self,
        subproblem: usize,
    ) -> Result<<P::Minorant as Minorant>::Primal, Error<M::Err, P::Err>> {
        let (tx, rx) = channel();
        self.tx.send(MasterTask::GetAggregatedPrimal { subproblem, tx })?;
        rx.recv()?.map_err(Error::Aggregation)
    }

    /// The main loop of the master process.
    fn master_main(
        master: M,
        master_config: MasterConfig,
        tx: &mut ClientSender<P, M>,
        rx: ToMasterReceiver<P, M>,
        subgradient_extender: &mut Option<Box<SubgradientExtender<P::Minorant, P::Err>>>,
    ) -> Result<(), Error<M::Err, P::Err>> {
        let mut master = master;

        // Initialize the master problem.
        master
            .set_num_subproblems(master_config.num_subproblems)
            .map_err(Error::Master)?;
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
                MasterTask::AddMinorant(i, mut m) => {
                    debug!("master: add minorant to subproblem {}", i);
                    // It may happen the number new minorant belongs to an earlier evaluation
                    // with less variables (i.e. new variables have been added
                    // after the start of the evaluation but before the new
                    // minorant is added, i.e. now). In this case we must add
                    // extend the minorant accordingly.
                    if m.linear.len() < master.num_variables() {
                        if let Some(ref mut sgext) = subgradient_extender {
                            let newinds = (m.linear.len()..master.num_variables()).collect::<Vec<_>>();
                            let new_subg =
                                sgext(i, &m.primal, &newinds).map_err(Error::<M::Err, P::Err>::SubgradientExtension)?;
                            m.linear.extend(new_subg);
                        }
                    }
                    master.add_minorant(i, m).map_err(Error::Master)?;
                }
                MasterTask::MoveCenter(alpha, d) => {
                    debug!("master: move center");
                    master.move_center(alpha, &d);







|

<
<
|
<







298
299
300
301
302
303
304
305
306


307

308
309
310
311
312
313
314
                MasterTask::AddMinorant(i, mut m) => {
                    debug!("master: add minorant to subproblem {}", i);
                    // It may happen the number new minorant belongs to an earlier evaluation
                    // with less variables (i.e. new variables have been added
                    // after the start of the evaluation but before the new
                    // minorant is added, i.e. now). In this case we must add
                    // extend the minorant accordingly.
                    if m.dim() < master.num_variables() {
                        if let Some(ref mut sgext) = subgradient_extender {


                            sgext(i, &mut m).map_err(Error::<M::Err, P::Err>::SubgradientExtension)?;

                        }
                    }
                    master.add_minorant(i, m).map_err(Error::Master)?;
                }
                MasterTask::MoveCenter(alpha, d) => {
                    debug!("master: move center");
                    master.move_center(alpha, &d);
Changes to src/solver/sync.rs.
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
#[cfg(feature = "crossbeam")]
use rs_crossbeam::channel::{unbounded as channel, RecvError};
#[cfg(not(feature = "crossbeam"))]
use std::sync::mpsc::{channel, RecvError};

use log::{debug, info, warn};
use num_cpus;
use num_traits::Float;
use std::sync::Arc;
use std::time::Instant;
use threadpool::ThreadPool;

use crate::{DVector, Real};

use super::channels::{
    ChannelResultSender, ChannelUpdateSender, ClientReceiver, ClientSender, EvalResult, Message, Update,
};
use super::masterprocess::{self, MasterConfig, MasterProcess, MasterResponse, Response};
use crate::master::{Builder as MasterBuilder, MasterProblem};
use crate::problem::{FirstOrderProblem, UpdateState};







|




|







20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
#[cfg(feature = "crossbeam")]
use rs_crossbeam::channel::{unbounded as channel, RecvError};
#[cfg(not(feature = "crossbeam"))]
use std::sync::mpsc::{channel, RecvError};

use log::{debug, info, warn};
use num_cpus;
use num_traits::{Float, Zero};
use std::sync::Arc;
use std::time::Instant;
use threadpool::ThreadPool;

use crate::{DVector, Minorant, Real};

use super::channels::{
    ChannelResultSender, ChannelUpdateSender, ClientReceiver, ClientSender, EvalResult, Message, Update,
};
use super::masterprocess::{self, MasterConfig, MasterProcess, MasterResponse, Response};
use crate::master::{Builder as MasterBuilder, MasterProblem};
use crate::problem::{FirstOrderProblem, UpdateState};
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
///
/// - `T` is the value type,
/// - `P` is the `FirstOrderProblem` associated with the solver,
/// - `M` is the `MasterBuilder` associated with the solver.
pub type Result<T, P, M> = std::result::Result<
    T,
    Error<
        <<M as MasterBuilder<<P as FirstOrderProblem>::Primal>>::MasterProblem as MasterProblem<
            <P as FirstOrderProblem>::Primal,
        >>::Err,
        <P as FirstOrderProblem>::Err,
    >,
>;

impl<MErr, PErr> std::fmt::Display for Error<MErr, PErr>
where







|
|







84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
///
/// - `T` is the value type,
/// - `P` is the `FirstOrderProblem` associated with the solver,
/// - `M` is the `MasterBuilder` associated with the solver.
pub type Result<T, P, M> = std::result::Result<
    T,
    Error<
        <<M as MasterBuilder<<P as FirstOrderProblem>::Minorant>>::MasterProblem as MasterProblem<
            <P as FirstOrderProblem>::Minorant,
        >>::Err,
        <P as FirstOrderProblem>::Err,
    >,
>;

impl<MErr, PErr> std::fmt::Display for Error<MErr, PErr>
where
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
    }
}

/// Implementation of a parallel bundle method.
pub struct Solver<P, T = StandardTerminator, W = HKWeighter, M = crate::master::FullMasterBuilder>
where
    P: FirstOrderProblem,
    M: MasterBuilder<P::Primal>,
    P::Err: 'static,
{
    /// Parameters for the solver.
    pub params: Parameters,

    /// Termination predicate.
    pub terminator: T,







|







393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
    }
}

/// Implementation of a parallel bundle method.
pub struct Solver<P, T = StandardTerminator, W = HKWeighter, M = crate::master::FullMasterBuilder>
where
    P: FirstOrderProblem,
    M: MasterBuilder<P::Minorant>,
    P::Err: 'static,
{
    /// Parameters for the solver.
    pub params: Parameters,

    /// Termination predicate.
    pub terminator: T,
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461

impl<P, T, W, M> Solver<P, T, W, M>
where
    P: FirstOrderProblem,
    P::Err: Send + 'static,
    T: Terminator<SolverData> + Default,
    W: Weighter<SolverData> + Default,
    M: MasterBuilder<P::Primal>,
{
    /// Create a new parallel bundle solver.
    pub fn new(problem: P) -> Self
    where
        M: Default,
    {
        Self::with_master(problem, M::default())







|







447
448
449
450
451
452
453
454
455
456
457
458
459
460
461

impl<P, T, W, M> Solver<P, T, W, M>
where
    P: FirstOrderProblem,
    P::Err: Send + 'static,
    T: Terminator<SolverData> + Default,
    W: Weighter<SolverData> + Default,
    M: MasterBuilder<P::Minorant>,
{
    /// Create a new parallel bundle solver.
    pub fn new(problem: P) -> Self
    where
        M: Default,
    {
        Self::with_master(problem, M::default())
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
            }
        }
    }

    /// Handle a response from a subproblem evaluation.
    ///
    /// The function returns `Ok(true)` if the final iteration count has been reached.
    fn handle_client_response(&mut self, msg: EvalResult<usize, P::Primal, P::Err>) -> Result<bool, P, M> {
        let master = self.master_proc.as_mut().ok_or(Error::NotInitialized)?;
        match msg {
            EvalResult::ObjectiveValue { index, value } => {
                debug!("Receive objective from subproblem {}: {}", index, value);
                if self.data.nxt_ubs[index].is_infinite() {
                    self.data.cnt_remaining_ubs -= 1;
                }
                self.data.nxt_ubs[index] = self.data.nxt_ubs[index].min(value);
            }
            EvalResult::Minorant { index, mut minorant } => {
                debug!("Receive minorant from subproblem {}", index);
                if self.data.nxt_cutvals[index].is_infinite() {
                    self.data.cnt_remaining_mins -= 1;
                }
                // move center of minorant to cur_y
                minorant.move_center(-1.0, &self.data.nxt_d);
                self.data.nxt_cutvals[index] = self.data.nxt_cutvals[index].max(minorant.constant);
                // add minorant to master problem
                master.add_minorant(index, minorant)?;
            }
            EvalResult::Done { .. } => return Ok(false), // nothing to do here
            EvalResult::Error { err, .. } => return Err(Error::Evaluation(err)),
        }








|
















|







650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
            }
        }
    }

    /// Handle a response from a subproblem evaluation.
    ///
    /// The function returns `Ok(true)` if the final iteration count has been reached.
    fn handle_client_response(&mut self, msg: EvalResult<usize, P::Minorant, P::Err>) -> Result<bool, P, M> {
        let master = self.master_proc.as_mut().ok_or(Error::NotInitialized)?;
        match msg {
            EvalResult::ObjectiveValue { index, value } => {
                debug!("Receive objective from subproblem {}: {}", index, value);
                if self.data.nxt_ubs[index].is_infinite() {
                    self.data.cnt_remaining_ubs -= 1;
                }
                self.data.nxt_ubs[index] = self.data.nxt_ubs[index].min(value);
            }
            EvalResult::Minorant { index, mut minorant } => {
                debug!("Receive minorant from subproblem {}", index);
                if self.data.nxt_cutvals[index].is_infinite() {
                    self.data.cnt_remaining_mins -= 1;
                }
                // move center of minorant to cur_y
                minorant.move_center(-1.0, &self.data.nxt_d);
                self.data.nxt_cutvals[index] = self.data.nxt_cutvals[index].max(minorant.constant());
                // add minorant to master problem
                master.add_minorant(index, minorant)?;
            }
            EvalResult::Done { .. } => return Ok(false), // nothing to do here
            EvalResult::Error { err, .. } => return Err(Error::Evaluation(err)),
        }

847
848
849
850
851
852
853
854

855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874

875
876
877
878
879
880
881
882
883
884
885
886
887

        let mut have_update = false;
        for msg in update_rx {
            if let Message::Update(update) = msg {
                match update {
                    Update::AddVariables { bounds, sgext, .. } => {
                        have_update = true;
                        let mut newvars = Vec::with_capacity(bounds.len());

                        for (lower, upper) in bounds {
                            if lower > upper {
                                return Err(Error::InvalidBounds { lower, upper });
                            }
                            let value = if lower > 0.0 {
                                lower
                            } else if upper < 0.0 {
                                upper
                            } else {
                                0.0
                            };
                            //self.bounds.push((lower, upper));
                            newvars.push((None, lower - value, upper - value, value));
                        }
                        if !newvars.is_empty() {
                            // modify moved variables
                            for (index, val) in newvars.iter().filter_map(|v| v.0.map(|i| (i, v.3))) {
                                self.data.cur_y[index] = val;
                            }


                            // add new variables
                            self.data
                                .cur_y
                                .extend(newvars.iter().filter(|v| v.0.is_none()).map(|v| v.3));

                            master_proc.add_vars(newvars.iter().map(|v| (v.0, v.1, v.2)).collect(), sgext)?;
                        }
                    }
                    Update::Done { .. } => (), // there's nothing to do
                    Update::Error { err, .. } => return Err(Error::Update(err)),
                }
            } else {
                unreachable!("Only update results allowed during update");







|
>
|
|
<
<
|
|
<
<
|
<
<
|
<
|
<
<
<
<
|
|
>

|
<
<

|







847
848
849
850
851
852
853
854
855
856
857


858
859


860


861

862




863
864
865
866
867


868
869
870
871
872
873
874
875
876

        let mut have_update = false;
        for msg in update_rx {
            if let Message::Update(update) = msg {
                match update {
                    Update::AddVariables { bounds, sgext, .. } => {
                        have_update = true;
                        let newvars = bounds
                            .into_iter()
                            .map(|(lower, upper)| {
                                if lower <= upper {


                                    let value = lower.max(Real::zero()).min(upper);
                                    Ok((lower - value, upper - value, value))


                                } else {


                                    Err(Error::InvalidBounds { lower, upper })

                                }




                            })
                            .collect::<Result<Vec<_>, P, M>>()?;
                        if !newvars.is_empty() {
                            // add new variables
                            self.data.cur_y.extend(newvars.iter().map(|v| v.2));



                            master_proc.add_vars(newvars.iter().map(|v| (v.0, v.1)).collect(), sgext)?;
                        }
                    }
                    Update::Done { .. } => (), // there's nothing to do
                    Update::Error { err, .. } => return Err(Error::Update(err)),
                }
            } else {
                unreachable!("Only update results allowed during update");
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
            self.data.nxt_mod,
            self.data.nxt_val,
            self.data.cur_val
        );
    }

    /// Return the aggregated primal of the given subproblem.
    pub fn aggregated_primal(&self, subproblem: usize) -> Result<P::Primal, P, M> {
        Ok(self
            .master_proc
            .as_ref()
            .ok_or(Error::NotSolved)?
            .get_aggregated_primal(subproblem)?)
    }
}







|







911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
            self.data.nxt_mod,
            self.data.nxt_val,
            self.data.cur_val
        );
    }

    /// Return the aggregated primal of the given subproblem.
    pub fn aggregated_primal(&self, subproblem: usize) -> Result<<P::Minorant as Minorant>::Primal, P, M> {
        Ok(self
            .master_proc
            .as_ref()
            .ok_or(Error::NotSolved)?
            .get_aggregated_primal(subproblem)?)
    }
}