00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "asl.h"
00012
00013 #include "nlp.h"
00014 #include "getstub.h"
00015 #include "opcode.hd"
00016
00017 #include "CoinHelperFunctions.hpp"
00018 #include "CoinTime.hpp"
00019
00020 #include "CouenneProblem.hpp"
00021
00022 #if defined HAVE_CSTDINT
00023 #include "cstdint"
00024 #elif defined HAVE_STDINT_H
00025 #include "stdint.h"
00026 #endif
00027
00028 #include "CouenneTypes.hpp"
00029 #include "CouenneExprSum.hpp"
00030 #include "CouenneExprOpp.hpp"
00031 #include "CouenneExprMul.hpp"
00032 #include "CouenneExprSub.hpp"
00033 #include "CouenneExprClone.hpp"
00034 #include "CouenneExprGroup.hpp"
00035
00036 #define OBJ_DE ((const ASL_fg *) asl) -> I.obj_de_
00037 #define VAR_E ((const ASL_fg *) asl) -> I.var_e_
00038 #define CON_DE ((const ASL_fg *) asl) -> I.con_de_
00039 #define OBJ_sense ((const ASL_fg *) asl) -> i.objtype_
00040
00041 #define THRESHOLD_OUTPUT_READNL 10000
00042
00043
00044
00045 using namespace Couenne;
00046
00047
00048 inline bool is_expr_zero (expr* e)
00049 {return ((e==NULL) || ((((Intcast (e->op)) == OPNUM) &&
00050 (fabs (((expr_n *)e) -> v) < COUENNE_EPS)
00051
00052
00053 )));}
00054
00055 void createCommonExpr (CouenneProblem *p, const ASL *asl, int i, int which);
00056
00057
00058 int CouenneProblem::readnl (const ASL *asl) {
00059
00060 problemName_ = filename;
00061
00062
00063 ndefined_ = como + comc + comb + como1 + comc1;
00064
00065
00066
00067
00068 if (nlvb >= 0) {
00069 for (int i = 0; i < nlvb - nlvbi; i++) addVariable (false, &domain_);
00070 for (int i = 0; i < nlvbi; i++) addVariable (true, &domain_);
00071 }
00072
00073
00074 if (nlvo > nlvc) {
00075 for (int i = 0; i < nlvc - (nlvb + nlvci); i++) addVariable (false, &domain_);
00076 for (int i = 0; i < nlvci; i++) addVariable (true, &domain_);
00077 for (int i = 0; i < nlvo - (nlvc + nlvoi); i++) addVariable (false, &domain_);
00078 for (int i = 0; i < nlvoi; i++) addVariable (true, &domain_);
00079 } else {
00080 for (int i = 0; i < nlvo - (nlvb + nlvoi); i++) addVariable (false, &domain_);
00081 for (int i = 0; i < nlvoi; i++) addVariable (true, &domain_);
00082 for (int i = 0; i < nlvc - (nlvo + nlvci); i++) addVariable (false, &domain_);
00083 for (int i = 0; i < nlvci; i++) addVariable (true, &domain_);
00084 }
00085
00086 for (int i = 0; i < nwv; i++) addVariable(false, &domain_);
00087 for (int i = n_var - (CoinMax (nlvc,nlvo) +niv+nbv+nwv); i--;) addVariable(false, &domain_);
00088 for (int i = 0; i < nbv; i++) addVariable(true, &domain_);
00089 for (int i = 0; i < niv; i++) addVariable(true, &domain_);
00090
00091
00092 for (int i = ndefined_; i--;) addVariable(false, &domain_);
00093
00094 double now = CoinCpuTime ();
00095
00096 if (nVars () > THRESHOLD_OUTPUT_READNL) {
00097 jnlst_ -> Printf (Ipopt::J_ERROR, J_COUENNE, "Reading problem: ");
00098 fflush(stdout);
00099 }
00100
00101
00102
00103 #ifdef DEBUG
00104 printf ("ndefd_ = %d\n", ndefined_);
00105 printf ("tot var = %d\n", variables_ . size ());
00106 printf ("c_vars_ = %d\n", ((const ASL_fg *) asl) -> i.c_vars_ );
00107 printf ("comb_ = %d\n", ((const ASL_fg *) asl) -> i.comb_ );
00108 printf ("combc_ = %d\n", ((const ASL_fg *) asl) -> i.combc_ );
00109 printf ("comc1_ = %d\n", ((const ASL_fg *) asl) -> i.comc1_ );
00110 printf ("comc_ = %d\n", ((const ASL_fg *) asl) -> i.comc_ );
00111 printf ("como1_ = %d\n", ((const ASL_fg *) asl) -> i.como1_ );
00112 printf ("como_ = %d\n", ((const ASL_fg *) asl) -> i.como_ );
00113 #endif
00114
00115
00116
00117
00118 for (int i = 0; i < como + comc + comb; i++)
00119 createCommonExpr (this, asl, i, 0);
00120
00121 for (int i = 0; i < como1 + comc1; i++)
00122 createCommonExpr (this, asl, i, 1);
00123
00124
00125
00126
00127
00128 if (n_obj == 0) {
00129
00130
00131
00132 jnlst_ -> Printf (Ipopt::J_ERROR, J_COUENNE, "Couenne: warning, no objective function found\nAdded fictitious function f(x)=0\n");
00133 addObjective (new exprConst (0.), "min");
00134 }
00135
00136 for (int i = 0; i < n_obj; i++) {
00137
00139 int nterms = 0;
00140
00141
00142
00143 for (ograd *objgrad = Ograd [i];
00144 objgrad;
00145 objgrad = objgrad -> next)
00146 if (fabs (objgrad -> coef) > COUENNE_EPS)
00147 nterms++;
00148
00149 expression
00150 *body,
00151 *nl = Simplified (nl2e (OBJ_DE [i] . e, asl));
00152
00153
00154
00155
00156
00157
00158
00159 if (nterms) {
00160
00161 int *indexL = new int [nterms+1];
00162 CouNumber *coeff = new CouNumber [nterms];
00163
00164 for (ograd *objgrad = Ograd [i]; objgrad; objgrad = objgrad -> next)
00165 if (fabs (objgrad -> coef) > COUENNE_EPS) {
00166
00167 *indexL++ = objgrad -> varno;
00168 *coeff++ = objgrad -> coef;
00169 }
00170
00171 *indexL = -1;
00172
00173 indexL -= nterms;
00174 coeff -= nterms;
00175
00176 std::vector <std::pair <exprVar *, CouNumber> > lcoeff;
00177 indcoe2vector (indexL, coeff, lcoeff);
00178
00179 if (nl -> code () == COU_EXPRSUM) {
00180 body = exprGroup::genExprGroup (0., lcoeff, nl -> ArgList (), nl -> nArgs ());
00181
00182 nl -> ArgList (NULL);
00183 delete nl;
00184 }
00185 else {
00186
00187 expression **nll = new expression * [1];
00188
00189 *nll = nl;
00190
00191
00192 body = exprGroup::genExprGroup (0., lcoeff, nll, 1);
00193
00194 }
00195
00196 delete [] indexL;
00197 delete [] coeff;
00198
00199 } else
00200
00201 body = nl;
00202
00203
00204
00205
00207
00208 expression *subst = Simplified (body);
00209
00210
00211
00212
00213
00214
00215
00216 addObjective (body, (OBJ_sense [i] == 0) ? "min" : "max");
00217 }
00218
00219
00220
00221 int *nterms = new int [n_con];
00222
00223
00224
00225
00226
00227 for (int i = n_con; i--;)
00228 *nterms++ = 0;
00229 nterms -= n_con;
00230
00231 cgrad *congrad;
00232
00233
00234 if (A_colstarts && A_vals)
00235 for (register int j = A_colstarts [n_var]; j--;) {
00236
00237 real coeff = A_vals [j];
00238
00239 if (fabs (coeff) > COUENNE_EPS)
00240 nterms [A_rownos [j]] ++;
00241 }
00242 else {
00243 for (register int i = 0; i < n_con; i++)
00244 for (congrad = Cgrad [i];
00245 congrad;
00246 congrad = congrad -> next)
00247 if (fabs (congrad -> coef) > COUENNE_EPS)
00248 nterms [i] ++;
00249 }
00250
00251
00252
00253 CouNumber **coeff = new CouNumber * [n_con];
00254 int **indexL = new int * [n_con];
00255
00256 for (register int i = n_con; i--;)
00257 *indexL++ = NULL;
00258
00259 indexL -= n_con;
00260
00261
00262
00263
00264 if (A_colstarts && A_vals)
00265 for (int j = 0; j < n_var; j++)
00266 for (register int i = A_colstarts [j], k = A_colstarts [j+1] - i; k--; i++) {
00267
00268 int rowno = A_rownos [i],
00269 nt = nterms [rowno] --;
00270
00271 CouNumber **cline = coeff + rowno;
00272 int **iline = indexL + rowno;
00273
00274 if (*iline==NULL) {
00275 *cline = new CouNumber [nt];
00276 *iline = new int [nt+1];
00277 (*iline) [nt] = -1;
00278 }
00279
00280 (*cline) [--nt] = A_vals [i];
00281 (*iline) [nt] = j;
00282
00283 }
00284 else {
00285 for (int i=0; i < n_con; i++) {
00286
00287 int nt = nterms [i];
00288
00289 CouNumber **cline = coeff + i;
00290 int **iline = indexL + i;
00291
00292 *cline = new CouNumber [nt];
00293 *iline = new int [nt+1];
00294 (*iline) [nt] = -1;
00295
00296 for (congrad = Cgrad [i]; congrad; congrad = congrad -> next)
00297 if (fabs (congrad -> coef) > COUENNE_EPS) {
00298 (*cline) [--nt] = congrad -> coef;
00299 (*iline) [nt] = congrad -> varno;
00300 }
00301 }
00302 }
00303
00304
00305
00306 for (int i = 0; i < n_con; i++) {
00307
00308 enum con_sign sign;
00309 double lb, ub;
00310
00311 if (Urhsx) {
00312 lb = LUrhs [i];
00313 ub = Urhsx [i];
00314 } else {
00315 int j = 2*i;
00316 lb = LUrhs [j];
00317 ub = LUrhs [j+1];
00318 }
00319
00320
00321 if (lb > negInfinity)
00322 if (ub < Infinity) sign = COUENNE_RNG;
00323 else sign = COUENNE_GE;
00324 else sign = COUENNE_LE;
00325
00326
00327 if (fabs (lb - ub) < COUENNE_EPS)
00328 sign = COUENNE_EQ;
00329
00330 expression
00331 *body,
00332 **nll = new expression * [1],
00333 *nls;
00334
00335 *nll = Simplified (nl2e (CON_DE [i] . e, asl));
00336
00337
00338
00339
00340
00341
00342
00343
00344 if (indexL [i] && (*(indexL [i]) >= 0)) {
00345
00346 int code = (*nll) -> code ();
00347
00348 std::vector <std::pair <exprVar *, CouNumber> > lcoeff;
00349 indcoe2vector (indexL [i], coeff [i], lcoeff);
00350
00351
00352
00353
00354
00355 if ((code == COU_EXPRSUM) ||
00356 (code == COU_EXPRGROUP)) {
00357
00358 body = exprGroup::genExprGroup (0., lcoeff, (*nll) -> ArgList (), (*nll) -> nArgs ());
00359
00360 (*nll) -> ArgList (NULL);
00361 delete *nll;
00362 delete [] nll;
00363 }
00364 else body = exprGroup::genExprGroup (0., lcoeff, nll, 1);
00365 }
00366 else {
00367 body = *nll;
00368 delete [] nll;
00369 }
00370
00371 expression *subst = Simplified (body);
00372
00373
00374
00375
00376
00377
00378
00379
00380 switch (sign) {
00381
00382 case COUENNE_EQ: addEQConstraint (body, new exprConst (ub)); break;
00383 case COUENNE_LE: addLEConstraint (body, new exprConst (ub)); break;
00384 case COUENNE_GE: addGEConstraint (body, new exprConst (lb)); break;
00385 case COUENNE_RNG: addRNGConstraint (body, new exprConst (lb),
00386 new exprConst (ub)); break;
00387 default: printf ("Could not recognize constraint\n"); return -1;
00388 }
00389
00390 delete [] indexL [i];
00391 delete [] coeff [i];
00392 }
00393
00394 delete [] indexL;
00395 delete [] coeff;
00396
00397
00398 CouNumber
00399 *x = (CouNumber *) malloc ((n_var + ndefined_) * sizeof (CouNumber)),
00400 *lb = (CouNumber *) malloc ((n_var + ndefined_) * sizeof (CouNumber)),
00401 *ub = (CouNumber *) malloc ((n_var + ndefined_) * sizeof (CouNumber));
00402
00403 for (int i = n_var + ndefined_; i--;) {
00404 x [i] = 0.;
00405 lb [i] = -COUENNE_INFINITY;
00406 ub [i] = COUENNE_INFINITY;
00407 }
00408
00409 domain_.push (n_var + ndefined_, x, lb, ub);
00410
00411 free (x); free (lb); free (ub);
00412
00413
00414
00415 if (LUv) {
00416
00417 real *Uvx_copy = Uvx;
00418
00419 if (!Uvx_copy)
00420 for (register int i=0; i<n_var; i++) {
00421
00422 register int j = 2*i;
00423
00424 Lb (i) = LUv [j];
00425 Ub (i) = LUv [j+1];
00426 }
00427 else
00428 for (register int i=n_var; i--;) {
00429 Lb (i) = LUv [i];
00430 Ub (i) = Uvx_copy [i];
00431 }
00432
00433 } else
00434 for (register int i=n_var; i--;) {
00435 Lb (i) = - COUENNE_INFINITY;
00436 Ub (i) = COUENNE_INFINITY;
00437 }
00438
00439
00440
00441 for (register int i=n_var; i--;)
00442
00443 if (X0 && havex0 [i]) X (i) = X0 [i];
00444
00445 else {
00446
00447 CouNumber x, l = Lb (i), u = Ub (i);
00448
00449 if (l < - COUENNE_INFINITY)
00450 if (u > COUENNE_INFINITY) x = 0.;
00451 else x = u;
00452 else if (u > COUENNE_INFINITY) x = l;
00453 else x = 0.5 * (l+u);
00454
00455 X (i) = x;
00456 }
00457
00458 for (register int i=n_var; i<ndefined_; i++) {
00459
00460 X (i) = 0.;
00461 Lb (i) = -COUENNE_INFINITY;
00462 Ub (i) = COUENNE_INFINITY;
00463 }
00464
00465 delete [] nterms;
00466
00467 if (nVars () > THRESHOLD_OUTPUT_READNL) {
00468 jnlst_ -> Printf (Ipopt::J_ERROR, J_COUENNE, "%.1f seconds\n", CoinCpuTime () - now);
00469 fflush(stdout);
00470 }
00471
00472 return 0;
00473 }
00474
00475
00476
00477
00478
00479
00480 void createCommonExpr (CouenneProblem *p, const ASL *asl, int i, int which) {
00481
00482 struct cexp *common = ((const ASL_fg *) asl) -> I.cexps_ + i;
00483 struct cexp1 *common1 = ((const ASL_fg *) asl) -> I.cexps1_ + i;
00484
00485 expression
00486 *nle = Simplified (p -> nl2e (which ? common1 -> e : common -> e, asl));
00487
00488
00489
00490
00491
00492
00493
00494
00495 #ifdef DEBUG
00496 printf ("cexp1 %d [%d]: ", i, p -> Variables () . size ()); nle -> print (); printf (" ||| ");
00497 #endif
00498
00499 int nlin = which ? common1 -> nlin : common -> nlin;
00500
00501 if (nlin > 0) {
00502
00503 linpart *L = which ? common1 -> L : common -> L;
00504
00505 std::vector <std::pair <exprVar *, CouNumber> > lcoeff;
00506
00507 for (int j = 0; j < nlin; j++) {
00508
00509
00510 int indVar = ((uintptr_t) (L [j].v.rp) - (uintptr_t) VAR_E) / sizeof (expr_v);
00511 CouNumber coeff = L [j]. fac;
00512
00513 lcoeff.push_back (std::pair <exprVar *, CouNumber> (p -> Var (indVar), coeff));
00514
00515 #ifdef DEBUG
00516 Printf( " %+g x_%d (%-3d)", L [j]. fac, indVar,
00517 (expr_v *) (L [j].v.rp) - VAR_E
00518
00519 );
00520 #endif
00521 }
00522
00523 expression **al = new expression * [1];
00524 *al = nle;
00525
00526 expression *eg;
00527
00528 if (lcoeff.size () == 1 &&
00529 nle -> Type () == CONST &&
00530 nle -> Value () == 0.) {
00531
00532 CouNumber coeff = lcoeff [0].second;
00533
00534 if (coeff == 1.) eg = new exprClone (lcoeff [0].first);
00535 else if (coeff == -1.) eg = new exprOpp (new exprClone (lcoeff [0].first));
00536 else eg = new exprMul (new exprConst (coeff), new exprClone (lcoeff [0].first));
00537 } else eg = exprGroup::genExprGroup (0, lcoeff, al, 1);
00538
00539 int indVar = p -> nVars () - p -> nDefVars () + p -> commonExprs () . size ();
00540
00541 if (eg -> Index () != indVar) {
00542
00543 expression *body = new exprSub (eg, new exprClone (p -> Var (indVar)));
00544 p -> addEQConstraint (body, new exprConst (0.));
00545 }
00546
00547 p -> commonExprs () . push_back (new exprClone (eg));
00548 }
00549 else {
00550
00551 int indVar = p -> nVars () - p -> nDefVars () + p -> commonExprs () . size ();
00552
00553 if (nle -> Index () != indVar) {
00554
00555 expression *body = new exprSub (nle, new exprClone (p -> Var (indVar)));
00556 p -> addEQConstraint (body, new exprConst (0.));
00557 }
00558
00559 p -> commonExprs () . push_back (new exprClone (nle));
00560 }
00561
00562 #ifdef DEBUG
00563 printf ("\n");
00564 #endif
00565
00566 }