RsBundle  Check-in [8dc408e351]

            obj += F[j] * y;
        }
        write!(s, "  objval = {:.2}", obj)?;
        info!("{}", s);
    }

    {
        let mut slv = ParallelSolver::<_>::new(CFLProblem::new());
        slv.terminator.termination_precision = 1e-9;
        slv.solve()?;
    }

    Ok(())
}

pub mod firstorderproblem;
pub use crate::firstorderproblem::{Evaluation, FirstOrderProblem, SimpleEvaluation, Update};

pub mod solver;
pub use crate::solver::{
    BundleState, DefaultSolver, IterationInfo, Solver, SolverParams, StandardTerminator, Step, Terminator, UpdateState,

};

pub mod parallel;

mod weighter;
pub use crate::weighter::{HKWeighter, Weighter};

mod master;

pub mod mcf;

use std::time::Instant;
use threadpool::ThreadPool;

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

use super::problem::{EvalResult, FirstOrderProblem};
use crate::master::{BoxedMasterProblem, CplexMaster, MasterProblem, UnconstrainedMasterProblem};
use crate::solver::{BundleState, SolverParams, StandardTerminator, Step, Terminator};
use crate::weighter::{HKWeightable, HKWeighter, Weighter};

/// The default iteration limit.
pub const DEFAULT_ITERATION_LIMIT: usize = 10_000;

type MasterProblemError = <BoxedMasterProblem<CplexMaster> as MasterProblem>::Err;

/// Error raised by the parallel bundle [`Solver`].

    sgnorm: Real,
}

type MasterSender = Sender<std::result::Result<MasterResponse, MasterProblemError>>;

type MasterReceiver<Pr> = Receiver<MasterTask<Pr>>;

type ClientSender<P> =
    Sender<std::result::Result<EvalResult<usize, <P as FirstOrderProblem>::Primal>, <P as FirstOrderProblem>::Err>>;

type ClientReceiver<P> =
    Receiver<std::result::Result<EvalResult<usize, <P as FirstOrderProblem>::Primal>, <P as FirstOrderProblem>::Err>>;

/// Parameters for tuning the solver.
pub type Parameters = SolverParams;

pub struct SolverData {
    /// Current center of stability.
    cur_y: DVector,

    /// Function value in the current point.
    cur_val: Real,

    /// Function value at the current candidate.
    nxt_val: Real,

    /// Model value at the current candidate.
    nxt_mod: Real,

    /// The value of the new minorant in the current center.
    new_cutval: Real,

    /// Norm of current aggregated subgradient.
    sgnorm: Real,

    /// The currently used master problem weight.
    cur_weight: Real,
}

impl HKWeightable for SolverData {
    fn current_weight(&self) -> Real {
        self.cur_weight
    }

    fn center(&self) -> &DVector {
        &self.cur_y
    }

    fn center_value(&self) -> Real {
        self.cur_val
    }

    fn candidate_value(&self) -> Real {
        self.nxt_val
    }

    fn candidate_model(&self) -> Real {
        self.nxt_mod
    }

    fn new_cutvalue(&self) -> Real {
        self.new_cutval
    }

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

/// Implementation of a parallel bundle method.
pub struct Solver<P, T = StandardTerminator, W = HKWeighter>
where
    P: FirstOrderProblem,
    T: Terminator,

{
    /// Parameters for the solver.
    pub params: Parameters,

    /// Termination predicate.
    pub terminator: T,

    /// Weighter heuristic.
    pub weighter: W,

    /// The first order problem.
    problem: P,



    /// The algorithm data.


    data: SolverData,






    /// The threadpool of the solver.
    threadpool: ThreadPool,

    /// The channel to transmit new tasks to the master problem.
    master_tx: Option<Sender<MasterTask<P::Primal>>>,

    /// The channel to receive solutions from the master problem.
    master_rx: Option<Receiver<std::result::Result<MasterResponse, MasterProblemError>>>,

    /// The channel to receive the evaluation results from subproblems.
    client_tx: Option<ClientSender<P>>,

    /// The channel to receive the evaluation results from subproblems.
    client_rx: Option<ClientReceiver<P>>,

    /// Number of descent steps.
    cnt_descent: usize,

    /// Number of null steps.
    cnt_null: usize,

impl<P, T, W> Solver<P, T, W>
where
    P: FirstOrderProblem,
    P::Primal: Send + Sync + 'static,
    P::Err: std::error::Error + Send + Sync + 'static,
    T: Terminator + Default,
    W: Weighter<SolverData> + Default,
{
    /// Create a new parallel bundle solver.
    pub fn new(problem: P) -> Self {
        Solver {
            params: Parameters::default(),
            terminator: Default::default(),
            weighter: Default::default(),
            problem,
            data: SolverData {
                cur_y: dvec![],
                cur_val: 0.0,
                nxt_val: 0.0,
                nxt_mod: 0.0,
                new_cutval: 0.0,
                sgnorm: 0.0,
                cur_weight: 1.0,
            },

            threadpool: ThreadPool::with_name("Parallel bundle solver".to_string(), num_cpus::get()),
            master_tx: None,
            master_rx: None,
            client_tx: None,
            client_rx: None,

            cnt_descent: 0,
            cnt_null: 0,
            cnt_evals: 0,

            start_time: Instant::now(),

        }
    }

    /// Return the underlying threadpool.
    ///
    /// In order to use the same threadpool for concurrent processes,
    /// just clone the returned `ThreadPool`.
    pub fn threadpool(&self) -> &ThreadPool {

    /// This function is automatically called by [`solve`].
    pub fn init(&mut self) -> Result<(), Error<P::Err>> {
        debug!("Initialize solver");

        let n = self.problem.num_variables();
        let m = self.problem.num_subproblems();

        self.data.cur_y.init0(n);
        self.cnt_descent = 0;
        self.cnt_null = 0;
        self.cnt_evals = 0;

        let (tx, rx) = channel();
        self.client_tx = Some(tx);
        self.client_rx = Some(rx);

                }
            }
            debug!("Master proces stopped");
        });

        debug!("Initial problem evaluation");
        // We need an initial evaluation of all oracles for the first center.
        let y = Arc::new(self.data.cur_y.clone());
        for i in 0..m {
            self.problem
                .evaluate(i, y.clone(), i, self.client_tx.clone().unwrap())
                .map_err(Error::Evaluation)?;
        }

        let mut have_minorants = vec![false; m];

                Err(err) => return Err(Error::Evaluation(err)),
            };
            if cnt_remaining_mins == 0 && cnt_remaining_objs == 0 {
                break;
            }
        }

        self.data.cur_weight = Real::infinity(); // gets initialized when the master problem is complete
        master_tx
            .send(MasterTask::SetWeight { weight: 1.0 })
            .map_err(|err| Error::Process(err.into()))?;
        master_tx
            .send(MasterTask::Solve {
                center_value: self.data.cur_val,
            })
            .map_err(|err| Error::Process(err.into()))?;

        debug!("Initialization complete");

        self.start_time = Instant::now();

        let client_rx = self.client_rx.as_ref().ok_or(Error::NotInitialized)?;
        let master_tx = self.master_tx.as_ref().ok_or(Error::NotInitialized)?;
        let master_rx = self.master_rx.as_ref().ok_or(Error::NotInitialized)?;

        let mut cnt_iter = 0;
        let mut nxt_ubs = vec![Real::infinity(); self.problem.num_subproblems()];
        let mut cnt_remaining_ubs = self.problem.num_subproblems();
        let mut nxt_cutvals = vec![-Real::infinity(); self.problem.num_subproblems()];
        let mut cnt_remaining_mins = self.problem.num_subproblems();
        let mut nxt_d = Arc::new(dvec![]);
        let mut nxt_y = Arc::new(dvec![]);

        let mut expected_progress = 0.0;

        loop {
            select! {
                recv(client_rx) -> msg => {
                    let msg = msg
                        .map_err(|err| Error::Process(err.into()))?
                        .map_err(Error::Evaluation)?;
                    match msg {
                        EvalResult::ObjectiveValue { index, value } => {
                                debug!("Receive objective from subproblem {}: {}", index, value);
                            if nxt_ubs[index].is_infinite() {
                                cnt_remaining_ubs -= 1;
                            }
                            nxt_ubs[index] = nxt_ubs[index].min(value);
                        }
                        EvalResult::Minorant { index, mut minorant, primal } => {
                            debug!("Receive minorant from subproblem {}", index);
                            if nxt_cutvals[index].is_infinite() {
                                cnt_remaining_mins -= 1;
                            }
                            // move center of minorant to cur_y
                            minorant.move_center(-1.0, &nxt_d);
                            nxt_cutvals[index] = nxt_cutvals[index].max(minorant.constant);
                            // add minorant to master problem
                            master_tx
                                .send(MasterTask::AddMinorant(index, minorant, primal))
                                .map_err(|err| Error::Process(err.into()))?;
                        }
                    }

                    if cnt_remaining_ubs == 0 && cnt_remaining_mins == 0 {
                        // All subproblems have been evaluated, do a step.

                        let nxt_ub = nxt_ubs.iter().sum::<Real>();
                        let descent_bnd = self.get_descent_bound();

                        self.data.nxt_val = nxt_ub;
                        self.data.new_cutval = nxt_cutvals.iter().sum::<Real>();

                        debug!("Step");
                        debug!("  cur_val    ={}", self.data.cur_val);
                        debug!("  nxt_mod    ={}", self.data.nxt_mod);
                        debug!("  nxt_ub     ={}", nxt_ub);
                        debug!("  descent_bnd={}", descent_bnd);


                        let step;
                        if nxt_ub <= descent_bnd {
                            step = Step::Descent;
                            self.cnt_descent += 1;

                            self.data.cur_y = nxt_y.as_ref().clone();
                            self.data.cur_val = nxt_ub;

                            master_tx
                                .send(MasterTask::MoveCenter(1.0, nxt_d.clone()))
                                .map_err(|err| Error::Process(err.into()))?;

                            debug!("Descent Step");
                            debug!("  dir ={}", nxt_d);
                            debug!("  newy={}", self.data.cur_y);
                            self.data.cur_weight = self.weighter.descent_weight(&self.data);
                        } else {
                            step = Step::Null;
                            self.cnt_null += 1;
                            self.data.cur_weight = self.weighter.null_weight(&self.data);
                        }















                        master_tx
                            .send(MasterTask::SetWeight { weight: self.data.cur_weight })
                            .map_err(|err| Error::Process(err.into()))?;

                        self.show_info(step, expected_progress, self.data.nxt_mod, nxt_ub, self.data.cur_val);
                        cnt_iter += 1;
                        if cnt_iter >= niter { break }

                        master_tx
                            .send(MasterTask::Solve { center_value: self.data.cur_val })
                            .map_err(|err| Error::Process(err.into()))?;
                    }
                },
                recv(master_rx) -> msg => {
                    debug!("Receive master response");
                    // Receive result (new candidate) from the master
                    let master_res = msg
                        .map_err(|err| Error::Process(err.into()))?
                        .map_err(Error::Master)?;

                    if self.data.cur_weight < Real::infinity() && self.terminator.terminate(&BundleState {
                        cur_y: &self.data.cur_y,
                        cur_val: self.data.cur_val,
                        nxt_y: &nxt_y,
                        nxt_mod: self.data.nxt_mod,
                        nxt_val: self.data.cur_val,    // hopefully does not matter
                        new_cutval: self.data.cur_val, // hopefully does not matter
                        sgnorm: self.data.sgnorm,
                        weight: self.data.cur_weight,
                        step: Step::Term,
                        expected_progress,
                    }, &self.params) {
                        info!("Termination criterion satisfied");
                        return Ok(true)
                    }

                    // Compress bundle
                    master_tx.send(MasterTask::Compress).map_err(|err| Error::Process(err.into()))?;

                    // Compute new candidate.
                    self.data.nxt_mod = master_res.nxt_mod;
                    expected_progress = self.data.cur_val - self.data.nxt_mod;
                    self.data.sgnorm = master_res.sgnorm;
                    let mut next_y = dvec![];
                    nxt_d = Arc::new(master_res.nxt_d);
                    next_y.add(&self.data.cur_y, &nxt_d);
                    nxt_y = Arc::new(next_y);

                    // Reset evaluation data.
                    nxt_ubs.clear();
                    nxt_ubs.resize(self.problem.num_subproblems(), Real::infinity());
                    cnt_remaining_ubs = self.problem.num_subproblems();
                    nxt_cutvals.clear();
                    nxt_cutvals.resize(self.problem.num_subproblems(), -Real::infinity());
                    cnt_remaining_mins = self.problem.num_subproblems();

                    // Start evaluation of all subproblems at the new candidate.
                    for i in 0..self.problem.num_subproblems() {
                        self.problem.evaluate(i, nxt_y.clone(), i, client_tx.clone()) .map_err(Error::Evaluation)?;
                    }

                    // Compute the real initial weight.
                    if self.data.cur_weight.is_infinite() {












                        let weight = self.weighter.initial_weight(&self.data);
                        self.data.cur_weight = weight;
                        master_tx
                            .send(MasterTask::SetWeight { weight })
                            .map_err(|err| Error::Process(err.into()))?;
                    }

                },
            }

    /// If the oracle guarantees that $f(\bar{y}) \le$ this bound, the
    /// bundle method will perform a descent step.
    ///
    /// This value is $f(\hat{y}) + \varrho \cdot \Delta$ where
    /// $\Delta = f(\hat{y}) - \hat{f}(\bar{y})$ is the expected
    /// progress and $\varrho$ is the `acceptance_factor`.
    fn get_descent_bound(&self) -> Real {
        self.data.cur_val - self.params.acceptance_factor * (self.data.cur_val - self.data.nxt_mod)
    }

    fn show_info(&self, step: Step, expected_progress: Real, nxt_mod: Real, nxt_val: Real, cur_val: Real) {
        let time = self.start_time.elapsed();
        info!(
            "{} {:0>2}:{:0>2}:{:0>2}.{:0>2} {:4} {:4} {:4}{:1}  {:9.4} {:9.4} \
             {:12.6e}({:12.6e}) {:12.6e}",
            if step == Step::Term { "_endit" } else { "endit " },
            time.as_secs() / 3600,
            (time.as_secs() / 60) % 60,
            time.as_secs() % 60,
            time.subsec_nanos() / 10_000_000,
            self.cnt_descent,
            self.cnt_descent + self.cnt_null,
            0, /*self.master.cnt_updates(),*/
            if step == Step::Descent { "*" } else { " " },
            self.data.cur_weight,
            expected_progress,
            nxt_mod,
            nxt_val,
            cur_val
        );
    }
}

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

//! The main bundle method solver.

use crate::{Aggregatable, DVector, Real};
use crate::{Evaluation, FirstOrderProblem, Update};

use crate::master::CplexMaster;
use crate::master::{BoxedMasterProblem, MasterProblem, MinimalMaster};

use crate::weighter::{HKWeightable, HKWeighter, Weighter};

use log::{debug, info, warn};

use std::error::Error;
use std::f64::{INFINITY, NEG_INFINITY};
use std::fmt;
use std::mem::swap;

            sgnorm: $slf.sgnorm,
            weight: $slf.master.weight(),
            step: $step,
            expected_progress: $slf.expected_progress,
        }
    };
}

impl<'a> HKWeightable for BundleState<'a> {
    fn current_weight(&self) -> Real {
        self.weight
    }

    fn center(&self) -> &DVector {
        self.cur_y
    }

    fn center_value(&self) -> Real {
        self.cur_val
    }

    fn candidate_value(&self) -> Real {
        self.nxt_val
    }

    fn candidate_model(&self) -> Real {
        self.nxt_mod
    }

    fn new_cutvalue(&self) -> Real {
        self.new_cutval
    }

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

/**
 * Termination predicate.
 *
 * Given the current state of the bundle method, this function returns
 * whether the solution process should be stopped.
 */

    fn default() -> StandardTerminator {
        StandardTerminator {
            termination_precision: 1e-3,
        }
    }
}












/// An invalid value for some parameter has been passes.
#[derive(Debug)]
pub struct ParameterError(String);

impl fmt::Display for ParameterError {
    fn fmt(&self, fmt: &mut fmt::Formatter) -> Result<(), fmt::Error> {
        write!(fmt, "{}", self.0)

    /// This is the last primal generated by the oracle.
    pub fn last_primal(&self, fidx: usize) -> Option<&Pr> {
        self.minorants[fidx].last().and_then(|m| m.primal.as_ref())
    }
}

/// The default bundle solver with general master problem.
pub type DefaultSolver<P> = Solver<P, HKWeighter, BoxedMasterProblem<CplexMaster>>;

/// A bundle solver with a minimal cutting plane model.
pub type NoBundleSolver<P> = Solver<P, HKWeighter, BoxedMasterProblem<MinimalMaster>>;

/**
 * Implementation of a bundle method.
 */
pub struct Solver<P, W, M = BoxedMasterProblem<CplexMaster>>
where
    P: FirstOrderProblem,
{
    /// The first order problem description.
    problem: P,

    /// The solver parameter.
    pub params: SolverParams,

    /// Termination predicate.
    pub terminator: Box<Terminator>,

    /// Weighter heuristic.
    pub weighter: W,

    /// Lower and upper bounds of all variables.
    bounds: Vec<(Real, Real)>,

    /// Current center of stability.
    cur_y: DVector,

    /// The active minorant indices for each subproblem.
    minorants: Vec<Vec<MinorantInfo<P::Primal>>>,

    /// Accumulated information about the last iteration.
    iterinfos: Vec<IterationInfo>,
}

impl<P, W, M> Solver<P, W, M>
where
    P: FirstOrderProblem,
    P::Err: Into<Box<std::error::Error + Send + Sync + 'static>>,
    W: for<'a> Weighter<BundleState<'a>> + Default,
    M: MasterProblem<MinorantIndex = usize>,
    M::Err: Into<Box<std::error::Error + Send + Sync + 'static>>,
{
    /**
     * Create a new solver for the given problem.
     *
     * Note that the solver owns the problem, so you cannot use the
     * same problem description elsewhere as long as it is assigned to
     * the solver. However, it is possible to get a reference to the
     * internally stored problem using `Solver::problem()`.
     */
    #[allow(clippy::type_complexity)]
    pub fn new_params(problem: P, params: SolverParams) -> Result<Solver<P, W, M>, SolverError<P::Err, M::Err>> {
        Ok(Solver {
            problem,
            params,
            terminator: Box::new(StandardTerminator::default()),
            weighter: W::default(),
            bounds: vec![],
            cur_y: dvec![],
            cur_val: 0.0,
            cur_mod: 0.0,
            cur_vals: dvec![],
            cur_mods: dvec![],
            cur_valid: false,

            master: M::new()?,
            minorants: vec![],
            iterinfos: vec![],
        })
    }

    /// A new solver with default parameter.
    #[allow(clippy::type_complexity)]
    pub fn new(problem: P) -> Result<Solver<P, W, M>, SolverError<P::Err, M::Err>> {
        Solver::new_params(problem, SolverParams::default())
    }

    /**
     * Set the first order problem description associated with this
     * solver.
     *

        // *not* be the initial weight for the first iteration.
        self.master.set_weight(1.0)?;
        self.master.solve(self.cur_val)?;
        self.sgnorm = self.master.get_dualoptnorm2().sqrt();

        // Compute the real initial weight.
        let state = current_state!(self, Step::Term);
        let new_weight = self.weighter.initial_weight(&state);
        self.master.set_weight(new_weight)?;

        debug!("Init master completed");

        Ok(())
    }

            }
        }
        Ok(())
    }

    /// Perform a descent step.
    fn descent_step(&mut self) -> Result<(), SolverError<P::Err, M::Err>> {
        let new_weight = self.weighter.descent_weight(&current_state!(self, Step::Descent));
        self.master.set_weight(new_weight)?;
        self.cnt_descent += 1;
        swap(&mut self.cur_y, &mut self.nxt_y);
        swap(&mut self.cur_val, &mut self.nxt_val);
        swap(&mut self.cur_mod, &mut self.nxt_mod);
        swap(&mut self.cur_vals, &mut self.nxt_vals);
        swap(&mut self.cur_mods, &mut self.nxt_mods);
        self.master.move_center(1.0, &self.nxt_d);
        debug!("Descent Step");
        debug!("  dir ={}", self.nxt_d);
        debug!("  newy={}", self.cur_y);
        Ok(())
    }

    /// Perform a null step.
    fn null_step(&mut self) -> Result<(), SolverError<P::Err, M::Err>> {
        let new_weight = self.weighter.null_weight(&current_state!(self, Step::Null));
        self.master.set_weight(new_weight)?;
        self.cnt_null += 1;
        debug!("Null Step");
        Ok(())
    }

    /// Perform one bundle iteration.

/*
 * 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/>
 */

//! Control strategies for the weight of the proximal term.

use crate::Real;

/// Bundle weight controller.
///
/// The weighter computes the new weight for the proximal term depending
/// on the current data provided by [`Weightable`].
pub trait Weighter<Weightable> {
    /// Return the initial weight of the proximal term.
    fn initial_weight(&mut self, w: &Weightable) -> Real;

    /// Return the new weight of the proximal term after a descent step.
    ///
    /// The weight must be within the given `bounds`.
    fn descent_weight(&mut self, w: &Weightable) -> Real;

    /// Return the new weight of the proximal term after a null step.
    ///
    /// The weight must be within the given `bounds`.
    fn null_weight(&mut self, w: &Weightable) -> Real;
}

mod hkweighter;
pub use self::hkweighter::{HKWeightable, HKWeighter};

// Copyright (c) 2016, 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

//!
//! The procedure is described in
//!
//! > Helmberg, C. and Kiwiel, K.C. (2002): A spectral bundle method
//! > with bounds, Math. Programming A 93, 173--194
//!

use super::Weighter;
use crate::{DVector, Real};


use log::debug;

use std::cmp::{max, min};
use std::f64::NEG_INFINITY;

const FACTOR: Real = 2.0;


/// Weight updating rule according to Helmberg and Kiwiel.

///
/// The procedure is described in

///
/// > Helmberg, C. and Kiwiel, K.C. (2002): A spectral bundle method
/// > with bounds, Math. Programming A 93, 173--194

pub struct HKWeighter {
    min_weight: Real,
    max_weight: Real,
    eps_weight: Real,
    m_r: Real,
    iter: isize,
    model_max: Real,
}

impl HKWeighter {
    /// Create a new HKWeighter with default weight $m_R = 0.5$.
    pub fn new() -> HKWeighter {
        HKWeighter::new_weight(0.5)
    }

    /// Create new HKWeighter with weight $m_R$.
    pub fn new_weight(m_r: Real) -> HKWeighter {
        assert!(m_r > 0.0);
        HKWeighter {
            min_weight: 1e-2,
            max_weight: 1e3,
            eps_weight: 1e30,
            m_r,
            iter: 0,
            model_max: NEG_INFINITY,
        }
    }

    /// Set weight bounds.
    pub fn set_weight_bounds(&mut self, min_weight: Real, max_weight: Real) {
        assert!(min_weight > 0.0);
        assert!(max_weight >= min_weight);
        self.min_weight = min_weight;
        self.max_weight = max_weight;
    }
}

impl Default for HKWeighter {
    fn default() -> Self {
        Self::new()
    }
}

/// Information required by [`HKWeighter`].
pub trait HKWeightable {
    /// Return the current weight.
    fn current_weight(&self) -> Real;

    /// Return the current center of stability.
    fn center(&self) -> &DVector;

    /// Return the function value in the center.
    fn center_value(&self) -> Real;

    /// Return the function value in the candidate.
    fn candidate_value(&self) -> Real;

    /// Return the model value in the candidate.
    fn candidate_model(&self) -> Real;

    /// The value of the new minorant in the current center.
    fn new_cutvalue(&self) -> Real;


    /// Return the norm of the current aggregated subgradient.
    fn sgnorm(&self) -> Real;
}

impl<Weightable> Weighter<Weightable> for HKWeighter
where
    Weightable: HKWeightable,
{
    fn initial_weight(&mut self, w: &Weightable) -> Real {
        debug!("Compute initial weight");

        self.eps_weight = 1e30;
        self.iter = 0;

        let cur_y = w.center();
        let sgnorm = w.sgnorm();
        if cur_y.len() == 0 || sgnorm < cur_y.len() as Real * 1e-10 {
            1.0
        } else {
            sgnorm.max(1e-4)
        }
        .max(self.min_weight)
        .min(self.max_weight)
    }

    fn null_weight(&mut self, w: &Weightable) -> Real {
        debug!("Compute new weight after null-step (iter: {})", self.iter);
        let cur_nxt = w.center_value() - w.candidate_value();
        let cur_mod = w.center_value() - w.candidate_model();
        let old_weight = w.current_weight();
        let new_weight = 2.0 * old_weight * (1.0 - cur_nxt / cur_mod);

        debug!("  cur_nxt={} cur_mod={} w={}", cur_nxt, cur_mod, new_weight);


        let sgnorm = w.sgnorm();
        let lin_err = w.center_value() - w.new_cutvalue();
        self.eps_weight = self.eps_weight.min(sgnorm + cur_mod - sgnorm * sgnorm / old_weight);
        let new_weight = if self.iter < -3 && lin_err > self.eps_weight.max(FACTOR * cur_mod) {
            new_weight
        } else {
            old_weight
        }
        .min(FACTOR * old_weight)
        .min(self.max_weight);

        if new_weight > old_weight {
            self.iter = -1
        } else {
            self.iter = min(self.iter - 1, -1);
        }

        debug!(
            "  sgnorm={} cur_val={} new_cutval={} lin_err={} eps_weight={}",
            sgnorm,
            w.center_value(),
            w.new_cutvalue(),
            lin_err,
            self.eps_weight
        );
        debug!("  new_weight={}", new_weight);

        new_weight
    }

    fn descent_weight(&mut self, w: &Weightable) -> Real {
        debug!("Compute new weight after null-step (iter: {})", self.iter);
        let cur_nxt = w.center_value() - w.candidate_value();
        let cur_mod = w.center_value() - w.candidate_model();
        let old_weight = w.current_weight();
        let new_weight = 2.0 * old_weight * (1.0 - cur_nxt / cur_mod);

        debug!("  cur_nxt={} cur_mod={} w={}", cur_nxt, cur_mod, new_weight);

        self.model_max = self.model_max.max(w.candidate_model());
        let new_weight = if self.iter > 0 && cur_nxt > self.m_r * cur_mod {
            new_weight
        } else if self.iter > 3 || w.candidate_value() < self.model_max {
            old_weight / 2.0
        } else {
            old_weight
        }
        .max(old_weight / FACTOR)
        .max(self.min_weight);
        self.eps_weight = self.eps_weight.max(2.0 * cur_mod);

        if new_weight < old_weight {
            self.iter = 1;
            self.model_max = NEG_INFINITY;
        } else {
            self.iter = max(self.iter + 1, 1);
        }

        debug!("  model_max={}", self.model_max);
        debug!("  new_weight={}", new_weight);

        new_weight

    }
}