VolVolume.hpp
Go to the documentation of this file.
1 // Copyright (C) 2000, International Business Machines
2 // Corporation and others. All Rights Reserved.
3 #ifndef __VOLUME_HPP__
4 #define __VOLUME_HPP__
5 
6 #include <cfloat>
7 #include <algorithm>
8 #include <cstdio>
9 #include <cmath>
10 
11 #ifndef VOL_DEBUG
12 // When VOL_DEBUG is 1, we check vector indices
13 #define VOL_DEBUG 0
14 #endif
15 
16 template <class T> static inline T
17 VolMax(register const T x, register const T y) {
18  return ((x) > (y)) ? (x) : (y);
19 }
20 
21 template <class T> static inline T
22 VolAbs(register const T x) {
23  return ((x) > 0) ? (x) : -(x);
24 }
25 
26 //############################################################################
27 
28 #if defined(VOL_DEBUG) && (VOL_DEBUG != 0)
29 #define VOL_TEST_INDEX(i, size) \
30 { \
31  if ((i) < 0 || (i) >= (size)) { \
32  printf("bad VOL_?vector index\n"); \
33  abort(); \
34  } \
35 }
36 #define VOL_TEST_SIZE(size) \
37 { \
38  if (s <= 0) { \
39  printf("bad VOL_?vector size\n"); \
40  abort(); \
41  } \
42 }
43 #else
44 #define VOL_TEST_INDEX(i, size)
45 #define VOL_TEST_SIZE(size)
46 #endif
47 
48 //############################################################################
49 
50 class VOL_dvector;
51 class VOL_ivector;
52 class VOL_primal;
53 class VOL_dual;
54 class VOL_swing;
55 class VOL_alpha_factor;
56 class VOL_vh;
57 class VOL_indc;
58 class VOL_user_hooks;
59 class VOL_problem;
60 
61 //############################################################################
62 
66 struct VOL_parms {
68  double lambdainit;
70  double alphainit;
72  double alphamin;
74  double alphafactor;
75 
77  double ubinit;
78 
86  double granularity;
87 
96 
99 
110  int printflag;
112  int printinvl;
114  int heurinvl;
115 
125 
127  int alphaint;
128 
131 };
132 
133 //############################################################################
134 
143 class VOL_dvector {
144 public:
146  double* v;
148  int sz;
149 
150 public:
152  VOL_dvector(const int s) {
153  VOL_TEST_SIZE(s);
154  v = new double[sz = s];
155  }
157  VOL_dvector() : v(0), sz(0) {}
159  VOL_dvector(const VOL_dvector& x) : v(0), sz(0) {
160  sz = x.sz;
161  if (sz > 0) {
162  v = new double[sz];
163  std::copy(x.v, x.v + sz, v);
164  }
165  }
167  ~VOL_dvector() { delete[] v; }
168 
170  inline int size() const {return sz;}
171 
173  inline double& operator[](const int i) {
174  VOL_TEST_INDEX(i, sz);
175  return v[i];
176  }
177 
179  inline double operator[](const int i) const {
180  VOL_TEST_INDEX(i, sz);
181  return v[i];
182  }
183 
186  inline void clear() {
187  delete[] v;
188  v = 0;
189  sz = 0;
190  }
193  inline void cc(const double gamma, const VOL_dvector& w) {
194  if (sz != w.sz) {
195  printf("bad VOL_dvector sizes\n");
196  abort();
197  }
198  double * p_v = v - 1;
199  const double * p_w = w.v - 1;
200  const double * const p_e = v + sz;
201  const double one_gamma = 1.0 - gamma;
202  while ( ++p_v != p_e ){
203  *p_v = one_gamma * (*p_v) + gamma * (*++p_w);
204  }
205  }
206 
209  inline void allocate(const int s) {
210  VOL_TEST_SIZE(s);
211  delete[] v;
212  v = new double[sz = s];
213  }
214 
216  inline void swap(VOL_dvector& w) {
217  std::swap(v, w.v);
218  std::swap(sz, w.sz);
219  }
220 
222  VOL_dvector& operator=(const VOL_dvector& w);
224  VOL_dvector& operator=(const double w);
225 };
226 
227 //-----------------------------------------------------------------------------
237 class VOL_ivector {
238 public:
240  int* v;
242  int sz;
243 public:
245  VOL_ivector(const int s) {
246  VOL_TEST_SIZE(s);
247  v = new int[sz = s];
248  }
250  VOL_ivector() : v(0), sz(0) {}
253  sz = x.sz;
254  if (sz > 0) {
255  v = new int[sz];
256  std::copy(x.v, x.v + sz, v);
257  }
258  }
261  delete [] v;
262  }
263 
265  inline int size() const { return sz; }
267  inline int& operator[](const int i) {
268  VOL_TEST_INDEX(i, sz);
269  return v[i];
270  }
271 
273  inline int operator[](const int i) const {
274  VOL_TEST_INDEX(i, sz);
275  return v[i];
276  }
277 
280  inline void clear() {
281  delete[] v;
282  v = 0;
283  sz = 0;
284  }
285 
288  inline void allocate(const int s) {
289  VOL_TEST_SIZE(s);
290  delete[] v;
291  v = new int[sz = s];
292  }
293 
295  inline void swap(VOL_ivector& w) {
296  std::swap(v, w.v);
297  std::swap(sz, w.sz);
298  }
299 
301  VOL_ivector& operator=(const VOL_ivector& v);
303  VOL_ivector& operator=(const int w);
304 };
305 
306 //############################################################################
307 // A class describing a primal solution. This class is used only internally
308 class VOL_primal {
309 public:
310  // objective value of this primal solution
311  double value;
312  // the largest of the v[i]'s
313  double viol;
314  // primal solution
316  // v=b-Ax, for the relaxed constraints
318 
319  VOL_primal(const int psize, const int dsize) : x(psize), v(dsize) {}
320  VOL_primal(const VOL_primal& primal) :
321  value(primal.value), viol(primal.viol), x(primal.x), v(primal.v) {}
323  inline VOL_primal& operator=(const VOL_primal& p) {
324  if (this == &p)
325  return *this;
326  value = p.value;
327  viol = p.viol;
328  x = p.x;
329  v = p.v;
330  return *this;
331  }
332 
333  // convex combination. data members in this will be overwritten
334  // convex combination between two primal solutions
335  // x <-- alpha x + (1 - alpha) p.x
336  // v <-- alpha v + (1 - alpha) p.v
337  inline void cc(const double alpha, const VOL_primal& p) {
338  value = alpha * p.value + (1.0 - alpha) * value;
339  x.cc(alpha, p.x);
340  v.cc(alpha, p.v);
341  }
342  // find maximum of v[i]
343  void find_max_viol(const VOL_dvector& dual_lb,
344  const VOL_dvector& dual_ub);
345 };
346 
347 //-----------------------------------------------------------------------------
348 // A class describing a dual solution. This class is used only internally
349 class VOL_dual {
350 public:
351  // lagrangian value
352  double lcost;
353  // reduced costs * (pstar-primal)
354  double xrc;
355  // this information is only printed
356  // dual vector
358 
359  VOL_dual(const int dsize) : u(dsize) { u = 0.0;}
360  VOL_dual(const VOL_dual& dual) :
361  lcost(dual.lcost), xrc(dual.xrc), u(dual.u) {}
363  inline VOL_dual& operator=(const VOL_dual& p) {
364  if (this == &p)
365  return *this;
366  lcost = p.lcost;
367  xrc = p.xrc;
368  u = p.u;
369  return *this;
370  }
371  // dual step
372  void step(const double target, const double lambda,
373  const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
374  const VOL_dvector& v);
375  double ascent(const VOL_dvector& v, const VOL_dvector& last_u) const;
376  void compute_xrc(const VOL_dvector& pstarx, const VOL_dvector& primalx,
377  const VOL_dvector& rc);
378 
379 };
380 
381 
382 //############################################################################
383 /* here we check whether an iteration is green, yellow or red. Also according
384  to this information we decide whether lambda should be changed */
385 class VOL_swing {
386 private:
387  VOL_swing(const VOL_swing&);
388  VOL_swing& operator=(const VOL_swing&);
389 public:
390  enum condition {green, yellow, red} lastswing;
391  int lastgreeniter, lastyellowiter, lastrediter;
392  int ngs, nrs, nys;
393  int rd;
394 
396  lastgreeniter = lastyellowiter = lastrediter = 0;
397  ngs = nrs = nys = 0;
398  }
400 
401  inline void cond(const VOL_dual& dual,
402  const double lcost, const double ascent, const int iter) {
403  double eps = 1.e-3;
404 
405  if (ascent > 0.0 && lcost > dual.lcost + eps) {
406  lastswing = green;
407  lastgreeniter = iter;
408  ++ngs;
409  rd = 0;
410  } else {
411  if (ascent <= 0 && lcost > dual.lcost) {
412  lastswing = yellow;
413  lastyellowiter = iter;
414  ++nys;
415  rd = 0;
416  } else {
417  lastswing = red;
418  lastrediter = iter;
419  ++nrs;
420  rd = 1;
421  }
422  }
423  }
424 
425  inline double
426  lfactor(const VOL_parms& parm, const double lambda, const int iter) {
427  double lambdafactor = 1.0;
428  double eps = 5.e-4;
429  int cons;
430 
431  switch (lastswing) {
432  case green:
433  cons = iter - VolMax(lastyellowiter, lastrediter);
434  if (parm.printflag & 4)
435  printf(" G: Consecutive Gs = %3d\n\n", cons);
436  if (cons >= parm.greentestinvl && lambda < 2.0) {
437  lastgreeniter = lastyellowiter = lastrediter = iter;
438  lambdafactor = 2.0;
439  if (parm.printflag & 2)
440  printf("\n ---- increasing lamda to %g ----\n\n",
441  lambda * lambdafactor);
442  }
443  break;
444 
445  case yellow:
446  cons = iter - VolMax(lastgreeniter, lastrediter);
447  if (parm.printflag & 4)
448  printf(" Y: Consecutive Ys = %3d\n\n", cons);
449  if (cons >= parm.yellowtestinvl) {
450  lastgreeniter = lastyellowiter = lastrediter = iter;
451  lambdafactor = 1.1;
452  if (parm.printflag & 2)
453  printf("\n **** increasing lamda to %g *****\n\n",
454  lambda * lambdafactor);
455  }
456  break;
457 
458  case red:
459  cons = iter - VolMax(lastgreeniter, lastyellowiter);
460  if (parm.printflag & 4)
461  printf(" R: Consecutive Rs = %3d\n\n", cons);
462  if (cons >= parm.redtestinvl && lambda > eps) {
463  lastgreeniter = lastyellowiter = lastrediter = iter;
464  lambdafactor = 0.67;
465  if (parm.printflag & 2)
466  printf("\n **** decreasing lamda to %g *****\n\n",
467  lambda * lambdafactor);
468  }
469  break;
470  }
471  return lambdafactor;
472  }
473 
474  inline void
475  print() {
476  printf("**** G= %i, Y= %i, R= %i ****\n", ngs, nys, nrs);
477  ngs = nrs = nys = 0;
478  }
479 };
480 
481 //############################################################################
482 /* alpha should be decreased if after some number of iterations the objective
483  has increased less that 1% */
485 private:
487  VOL_alpha_factor& operator=(const VOL_alpha_factor&);
488 public:
489  double lastvalue;
490 
491  VOL_alpha_factor() {lastvalue = -DBL_MAX;}
493 
494  inline double factor(const VOL_parms& parm, const double lcost,
495  const double alpha) {
496  if (alpha < parm.alphamin)
497  return 1.0;
498  register const double ll = VolAbs(lcost);
499  const double x = ll > 10 ? (lcost-lastvalue)/ll : (lcost-lastvalue);
500  lastvalue = lcost;
501  return (x <= 0.01) ? parm.alphafactor : 1.0;
502  }
503 };
504 
505 //############################################################################
506 /* here we compute the norm of the conjugate direction -hh-, the norm of the
507  subgradient -norm-, the inner product between the subgradient and the
508  last conjugate direction -vh-, and the inner product between the new
509  conjugate direction and the subgradient */
510 class VOL_vh {
511 private:
512  VOL_vh(const VOL_vh&);
513  VOL_vh& operator=(const VOL_vh&);
514 public:
515  double hh;
516  double norm;
517  double vh;
518  double asc;
519 
520  VOL_vh(const double alpha,
521  const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
522  const VOL_dvector& v, const VOL_dvector& vstar,
523  const VOL_dvector& u);
525 };
526 
527 //############################################################################
528 /* here we compute different parameter to be printed. v2 is the square of
529  the norm of the subgradient. vu is the inner product between the dual
530  variables and the subgradient. vabs is the maximum absolute value of
531  the violations of pstar. asc is the inner product between the conjugate
532  direction and the subgradient */
533 class VOL_indc {
534 private:
535  VOL_indc(const VOL_indc&);
536  VOL_indc& operator=(const VOL_indc&);
537 public:
538  double v2;
539  double vu;
540  double vabs;
541  double asc;
542 
543 public:
544  VOL_indc(const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
545  const VOL_primal& primal, const VOL_primal& pstar,
546  const VOL_dual& dual);
548 };
549 
550 //#############################################################################
551 
559 public:
560  virtual ~VOL_user_hooks() {}
561 public:
562  // for all hooks: return value of -1 means that volume should quit
567  virtual int compute_rc(const VOL_dvector& u, VOL_dvector& rc) = 0;
568 
577  virtual int solve_subproblem(const VOL_dvector& dual, const VOL_dvector& rc,
578  double& lcost, VOL_dvector& x, VOL_dvector& v,
579  double& pcost) = 0;
586  virtual int heuristics(const VOL_problem& p,
587  const VOL_dvector& x, double& heur_val) = 0;
588 };
589 
590 //#############################################################################
591 
600 class VOL_problem {
601 private:
602  VOL_problem(const VOL_problem&);
603  VOL_problem& operator=(const VOL_problem&);
604  void set_default_parm();
605  // ############ INPUT fields ########################
606 public:
610  VOL_problem();
613  VOL_problem(const char *filename);
615  ~VOL_problem();
617 
623  int solve(VOL_user_hooks& hooks, const bool use_preset_dual = false);
625 
626 private:
630  double alpha_;
632  double lambda_;
633  // This union is here for padding (so that data members would be
634  // double-aligned on x86 CPU
635  union {
637  int iter_;
638  double __pad0;
639  };
641 
642 public:
643 
647  double value;
655 
661  int psize;
663  int dsize;
671 
672 public:
676  int iter() const { return iter_; }
678  double alpha() const { return alpha_; }
680  double lambda() const { return lambda_; }
682 
683 private:
687  void read_params(const char* filename);
688 
690  int initialize(const bool use_preset_dual);
691 
693  void print_info(const int iter,
694  const VOL_primal& primal, const VOL_primal& pstar,
695  const VOL_dual& dual);
696 
699  double readjust_target(const double oldtarget, const double lcost) const;
700 
708  double power_heur(const VOL_primal& primal, const VOL_primal& pstar,
709  const VOL_dual& dual) const;
711 };
712 
713 #endif
double value
Definition: VolVolume.hpp:311
int size() const
Return the size of the vector.
Definition: VolVolume.hpp:170
void swap(VOL_dvector &w)
swaps the vector with w.
Definition: VolVolume.hpp:216
#define VOL_TEST_SIZE(size)
Definition: VolVolume.hpp:45
double & operator[](const int i)
Return a reference to the i-th entry.
Definition: VolVolume.hpp:173
char * temp_dualfile
name of file for saving dual solution
Definition: VolVolume.hpp:130
VOL_dvector viol
violations (b-Ax) for the relaxed constraints
Definition: VolVolume.hpp:653
VOL_dvector dual_ub
upper bounds for the duals (if 0 length, then filled with +inf) (INPUT)
Definition: VolVolume.hpp:669
int iter_
iteration number
Definition: VolVolume.hpp:637
vector of ints.
Definition: VolVolume.hpp:237
int * v
The array holding the vector.
Definition: VolVolume.hpp:240
VOL_ivector(const int s)
Construct a vector of size s.
Definition: VolVolume.hpp:245
VOL_dual(const int dsize)
Definition: VolVolume.hpp:359
void cond(const VOL_dual &dual, const double lcost, const double ascent, const int iter)
Definition: VolVolume.hpp:401
double alpha_
value of alpha
Definition: VolVolume.hpp:630
double v2
Definition: VolVolume.hpp:538
int sz
The size of the vector.
Definition: VolVolume.hpp:242
VOL_dvector(const VOL_dvector &x)
Copy constructor makes a replica of x.
Definition: VolVolume.hpp:159
double hh
Definition: VolVolume.hpp:515
double alphafactor
when little progress is being done, we multiply alpha by alphafactor
Definition: VolVolume.hpp:74
double * v
The array holding the vector.
Definition: VolVolume.hpp:146
double lfactor(const VOL_parms &parm, const double lambda, const int iter)
Definition: VolVolume.hpp:426
VOL_dvector psol
final primal solution (OUTPUT)
Definition: VolVolume.hpp:651
virtual ~VOL_user_hooks()
Definition: VolVolume.hpp:560
int & operator[](const int i)
Return a reference to the i-th entry.
Definition: VolVolume.hpp:267
VOL_dual(const VOL_dual &dual)
Definition: VolVolume.hpp:360
double gap_rel_precision
accept if rel gap is less than this
Definition: VolVolume.hpp:84
VOL_dual & operator=(const VOL_dual &p)
Definition: VolVolume.hpp:363
double alpha() const
returns the value of alpha
Definition: VolVolume.hpp:678
VOL_dvector dual_lb
lower bounds for the duals (if 0 length, then filled with -inf) (INPUT)
Definition: VolVolume.hpp:666
VOL_ivector()
Default constructor creates a vector of size 0.
Definition: VolVolume.hpp:250
int size() const
Return the size of the vector.
Definition: VolVolume.hpp:265
void print()
Definition: VolVolume.hpp:475
VOL_ivector(const VOL_ivector &x)
Copy constructor makes a replica of x.
Definition: VolVolume.hpp:252
~VOL_dvector()
The destructor deletes the data array.
Definition: VolVolume.hpp:167
double lambda() const
returns the value of lambda
Definition: VolVolume.hpp:680
int ascent_first_check
when to check for sufficient relative ascent the first time
Definition: VolVolume.hpp:92
VOL_primal(const VOL_primal &primal)
Definition: VolVolume.hpp:320
double value
final lagrangian value (OUTPUT)
Definition: VolVolume.hpp:647
double alphainit
initial value of alpha
Definition: VolVolume.hpp:70
vector of doubles.
Definition: VolVolume.hpp:143
VOL_dvector u
Definition: VolVolume.hpp:357
void clear()
Delete the content of the vector and replace it with a vector of length 0.
Definition: VolVolume.hpp:280
double gap_abs_precision
accept if abs gap is less than this
Definition: VolVolume.hpp:82
double primal_abs_precision
accept if max abs viol is less than this
Definition: VolVolume.hpp:80
static T VolAbs(register const T x)
Definition: VolVolume.hpp:22
double lambdainit
initial value of lambda
Definition: VolVolume.hpp:68
double vh
Definition: VolVolume.hpp:517
void clear()
Delete the content of the vector and replace it with a vector of length 0.
Definition: VolVolume.hpp:186
This class holds every data for the Volume Algorithm and its solve method must be invoked to solve th...
Definition: VolVolume.hpp:600
double minimum_rel_ascent
terminate if the relative increase in lcost through ascent_check_invl steps is less than this ...
Definition: VolVolume.hpp:90
void cc(const double alpha, const VOL_primal &p)
Definition: VolVolume.hpp:337
double asc
Definition: VolVolume.hpp:518
int psize
length of primal solution (INPUT)
Definition: VolVolume.hpp:661
int redtestinvl
how many consecutive red iterations are allowed before changing lambda
Definition: VolVolume.hpp:124
double norm
Definition: VolVolume.hpp:516
void allocate(const int s)
delete the current vector and allocate space for a vector of size s.
Definition: VolVolume.hpp:288
double viol
Definition: VolVolume.hpp:313
VOL_primal(const int psize, const int dsize)
Definition: VolVolume.hpp:319
double vu
Definition: VolVolume.hpp:539
int lastyellowiter
Definition: VolVolume.hpp:391
int yellowtestinvl
how many consecutive yellow iterations are allowed before changing lambda
Definition: VolVolume.hpp:121
int dsize
length of dual solution (INPUT)
Definition: VolVolume.hpp:663
double lambda_
value of lambda
Definition: VolVolume.hpp:632
double ubinit
initial upper bound of the value of an integer solution
Definition: VolVolume.hpp:77
This class contains the parameters controlling the Volume Algorithm.
Definition: VolVolume.hpp:66
int greentestinvl
how many consecutive green iterations are allowed before changing lambda
Definition: VolVolume.hpp:118
VOL_primal & operator=(const VOL_primal &p)
Definition: VolVolume.hpp:323
double granularity
terminate if best_ub - lcost < granularity
Definition: VolVolume.hpp:86
int ascent_check_invl
through how many iterations does the relative ascent have to reach a minimum
Definition: VolVolume.hpp:95
int heurinvl
controls how often we run the primal heuristic
Definition: VolVolume.hpp:114
VOL_dvector(const int s)
Construct a vector of size s.
Definition: VolVolume.hpp:152
#define VOL_TEST_INDEX(i, size)
Definition: VolVolume.hpp:44
int iter() const
returns the iteration number
Definition: VolVolume.hpp:676
VOL_parms parm
The parameters controlling the Volume Algorithm (INPUT)
Definition: VolVolume.hpp:659
double lcost
Definition: VolVolume.hpp:352
double factor(const VOL_parms &parm, const double lcost, const double alpha)
Definition: VolVolume.hpp:494
~VOL_ivector()
The destructor deletes the data array.
Definition: VolVolume.hpp:260
void swap(VOL_ivector &w)
swaps the vector with w.
Definition: VolVolume.hpp:295
double __pad0
Definition: VolVolume.hpp:638
void cc(const double gamma, const VOL_dvector &w)
Convex combination.
Definition: VolVolume.hpp:193
VOL_dvector dsol
final dual solution (INPUT/OUTPUT)
Definition: VolVolume.hpp:649
int printinvl
controls how often do we print
Definition: VolVolume.hpp:112
double vabs
Definition: VolVolume.hpp:540
The user hooks should be overridden by the user to provide the problem specific routines for the volu...
Definition: VolVolume.hpp:558
double operator[](const int i) const
Return the i-th entry.
Definition: VolVolume.hpp:179
int sz
The size of the vector.
Definition: VolVolume.hpp:148
double asc
Definition: VolVolume.hpp:541
VOL_dvector x
Definition: VolVolume.hpp:315
int operator[](const int i) const
Return the i-th entry.
Definition: VolVolume.hpp:273
int alphaint
number of iterations before we check if alpha should be decreased
Definition: VolVolume.hpp:127
static T VolMax(register const T x, register const T y)
Definition: VolVolume.hpp:17
double xrc
Definition: VolVolume.hpp:354
double alphamin
minimum value for alpha
Definition: VolVolume.hpp:72
int printflag
controls the level of printing.
Definition: VolVolume.hpp:110
VOL_dvector v
Definition: VolVolume.hpp:317
int maxsgriters
maximum number of iterations
Definition: VolVolume.hpp:98
void allocate(const int s)
delete the current vector and allocate space for a vector of size s.
Definition: VolVolume.hpp:209
VOL_dvector()
Default constructor creates a vector of size 0.
Definition: VolVolume.hpp:157