source: stable/0.2/Couenne/src/standardize/linStandardize.cpp @ 159

Last change on this file since 159 was 159, checked in by pbelotti, 11 years ago

created new stable branch 0.2 from trunk (rev. 157)

File size: 3.6 KB
Line 
1/* $Id: linStandardize.cpp 154 2009-06-16 18:52:53Z pbelotti $ */
2/*
3 * Name:    linStandardize.cpp
4 * Author:  Pietro Belotti
5 * Purpose: standardize sum expressions (expr{Sum,Sub,Quad,Group,Opp})
6 *
7 * (C) Carnegie-Mellon University, 2007.
8 * This file is licensed under the Common Public License (CPL)
9 */
10
11#include "CouenneProblem.hpp"
12#include "exprSum.hpp"
13#include "exprSub.hpp"
14#include "exprOpp.hpp"
15#include "exprMul.hpp"
16#include "exprPow.hpp"
17#include "exprGroup.hpp"
18#include "exprQuad.hpp"
19#include "lqelems.hpp"
20
21//#define DEBUG
22
23/// standardization of linear exprOp's
24exprAux *CouenneProblem::linStandardize (bool addAux, 
25                                         CouNumber c0, 
26                                         LinMap  &lmap,
27                                         QuadMap &qmap) {
28
29  // check if quadratic forms are dense enough ///////////////////////////////////////////
30
31  analyzeSparsity (c0, lmap, qmap);
32
33  ////////////////////////////////////////////////////////////////////////////////////////
34
35  int  nq = qmap.Map().size (),     /// data for exprQuad
36      *qi = new int [nq+1], 
37      *qj = new int [nq+1];
38
39  CouNumber *qc = new CouNumber [nq];
40
41  int 
42     nl = lmap.Map().size(),      /// data for exprGroup
43    *li = new int [nl+1];
44
45  CouNumber *lc = new CouNumber [nl];
46
47  // terminate arrays with negative index
48  qi [nq] = li [nl] = -1; 
49
50  std::map <int, CouNumber>::iterator lit = lmap.Map().begin (); 
51
52  // fill in arrays for linear part
53  for (int i=0; i<nl; i++, lit++) {
54
55    li [i] = lit -> first;
56    lc [i] = lit -> second;
57  }
58
59  std::map <std::pair <int, int>, CouNumber>::iterator qit = qmap.Map().begin (); 
60
61  // fill in arrays for quadratic part
62  for (int i=0; i < nq; i++, qit++) {
63    qi [i] = qit -> first. first;
64    qj [i] = qit -> first. second;
65    qc [i] = qit -> second;
66  }
67
68  nl = lmap.Map().size ();
69  nq = qmap.Map().size ();
70
71  // particular cases ///////////////////////////////////////////////////////////
72
73  expression *ret;
74
75  // a constant
76  if ((nq==0) && (nl==0)) 
77
78    ret = new exprConst (c0); // a constant auxiliary?
79
80  else if ((nq==0) && (fabs (c0) < COUENNE_EPS) && (nl==1)) { // a linear monomial, cx
81
82    if      (fabs (*lc - 1.) < COUENNE_EPS) ret = new exprClone (Var (*li));
83    else if (fabs (*lc + 1.) < COUENNE_EPS) ret = new exprOpp (new exprClone (Var (*li)));
84    else                                    ret = new exprMul (new exprConst (*lc), 
85                                                               new exprClone (Var (*li)));
86
87  } else if ((nl==0) && (fabs (c0) < COUENNE_EPS) && (nq==1)) { 
88
89    // a bilinear/quadratic monomial, cx^2 or cxy
90
91    expression *quad;
92
93    if (*qi == *qj) quad = new exprPow (new exprClone (Var (*qi)), new exprConst (2.));
94    else            quad = new exprMul (new exprClone (Var (*qi)), 
95                                        new exprClone (Var (*qj)));
96
97    if (fabs (*qc - 1) < COUENNE_EPS) 
98      ret    = quad;
99    else {
100      quad = addAuxiliary (quad);
101      ret  = new exprMul (new exprConst (*qc), new exprClone (quad));
102    }
103
104  } else {
105
106    // general case ///////////////////////////////////////////////////////////////
107
108    std::vector <std::pair <exprVar *, CouNumber> > lcoeff;
109    indcoe2vector (li, lc, lcoeff);
110
111    std::vector <quadElem> qcoeff;
112    indcoe2vector (qi, qj, qc, qcoeff);
113
114    if (!nq) ret = new exprGroup (c0, lcoeff); // create exprGroup
115    else     ret = new exprQuad  (c0, lcoeff, qcoeff); // create exprQuad
116  }
117
118  delete [] li;
119  delete [] lc;
120  delete [] qi;
121  delete [] qj;
122  delete [] qc;
123
124#ifdef DEBUG
125  printf ("\nlinstand ==> "); 
126  ret -> print (); printf ("\n"); 
127  //  ret -> Image () -> print (); printf ("\n");
128#endif
129
130  //if (ret -> Type () == AUX)
131  //return ret;
132
133  return (addAux ? addAuxiliary (ret) : new exprAux (ret, &domain_));
134}
Note: See TracBrowser for help on using the repository browser.