QPSolver.hxx

Go to the documentation of this file.

#ifndef CONICBUNDLE_QPSOLVER_HXX
#define CONICBUNDLE_QPSOLVER_HXX

#include "QPSolverBasicStructures.hxx"
#include "Groundset.hxx"
#include "LPGroundsetModification.hxx"

namespace ConicBundle {




class GroundsetModification;

  class QPSolver:public virtual QPSolverBasicStructures
{
private:
  class QPProblemData {
  public:
    CH_Matrix_Classes::Integer dim;  
    CH_Matrix_Classes::Matrix lby;    
    CH_Matrix_Classes::Matrix uby;    
    CH_Matrix_Classes::Indexmatrix lbindex;  
    CH_Matrix_Classes::Indexmatrix ubindex;  
    CH_Matrix_Classes::Sparsemat A;  
    CH_Matrix_Classes::Matrix rhslb; 
    CH_Matrix_Classes::Matrix rhsub; 
    CH_Matrix_Classes::Indexmatrix rhslbindex;  
    CH_Matrix_Classes::Indexmatrix rhsubindex;  
    CH_Matrix_Classes::Indexmatrix rhseqindex;
    
    
  //data for describing the quadratic cost matrix Q
    BundleProxObject* Hp;  
    CH_Matrix_Classes::Matrix c; 
    CH_Matrix_Classes::Real gamma; 
    
    CH_Matrix_Classes::Matrix groundset_c; 
    CH_Matrix_Classes::Real groundset_gamma; 

  };
  
  QPProblemData original_data;

  QPProblemData preproc_data;
  
  CH_Matrix_Classes::Indexmatrix preproc_fixed; 
  CH_Matrix_Classes::Matrix preproc_fixedval; 
  CH_Matrix_Classes::Indexmatrix preproc_indices; 
  AffineFunctionTransformation preproc_aft; 
  std::map<MinorantPointer,MinorantPointer> preproc_bundle_projection; 
  
  QPProblemData* qp_data; 


  CH_Matrix_Classes::Real sol_val_lb;
  CH_Matrix_Classes::Real sol_val_ub;
  CH_Matrix_Classes::Matrix sol_point;
  CH_Matrix_Classes::Real gs_aggr_offset;
  CH_Matrix_Classes::Matrix gs_aggr_gradient;
  
  
  
  CH_Matrix_Classes::Real primal_infeasibility_eps; 
  CH_Matrix_Classes::Real dual_infeasibility_eps; 
  CH_Matrix_Classes::Real lower_bound_gap_eps; 
  CH_Matrix_Classes::Real upper_bound_gap_eps; 
  CH_Matrix_Classes::Real objective_gap_eps; 


  
  int determine_indices(QPProblemData& qpd);

  int preprocess_data(const CH_Matrix_Classes::Matrix& center_y,
                      CH_Matrix_Classes::Indexmatrix* yfixed,
                      bool& no_changes);

  int postprocess_data(bool round_to_active_bounds);

  
protected:

  CH_Matrix_Classes::Integer QPget_xdim() const {return qp_data->dim;}
  CH_Matrix_Classes::Integer QPget_ydim() const {return qp_data->A.rowdim();}
 
  CH_Matrix_Classes::Matrix&  QPadd_Qx(const CH_Matrix_Classes::Matrix& xin,CH_Matrix_Classes::Matrix& outplusQx) const 
  {return qp_data->Hp->add_Hx(xin,outplusQx);}
 
  const CH_Matrix_Classes::Matrix& QPget_c() const {return qp_data->c;}

  CH_Matrix_Classes::Real QPget_gamma() const {return qp_data->gamma;}

  const CH_Matrix_Classes::Matrix& QPget_rhslb() const {return qp_data->rhslb;}

  const CH_Matrix_Classes::Matrix& QPget_rhsub() const {return qp_data->rhsub;}

  const CH_Matrix_Classes::Indexmatrix& QPget_rhslbind() const {return qp_data->rhslbindex;}

  const CH_Matrix_Classes::Indexmatrix& QPget_rhsubind() const {return qp_data->rhsubindex;}

  const CH_Matrix_Classes::Matrix& QPget_lb() const {return qp_data->lby;}

  const CH_Matrix_Classes::Matrix& QPget_ub() const {return qp_data->uby;}

  const CH_Matrix_Classes::Indexmatrix& QPget_lbind() const {return qp_data->lbindex;}

  const CH_Matrix_Classes::Indexmatrix& QPget_ubind() const {return qp_data->ubindex;}

  CH_Matrix_Classes::Matrix& QPadd_Ax(const CH_Matrix_Classes::Matrix& xin,CH_Matrix_Classes::Matrix& outplusAx) const 
  {return CH_Matrix_Classes::genmult(qp_data->A,xin,outplusAx,1.,1.);}
  
  CH_Matrix_Classes::Matrix& QPadd_Aty(const CH_Matrix_Classes::Matrix& yin,CH_Matrix_Classes::Matrix& outplusAty) const 
  {return CH_Matrix_Classes::genmult(qp_data->A,yin,outplusAty,1.,1.,1);}
  

public:
  
  void QPclear();

  QPSolver(CBout* cb=0,int cbinc=-1):
    CBout(cb,cbinc),QPSolverBasicStructures(0,cb)
  { original_data.Hp=0;preproc_data.Hp=0; QPclear();}
           
  ~QPSolver(){}

  //------------  QPSolverOject routines

  int QPset_parameters(QPSolverParametersObject* params)
  {
    QPSolverParameters* qpp=dynamic_cast<QPSolverParameters*>(params);
    if (qpp==0)
      return 1;
    return QPset_solver_parameters(qpp);
  }

  bool QPsupports_yfixing()
  { return false;}

  virtual bool QPsupports_updates()
  { return false;}
  
  CH_Matrix_Classes::Real QPget_lower_bound()
  { return QPget_dualval(); }

  int QPsolve(const CH_Matrix_Classes::Matrix& center_y,
              CH_Matrix_Classes::Real lower_bound,
              CH_Matrix_Classes::Real upper_bound,
              CH_Matrix_Classes::Real relprec,
              QPSolverProxObject* Hp,
              const MinorantPointer& gs_aggr,
              CH_Matrix_Classes::Indexmatrix* yfixed);

  int QPupdate(const CH_Matrix_Classes::Matrix& center_y,
               CH_Matrix_Classes::Real lower_bound,
               CH_Matrix_Classes::Real upper_bound,
               CH_Matrix_Classes::Real relprec,
               QPSolverProxObject* Hp,
               const MinorantPointer& gs_aggr,
               CH_Matrix_Classes::Indexmatrix* yfixed,
               const MinorantPointer& delta_gs_aggr,
               const CH_Matrix_Classes::Indexmatrix& delta_index);

  int QPresolve(CH_Matrix_Classes::Real lower_bound,
                CH_Matrix_Classes::Real upper_bound,
                CH_Matrix_Classes::Real relprec);

  int QPget_solution(CH_Matrix_Classes::Real& augval_lb,
                     CH_Matrix_Classes::Real& augval_ub,
                     CH_Matrix_Classes::Matrix& new_point,
                     CH_Matrix_Classes::Real& gsaggr_offset,
                     CH_Matrix_Classes::Matrix& gsaggr_gradient);

  std::ostream& QPprint_statistics(std::ostream& out,int /* printlevel*/ =0)
  {return out;}
  
  GroundsetModification* QPstart_modification()
  {return new LPGroundsetModification(original_data.dim,original_data.A.rowdim(),this);}

  int QPapply_modification(const GroundsetModification& mdf);

  bool QPis_feasible(const CH_Matrix_Classes::Matrix& y,
                   CH_Matrix_Classes::Real relprec=1e-10);

  int QPensure_feasibility(CH_Matrix_Classes::Matrix& y,
                           bool& ychanged,
                           QPSolverProxObject* inHp,
                           CH_Matrix_Classes::Real relprec=1e-10);

  int solve(BundleProxObject* Hp,
            const CH_Matrix_Classes::Matrix& c,
            CH_Matrix_Classes::Real gamma,
            CH_Matrix_Classes::Real lowerbound,
            CH_Matrix_Classes::Real upperbound,
            CH_Matrix_Classes::Real relprec,
            CH_Matrix_Classes::Real skip_factor);
           
  bool QPprefer_UQPSolver(QPSolverProxObject*) const;

  bool QPconstrained() const
  {return (original_data.lbindex.rowdim()+original_data.ubindex.rowdim()+original_data.A.rowdim()==0)?false:true;}

  CH_Matrix_Classes::Integer rowdim() const {return original_data.A.rowdim();}

  const CH_Matrix_Classes::Matrix& get_lby() const { return original_data.lby;}

  const CH_Matrix_Classes::Matrix& get_uby() const { return original_data.uby;}

  const CH_Matrix_Classes::Indexmatrix& get_lbindex() const { return original_data.lbindex;}

  const CH_Matrix_Classes::Indexmatrix& get_ubindex() const { return original_data.ubindex;}

  const CH_Matrix_Classes::Sparsemat& get_A() const { return original_data.A;}

  const CH_Matrix_Classes::Matrix& get_rhslb() const {return original_data.rhslb;}

  const CH_Matrix_Classes::Matrix& get_rhsub() const {return original_data.rhsub;}

  const CH_Matrix_Classes::Indexmatrix& get_rhslbind() const {return original_data.rhslbindex;}

  const CH_Matrix_Classes::Indexmatrix& get_rhsubind() const {return original_data.rhsubindex;}

  bool QPboxconstrained(const CH_Matrix_Classes::Matrix*& lb, 
                        const CH_Matrix_Classes::Matrix*& ub, 
                        const CH_Matrix_Classes::Indexmatrix*& lbind, 
                        const CH_Matrix_Classes::Indexmatrix*& ubind 
                        ) const
  {
    lb=&get_lby();
    ub=&get_uby();
    lbind=&get_lbindex();
    ubind=&get_ubindex();
    if ((!QPconstrained())||(original_data.A.rowdim()!=0))
      return false;
    return true;
  }

  int mfile_data(std::ostream& out) const;

  void set_cbout(const CBout* cb, int incr=-1)
  {
    CBout::set_cbout(cb,incr);
    if (model_block)
      model_block->set_cbout(cb,incr);
  }

};





}

#endif