QPKKTSolverComparison.hxx

Go to the documentation of this file.

#ifndef CONICBUNDLE_QPKKTSOLVERCOMPARISON_HXX
#define CONICBUNDLE_QPKKTSOLVERCOMPARISON_HXX

#include <string>
#include "QPKKTSolverObject.hxx"
#include "clock.hxx"

namespace ConicBundle {




  class QPKKT_SolverStats
  {
  public:
    CH_Tools::Microseconds preptime;  
    CH_Tools::Microseconds predtime;  
    CH_Tools::Microseconds corrtime;  
    CH_Matrix_Classes::Integer prepnmult; 
    CH_Matrix_Classes::Integer prednmult; 
    CH_Matrix_Classes::Integer corrnmult; 
    CH_Matrix_Classes::Real cond;   
    CH_Matrix_Classes::Integer rank;  
    CH_Matrix_Classes::Real Hviol;   
    CH_Matrix_Classes::Real Aviol;   
    CH_Matrix_Classes::Real Bviol;   
    CH_Matrix_Classes::Real Cviol;   
    CH_Matrix_Classes::Real sysviol; 

    QPKKT_SolverStats()
    {
      preptime=predtime=corrtime=0;
      prepnmult=prednmult=corrnmult=0;
      cond=-1.;
      rank=-1;
      Hviol=Aviol=Bviol=Cviol=sysviol=0.;
    }

    ~QPKKT_SolverStats(){}

    friend std::ostream& operator<<(std::ostream& out,const QPKKT_SolverStats& s);

    friend std::istream& operator>>(std::ostream& in,QPKKT_SolverStats& s);
  };

  std::ostream& operator<<(std::ostream& out,const QPKKT_SolverStats& s);

  std::istream& operator>>(std::ostream& in,QPKKT_SolverStats& s);



  
  class QPKKT_KKTStats
  {
  public:
    CH_Matrix_Classes::Real prec;       
    CH_Matrix_Classes::Real mu;         
    std::vector<QPKKT_SolverStats> sdata; 

    QPKKT_KKTStats()
    { prec=mu=0.; }

    ~QPKKT_KKTStats(){}

    int append_col(CH_Matrix_Classes::Matrix& muvals, 
                   CH_Matrix_Classes::Matrix& prepsecs, 
                   CH_Matrix_Classes::Matrix& predsecs, 
                   CH_Matrix_Classes::Matrix& corrsecs, 
                   CH_Matrix_Classes::Indexmatrix& predcalls, 
                   CH_Matrix_Classes::Indexmatrix& corrcalls, 
                   CH_Matrix_Classes::Matrix& cond, 
                   CH_Matrix_Classes::Indexmatrix& pccols, 
                   CH_Matrix_Classes::Matrix& sysviol) 
    {
      CH_Matrix_Classes::Integer ind=muvals.coldim();
      muvals.concat_right(mu);
      assert(prepsecs.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(predsecs.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(corrsecs.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(predcalls.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(corrcalls.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(cond.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(pccols.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      assert(sysviol.rowdim()==CH_Matrix_Classes::Integer(sdata.size()));
      prepsecs.enlarge_right(1,0.);
      predsecs.enlarge_right(1,0.);
      corrsecs.enlarge_right(1,0.);
      predcalls.enlarge_right(1,CH_Matrix_Classes::Integer(0));
      corrcalls.enlarge_right(1,CH_Matrix_Classes::Integer(0));
      cond.enlarge_right(1,CH_Matrix_Classes::Integer(0));
      pccols.enlarge_right(1,CH_Matrix_Classes::Integer(0));
      sysviol.enlarge_right(1,CH_Matrix_Classes::Integer(0));
      for(CH_Matrix_Classes::Integer i=0;i<CH_Matrix_Classes::Integer(sdata.size());i++){
        prepsecs(i,ind)+=sdata[unsigned(i)].preptime;
        predsecs(i,ind)+=sdata[unsigned(i)].predtime;
        corrsecs(i,ind)=sdata[unsigned(i)].corrtime;
        predcalls(i,ind)=sdata[unsigned(i)].prednmult;
        corrcalls(i,ind)=sdata[unsigned(i)].corrnmult;
        cond(i,ind)=sdata[unsigned(i)].cond;
        pccols(i,ind)=sdata[unsigned(i)].rank;
        sysviol(i,ind)=sdata[unsigned(i)].sysviol;
      }
      return 0;
    }

    int add_col(CH_Matrix_Classes::Matrix& prepsecs, 
                CH_Matrix_Classes::Matrix& predsecs, 
                CH_Matrix_Classes::Matrix& corrsecs, 
                CH_Matrix_Classes::Indexmatrix& predcalls, 
                CH_Matrix_Classes::Indexmatrix& corrcalls 
                )
    {
      CH_Matrix_Classes::Integer ind=prepsecs.coldim()-1;
      for(CH_Matrix_Classes::Integer i=0;i<CH_Matrix_Classes::Integer(sdata.size());i++){
        prepsecs(i,ind)+=sdata[unsigned(i)].preptime;
        predsecs(i,ind)+=sdata[unsigned(i)].predtime;
        corrsecs(i,ind)+=sdata[unsigned(i)].corrtime;
        predcalls(i,ind)+=sdata[unsigned(i)].prednmult;
        corrcalls(i,ind)+=sdata[unsigned(i)].corrnmult;
      }
      return 0;
    }
    
    friend std::ostream& operator<<(std::ostream& out,const QPKKT_KKTStats& k);
    
    friend std::istream& operator>>(std::ostream& in,QPKKT_KKTStats& k);
    
  };
    
  std::ostream& operator<<(std::ostream& out,const QPKKT_KKTStats& k);
    
  std::istream& operator>>(std::istream& in,QPKKT_KKTStats& k);



  
  class QPKKT_ProbStats
  {
  public:
    CH_Matrix_Classes::Integer Qdim;       
    CH_Matrix_Classes::Integer Vdim;       
    CH_Matrix_Classes::Integer Arowdim;    
    CH_Matrix_Classes::Integer Aeqdim;     
    CH_Matrix_Classes::Integer Bdim;       
    CH_Matrix_Classes::Integer Cdim;       
    std::vector<CH_Tools::Microseconds> inittime; 
    std::vector<QPKKT_KKTStats> kktdata; 

    QPKKT_ProbStats()
    {Qdim=Vdim=Arowdim=Aeqdim=Bdim=Cdim=0; }

    ~QPKKT_ProbStats(){}

    CH_Matrix_Classes::Integer cnt_mu_cols(CH_Matrix_Classes::Real lbmu,
                                           CH_Matrix_Classes::Real ubmu)
    {
      CH_Matrix_Classes::Integer sumval=0;
      for (unsigned int i=0;i<kktdata.size();i++)
        if ((lbmu<=kktdata[i].mu)&&(kktdata[i].mu<ubmu))
          sumval++;
      return sumval;
    }
    
    int get_mu_stats(CH_Matrix_Classes::Real lbmu,
                     CH_Matrix_Classes::Real ubmu,
                     CH_Matrix_Classes::Indexmatrix& dims,
                     CH_Matrix_Classes::Matrix& mu,
                     CH_Matrix_Classes::Matrix& prepsecs,
                     CH_Matrix_Classes::Matrix& predsecs,
                     CH_Matrix_Classes::Matrix& corrsecs,
                     CH_Matrix_Classes::Indexmatrix& predcalls,
                     CH_Matrix_Classes::Indexmatrix& corrcalls,
                     CH_Matrix_Classes::Matrix& cond,
                     CH_Matrix_Classes::Indexmatrix& pccols,
                     CH_Matrix_Classes::Matrix& sysviol)
    {
      CH_Matrix_Classes::Integer cnt=0;
      for (unsigned int i=0;i<kktdata.size();i++){
        if ((lbmu<=kktdata[i].mu)&&(kktdata[i].mu<ubmu)){
          kktdata[i].append_col(mu,prepsecs,predsecs,corrsecs,predcalls,corrcalls,cond,pccols,sysviol);
          cnt++;
        }
      }
      assert(dims.coldim()+cnt==mu.coldim());
      CH_Matrix_Classes::Integer startcol=dims.coldim();
      if (dims.coldim()==0)
        dims.init(6,cnt,CH_Matrix_Classes::Integer(0));
      else {
        assert(dims.rowdim()==6);
        dims.enlarge_right(cnt,CH_Matrix_Classes::Integer(0));
      }
      for (CH_Matrix_Classes::Integer i=0;i<cnt;i++){
        dims(0,startcol+i)=Qdim;
        dims(1,startcol+i)=Vdim;
        dims(2,startcol+i)=Arowdim;
        dims(3,startcol+i)=Aeqdim;
        dims(4,startcol+i)=Bdim;
        dims(5,startcol+i)=Cdim;
      }
      return 0;
    }

    int get_prob_stats(CH_Matrix_Classes::Indexmatrix& dims, 
                       CH_Matrix_Classes::Indexmatrix& iterations, 
                       CH_Matrix_Classes::Matrix& lastmu, 
                       CH_Matrix_Classes::Matrix& prepsecs, 
                       CH_Matrix_Classes::Matrix& predsecs, 
                       CH_Matrix_Classes::Matrix& corrsecs, 
                       CH_Matrix_Classes::Indexmatrix& predcalls, 
                       CH_Matrix_Classes::Indexmatrix& corrcalls 
                       )
    {
      assert(dims.coldim()==iterations.coldim());
      assert(dims.coldim()==lastmu.coldim());
      assert(dims.coldim()==prepsecs.coldim());
      assert(dims.coldim()==predsecs.coldim());
      assert(dims.coldim()==corrsecs.coldim());
      assert(dims.coldim()==predcalls.coldim());
      assert(dims.coldim()==corrcalls.coldim());
      CH_Matrix_Classes::Integer ind=dims.coldim();
      if (dims.coldim()==0)
        dims.init(6,1,CH_Matrix_Classes::Integer(0));
      else {
        assert(dims.rowdim()==6);
        dims.enlarge_right(1,CH_Matrix_Classes::Integer(0));
      }
      dims(0,ind)=Qdim;
      dims(1,ind)=Vdim;
      dims(2,ind)=Arowdim;
      dims(3,ind)=Aeqdim;
      dims(4,ind)=Bdim;
      dims(5,ind)=Cdim;
      iterations.concat_right(CH_Matrix_Classes::Integer(kktdata.size()));
      lastmu.concat_right(kktdata[kktdata.size()-1].mu);
      assert(inittime.size()==unsigned(prepsecs.rowdim()));
      prepsecs.enlarge_right(1,0.);
      for(unsigned int i=0;i<inittime.size();i++)
        prepsecs(CH_Matrix_Classes::Integer(i),ind)+=inittime[i];
      predsecs.enlarge_right(1,0.);
      corrsecs.enlarge_right(1,0.);
      predcalls.enlarge_right(1,0.);
      corrcalls.enlarge_right(1,0.);
      for (unsigned int i=0;i<kktdata.size();i++){
        kktdata[i].add_col(prepsecs,predsecs,corrsecs,predcalls,corrcalls);
      }
      return 0;
    }
 

    friend std::ostream& operator<<(std::ostream& out,const QPKKT_ProbStats& p);

    friend std::istream& operator>>(std::istream& in,QPKKT_ProbStats& p);
};

  std::ostream& operator<<(std::ostream& out,const QPKKT_ProbStats& p);

  std::istream& operator>>(std::istream& in,QPKKT_ProbStats& p);





class QPKKTSolverComparison: public QPKKTSolverObject
{
private:
  std::vector<QPKKTSolverObject*> solver;  
  std::vector<std::string> solvername;    
  std::vector<QPKKT_ProbStats> probdata;  
  std::vector<QPModelBlockObject*> model; 
  
  QPSolverProxObject* testHp; 
  const CH_Matrix_Classes::Sparsemat* testA; 
  CH_Matrix_Classes::Matrix testKKTdiagx; 
  CH_Matrix_Classes::Matrix testKKTdiagy; 
  CH_Matrix_Classes::Real testHfactor; 
  CH_Matrix_Classes::Matrix testprimalrhs; 
  CH_Matrix_Classes::Matrix testdualrhs; 
  
  CH_Tools::Clock clock; 
 
  CH_Matrix_Classes::Matrix insolx;   
  CH_Matrix_Classes::Matrix insoly;  
  CH_Matrix_Classes::Matrix outsolx; 
  CH_Matrix_Classes::Matrix outsoly; 

  int violation(CH_Matrix_Classes::Matrix& violvec,const CH_Matrix_Classes::Matrix& solx,const CH_Matrix_Classes::Matrix& soly,QPModelBlockObject* inmodel,QPKKT_SolverStats& stats);
  
public:
  virtual void clear();
  

  QPKKTSolverComparison(CBout* cb=0,int cbinc=-1);

  virtual ~QPKKTSolverComparison();
  
  virtual int add_solver(QPKKTSolverObject* solver,const char* name);
  
  virtual int QPinit_KKTdata(QPSolverProxObject* Hp, 
                             QPModelBlockObject* model, 
                             const CH_Matrix_Classes::Sparsemat* A, 
                             const CH_Matrix_Classes::Indexmatrix* eq_indices 
                             );

  virtual int QPinit_KKTsystem(const CH_Matrix_Classes::Matrix& KKTdiagx,
                               const CH_Matrix_Classes::Matrix& KKTdiagy,
                               CH_Matrix_Classes::Real Hfactor,
                               CH_Matrix_Classes::Real prec,
                               QPSolverParameters* params);

  virtual int QPsolve_KKTsystem(CH_Matrix_Classes::Matrix& solx,
                                 CH_Matrix_Classes::Matrix& soly,
                                 const CH_Matrix_Classes::Matrix& primalrhs,
                                 const CH_Matrix_Classes::Matrix& dualrhs,
                                 CH_Matrix_Classes::Real rhsmu,
                                 CH_Matrix_Classes::Real rhscorr,
                                 CH_Matrix_Classes::Real prec,
                                 QPSolverParameters* params);
  
  virtual CH_Matrix_Classes::Real QPget_blockH_norm();

  virtual CH_Matrix_Classes::Real QPget_blockA_norm();

  //------- for later statistical evaluation purposes

  std::vector<std::string>& get_solvernames()
  {return solvername;} 

  std::vector<QPKKT_ProbStats>& get_probdata()
  {return probdata;}

  int get_mu_stats(CH_Matrix_Classes::Real lbmu, 
                   CH_Matrix_Classes::Real ubmu, 
                   CH_Matrix_Classes::Indexmatrix& dims, 
                   CH_Matrix_Classes::Matrix& mu, 
                   CH_Matrix_Classes::Matrix& prepsecs, 
                   CH_Matrix_Classes::Matrix& predsecs, 
                   CH_Matrix_Classes::Matrix& corrsecs, 
                   CH_Matrix_Classes::Indexmatrix& predcalls, 
                   CH_Matrix_Classes::Indexmatrix& corrcalls, 
                   CH_Matrix_Classes::Matrix& cond, 
                   CH_Matrix_Classes::Indexmatrix& pccols, 
                     CH_Matrix_Classes::Matrix& sysviol 
                   );

  int get_prob_stats(CH_Matrix_Classes::Indexmatrix& dims,
                     CH_Matrix_Classes::Indexmatrix& iterations, 
                     CH_Matrix_Classes::Matrix& lastmu, 
                     CH_Matrix_Classes::Matrix& prepsecs, 
                     CH_Matrix_Classes::Matrix& predsecs, 
                     CH_Matrix_Classes::Matrix& corrsecs, 
                     CH_Matrix_Classes::Indexmatrix& predcalls, 
                     CH_Matrix_Classes::Indexmatrix& corrcalls 
                     );
  
  friend std::ostream& operator<<(std::ostream& out,const QPKKTSolverComparison& q);

  friend std::istream& operator>>(std::istream& in,QPKKTSolverComparison& q);
};


  std::ostream& operator<<(std::ostream& out,const QPKKTSolverComparison& q);

  std::istream& operator>>(std::istream& in,QPKKTSolverComparison& q);




}

#endif