opt.h

Go to the documentation of this file.
00001 // Copyright (C) 2006 Ivo Nowak and Stefan Vigerske
00002 // All Rights Reserved.
00003 // This code is published under the Common Public License.
00004 //
00005 // Author: Ivo Nowak, Stefan Vigerske
00006 
00007 #ifndef OPT_H
00008 #define OPT_H
00009 
00010 #include "standard.h"
00011 #include "usermatrix.h"
00012 #include "func.h"
00013 #include "param.h"
00014 #include "dualqqp.h"
00015 
00020 class Solver {
00027       friend ostream& operator << (ostream& out, const Solver& a)
00028       {  a.print(out); return out; };
00029 
00030    protected:
00034       int dim_;
00035 
00039       double time_;
00040 
00044       double opt_val_;
00045 
00049       int iter_;
00050    public:
00054       Pointer<ostream> out_solver_p;
00058       Pointer<ostream> out_solver_log_p;
00064 #define out_solver if (out_solver_p) (*out_solver_p)
00065 
00070 #define out_solver_log if (out_solver_log_p) (*out_solver_log_p)
00071 
00075       double tol;
00076 
00080       int iter_max;
00081 
00084       dvector sol_point;
00085 
00092       Solver(int n, Pointer<ostream> out_solver_p_=out_out_p, Pointer<ostream> out_solver_log_p_=out_log_p)
00093       : dim_(n), sol_point(n), iter_(0), time_(-1), tol(-1.), opt_val_(INFINITY), iter_max(INF),
00094         out_solver_p(out_solver_p_), out_solver_log_p(out_solver_log_p_)
00095       {  };
00096 
00099       virtual ~Solver() { };
00100 
00104       int dim() const { return dim_;  };
00105 
00109       double opt_val() { return opt_val_;  };
00110 
00113       double time() { return time_;  };
00114 
00117       int iter() { return iter_;  };
00118 
00124       virtual int solve()=0;
00125 
00132       virtual int solve(dvector &x) { sol_point=x; return solve(); };
00133 
00138       virtual void print(ostream &out) const {
00139                                 out << "dim: " << dim() << " tol: " << tol << endl;
00140                                 out << "solution time: " << time_ << endl;
00141                                 out << "iterations: " << iter_ << endl;
00142                         };
00143 };
00144 
00145 
00148 class DualSolver : public Solver {
00149   private:
00152     bool threshold_cntrl;
00153 
00156     bool conv_rate_cntrl;
00159     double stopping_rho;
00162     int minor_iter;
00165     double last_major_val;
00168     double last_val;
00171     double first_major_val;
00174     double max_rel_improvement;
00177     int improve_iter;
00178 
00179   protected:
00182     Param& param;
00183 
00186     DualFunc& obj;
00187 
00195     void do_log();
00196 
00205     int check(double val);
00206 
00207   public:
00210     double last_rel_improvement;
00211 
00214     double threshold;
00215 
00218     dvector lower_bound;
00219 
00222     Pointer<vector<double> > dual_vals;
00225     Pointer<vector<dvector> > dual_points;
00228     Pointer<vector<vector<dvector> > > orig_points;
00233     int log_frequency;
00237     bool store_primals;
00238 
00247     DualSolver(DualFunc& f, vector<double> lower_bound_, Param& param_, Pointer<ostream> out_solver_p_=out_out_p, Pointer<ostream> out_solver_log_p_=out_log_p)
00248     : Solver(f.dim(), out_solver_p_, out_solver_log_p_), threshold(INFINITY), lower_bound(f.dim()), obj(f),
00249       log_frequency(0), store_primals(true), param(param_), last_major_val(-INFINITY), last_val(-INFINITY), improve_iter(-1)
00250     { for (int i=0; i<dim(); i++) lower_bound[i]=lower_bound_[i];
00251       param.read();
00252       threshold_cntrl=param.get_i("control threshold", 1);
00253       conv_rate_cntrl=param.get_i("control convergence rate", 1);
00254                         stopping_rho=param.get_d("stopping rho", 0.1);
00255                         minor_iter=param_.get_i("minor iterations", 5);
00256                 }
00257 
00258                 virtual ~DualSolver() { }
00259 
00260                 double rel_improvement1() {
00261                   return (last_major_val-first_major_val)/fabs(first_major_val);
00262                 }
00263 
00264                 double rel_improvement2() {
00265                   return last_rel_improvement;
00266                 }
00267 
00268                 double rel_improvement_max() {
00269                   return max_rel_improvement;
00270                 }
00271 
00272                 int serious_steps() {
00273                   return improve_iter;
00274                 }
00275 
00276 };
00277 
00281 class LocOpt : public Solver {
00282         protected:
00283                 Timer timer;
00284         public:
00285                 static bool nlp_solver_available();
00288                 static Pointer<LocOpt> get_solver(const Pointer<MinlpProblem> prob, Pointer<Param> param, const char* param_prefix=NULL, Pointer<ostream> out_solver_p_=out_out_p, Pointer<ostream> out_solver_log_p_=out_log_p);
00289 
00292                 static Pointer<LocOpt> get_lp_solver(const Pointer<MinlpProblem> prob, Pointer<Param> param, const char* param_prefix=NULL, Pointer<ostream> out_solver_p_=out_out_p, Pointer<ostream> out_solver_log_p_=out_log_p);
00293 
00297                 static Pointer<LocOpt> get_solver_origprob(const Pointer<MinlpProblem> prob, Pointer<Param> param, const char* param_prefix=NULL, Pointer<ostream> out_solver_p_=out_out_p, Pointer<ostream> out_solver_log_p_=out_log_p);
00298 
00299                 LocOpt(int n, Pointer<ostream> out_solver_p_=out_out_p, Pointer<ostream> out_solver_log_p_=out_log_p)
00300                 : Solver(n, out_solver_p_, out_solver_log_p_)
00301                 { }
00302 
00303                 virtual ~LocOpt() { }
00304 
00305                 virtual double time() { return timer; }
00306 
00307                 virtual void reinit() { };
00308 
00309                 virtual dvector get_lag_multipliers()=0;
00310 };
00311 
00312 class SimpleCut;
00313 
00314 
00317 class MIPSolver {
00318         public:
00327                 typedef enum { SOLVED, FEASIBLE, UNBOUNDED, INFEASIBLE, ITERATIONLIMITEXCEEDED, ABORTED, UNKNOWN } SolutionStatus;
00328 
00329                 friend ostream& operator<<(ostream& out, const SolutionStatus status) {
00330                         switch(status) {
00331                                 case SOLVED: out << "Solved"; break;
00332                                 case FEASIBLE: out << "Feasible"; break;
00333                                 case UNBOUNDED: out << "Unbounded"; break;
00334                                 case INFEASIBLE: out << "Infeasible"; break;
00335                                 case ITERATIONLIMITEXCEEDED: out << "Iterationlimit exceeded"; break;
00336                                 case ABORTED: out << "Aborted"; break;
00337                                 case UNKNOWN: out << "Unknown"; break;
00338                                 default: out << "SolutionStatus not known!!";
00339                         }
00340                         return out;
00341                 }
00342 
00347                 class RowItem {
00348                         public:
00349                                 RowItem() { }
00350                                 virtual ~RowItem() { }
00351                 };
00352 
00357                 class ColItem {
00358                         public:
00359                                 ColItem() { }
00360                                 virtual ~ColItem() { }
00361                 };
00362 
00363                 static Pointer<MIPSolver> get_solver(const MipProblem& mip, Pointer<Param> param);
00364 
00365                 virtual ~MIPSolver() { }
00366 
00367                 virtual int nr_col()=0;
00368                 virtual int nr_row()=0;
00369 
00370                 virtual void set_tol(double tol)=0;
00371                 virtual void set_maxiter(int maxiter)=0;
00372 
00376                 virtual void reset()=0;
00377 
00380                 virtual SolutionStatus solve()=0;
00381 
00384                 virtual SolutionStatus solve(const UserVector<double>& x) { return solve(); };
00385 
00386                 virtual SolutionStatus solveMIP()=0;
00387 
00391                 virtual SolutionStatus feasible()=0;
00392 
00397                 virtual void get_primal(UserVector<double>& x)=0;
00398                 virtual dvector get_primal() { dvector x(nr_col()); get_primal(x); return x; }
00401                 virtual double get_primal(const ColItem& colitem)=0;
00402 
00403                 virtual int get_colindex(const ColItem& colitem)=0;
00404 
00409                 virtual void get_dual(UserVector<double>& mu)=0;
00412                 virtual double get_dual(const RowItem& rowitem)=0;
00413 
00418                 virtual void get_reducedcosts(UserVector<double>& rc)=0;
00421                 virtual double get_reducedcosts(const ColItem& colitem)=0;
00422                 
00426                 virtual void get_rowactivity(UserVector<double>& rowact)=0;
00429                 virtual double get_rowactivity(const RowItem& rowitem)=0;
00430                 
00431 
00434                 virtual double get_optval()=0;
00435 
00436                 virtual int get_iter()=0;
00437 
00438                 virtual void set_obj(const UserVector<double>& obj, double obj_const=0.)=0;
00439                 virtual void modify_obj(int i, double coeff)=0;
00440 
00447                 virtual const RowItem* add_row(const UserVector<double>& row, double low, double up)=0;
00448                 virtual const RowItem* add_row(const UserVector<double>& row, const ivector& indices, double low, double up) {
00449                         dvector bigrow(nr_col()); bigrow.set_block(row, indices);
00450                         return add_row(bigrow, low, up);
00451                 };
00452                 virtual void add_rows(const vector<Pointer<UserVector<double> > >& rows, const dvector& low, const dvector& up) {
00453                         for (int i=0; i<rows.size(); i++) add_row(*rows[i], low[i], up[i]);
00454                 }
00455                 virtual void add_rows(const vector<pair<dvector, ivector> >& rows, const dvector& low, const dvector& up) {
00456                         for (int i=0; i<rows.size(); i++) add_row(rows[i].first, rows[i].second, low[i], up[i]);
00457                 } 
00458                 virtual void add_rows(list<const RowItem*>& rowitems, const vector<pair<dvector, ivector> >& rows, const dvector& low, const dvector& up) {
00459                         for (int i=0; i<rows.size(); i++) rowitems.push_back(add_row(rows[i].first, rows[i].second, low[i], up[i]));
00460                 }
00461                 virtual void delete_row(const RowItem& rowitem)=0;
00462                 virtual void delete_rows(const list<const RowItem*>& rowitems) {
00463                         for (list<const RowItem*>::const_iterator it(rowitems.begin()); it!=rowitems.end(); ++it) delete_row(**it);
00464                 }
00465                 virtual void modify_row(const RowItem& rowitem, double low, double up)=0;
00466                 virtual void modify_row(int index, double low, double up)=0;
00467 
00468                 virtual const ColItem* add_col(double low, double up, MipProblem::VarType=MipProblem::CONTINUOUS)=0;
00469                 virtual const ColItem* add_col(const UserVector<double>& col, double obj_coeff, double low, double up, MipProblem::VarType=MipProblem::CONTINUOUS)=0;
00475                 virtual void add_cols(list<const ColItem*>& colitems, const dvector& low, const dvector& up)=0;
00476                 virtual void add_cols(list<const ColItem*>& colitems, vector<Pointer<UserVector<double> > >& cols, const vector<double>& obj_coeff, const dvector& low, const dvector& up)=0;
00477                 virtual void delete_col(const ColItem& colitem)=0;
00478                 virtual void delete_cols(const list<const ColItem*>& colitems) {
00479                         for (list<const ColItem*>::const_iterator it(colitems.begin()); it!=colitems.end(); ++it) delete_col(**it);
00480                 }
00481                 virtual void modify_col(const ColItem& colitem, double low, double up, MipProblem::VarType vartype)=0;
00482                 virtual void modify_col(const ColItem& colitem, const UserVector<double>& col, double obj_coeff, double low, double up, MipProblem::VarType vartype)=0;
00485                 virtual void modify_col(int index, double low, double up, MipProblem::VarType type)=0;
00486                 
00487                 virtual double get_collow(int index)=0;
00488                 virtual double get_colup(int index)=0;
00489                 
00490                 virtual int generate_cuts(list<Pointer<SimpleCut> >& rowcuts)=0;
00491 };
00492 
00495 class LPSolver : public LocOpt {
00496         private:
00497                 Pointer<MIPSolver> mipsolver;
00498                 const Pointer<MinlpProblem> prob;
00499 //      Pointer<Param> param;
00500 //              char* param_prefix;
00501 
00502                 void initmip();
00503 
00504                 // called from solve or solve(x)
00505                 int solved(MIPSolver::SolutionStatus);
00506 
00507         public:
00508                 LPSolver(const Pointer<MinlpProblem> prob_/*, Pointer<Param> param_, char* param_prefix_=NULL*/);
00509 
00510                 dvector get_lag_multipliers();
00511 
00512                 void reinit() { initmip(); }
00513 
00514                 int solve() { return solved(mipsolver->solve()); }
00515                 int solve(dvector& start) { return solved(mipsolver->solve(start)); }
00516 };
00517 
00521 class BoxMinimizer : public Solver {
00522         protected:
00525           const Func& f;
00526 
00529           Pointer<dvector> lower, upper;
00530 
00531         public:
00532           BoxMinimizer(const Func& f_, Pointer<dvector> lower_, Pointer<dvector> upper_)
00533           : Solver(f_.dim()), f(f_), lower(lower_), upper(upper_)
00534           { }
00535 
00536           virtual void set_box(Pointer<dvector> lower_, Pointer<dvector> upper_) {
00537                 lower=lower_;
00538                 upper=upper_;
00539           }
00540 
00541 
00542 };
00543 
00544 class BoxLocOpt : public BoxMinimizer {
00545         private:
00546                 Pointer<LocOpt> locopt;
00547                 Pointer<MinlpProblem> prob;
00548                 Pointer<Param> param;
00549 
00550         public:
00551                 BoxLocOpt(Func &f_, Pointer<dvector> lower_, Pointer<dvector> upper_, Pointer<Param> param_=NULL, char* param_prefix=NULL)
00552                 : BoxMinimizer(f_, lower_, upper_), param(param_),
00553                   prob(new MinlpProblem(new SepQcFunc(NULL, NULL, Pointer<Func>(&f_, false)), lower_ ? *lower_ : dvector(f_.dim()), upper_ ? *upper_ : dvector(f_.dim())))
00554                 { locopt=LocOpt::get_solver(prob, param, param_prefix, NULL, NULL);
00555                 }
00556 
00557                 BoxLocOpt(const SepQcFunc& f_, Pointer<dvector> lower_, Pointer<dvector> upper_, Pointer<Param> param_=NULL, char* param_prefix=NULL)
00558                 : BoxMinimizer(f_, lower_, upper_), param(param_),
00559                         prob(new MinlpProblem(Pointer<SepQcFunc>((SepQcFunc*)&f_, false), lower_ ? *lower_ : dvector(f_.dim()), upper_ ? *upper_ : dvector(f_.dim())))
00560                 { locopt=LocOpt::get_solver(prob, param, param_prefix, NULL, NULL);
00561                 }
00562 
00563                 BoxLocOpt(SepQcFunc& f_, Pointer<dvector> lower_, Pointer<dvector> upper_, Pointer<Param> param_=NULL, char* param_prefix=NULL)
00564                 : BoxMinimizer(f_, lower_, upper_), param(param_),
00565                         prob(new MinlpProblem(Pointer<SepQcFunc>(&f_, false), lower_ ? *lower_ : dvector(f_.dim()), upper_ ? *upper_ : dvector(f_.dim())))
00566                 { locopt=LocOpt::get_solver(prob, param, param_prefix, NULL, NULL);
00567                 }
00568 
00569                 void set_box(Pointer<dvector> lower_, Pointer<dvector> upper_) {
00570                         if (lower==lower_ && upper==upper_) return; // same bounds, so no change needed
00571       BoxMinimizer::set_box(lower_, upper_);
00572       prob->lower=*lower;
00573       prob->upper=*upper;
00574       locopt->iter_max=iter_max;
00575       locopt->reinit();
00576                 }
00577 
00578                 int solve() {
00579                         locopt->iter_max=iter_max;
00580                         int ret=locopt->solve();
00581                         sol_point=locopt->sol_point;
00582                         opt_val_=locopt->opt_val();
00583                         return ret;
00584                 }
00585 
00586                 int solve(dvector& x) {
00587                         locopt->iter_max=iter_max;
00588                         int ret=locopt->solve(x);
00589                         sol_point=locopt->sol_point;
00590                         opt_val_=locopt->opt_val();
00591                         return ret;
00592                 }
00593 
00594 };
00595 
00596 #endif // OPT_H

Generated on Tue Oct 21 03:12:10 2008 for LaGO by  doxygen 1.4.7