/home/coin/SVN-release/OS-2.4.0/Couenne/src/expression/operators/exprMul.cpp

Go to the documentation of this file.
00001 /* $Id: exprMul.cpp 615 2011-06-08 20:36:24Z pbelotti $
00002  *
00003  * Name:    exprMul.cpp
00004  * Author:  Pietro Belotti
00005  * Purpose: definition of multiplications
00006  *
00007  * (C) Carnegie-Mellon University, 2006. 
00008  * This file is licensed under the Eclipse Public License (EPL)
00009  */
00010 
00011 #include <stdlib.h>
00012 #include <assert.h>
00013 
00014 #include "CouenneExprMul.hpp"
00015 #include "CouenneExprSum.hpp"
00016 #include "CouenneExprConst.hpp"
00017 #include "CouenneExprClone.hpp"
00018 #include "CouennePrecisions.hpp"
00019 
00020 #include "CouenneConfig.h"
00021 #include "CoinHelperFunctions.hpp"
00022 #include "CoinFinite.hpp"
00023 
00024 using namespace Couenne;
00025 
00027 exprMul::exprMul  (expression **al, int n): 
00028   exprOp (al, n) { //< non-leaf expression, with argument list
00029 
00030   // commutative operator, sort elements
00031   qsort (arglist_, nargs_, sizeof (expression*), compareExpr);
00032 }
00033 
00035 exprMul::exprMul (expression *arg0, expression *arg1):
00036   exprOp (arg0, arg1) {
00037 
00038   if (compareExpr (arglist_, arglist_ + 1) > 0) {
00039 
00040     register expression
00041            *swap = arglist_ [0];
00042     arglist_ [0] = arglist_ [1];
00043     arglist_ [1] = swap;
00044   }
00045 }
00046 
00047 
00049 
00050 expression *exprMul::simplify () {
00051 
00052   exprOp:: simplify ();
00053 
00054   if (nargs_ == 1) {
00055 
00056     expression *ret = arglist_ [0];
00057     arglist_ [0] = NULL;
00058     return ret;
00059   }
00060 
00061   CouNumber prod = 1;
00062 
00063   bool found_one = false;
00064 
00065   for (register int i=0; i<nargs_; i++) {
00066 
00067     // check for null operands in multiplications
00068 
00069     if (arglist_ [i] -> Type () == CONST) {
00070 
00071       found_one = true;
00072 
00073       CouNumber c = arglist_ [i] -> Value ();
00074       prod *= c;
00075 
00076       if (fabs (c) < COUENNE_EPS_SIMPL) {
00077 
00078         for (int j=0; j<nargs_; j++)
00079           if (arglist_ [j]) {
00080             delete arglist_ [j];
00081             arglist_ [j] = NULL;
00082           }
00083 
00084         return new exprConst (0.);
00085       }
00086 
00087       // check for nonzero constants in multiplications
00088 
00089       delete arglist_ [i];
00090       arglist_ [i] = NULL;
00091     }
00092   }
00093   /*
00094   if (found_one && shrink_arglist (prod, 1))
00095     return new exprConst (arglist_ [0] -> Value ());
00096   */
00097   if (found_one && shrink_arglist (prod, 1)) {
00098     expression *ret = arglist_ [0];
00099     arglist_ [0] = NULL;
00100     return ret;
00101   }
00102   else return NULL;
00103 }
00104 
00105 
00106 // differentiate product of expressions
00107 
00108 expression *exprMul::differentiate (int index) {
00109 
00110   expression **als = new expression * [nargs_];
00111   int nonconst = 0;
00112 
00113   for (int i = 0; i < nargs_; i++) 
00114 
00115     if (arglist_ [i] -> dependsOn (index)) {
00116 
00117       expression **alm = new expression * [nargs_];
00118 
00119       alm [i] = arglist_ [i] -> differentiate (index);
00120 
00121       for (int j = 0; j < nargs_; j++) 
00122         if (i!=j)
00123           alm [j] = arglist_ [j] -> clone (); //new exprClone (arglist_ [j]);
00124 
00125       als [nonconst++] = new exprMul (alm, nargs_);
00126     }
00127 
00128   if (nonconst) 
00129     return new exprSum (als, nonconst);
00130   else {
00131     delete [] als;
00132     return new exprConst (0.);
00133   }
00134 }
00135 
00136 
00137 // get a measure of "how linear" the expression is:
00138 //
00139 // ZERO      = 0: constant 0
00140 // CONSTANT  = 1: a constant
00141 // LINEAR    = 2: linear
00142 // QUADRATIC = 3: quadratic
00143 // NONLINER  = 4: nonlinear non-quadratic
00144 
00145 int exprMul::Linearity () {
00146 
00147   int lin0 = arglist_ [0] -> Linearity ();
00148 
00149   if (lin0 >= NONLINEAR) return NONLINEAR;
00150   if (lin0 == ZERO)      return ZERO;
00151 
00152   for (register int i=1; i<nargs_; i++) {
00153 
00154     int lin = arglist_ [i] -> Linearity ();
00155 
00156     switch (lin) {
00157     case NONLINEAR: return NONLINEAR;
00158     case ZERO:      return ZERO;
00159     case LINEAR:    lin0++; break;
00160     case QUADRATIC: lin0 += 2; break;
00161     default: break;
00162     }
00163 
00164     if (lin0 >= NONLINEAR) return NONLINEAR;
00165   }
00166   return lin0;
00167 }
00168 
00169 
00171 void exprMul::closestFeasible (expression *varind,
00172                                expression *vardep,
00173                                CouNumber &left,
00174                                CouNumber &right) const {
00175 
00176   expression *varoth = arglist_ [0]; // suppose $w = cy$;
00177 
00178   if (varoth -> Index () == varind -> Index ())
00179     varoth = arglist_ [1]; // actually no, it's $w = x*c$
00180 
00181   assert (varoth -> Index () >= 0);
00182 
00183   CouNumber
00184     x = (*varind) (),
00185     y = (*vardep) (),
00186     c = (*varoth) ();
00187 
00188   if (c < 0.)
00189     if (y < c*x) {assert (y/c >= right); right = y/c;}
00190     else         {assert (y/c <= left);  left  = y/c;}
00191   else if (c > 0.)
00192     if (y < c*x) {assert (y/c <= left);  left  = y/c;}
00193     else         {assert (y/c >= right); right = y/c;}
00194   else left = - (right = COIN_DBL_MAX);
00195 }
00196 
00197 
00199 CouNumber exprMul::gradientNorm (const double *x) {
00200 
00201   int 
00202     ind0 = arglist_ [0] -> Index (),
00203     ind1 = arglist_ [1] -> Index ();
00204 
00205   CouNumber 
00206     x0 = (ind0 < 0) ? arglist_ [0] -> Value () : x [ind0],
00207     x1 = (ind1 < 0) ? arglist_ [1] -> Value () : x [ind1];
00208 
00209   if (ind0 < 0)
00210     if (ind1 < 0) return 0.;                   // c*d
00211     else          return fabs (x0);            // c*y
00212   else 
00213     if (ind1 < 0) return fabs (x1);            // x*d
00214     else          return sqrt (x0*x0 + x1*x1); // x*y
00215 }
00216 

Generated on Thu Sep 22 03:05:58 2011 by  doxygen 1.4.7