Vol  1.5.4
VolVolume.hpp
Go to the documentation of this file.
1 // $Id: VolVolume.hpp 375 2013-11-22 15:46:16Z tkr $
2 // Copyright (C) 2000, International Business Machines
3 // Corporation and others. All Rights Reserved.
4 // This code is licensed under the terms of the Eclipse Public License (EPL).
5 
6 #ifndef __VOLUME_HPP__
7 #define __VOLUME_HPP__
8 
9 #include <cfloat>
10 #include <algorithm>
11 #include <cstdio>
12 #include <cmath>
13 #include <cstring>
14 #include <cstdlib>
15 
16 #ifndef VOL_DEBUG
17 // When VOL_DEBUG is 1, we check vector indices
18 #define VOL_DEBUG 0
19 #endif
20 
21 template <class T> static inline T
22 VolMax(register const T x, register const T y) {
23  return ((x) > (y)) ? (x) : (y);
24 }
25 
26 template <class T> static inline T
27 VolAbs(register const T x) {
28  return ((x) > 0) ? (x) : -(x);
29 }
30 
31 //############################################################################
32 
33 #if defined(VOL_DEBUG) && (VOL_DEBUG != 0)
34 #define VOL_TEST_INDEX(i, size) \
35 { \
36  if ((i) < 0 || (i) >= (size)) { \
37  printf("bad VOL_?vector index\n"); \
38  abort(); \
39  } \
40 }
41 #define VOL_TEST_SIZE(size) \
42 { \
43  if (s <= 0) { \
44  printf("bad VOL_?vector size\n"); \
45  abort(); \
46  } \
47 }
48 #else
49 #define VOL_TEST_INDEX(i, size)
50 #define VOL_TEST_SIZE(size)
51 #endif
52 
53 //############################################################################
54 
55 class VOL_dvector;
56 class VOL_ivector;
57 class VOL_primal;
58 class VOL_dual;
59 class VOL_swing;
60 class VOL_alpha_factor;
61 class VOL_vh;
62 class VOL_indc;
63 class VOL_user_hooks;
64 class VOL_problem;
65 
66 //############################################################################
67 
71 struct VOL_parms {
73  double lambdainit;
75  double alphainit;
77  double alphamin;
79  double alphafactor;
80 
82  double ubinit;
83 
85  double primal_abs_precision;
87  double gap_abs_precision;
89  double gap_rel_precision;
91  double granularity;
92 
95  double minimum_rel_ascent;
100  int ascent_check_invl;
101 
103  int maxsgriters;
104 
115  int printflag;
117  int printinvl;
119  int heurinvl;
120 
123  int greentestinvl;
126  int yellowtestinvl;
129  int redtestinvl;
130 
132  int alphaint;
133 
135  char* temp_dualfile;
136 };
137 
138 //############################################################################
139 
148 class VOL_dvector {
149 public:
151  double* v;
153  int sz;
154 
155 public:
157  VOL_dvector(const int s) {
158  VOL_TEST_SIZE(s);
159  v = new double[sz = s];
160  }
162  VOL_dvector() : v(0), sz(0) {}
164  VOL_dvector(const VOL_dvector& x) : v(0), sz(0) {
165  sz = x.sz;
166  if (sz > 0) {
167  v = new double[sz];
168  std::memcpy(v, x.v, sz * sizeof(double));
169  }
170  }
172  ~VOL_dvector() { delete[] v; }
173 
175  inline int size() const {return sz;}
176 
178  inline double& operator[](const int i) {
179  VOL_TEST_INDEX(i, sz);
180  return v[i];
181  }
182 
184  inline double operator[](const int i) const {
185  VOL_TEST_INDEX(i, sz);
186  return v[i];
187  }
188 
191  inline void clear() {
192  delete[] v;
193  v = 0;
194  sz = 0;
195  }
198  inline void cc(const double gamma, const VOL_dvector& w) {
199  if (sz != w.sz) {
200  printf("bad VOL_dvector sizes\n");
201  abort();
202  }
203  double * p_v = v - 1;
204  const double * p_w = w.v - 1;
205  const double * const p_e = v + sz;
206  const double one_gamma = 1.0 - gamma;
207  while ( ++p_v != p_e ){
208  *p_v = one_gamma * (*p_v) + gamma * (*++p_w);
209  }
210  }
211 
214  inline void allocate(const int s) {
215  VOL_TEST_SIZE(s);
216  delete[] v;
217  v = new double[sz = s];
218  }
219 
221  inline void swap(VOL_dvector& w) {
222  std::swap(v, w.v);
223  std::swap(sz, w.sz);
224  }
225 
227  VOL_dvector& operator=(const VOL_dvector& w);
229  VOL_dvector& operator=(const double w);
230 };
231 
232 //-----------------------------------------------------------------------------
242 class VOL_ivector {
243 public:
245  int* v;
247  int sz;
248 public:
250  VOL_ivector(const int s) {
251  VOL_TEST_SIZE(s);
252  v = new int[sz = s];
253  }
255  VOL_ivector() : v(0), sz(0) {}
257  VOL_ivector(const VOL_ivector& x) {
258  sz = x.sz;
259  if (sz > 0) {
260  v = new int[sz];
261  std::memcpy(v, x.v, sz * sizeof(int));
262  }
263  }
265  ~VOL_ivector(){
266  delete [] v;
267  }
268 
270  inline int size() const { return sz; }
272  inline int& operator[](const int i) {
273  VOL_TEST_INDEX(i, sz);
274  return v[i];
275  }
276 
278  inline int operator[](const int i) const {
279  VOL_TEST_INDEX(i, sz);
280  return v[i];
281  }
282 
285  inline void clear() {
286  delete[] v;
287  v = 0;
288  sz = 0;
289  }
290 
293  inline void allocate(const int s) {
294  VOL_TEST_SIZE(s);
295  delete[] v;
296  v = new int[sz = s];
297  }
298 
300  inline void swap(VOL_ivector& w) {
301  std::swap(v, w.v);
302  std::swap(sz, w.sz);
303  }
304 
308  VOL_ivector& operator=(const int w);
309 };
310 
311 //############################################################################
312 // A class describing a primal solution. This class is used only internally
313 class VOL_primal {
314 public:
315  // objective value of this primal solution
316  double value;
317  // the largest of the v[i]'s
318  double viol;
319  // primal solution
320  VOL_dvector x;
321  // v=b-Ax, for the relaxed constraints
322  VOL_dvector v;
323 
324  VOL_primal(const int psize, const int dsize) : x(psize), v(dsize) {}
325  VOL_primal(const VOL_primal& primal) :
326  value(primal.value), viol(primal.viol), x(primal.x), v(primal.v) {}
327  ~VOL_primal() {}
328  inline VOL_primal& operator=(const VOL_primal& p) {
329  if (this == &p)
330  return *this;
331  value = p.value;
332  viol = p.viol;
333  x = p.x;
334  v = p.v;
335  return *this;
336  }
337 
338  // convex combination. data members in this will be overwritten
339  // convex combination between two primal solutions
340  // x <-- alpha x + (1 - alpha) p.x
341  // v <-- alpha v + (1 - alpha) p.v
342  inline void cc(const double alpha, const VOL_primal& p) {
343  value = alpha * p.value + (1.0 - alpha) * value;
344  x.cc(alpha, p.x);
345  v.cc(alpha, p.v);
346  }
347  // find maximum of v[i]
348  void find_max_viol(const VOL_dvector& dual_lb,
349  const VOL_dvector& dual_ub);
350 };
351 
352 //-----------------------------------------------------------------------------
353 // A class describing a dual solution. This class is used only internally
354 class VOL_dual {
355 public:
356  // lagrangian value
357  double lcost;
358  // reduced costs * (pstar-primal)
359  double xrc;
360  // this information is only printed
361  // dual vector
362  VOL_dvector u;
363 
364  VOL_dual(const int dsize) : u(dsize) { u = 0.0;}
365  VOL_dual(const VOL_dual& dual) :
366  lcost(dual.lcost), xrc(dual.xrc), u(dual.u) {}
367  ~VOL_dual() {}
368  inline VOL_dual& operator=(const VOL_dual& p) {
369  if (this == &p)
370  return *this;
371  lcost = p.lcost;
372  xrc = p.xrc;
373  u = p.u;
374  return *this;
375  }
376  // dual step
377  void step(const double target, const double lambda,
378  const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
379  const VOL_dvector& v);
380  double ascent(const VOL_dvector& v, const VOL_dvector& last_u) const;
381  void compute_xrc(const VOL_dvector& pstarx, const VOL_dvector& primalx,
382  const VOL_dvector& rc);
383 
384 };
385 
386 
387 //############################################################################
388 /* here we check whether an iteration is green, yellow or red. Also according
389  to this information we decide whether lambda should be changed */
390 class VOL_swing {
391 private:
392  VOL_swing(const VOL_swing&);
393  VOL_swing& operator=(const VOL_swing&);
394 public:
397  int ngs, nrs, nys;
398  int rd;
399 
400  VOL_swing() {
402  ngs = nrs = nys = 0;
403  }
404  ~VOL_swing(){}
405 
406  inline void cond(const VOL_dual& dual,
407  const double lcost, const double ascent, const int iter) {
408  double eps = 1.e-3;
409 
410  if (ascent > 0.0 && lcost > dual.lcost + eps) {
411  lastswing = green;
412  lastgreeniter = iter;
413  ++ngs;
414  rd = 0;
415  } else {
416  if (ascent <= 0 && lcost > dual.lcost) {
417  lastswing = yellow;
418  lastyellowiter = iter;
419  ++nys;
420  rd = 0;
421  } else {
422  lastswing = red;
423  lastrediter = iter;
424  ++nrs;
425  rd = 1;
426  }
427  }
428  }
429 
430  inline double
431  lfactor(const VOL_parms& parm, const double lambda, const int iter) {
432  double lambdafactor = 1.0;
433  double eps = 5.e-4;
434  int cons;
435 
436  switch (lastswing) {
437  case green:
438  cons = iter - VolMax(lastyellowiter, lastrediter);
439  if (parm.printflag & 4)
440  printf(" G: Consecutive Gs = %3d\n\n", cons);
441  if (cons >= parm.greentestinvl && lambda < 2.0) {
443  lambdafactor = 2.0;
444  if (parm.printflag & 2)
445  printf("\n ---- increasing lamda to %g ----\n\n",
446  lambda * lambdafactor);
447  }
448  break;
449 
450  case yellow:
451  cons = iter - VolMax(lastgreeniter, lastrediter);
452  if (parm.printflag & 4)
453  printf(" Y: Consecutive Ys = %3d\n\n", cons);
454  if (cons >= parm.yellowtestinvl) {
456  lambdafactor = 1.1;
457  if (parm.printflag & 2)
458  printf("\n **** increasing lamda to %g *****\n\n",
459  lambda * lambdafactor);
460  }
461  break;
462 
463  case red:
464  cons = iter - VolMax(lastgreeniter, lastyellowiter);
465  if (parm.printflag & 4)
466  printf(" R: Consecutive Rs = %3d\n\n", cons);
467  if (cons >= parm.redtestinvl && lambda > eps) {
469  lambdafactor = 0.67;
470  if (parm.printflag & 2)
471  printf("\n **** decreasing lamda to %g *****\n\n",
472  lambda * lambdafactor);
473  }
474  break;
475  }
476  return lambdafactor;
477  }
478 
479  inline void
480  print() {
481  printf("**** G= %i, Y= %i, R= %i ****\n", ngs, nys, nrs);
482  ngs = nrs = nys = 0;
483  }
484 };
485 
486 //############################################################################
487 /* alpha should be decreased if after some number of iterations the objective
488  has increased less that 1% */
489 class VOL_alpha_factor {
490 private:
493 public:
494  double lastvalue;
495 
496  VOL_alpha_factor() {lastvalue = -DBL_MAX;}
497  ~VOL_alpha_factor() {}
498 
499  inline double factor(const VOL_parms& parm, const double lcost,
500  const double alpha) {
501  if (alpha < parm.alphamin)
502  return 1.0;
503  register const double ll = VolAbs(lcost);
504  const double x = ll > 10 ? (lcost-lastvalue)/ll : (lcost-lastvalue);
505  lastvalue = lcost;
506  return (x <= 0.01) ? parm.alphafactor : 1.0;
507  }
508 };
509 
510 //############################################################################
511 /* here we compute the norm of the conjugate direction -hh-, the norm of the
512  subgradient -norm-, the inner product between the subgradient and the
513  last conjugate direction -vh-, and the inner product between the new
514  conjugate direction and the subgradient */
515 class VOL_vh {
516 private:
517  VOL_vh(const VOL_vh&);
518  VOL_vh& operator=(const VOL_vh&);
519 public:
520  double hh;
521  double norm;
522  double vh;
523  double asc;
524 
525  VOL_vh(const double alpha,
526  const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
527  const VOL_dvector& v, const VOL_dvector& vstar,
528  const VOL_dvector& u);
529  ~VOL_vh(){}
530 };
531 
532 //############################################################################
533 /* here we compute different parameter to be printed. v2 is the square of
534  the norm of the subgradient. vu is the inner product between the dual
535  variables and the subgradient. vabs is the maximum absolute value of
536  the violations of pstar. asc is the inner product between the conjugate
537  direction and the subgradient */
538 class VOL_indc {
539 private:
540  VOL_indc(const VOL_indc&);
541  VOL_indc& operator=(const VOL_indc&);
542 public:
543  double v2;
544  double vu;
545  double vabs;
546  double asc;
547 
548 public:
549  VOL_indc(const VOL_dvector& dual_lb, const VOL_dvector& dual_ub,
550  const VOL_primal& primal, const VOL_primal& pstar,
551  const VOL_dual& dual);
552  ~VOL_indc() {}
553 };
554 
555 //#############################################################################
556 
563 class VOL_user_hooks {
564 public:
565  virtual ~VOL_user_hooks() {}
566 public:
567  // for all hooks: return value of -1 means that volume should quit
572  virtual int compute_rc(const VOL_dvector& u, VOL_dvector& rc) = 0;
573 
582  virtual int solve_subproblem(const VOL_dvector& dual, const VOL_dvector& rc,
583  double& lcost, VOL_dvector& x, VOL_dvector& v,
584  double& pcost) = 0;
591  virtual int heuristics(const VOL_problem& p,
592  const VOL_dvector& x, double& heur_val) = 0;
593 };
594 
595 //#############################################################################
596 
605 class VOL_problem {
606 private:
607  VOL_problem(const VOL_problem&);
609  void set_default_parm();
610  // ############ INPUT fields ########################
611 public:
615  VOL_problem();
618  VOL_problem(const char *filename);
620  ~VOL_problem();
622 
628  int solve(VOL_user_hooks& hooks, const bool use_preset_dual = false);
630 
631 private:
635  double alpha_;
637  double lambda_;
638  // This union is here for padding (so that data members would be
639  // double-aligned on x86 CPU
640  union {
642  int iter_;
643  double __pad0;
644  };
646 
647 public:
648 
652  double value;
660 
664  VOL_parms parm;
666  int psize;
668  int dsize;
676 
677 public:
681  int iter() const { return iter_; }
683  double alpha() const { return alpha_; }
685  double lambda() const { return lambda_; }
687 
688 private:
692  void read_params(const char* filename);
693 
695  int initialize(const bool use_preset_dual);
696 
698  void print_info(const int iter,
699  const VOL_primal& primal, const VOL_primal& pstar,
700  const VOL_dual& dual);
701 
704  double readjust_target(const double oldtarget, const double lcost) const;
705 
713  double power_heur(const VOL_primal& primal, const VOL_primal& pstar,
714  const VOL_dual& dual) const;
716 };
717 
718 #endif
#define VOL_TEST_INDEX(i, size)
Definition: VolVolume.hpp:49
static T VolMax(register const T x, register const T y)
Definition: VolVolume.hpp:22
static T VolAbs(register const T x)
Definition: VolVolume.hpp:27
#define VOL_TEST_SIZE(size)
Definition: VolVolume.hpp:50
double factor(const VOL_parms &parm, const double lcost, const double alpha)
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
Definition: VolVolume.hpp:362
double xrc
Definition: VolVolume.hpp:359
VOL_dual & operator=(const VOL_dual &p)
void step(const double target, const double lambda, const VOL_dvector &dual_lb, const VOL_dvector &dual_ub, const VOL_dvector &v)
double ascent(const VOL_dvector &v, const VOL_dvector &last_u) const
VOL_dual(const int dsize)
double lcost
Definition: VolVolume.hpp:357
void cc(const double gamma, const VOL_dvector &w)
void allocate(const int s)
double * v
Definition: VolVolume.hpp:151
double & operator[](const int i)
void clear()
void swap(VOL_dvector &w)
VOL_dvector & operator=(const VOL_dvector &w)
int size() const
double vabs
Definition: VolVolume.hpp:545
double asc
Definition: VolVolume.hpp:546
VOL_indc & operator=(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)
double vu
Definition: VolVolume.hpp:544
double v2
Definition: VolVolume.hpp:543
void swap(VOL_ivector &w)
void clear()
VOL_ivector & operator=(const VOL_ivector &v)
int size() const
int & operator[](const int i)
void allocate(const int s)
VOL_primal(const int psize, const int dsize)
double value
Definition: VolVolume.hpp:316
VOL_primal & operator=(const VOL_primal &p)
VOL_dvector v
Definition: VolVolume.hpp:322
void find_max_viol(const VOL_dvector &dual_lb, const VOL_dvector &dual_ub)
VOL_dvector x
Definition: VolVolume.hpp:320
double viol
Definition: VolVolume.hpp:318
void cc(const double alpha, const VOL_primal &p)
double alpha_
Definition: VolVolume.hpp:635
double value
Definition: VolVolume.hpp:652
void print_info(const int iter, const VOL_primal &primal, const VOL_primal &pstar, const VOL_dual &dual)
VOL_dvector psol
Definition: VolVolume.hpp:656
void read_params(const char *filename)
VOL_dvector dual_lb
Definition: VolVolume.hpp:671
double lambda() const
VOL_parms parm
Definition: VolVolume.hpp:664
int iter() const
double alpha() const
double readjust_target(const double oldtarget, const double lcost) const
int initialize(const bool use_preset_dual)
VOL_problem & operator=(const VOL_problem &)
VOL_dvector viol
Definition: VolVolume.hpp:658
double power_heur(const VOL_primal &primal, const VOL_primal &pstar, const VOL_dual &dual) const
double __pad0
Definition: VolVolume.hpp:643
VOL_dvector dual_ub
Definition: VolVolume.hpp:674
VOL_dvector dsol
Definition: VolVolume.hpp:654
int solve(VOL_user_hooks &hooks, const bool use_preset_dual=false)
void set_default_parm()
double lambda_
Definition: VolVolume.hpp:637
int lastrediter
Definition: VolVolume.hpp:396
double lfactor(const VOL_parms &parm, const double lambda, const int iter)
enum VOL_swing::condition lastswing
VOL_swing & operator=(const VOL_swing &)
void print()
int lastgreeniter
Definition: VolVolume.hpp:396
void cond(const VOL_dual &dual, const double lcost, const double ascent, const int iter)
int lastyellowiter
Definition: VolVolume.hpp:396
virtual int compute_rc(const VOL_dvector &u, VOL_dvector &rc)=0
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
virtual int heuristics(const VOL_problem &p, const VOL_dvector &x, double &heur_val)=0
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)
double hh
Definition: VolVolume.hpp:520
double asc
Definition: VolVolume.hpp:523
VOL_vh & operator=(const VOL_vh &)
double vh
Definition: VolVolume.hpp:522
double norm
Definition: VolVolume.hpp:521
int printflag
Definition: VolVolume.hpp:115
double lambdainit
Definition: VolVolume.hpp:73
double gap_abs_precision
Definition: VolVolume.hpp:87
char * temp_dualfile
Definition: VolVolume.hpp:135
double alphamin
Definition: VolVolume.hpp:77
int greentestinvl
Definition: VolVolume.hpp:123
int maxsgriters
Definition: VolVolume.hpp:103
double minimum_rel_ascent
Definition: VolVolume.hpp:95
int heurinvl
Definition: VolVolume.hpp:119
double granularity
Definition: VolVolume.hpp:91
int alphaint
Definition: VolVolume.hpp:132
int redtestinvl
Definition: VolVolume.hpp:129
int printinvl
Definition: VolVolume.hpp:117
double alphafactor
Definition: VolVolume.hpp:79
double gap_rel_precision
Definition: VolVolume.hpp:89
int ascent_first_check
Definition: VolVolume.hpp:97
double alphainit
Definition: VolVolume.hpp:75
int ascent_check_invl
Definition: VolVolume.hpp:100
double primal_abs_precision
Definition: VolVolume.hpp:85
int yellowtestinvl
Definition: VolVolume.hpp:126
double ubinit
Definition: VolVolume.hpp:82