Osi 0.108.9
Loading...
Searching...
No Matches
OsiTestSolver.hpp
Go to the documentation of this file.
1// Copyright (C) 2000, International Business Machines
2// Corporation and others. All Rights Reserved.
3// This file is licensed under the terms of Eclipse Public License (EPL).
4
5// this is a copy of VolVolume (stable/1.1 rev. 233)
6
7#ifndef __OSITESTSOLVER_HPP__
8#define __OSITESTSOLVER_HPP__
9
10#include <cfloat>
11#include <algorithm>
12#include <cstdio>
13#include <cmath>
14#include <cstdlib>
15
16#include "CoinFinite.hpp"
17
18#ifndef VOL_DEBUG
19// When VOL_DEBUG is 1, we check vector indices
20#define VOL_DEBUG 0
21#endif
22
23template <class T> static inline T
24VolMax(const T x, const T y) {
25 return ((x) > (y)) ? (x) : (y);
26}
27
28template <class T> static inline T
29VolAbs(const T x) {
30 return ((x) > 0) ? (x) : -(x);
31}
32
33//############################################################################
34
35#if defined(VOL_DEBUG) && (VOL_DEBUG != 0)
36#define VOL_TEST_INDEX(i, size) \
37{ \
38 if ((i) < 0 || (i) >= (size)) { \
39 printf("bad VOL_?vector index\n"); \
40 abort(); \
41 } \
42}
43#define VOL_TEST_SIZE(size) \
44{ \
45 if (s <= 0) { \
46 printf("bad VOL_?vector size\n"); \
47 abort(); \
48 } \
49}
50#else
51#define VOL_TEST_INDEX(i, size)
52#define VOL_TEST_SIZE(size)
53#endif
54
55//############################################################################
56
57class VOL_dvector;
58class VOL_ivector;
59class VOL_primal;
60class VOL_dual;
61class VOL_swing;
63class VOL_vh;
64class VOL_indc;
65class VOL_user_hooks;
66class VOL_problem;
67
68//############################################################################
69
139
140//############################################################################
141
151public:
153 double* v;
155 int sz;
156
157public:
159 VOL_dvector(const int s) {
160 VOL_TEST_SIZE(s);
161 v = new double[sz = s];
162 }
164 VOL_dvector() : v(0), sz(0) {}
166 VOL_dvector(const VOL_dvector& x) : v(0), sz(0) {
167 sz = x.sz;
168 if (sz > 0) {
169 v = new double[sz];
170 std::copy(x.v, x.v + sz, v);
171 }
172 }
174 ~VOL_dvector() { delete[] v; }
175
177 inline int size() const {return sz;}
178
180 inline double& operator[](const int i) {
181 VOL_TEST_INDEX(i, sz);
182 return v[i];
183 }
184
186 inline double operator[](const int i) const {
187 VOL_TEST_INDEX(i, sz);
188 return v[i];
189 }
190
193 inline void clear() {
194 delete[] v;
195 v = 0;
196 sz = 0;
197 }
200 inline void cc(const double gamma, const VOL_dvector& w) {
201 if (sz != w.sz) {
202 printf("bad VOL_dvector sizes\n");
203 abort();
204 }
205 double * p_v = v - 1;
206 const double * p_w = w.v - 1;
207 const double * const p_e = v + sz;
208 const double one_gamma = 1.0 - gamma;
209 while ( ++p_v != p_e ){
210 *p_v = one_gamma * (*p_v) + gamma * (*++p_w);
211 }
212 }
213
216 inline void allocate(const int s) {
217 VOL_TEST_SIZE(s);
218 delete[] v;
219 v = new double[sz = s];
220 }
221
223 inline void swap(VOL_dvector& w) {
224 std::swap(v, w.v);
225 std::swap(sz, w.sz);
226 }
227
231 VOL_dvector& operator=(const double w);
232};
233
234//-----------------------------------------------------------------------------
245public:
247 int* v;
249 int sz;
250public:
252 VOL_ivector(const int s) {
253 VOL_TEST_SIZE(s);
254 v = new int[sz = s];
255 }
257 VOL_ivector() : v(0), sz(0) {}
260 sz = x.sz;
261 if (sz > 0) {
262 v = new int[sz];
263 std::copy(x.v, x.v + sz, v);
264 }
265 }
268 delete [] v;
269 }
270
272 inline int size() const { return sz; }
274 inline int& operator[](const int i) {
275 VOL_TEST_INDEX(i, sz);
276 return v[i];
277 }
278
280 inline int operator[](const int i) const {
281 VOL_TEST_INDEX(i, sz);
282 return v[i];
283 }
284
287 inline void clear() {
288 delete[] v;
289 v = 0;
290 sz = 0;
291 }
292
295 inline void allocate(const int s) {
296 VOL_TEST_SIZE(s);
297 delete[] v;
298 v = new int[sz = s];
299 }
300
302 inline void swap(VOL_ivector& w) {
303 std::swap(v, w.v);
304 std::swap(sz, w.sz);
305 }
306
310 VOL_ivector& operator=(const int w);
311};
312
313//############################################################################
314// A class describing a primal solution. This class is used only internally
316public:
317 // objective value of this primal solution
318 double value;
319 // the largest of the v[i]'s
320 double viol;
321 // primal solution
323 // v=b-Ax, for the relaxed constraints
325
326 VOL_primal(const int psize, const int dsize) : x(psize), v(dsize) {}
327 VOL_primal(const VOL_primal& primal) :
328 value(primal.value), viol(primal.viol), x(primal.x), v(primal.v) {}
330 inline VOL_primal& operator=(const VOL_primal& p) {
331 if (this == &p)
332 return *this;
333 value = p.value;
334 viol = p.viol;
335 x = p.x;
336 v = p.v;
337 return *this;
338 }
339
340 // convex combination. data members in this will be overwritten
341 // convex combination between two primal solutions
342 // x <-- alpha x + (1 - alpha) p.x
343 // v <-- alpha v + (1 - alpha) p.v
344 inline void cc(const double alpha, const VOL_primal& p) {
345 value = alpha * p.value + (1.0 - alpha) * value;
346 x.cc(alpha, p.x);
347 v.cc(alpha, p.v);
348 }
349 // find maximum of v[i]
350 void find_max_viol(const VOL_dvector& dual_lb,
351 const VOL_dvector& dual_ub);
352};
353
354//-----------------------------------------------------------------------------
355// A class describing a dual solution. This class is used only internally
356class VOL_dual {
357public:
358 // lagrangian value
359 double lcost;
360 // reduced costs * (pstar-primal)
361 double xrc;
362 // this information is only printed
363 // dual vector
365
366 VOL_dual(const int dsize) : u(dsize) { u = 0.0;}
367 VOL_dual(const VOL_dual& dual) :
368 lcost(dual.lcost), xrc(dual.xrc), u(dual.u) {}
370 inline VOL_dual& operator=(const VOL_dual& p) {
371 if (this == &p)
372 return *this;
373 lcost = p.lcost;
374 xrc = p.xrc;
375 u = p.u;
376 return *this;
377 }
378 // dual step
379 void step(const double target, const double lambda,
380 const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
381 const VOL_dvector& v);
382 double ascent(const VOL_dvector& v, const VOL_dvector& last_u) const;
383 void compute_xrc(const VOL_dvector& pstarx, const VOL_dvector& primalx,
384 const VOL_dvector& rc);
385
386};
387
388
389//############################################################################
390/* here we check whether an iteration is green, yellow or red. Also according
391 to this information we decide whether lambda should be changed */
393private:
396public:
399 int ngs, nrs, nys;
400 int rd;
401
404 ngs = nrs = nys = 0;
405 }
407
408 inline void cond(const VOL_dual& dual,
409 const double lcost, const double ascent, const int iter) {
410 double eps = 1.e-3;
411
412 if (ascent > 0.0 && lcost > dual.lcost + eps) {
414 lastgreeniter = iter;
415 ++ngs;
416 rd = 0;
417 } else {
418 if (ascent <= 0 && lcost > dual.lcost) {
420 lastyellowiter = iter;
421 ++nys;
422 rd = 0;
423 } else {
424 lastswing = red;
425 lastrediter = iter;
426 ++nrs;
427 rd = 1;
428 }
429 }
430 }
431
432 inline double
433 lfactor(const VOL_parms& parm, const double lambda, const int iter) {
434 double lambdafactor = 1.0;
435 double eps = 5.e-4;
436 int cons;
437
438 switch (lastswing) {
439 case green:
440 cons = iter - VolMax(lastyellowiter, lastrediter);
441 if (parm.printflag & 4)
442 printf(" G: Consecutive Gs = %3d\n\n", cons);
443 if (cons >= parm.greentestinvl && lambda < 2.0) {
445 lambdafactor = 2.0;
446 if (parm.printflag & 2)
447 printf("\n ---- increasing lamda to %g ----\n\n",
448 lambda * lambdafactor);
449 }
450 break;
451
452 case yellow:
453 cons = iter - VolMax(lastgreeniter, lastrediter);
454 if (parm.printflag & 4)
455 printf(" Y: Consecutive Ys = %3d\n\n", cons);
456 if (cons >= parm.yellowtestinvl) {
458 lambdafactor = 1.1;
459 if (parm.printflag & 2)
460 printf("\n **** increasing lamda to %g *****\n\n",
461 lambda * lambdafactor);
462 }
463 break;
464
465 case red:
466 cons = iter - VolMax(lastgreeniter, lastyellowiter);
467 if (parm.printflag & 4)
468 printf(" R: Consecutive Rs = %3d\n\n", cons);
469 if (cons >= parm.redtestinvl && lambda > eps) {
471 lambdafactor = 0.67;
472 if (parm.printflag & 2)
473 printf("\n **** decreasing lamda to %g *****\n\n",
474 lambda * lambdafactor);
475 }
476 break;
477 }
478 return lambdafactor;
479 }
480
481 inline void
483 printf("**** G= %i, Y= %i, R= %i ****\n", ngs, nys, nrs);
484 ngs = nrs = nys = 0;
485 }
486};
487
488//############################################################################
489/* alpha should be decreased if after some number of iterations the objective
490 has increased less that 1% */
492private:
495public:
496 double lastvalue;
497
498 VOL_alpha_factor() {lastvalue = -COIN_DBL_MAX;}
500
501 inline double factor(const VOL_parms& parm, const double lcost,
502 const double alpha) {
503 if (alpha < parm.alphamin)
504 return 1.0;
505 const double ll = VolAbs(lcost);
506 const double x = ll > 10 ? (lcost-lastvalue)/ll : (lcost-lastvalue);
507 lastvalue = lcost;
508 return (x <= 0.01) ? parm.alphafactor : 1.0;
509 }
510};
511
512//############################################################################
513/* here we compute the norm of the conjugate direction -hh-, the norm of the
514 subgradient -norm-, the inner product between the subgradient and the
515 last conjugate direction -vh-, and the inner product between the new
516 conjugate direction and the subgradient */
517class VOL_vh {
518private:
519 VOL_vh(const VOL_vh&);
521public:
522 double hh;
523 double norm;
524 double vh;
525 double asc;
526
527 VOL_vh(const double alpha,
528 const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
529 const VOL_dvector& v, const VOL_dvector& vstar,
530 const VOL_dvector& u);
532};
533
534//############################################################################
535/* here we compute different parameter to be printed. v2 is the square of
536 the norm of the subgradient. vu is the inner product between the dual
537 variables and the subgradient. vabs is the maximum absolute value of
538 the violations of pstar. asc is the inner product between the conjugate
539 direction and the subgradient */
540class VOL_indc {
541private:
544public:
545 double v2;
546 double vu;
547 double vabs;
548 double asc;
549
550public:
551 VOL_indc(const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
552 const VOL_primal& primal, const VOL_primal& pstar,
553 const VOL_dual& dual);
555};
556
557//#############################################################################
558
566public:
567 virtual ~VOL_user_hooks() {}
568public:
569 // for all hooks: return value of -1 means that volume should quit
574 virtual int compute_rc(const VOL_dvector& u, VOL_dvector& rc) = 0;
575
584 virtual int solve_subproblem(const VOL_dvector& dual, const VOL_dvector& rc,
585 double& lcost, VOL_dvector& x, VOL_dvector& v,
586 double& pcost) = 0;
593 virtual int heuristics(const VOL_problem& p,
594 const VOL_dvector& x, double& heur_val) = 0;
595};
596
597//#############################################################################
598
608private:
612 // ############ INPUT fields ########################
613public:
620 VOL_problem(const char *filename);
624
630 int solve(VOL_user_hooks& hooks, const bool use_preset_dual = false);
632
633private:
637 double alpha_;
639 double lambda_;
640 // This union is here for padding (so that data members would be
641 // double-aligned on x86 CPU
642 union {
644 int iter_;
645 double __pad0;
646 };
648
649public:
650
654 double value;
662
668 int psize;
670 int dsize;
678
679public:
683 int iter() const { return iter_; }
685 double alpha() const { return alpha_; }
687 double lambda() const { return lambda_; }
689
690private:
694 void read_params(const char* filename);
695
697 int initialize(const bool use_preset_dual);
698
700 void print_info(const int iter,
701 const VOL_primal& primal, const VOL_primal& pstar,
702 const VOL_dual& dual);
703
706 double readjust_target(const double oldtarget, const double lcost) const;
707
715 double power_heur(const VOL_primal& primal, const VOL_primal& pstar,
716 const VOL_dual& dual) const;
718};
719
720#endif
static T VolAbs(const T x)
#define VOL_TEST_INDEX(i, size)
static T VolMax(const T x, const T y)
#define VOL_TEST_SIZE(size)
double factor(const VOL_parms &parm, const double lcost, const double alpha)
VOL_alpha_factor(const VOL_alpha_factor &)
VOL_alpha_factor & operator=(const VOL_alpha_factor &)
void compute_xrc(const VOL_dvector &pstarx, const VOL_dvector &primalx, const VOL_dvector &rc)
VOL_dvector u
void step(const double target, const double lambda, const VOL_dvector &dual_lb, const VOL_dvector &dual_ub, const VOL_dvector &v)
VOL_dual(const VOL_dual &dual)
double ascent(const VOL_dvector &v, const VOL_dvector &last_u) const
VOL_dual & operator=(const VOL_dual &p)
VOL_dual(const int dsize)
vector of doubles.
void cc(const double gamma, const VOL_dvector &w)
Convex combination.
void allocate(const int s)
delete the current vector and allocate space for a vector of size s.
double * v
The array holding the vector.
double & operator[](const int i)
Return a reference to the i-th entry.
VOL_dvector(const int s)
Construct a vector of size s.
double operator[](const int i) const
Return the i-th entry.
void clear()
Delete the content of the vector and replace it with a vector of length 0.
~VOL_dvector()
The destructor deletes the data array.
VOL_dvector(const VOL_dvector &x)
Copy constructor makes a replica of x.
VOL_dvector & operator=(const VOL_dvector &w)
Copy w into the vector.
VOL_dvector()
Default constructor creates a vector of size 0.
VOL_dvector & operator=(const double w)
Replace every entry in the vector with w.
void swap(VOL_dvector &w)
swaps the vector with w.
int sz
The size of the vector.
int size() const
Return the size of the vector.
VOL_indc(const VOL_indc &)
VOL_indc(const VOL_dvector &dual_lb, const VOL_dvector &dual_ub, const VOL_primal &primal, const VOL_primal &pstar, const VOL_dual &dual)
VOL_indc & operator=(const VOL_indc &)
vector of ints.
VOL_ivector(const VOL_ivector &x)
Copy constructor makes a replica of x.
void swap(VOL_ivector &w)
swaps the vector with w.
~VOL_ivector()
The destructor deletes the data array.
VOL_ivector()
Default constructor creates a vector of size 0.
VOL_ivector & operator=(const VOL_ivector &v)
Copy w into the vector.
void clear()
Delete the content of the vector and replace it with a vector of length 0.
VOL_ivector & operator=(const int w)
Replace every entry in the vector with w.
int operator[](const int i) const
Return the i-th entry.
int * v
The array holding the vector.
VOL_ivector(const int s)
Construct a vector of size s.
int size() const
Return the size of the vector.
int sz
The size of the vector.
int & operator[](const int i)
Return a reference to the i-th entry.
void allocate(const int s)
delete the current vector and allocate space for a vector of size s.
VOL_primal(const int psize, const int dsize)
VOL_dvector v
void find_max_viol(const VOL_dvector &dual_lb, const VOL_dvector &dual_ub)
VOL_dvector x
VOL_primal & operator=(const VOL_primal &p)
VOL_primal(const VOL_primal &primal)
void cc(const double alpha, const VOL_primal &p)
This class holds every data for the Volume Algorithm and its solve method must be invoked to solve th...
double alpha_
value of alpha
double value
final lagrangian value (OUTPUT)
void print_info(const int iter, const VOL_primal &primal, const VOL_primal &pstar, const VOL_dual &dual)
print volume info every parm.printinvl iterations
VOL_dvector psol
final primal solution (OUTPUT)
void read_params(const char *filename)
Read in the parameters from the file filename.
int dsize
length of dual solution (INPUT)
VOL_dvector dual_lb
lower bounds for the duals (if 0 length, then filled with -inf) (INPUT)
double lambda() const
returns the value of lambda
VOL_parms parm
The parameters controlling the Volume Algorithm (INPUT)
int iter() const
returns the iteration number
double alpha() const
returns the value of alpha
VOL_problem & operator=(const VOL_problem &)
double readjust_target(const double oldtarget, const double lcost) const
Checks if lcost is close to the target, if so it increases the target.
int initialize(const bool use_preset_dual)
initializes duals, bounds for the duals, alpha, lambda
VOL_problem(const VOL_problem &)
VOL_dvector viol
violations (b-Ax) for the relaxed constraints
VOL_problem(const char *filename)
Create a a VOL_problem object and read in the parameters from filename.
double power_heur(const VOL_primal &primal, const VOL_primal &pstar, const VOL_dual &dual) const
Here we decide the value of alpha1 to be used in the convex combination.
VOL_problem()
Default constructor.
VOL_dvector dual_ub
upper bounds for the duals (if 0 length, then filled with +inf) (INPUT)
int psize
length of primal solution (INPUT)
VOL_dvector dsol
final dual solution (INPUT/OUTPUT)
int solve(VOL_user_hooks &hooks, const bool use_preset_dual=false)
Solve the problem using the hooks.
void set_default_parm()
double lambda_
value of lambda
int iter_
iteration number
~VOL_problem()
Destruct the object.
double lfactor(const VOL_parms &parm, const double lambda, const int iter)
enum VOL_swing::condition lastswing
VOL_swing & operator=(const VOL_swing &)
VOL_swing(const VOL_swing &)
void cond(const VOL_dual &dual, const double lcost, const double ascent, const int iter)
The user hooks should be overridden by the user to provide the problem specific routines for the volu...
virtual int compute_rc(const VOL_dvector &u, VOL_dvector &rc)=0
compute reduced costs
virtual ~VOL_user_hooks()
virtual int solve_subproblem(const VOL_dvector &dual, const VOL_dvector &rc, double &lcost, VOL_dvector &x, VOL_dvector &v, double &pcost)=0
Solve the subproblem for the subgradient step.
virtual int heuristics(const VOL_problem &p, const VOL_dvector &x, double &heur_val)=0
Starting from the primal vector x, run a heuristic to produce an integer solution
VOL_vh(const double alpha, const VOL_dvector &dual_lb, const VOL_dvector &dual_ub, const VOL_dvector &v, const VOL_dvector &vstar, const VOL_dvector &u)
VOL_vh & operator=(const VOL_vh &)
double norm
VOL_vh(const VOL_vh &)
This class contains the parameters controlling the Volume Algorithm.
int printflag
controls the level of printing.
double lambdainit
initial value of lambda
double gap_abs_precision
accept if abs gap is less than this
char * temp_dualfile
name of file for saving dual solution
double alphamin
minimum value for alpha
int greentestinvl
how many consecutive green iterations are allowed before changing lambda
int maxsgriters
maximum number of iterations
double minimum_rel_ascent
terminate if the relative increase in lcost through ascent_check_invl steps is less than this
int heurinvl
controls how often we run the primal heuristic
double granularity
terminate if best_ub - lcost < granularity
int alphaint
number of iterations before we check if alpha should be decreased
int redtestinvl
how many consecutive red iterations are allowed before changing lambda
int printinvl
controls how often do we print
double alphafactor
when little progress is being done, we multiply alpha by alphafactor
double gap_rel_precision
accept if rel gap is less than this
int ascent_first_check
when to check for sufficient relative ascent the first time
double alphainit
initial value of alpha
int ascent_check_invl
through how many iterations does the relative ascent have to reach a minimum
double primal_abs_precision
accept if max abs viol is less than this
int yellowtestinvl
how many consecutive yellow iterations are allowed before changing lambda
double ubinit
initial upper bound of the value of an integer solution