1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ 2 /* */ 3 /* This file is part of the program and library */ 4 /* SCIP --- Solving Constraint Integer Programs */ 5 /* */ 6 /* Copyright (c) 2002-2023 Zuse Institute Berlin (ZIB) */ 7 /* */ 8 /* Licensed under the Apache License, Version 2.0 (the "License"); */ 9 /* you may not use this file except in compliance with the License. */ 10 /* You may obtain a copy of the License at */ 11 /* */ 12 /* http://www.apache.org/licenses/LICENSE-2.0 */ 13 /* */ 14 /* Unless required by applicable law or agreed to in writing, software */ 15 /* distributed under the License is distributed on an "AS IS" BASIS, */ 16 /* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */ 17 /* See the License for the specific language governing permissions and */ 18 /* limitations under the License. */ 19 /* */ 20 /* You should have received a copy of the Apache-2.0 license */ 21 /* along with SCIP; see the file LICENSE. If not visit scipopt.org. */ 22 /* */ 23 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ 24 25 /**@file expr.c 26 * @ingroup OTHER_CFILES 27 * @brief functions for algebraic expressions 28 * @author Ksenia Bestuzheva 29 * @author Benjamin Mueller 30 * @author Felipe Serrano 31 * @author Stefan Vigerske 32 */ 33 34 #include <assert.h> 35 #include <ctype.h> 36 37 #include "scip/expr.h" 38 #include "scip/struct_expr.h" 39 #include "scip/pub_misc.h" 40 #include "scip/clock.h" 41 #include "scip/set.h" 42 #include "scip/pub_var.h" 43 #include "scip/sol.h" 44 #include "scip/tree.h" 45 #include "scip/struct_set.h" 46 #include "scip/struct_stat.h" 47 #include "scip/lapack_calls.h" 48 49 /*lint -e440*/ 50 /*lint -e441*/ 51 /*lint -e777*/ 52 53 /* 54 * Data structures 55 */ 56 57 /** printing to file data */ 58 struct SCIP_ExprPrintData 59 { 60 FILE* file; /**< file to print to */ 61 SCIP_EXPRITER* iterator; /**< iterator to use */ 62 SCIP_Bool closefile; /**< whether file need to be closed when finished printing */ 63 SCIP_HASHMAP* leaveexprs; /**< hashmap storing leave (no children) expressions */ 64 SCIP_EXPRPRINT_WHAT whattoprint; /**< flags that indicate what to print for each expression */ 65 }; 66 67 /* 68 * Local methods 69 */ 70 71 /** frees an expression */ 72 static 73 SCIP_RETCODE freeExpr( 74 BMS_BLKMEM* blkmem, /**< block memory */ 75 SCIP_EXPR** expr /**< pointer to free the expression */ 76 ) 77 { 78 assert(expr != NULL); 79 assert(*expr != NULL); 80 assert((*expr)->nuses == 1); 81 assert((*expr)->quaddata == NULL); 82 assert((*expr)->ownerdata == NULL); 83 84 /* free children array, if any */ 85 BMSfreeBlockMemoryArrayNull(blkmem, &(*expr)->children, (*expr)->childrensize); 86 87 BMSfreeBlockMemory(blkmem, expr); 88 assert(*expr == NULL); 89 90 return SCIP_OKAY; 91 } 92 93 /* 94 * quadratic representation of expression 95 */ 96 97 /** first time seen quadratically and 98 * seen before linearly --> --nlinterms; assign 2; ++nquadterms 99 * not seen before linearly --> assing 1; ++nquadterms 100 * 101 * seen before --> assign += 1 102 */ 103 static 104 SCIP_RETCODE quadDetectProcessExpr( 105 SCIP_EXPR* expr, /**< the expression */ 106 SCIP_HASHMAP* seenexpr, /**< hash map */ 107 int* nquadterms, /**< number of quadratic terms */ 108 int* nlinterms /**< number of linear terms */ 109 ) 110 { 111 if( SCIPhashmapExists(seenexpr, (void*)expr) ) 112 { 113 int nseen = SCIPhashmapGetImageInt(seenexpr, (void*)expr); 114 115 if( nseen < 0 ) 116 { 117 /* only seen linearly before */ 118 assert(nseen == -1); 119 120 --*nlinterms; 121 ++*nquadterms; 122 SCIP_CALL( SCIPhashmapSetImageInt(seenexpr, (void*)expr, 2) ); 123 } 124 else 125 { 126 assert(nseen > 0); 127 SCIP_CALL( SCIPhashmapSetImageInt(seenexpr, (void*)expr, nseen + 1) ); 128 } 129 } 130 else 131 { 132 ++*nquadterms; 133 SCIP_CALL( SCIPhashmapInsertInt(seenexpr, (void*)expr, 1) ); 134 } 135 136 return SCIP_OKAY; 137 } 138 139 /** returns a quadexprterm that contains the expr 140 * 141 * it either finds one that already exists or creates a new one 142 */ 143 static 144 SCIP_RETCODE quadDetectGetQuadexprterm( 145 BMS_BLKMEM* blkmem, /**< block memory */ 146 SCIP_EXPR* expr, /**< the expression */ 147 SCIP_HASHMAP* expr2idx, /**< map: expr to index in quadexpr->quadexprterms */ 148 SCIP_HASHMAP* seenexpr, /**< map: expr to number of times it was seen */ 149 SCIP_QUADEXPR* quadexpr, /**< data of quadratic representation of expression */ 150 SCIP_QUADEXPR_QUADTERM** quadexprterm /**< buffer to store quadexprterm */ 151 ) 152 { 153 assert(expr != NULL); 154 assert(expr2idx != NULL); 155 assert(quadexpr != NULL); 156 assert(quadexprterm != NULL); 157 158 if( SCIPhashmapExists(expr2idx, (void*)expr) ) 159 { 160 *quadexprterm = &quadexpr->quadexprterms[SCIPhashmapGetImageInt(expr2idx, (void*)expr)]; 161 assert((*quadexprterm)->expr == expr); 162 } 163 else 164 { 165 SCIP_CALL( SCIPhashmapInsertInt(expr2idx, expr, quadexpr->nquadexprs) ); 166 *quadexprterm = &quadexpr->quadexprterms[quadexpr->nquadexprs]; 167 ++quadexpr->nquadexprs; 168 169 (*quadexprterm)->expr = expr; 170 (*quadexprterm)->sqrcoef = 0.0; 171 (*quadexprterm)->sqrexpr = NULL; 172 (*quadexprterm)->lincoef = 0.0; 173 (*quadexprterm)->nadjbilin = 0; 174 (*quadexprterm)->adjbilinsize = SCIPhashmapGetImageInt(seenexpr, (void*)expr); 175 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &(*quadexprterm)->adjbilin, (*quadexprterm)->adjbilinsize) ); 176 } 177 178 return SCIP_OKAY; 179 } 180 181 182 /** evaluate and forward-differentiate expression 183 * 184 * also initializes derivative and bardot to 0.0 185 */ 186 static 187 SCIP_RETCODE evalAndDiff( 188 SCIP_SET* set, /**< global SCIP settings */ 189 SCIP_STAT* stat, /**< dynamic problem statistics */ 190 BMS_BLKMEM* blkmem, /**< block memory */ 191 SCIP_EXPR* expr, /**< expression to be evaluated */ 192 SCIP_SOL* sol, /**< solution to be evaluated */ 193 SCIP_Longint soltag, /**< tag that uniquely identifies the solution (with its values), or 0. */ 194 SCIP_SOL* direction /**< direction for directional derivative */ 195 ) 196 { 197 SCIP_EXPRITER* it; 198 199 assert(set != NULL); 200 assert(stat != NULL); 201 assert(blkmem != NULL); 202 assert(expr != NULL); 203 204 /* assume we'll get a domain error, so we don't have to get this expr back if we abort the iteration 205 * if there is no domain error, then we will overwrite the evalvalue in the last leaveexpr stage 206 */ 207 expr->evalvalue = SCIP_INVALID; 208 expr->evaltag = soltag; 209 expr->dot = SCIP_INVALID; 210 211 /* start a new difftag */ 212 ++stat->exprlastdifftag; 213 214 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 215 SCIP_CALL( SCIPexpriterInit(it, expr, SCIP_EXPRITER_DFS, TRUE) ); 216 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_LEAVEEXPR); 217 218 for( expr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it); expr = SCIPexpriterGetNext(it) ) 219 { 220 /* evaluate expression only if necessary */ 221 if( soltag == 0 || expr->evaltag != soltag ) 222 { 223 SCIP_CALL( SCIPexprhdlrEvalExpr(expr->exprhdlr, set, NULL, expr, &expr->evalvalue, NULL, sol) ); 224 225 expr->evaltag = soltag; 226 } 227 228 if( expr->evalvalue == SCIP_INVALID ) 229 break; 230 231 if( expr->difftag != stat->exprlastdifftag ) 232 { 233 /* compute forward diff */ 234 SCIP_CALL( SCIPexprhdlrFwDiffExpr(expr->exprhdlr, set, expr, &expr->dot, direction) ); 235 236 if( expr->dot == SCIP_INVALID ) 237 break; 238 239 expr->derivative = 0.0; 240 expr->bardot = 0.0; 241 expr->difftag = stat->exprlastdifftag; 242 } 243 } 244 245 SCIPexpriterFree(&it); 246 247 return SCIP_OKAY; 248 } 249 250 251 /* 252 * Public methods 253 */ 254 255 /* Undo the defines from pub_expr.h, which exist if NDEBUG is defined. */ 256 #ifdef NDEBUG 257 #undef SCIPexprhdlrSetCopyFreeHdlr 258 #undef SCIPexprhdlrSetCopyFreeData 259 #undef SCIPexprhdlrSetPrint 260 #undef SCIPexprhdlrSetParse 261 #undef SCIPexprhdlrSetCurvature 262 #undef SCIPexprhdlrSetMonotonicity 263 #undef SCIPexprhdlrSetIntegrality 264 #undef SCIPexprhdlrSetHash 265 #undef SCIPexprhdlrSetCompare 266 #undef SCIPexprhdlrSetDiff 267 #undef SCIPexprhdlrSetIntEval 268 #undef SCIPexprhdlrSetSimplify 269 #undef SCIPexprhdlrSetReverseProp 270 #undef SCIPexprhdlrSetEstimate 271 #undef SCIPexprhdlrSetGetSymdata 272 #undef SCIPexprhdlrGetName 273 #undef SCIPexprhdlrGetDescription 274 #undef SCIPexprhdlrGetPrecedence 275 #undef SCIPexprhdlrGetData 276 #undef SCIPexprhdlrHasPrint 277 #undef SCIPexprhdlrHasBwdiff 278 #undef SCIPexprhdlrHasFwdiff 279 #undef SCIPexprhdlrHasIntEval 280 #undef SCIPexprhdlrHasEstimate 281 #undef SCIPexprhdlrHasInitEstimates 282 #undef SCIPexprhdlrHasSimplify 283 #undef SCIPexprhdlrHasCurvature 284 #undef SCIPexprhdlrHasMonotonicity 285 #undef SCIPexprhdlrHasReverseProp 286 #undef SCIPexprhdlrGetNCreated 287 #undef SCIPexprhdlrGetNIntevalCalls 288 #undef SCIPexprhdlrGetIntevalTime 289 #undef SCIPexprhdlrGetNReversepropCalls 290 #undef SCIPexprhdlrGetReversepropTime 291 #undef SCIPexprhdlrGetNCutoffs 292 #undef SCIPexprhdlrGetNDomainReductions 293 #undef SCIPexprhdlrIncrementNDomainReductions 294 #undef SCIPexprhdlrGetNEstimateCalls 295 #undef SCIPexprhdlrGetEstimateTime 296 #undef SCIPexprhdlrGetNBranchings 297 #undef SCIPexprhdlrIncrementNBranchings 298 #undef SCIPexprhdlrGetNSimplifyCalls 299 #undef SCIPexprhdlrGetSimplifyTime 300 #undef SCIPexprhdlrGetNSimplifications 301 #endif 302 303 /** create expression handler */ 304 SCIP_RETCODE SCIPexprhdlrCreate( 305 BMS_BLKMEM* blkmem, /**< block memory */ 306 SCIP_EXPRHDLR** exprhdlr, /**< buffer where to store created expression handler */ 307 const char* name, /**< name of expression handler (must not be NULL) */ 308 const char* desc, /**< description of expression handler (can be NULL) */ 309 unsigned int precedence, /**< precedence of expression operation (used for printing) */ 310 SCIP_DECL_EXPREVAL((*eval)), /**< point evaluation callback (must not be NULL) */ 311 SCIP_EXPRHDLRDATA* data /**< data of expression handler (can be NULL) */ 312 ) 313 { 314 assert(exprhdlr != NULL); 315 assert(name != NULL); 316 317 SCIP_ALLOC( BMSallocClearBlockMemory(blkmem, exprhdlr) ); 318 319 SCIP_ALLOC( BMSduplicateMemoryArray(&(*exprhdlr)->name, name, strlen(name)+1) ); 320 if( desc != NULL ) 321 { 322 SCIP_ALLOC( BMSduplicateMemoryArray(&(*exprhdlr)->desc, desc, strlen(desc)+1) ); 323 } 324 325 (*exprhdlr)->precedence = precedence; 326 (*exprhdlr)->eval = eval; 327 (*exprhdlr)->data = data; 328 329 /* create clocks */ 330 SCIP_CALL( SCIPclockCreate(&(*exprhdlr)->estimatetime, SCIP_CLOCKTYPE_DEFAULT) ); 331 SCIP_CALL( SCIPclockCreate(&(*exprhdlr)->intevaltime, SCIP_CLOCKTYPE_DEFAULT) ); 332 SCIP_CALL( SCIPclockCreate(&(*exprhdlr)->proptime, SCIP_CLOCKTYPE_DEFAULT) ); 333 SCIP_CALL( SCIPclockCreate(&(*exprhdlr)->simplifytime, SCIP_CLOCKTYPE_DEFAULT) ); 334 335 return SCIP_OKAY; 336 } 337 338 /** frees expression handler */ 339 SCIP_RETCODE SCIPexprhdlrFree( 340 SCIP_EXPRHDLR** exprhdlr, /**< pointer to expression handler to be freed */ 341 SCIP_SET* set, /**< global SCIP settings */ 342 BMS_BLKMEM* blkmem /**< block memory */ 343 ) 344 { 345 if( (*exprhdlr)->freehdlr != NULL ) 346 { 347 SCIP_CALL( (*exprhdlr)->freehdlr(set->scip, *exprhdlr, &(*exprhdlr)->data) ); 348 } 349 350 /* free clocks */ 351 SCIPclockFree(&(*exprhdlr)->simplifytime); 352 SCIPclockFree(&(*exprhdlr)->intevaltime); 353 SCIPclockFree(&(*exprhdlr)->proptime); 354 SCIPclockFree(&(*exprhdlr)->estimatetime); 355 356 BMSfreeMemoryArrayNull(&(*exprhdlr)->desc); 357 BMSfreeMemoryArray(&(*exprhdlr)->name); 358 359 BMSfreeBlockMemory(blkmem, exprhdlr); 360 361 return SCIP_OKAY; 362 } 363 364 /**@addtogroup PublicExprHandlerMethods 365 * @{ 366 */ 367 368 /** set the expression handler callbacks to copy and free an expression handler */ 369 void SCIPexprhdlrSetCopyFreeHdlr( 370 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 371 SCIP_DECL_EXPRCOPYHDLR((*copyhdlr)), /**< handler copy callback (can be NULL) */ 372 SCIP_DECL_EXPRFREEHDLR((*freehdlr)) /**< handler free callback (can be NULL) */ 373 ) 374 { 375 assert(exprhdlr != NULL); 376 377 exprhdlr->copyhdlr = copyhdlr; 378 exprhdlr->freehdlr = freehdlr; 379 } 380 381 /** set the expression handler callbacks to copy and free expression data */ 382 void SCIPexprhdlrSetCopyFreeData( 383 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 384 SCIP_DECL_EXPRCOPYDATA((*copydata)), /**< expression data copy callback (can be NULL for expressions without data) */ 385 SCIP_DECL_EXPRFREEDATA((*freedata)) /**< expression data free callback (can be NULL if data does not need to be freed) */ 386 ) 387 { /*lint --e{715}*/ 388 assert(exprhdlr != NULL); 389 390 exprhdlr->copydata = copydata; 391 exprhdlr->freedata = freedata; 392 } 393 394 /** set the print callback of an expression handler */ 395 void SCIPexprhdlrSetPrint( 396 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 397 SCIP_DECL_EXPRPRINT((*print)) /**< print callback (can be NULL) */ 398 ) 399 { 400 assert(exprhdlr != NULL); 401 402 exprhdlr->print = print; 403 } 404 405 /** set the parse callback of an expression handler */ 406 void SCIPexprhdlrSetParse( 407 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 408 SCIP_DECL_EXPRPARSE((*parse)) /**< parse callback (can be NULL) */ 409 ) 410 { 411 assert(exprhdlr != NULL); 412 413 exprhdlr->parse = parse; 414 } 415 416 /** set the curvature detection callback of an expression handler */ 417 void SCIPexprhdlrSetCurvature( 418 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 419 SCIP_DECL_EXPRCURVATURE((*curvature)) /**< curvature detection callback (can be NULL) */ 420 ) 421 { 422 assert(exprhdlr != NULL); 423 424 exprhdlr->curvature = curvature; 425 } 426 427 /** set the monotonicity detection callback of an expression handler */ 428 void SCIPexprhdlrSetMonotonicity( 429 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 430 SCIP_DECL_EXPRMONOTONICITY((*monotonicity)) /**< monotonicity detection callback (can be NULL) */ 431 ) 432 { 433 assert(exprhdlr != NULL); 434 435 exprhdlr->monotonicity = monotonicity; 436 } 437 438 /** set the integrality detection callback of an expression handler */ 439 void SCIPexprhdlrSetIntegrality( 440 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 441 SCIP_DECL_EXPRINTEGRALITY((*integrality)) /**< integrality detection callback (can be NULL) */ 442 ) 443 { 444 assert(exprhdlr != NULL); 445 446 exprhdlr->integrality = integrality; 447 } 448 449 /** set the hash callback of an expression handler */ 450 void SCIPexprhdlrSetHash( 451 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 452 SCIP_DECL_EXPRHASH((*hash)) /**< hash callback (can be NULL) */ 453 ) 454 { 455 assert(exprhdlr != NULL); 456 457 exprhdlr->hash = hash; 458 } 459 460 /** set the compare callback of an expression handler */ 461 void SCIPexprhdlrSetCompare( 462 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 463 SCIP_DECL_EXPRCOMPARE((*compare)) /**< compare callback (can be NULL) */ 464 ) 465 { 466 assert(exprhdlr != NULL); 467 468 exprhdlr->compare = compare; 469 } 470 471 /** set differentiation callbacks of an expression handler */ 472 void SCIPexprhdlrSetDiff( 473 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 474 SCIP_DECL_EXPRBWDIFF((*bwdiff)), /**< backward derivative evaluation callback (can be NULL) */ 475 SCIP_DECL_EXPRFWDIFF((*fwdiff)), /**< forward derivative evaluation callback (can be NULL) */ 476 SCIP_DECL_EXPRBWFWDIFF((*bwfwdiff)) /**< backward-forward derivative evaluation callback (can be NULL) */ 477 ) 478 { 479 assert(exprhdlr != NULL); 480 481 exprhdlr->bwdiff = bwdiff; 482 exprhdlr->fwdiff = fwdiff; 483 exprhdlr->bwfwdiff = bwfwdiff; 484 } 485 486 /** set the interval evaluation callback of an expression handler */ 487 void SCIPexprhdlrSetIntEval( 488 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 489 SCIP_DECL_EXPRINTEVAL((*inteval)) /**< interval evaluation callback (can be NULL) */ 490 ) 491 { 492 assert(exprhdlr != NULL); 493 494 exprhdlr->inteval = inteval; 495 } 496 497 /** set the simplify callback of an expression handler */ 498 void SCIPexprhdlrSetSimplify( 499 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 500 SCIP_DECL_EXPRSIMPLIFY((*simplify)) /**< simplify callback (can be NULL) */ 501 ) 502 { 503 assert(exprhdlr != NULL); 504 505 exprhdlr->simplify = simplify; 506 } 507 508 /** set the reverse propagation callback of an expression handler */ 509 void SCIPexprhdlrSetReverseProp( 510 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 511 SCIP_DECL_EXPRREVERSEPROP((*reverseprop)) /**< reverse propagation callback (can be NULL) */ 512 ) 513 { 514 assert(exprhdlr != NULL); 515 516 exprhdlr->reverseprop = reverseprop; 517 } 518 519 /** set the symmetry information callback of an expression handler */ 520 void SCIPexprhdlrSetGetSymdata( 521 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 522 SCIP_DECL_EXPRGETSYMDATA((*getsymdata)) /**< symmetry information callback (can be NULL) */ 523 ) 524 { 525 assert(exprhdlr != NULL); 526 527 exprhdlr->getsymdata = getsymdata; 528 } 529 530 /** set the estimation callbacks of an expression handler */ 531 void SCIPexprhdlrSetEstimate( 532 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 533 SCIP_DECL_EXPRINITESTIMATES((*initestimates)), /**< initial estimators callback (can be NULL) */ 534 SCIP_DECL_EXPRESTIMATE((*estimate)) /**< estimator callback (can be NULL) */ 535 ) 536 { 537 assert(exprhdlr != NULL); 538 539 exprhdlr->initestimates = initestimates; 540 exprhdlr->estimate = estimate; 541 } 542 543 /** gives the name of an expression handler */ 544 const char* SCIPexprhdlrGetName( 545 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 546 ) 547 { 548 assert(exprhdlr != NULL); 549 550 return exprhdlr->name; 551 } 552 553 /** gives the description of an expression handler (can be NULL) */ 554 const char* SCIPexprhdlrGetDescription( 555 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 556 ) 557 { 558 assert(exprhdlr != NULL); 559 560 return exprhdlr->desc; 561 } 562 563 /** gives the precedence of an expression handler */ 564 unsigned int SCIPexprhdlrGetPrecedence( 565 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 566 ) 567 { 568 assert(exprhdlr != NULL); 569 570 return exprhdlr->precedence; 571 } 572 573 /** gives the data of an expression handler */ 574 SCIP_EXPRHDLRDATA* SCIPexprhdlrGetData( 575 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 576 ) 577 { 578 assert(exprhdlr != NULL); 579 580 return exprhdlr->data; 581 } 582 583 /** returns whether expression handler implements the print callback */ 584 SCIP_Bool SCIPexprhdlrHasPrint( 585 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 586 ) 587 { 588 assert(exprhdlr != NULL); 589 590 return exprhdlr->print != NULL; 591 } 592 593 /** returns whether expression handler implements the backward differentiation callback */ 594 SCIP_Bool SCIPexprhdlrHasBwdiff( 595 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 596 ) 597 { 598 assert(exprhdlr != NULL); 599 600 return exprhdlr->bwdiff != NULL; 601 } 602 603 /** returns whether expression handler implements the forward differentiation callback */ 604 SCIP_Bool SCIPexprhdlrHasFwdiff( 605 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 606 ) 607 { 608 assert(exprhdlr != NULL); 609 610 return exprhdlr->fwdiff != NULL; 611 } 612 613 /** returns whether expression handler implements the interval evaluation callback */ 614 SCIP_Bool SCIPexprhdlrHasIntEval( 615 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 616 ) 617 { 618 assert(exprhdlr != NULL); 619 620 return exprhdlr->inteval != NULL; 621 } 622 623 /** returns whether expression handler implements the estimator callback */ 624 SCIP_Bool SCIPexprhdlrHasEstimate( 625 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 626 ) 627 { 628 assert(exprhdlr != NULL); 629 630 return exprhdlr->estimate != NULL; 631 } 632 633 /** returns whether expression handler implements the initial estimators callback */ 634 SCIP_Bool SCIPexprhdlrHasInitEstimates( 635 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 636 ) 637 { 638 assert(exprhdlr != NULL); 639 640 return exprhdlr->initestimates != NULL; 641 } 642 643 /** returns whether expression handler implements the simplification callback */ 644 SCIP_Bool SCIPexprhdlrHasSimplify( 645 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 646 ) 647 { 648 assert(exprhdlr != NULL); 649 650 return exprhdlr->simplify != NULL; 651 } 652 653 /** returns whether expression handler implements the curvature callback */ 654 SCIP_Bool SCIPexprhdlrHasCurvature( 655 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 656 ) 657 { 658 assert(exprhdlr != NULL); 659 660 return exprhdlr->curvature != NULL; 661 } 662 663 /** returns whether expression handler implements the monotonicity callback */ 664 SCIP_Bool SCIPexprhdlrHasMonotonicity( 665 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 666 ) 667 { 668 assert(exprhdlr != NULL); 669 670 return exprhdlr->monotonicity != NULL; 671 } 672 673 /** returns whether expression handler implements the reverse propagation callback */ 674 SCIP_Bool SCIPexprhdlrHasReverseProp( 675 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 676 ) 677 { 678 assert(exprhdlr != NULL); 679 680 return exprhdlr->reverseprop != NULL; 681 } 682 683 /** returns whether expression handler implements the symmetry information callback */ 684 SCIP_Bool SCIPexprhdlrHasGetSymData( 685 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 686 ) 687 { 688 assert(exprhdlr != NULL); 689 690 return exprhdlr->getsymdata != NULL; 691 } 692 693 /** compares two expression handler w.r.t. their name */ 694 SCIP_DECL_SORTPTRCOMP(SCIPexprhdlrComp) 695 { 696 return strcmp(((SCIP_EXPRHDLR*)elem1)->name, ((SCIP_EXPRHDLR*)elem2)->name); 697 } 698 699 /** gets number of times an expression has been created with given expression handler */ 700 unsigned int SCIPexprhdlrGetNCreated( 701 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 702 ) 703 { 704 assert(exprhdlr != NULL); 705 706 return exprhdlr->ncreated; 707 } 708 709 /** gets number of times the interval evaluation callback was called */ 710 SCIP_Longint SCIPexprhdlrGetNIntevalCalls( 711 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 712 ) 713 { 714 assert(exprhdlr != NULL); 715 716 return exprhdlr->nintevalcalls; 717 } 718 719 /** gets time spend in interval evaluation callback */ 720 SCIP_Real SCIPexprhdlrGetIntevalTime( 721 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 722 ) 723 { 724 assert(exprhdlr != NULL); 725 726 return SCIPclockGetTime(exprhdlr->intevaltime); 727 } 728 729 /** gets number of times the reverse propagation callback was called */ 730 SCIP_Longint SCIPexprhdlrGetNReversepropCalls( 731 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 732 ) 733 { 734 assert(exprhdlr != NULL); 735 736 return exprhdlr->npropcalls; 737 } 738 739 /** gets time spend in reverse propagation callback */ 740 SCIP_Real SCIPexprhdlrGetReversepropTime( 741 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 742 ) 743 { 744 assert(exprhdlr != NULL); 745 746 return SCIPclockGetTime(exprhdlr->proptime); 747 } 748 749 /** gets number of times an empty interval was found in reverse propagation */ 750 SCIP_Longint SCIPexprhdlrGetNCutoffs( 751 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 752 ) 753 { 754 assert(exprhdlr != NULL); 755 756 return exprhdlr->ncutoffs; 757 } 758 759 /** gets number of times a bound reduction was found in reverse propagation (and accepted by caller) */ 760 SCIP_Longint SCIPexprhdlrGetNDomainReductions( 761 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 762 ) 763 { 764 assert(exprhdlr != NULL); 765 766 return exprhdlr->ndomreds; 767 } 768 769 /** increments the domain reductions count of an expression handler */ 770 void SCIPexprhdlrIncrementNDomainReductions( 771 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 772 int nreductions /**< number of reductions to add to counter */ 773 ) 774 { 775 assert(exprhdlr != NULL); 776 assert(nreductions >= 0); 777 778 exprhdlr->ndomreds += nreductions; 779 } 780 781 /** gets number of times the estimation callback was called */ 782 SCIP_Longint SCIPexprhdlrGetNEstimateCalls( 783 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 784 ) 785 { 786 assert(exprhdlr != NULL); 787 788 return exprhdlr->nestimatecalls; 789 } 790 791 /** gets time spend in estimation callback */ 792 SCIP_Real SCIPexprhdlrGetEstimateTime( 793 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 794 ) 795 { 796 assert(exprhdlr != NULL); 797 798 return SCIPclockGetTime(exprhdlr->estimatetime); 799 } 800 801 /** gets number of times branching candidates reported by of this expression handler were used to 802 * assemble branching candidates 803 * 804 * that is, how often did we consider branching on a child of this expression 805 */ 806 SCIP_Longint SCIPexprhdlrGetNBranchings( 807 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 808 ) 809 { 810 assert(exprhdlr != NULL); 811 812 return exprhdlr->nbranchscores; 813 } 814 815 /** increments the branching candidates count of an expression handler */ 816 void SCIPexprhdlrIncrementNBranchings( 817 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 818 ) 819 { 820 assert(exprhdlr != NULL); 821 822 ++exprhdlr->nbranchscores; 823 } 824 825 /** gets number of times the simplify callback was called */ 826 SCIP_Longint SCIPexprhdlrGetNSimplifyCalls( 827 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 828 ) 829 { 830 assert(exprhdlr != NULL); 831 832 return exprhdlr->nsimplifycalls; 833 } 834 835 /** gets time spend in simplify callback */ 836 SCIP_Real SCIPexprhdlrGetSimplifyTime( 837 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 838 ) 839 { 840 assert(exprhdlr != NULL); 841 842 return SCIPclockGetTime(exprhdlr->simplifytime); 843 } 844 845 /** gets number of times the simplify callback found a simplification */ 846 SCIP_Longint SCIPexprhdlrGetNSimplifications( 847 SCIP_EXPRHDLR* exprhdlr /**< expression handler */ 848 ) 849 { 850 assert(exprhdlr != NULL); 851 852 return exprhdlr->nsimplified; 853 } 854 855 /** @} */ 856 857 /** copies the given expression handler to a new scip */ 858 SCIP_RETCODE SCIPexprhdlrCopyInclude( 859 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 860 SCIP_SET* targetset /**< SCIP_SET of SCIP to copy to */ 861 ) 862 { 863 assert(exprhdlr != NULL); 864 assert(targetset != NULL); 865 assert(targetset->scip != NULL); 866 867 if( exprhdlr->copyhdlr != NULL ) 868 { 869 SCIPsetDebugMsg(targetset, "including expression handler <%s> in subscip %p\n", 870 SCIPexprhdlrGetName(exprhdlr), (void*)targetset->scip); 871 SCIP_CALL( exprhdlr->copyhdlr(targetset->scip, exprhdlr) ); 872 } 873 else 874 { 875 SCIPsetDebugMsg(targetset, "expression handler <%s> cannot be copied to subscip %p due " 876 "to missing copyhdlr callback\n", SCIPexprhdlrGetName(exprhdlr), (void*)targetset->scip); 877 } 878 879 return SCIP_OKAY; 880 } 881 882 /** initialization of expression handler (resets statistics) */ 883 void SCIPexprhdlrInit( 884 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 885 SCIP_SET* set /**< global SCIP settings */ 886 ) 887 { 888 assert(exprhdlr != NULL); 889 890 if( set->misc_resetstat ) 891 { 892 exprhdlr->ncreated = 0; 893 exprhdlr->nestimatecalls = 0; 894 exprhdlr->nintevalcalls = 0; 895 exprhdlr->npropcalls = 0; 896 exprhdlr->ncutoffs = 0; 897 exprhdlr->ndomreds = 0; 898 exprhdlr->nbranchscores = 0; 899 exprhdlr->nsimplifycalls = 0; 900 exprhdlr->nsimplified = 0; 901 902 SCIPclockReset(exprhdlr->estimatetime); 903 SCIPclockReset(exprhdlr->intevaltime); 904 SCIPclockReset(exprhdlr->proptime); 905 SCIPclockReset(exprhdlr->simplifytime); 906 } 907 } 908 909 /** calls the print callback of an expression handler 910 * 911 * The method prints an expression. 912 * It is called while iterating over the expression graph at different stages. 913 * 914 * @see SCIP_DECL_EXPRPRINT 915 */ 916 SCIP_RETCODE SCIPexprhdlrPrintExpr( 917 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 918 SCIP_SET* set, /**< global SCIP settings */ 919 SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */ 920 SCIP_EXPR* expr, /**< expression */ 921 SCIP_EXPRITER_STAGE stage, /**< stage of expression iteration */ 922 int currentchild, /**< index of current child if in stage visitingchild or visitedchild */ 923 unsigned int parentprecedence, /**< precedence of parent */ 924 FILE* file /**< the file to print to */ 925 ) 926 { 927 assert(exprhdlr != NULL); 928 assert(set != NULL); 929 assert(expr != NULL); 930 assert(expr->exprhdlr == exprhdlr); 931 assert(messagehdlr != NULL); 932 933 if( SCIPexprhdlrHasPrint(exprhdlr) ) 934 { 935 SCIP_CALL( exprhdlr->print(set->scip, expr, stage, currentchild, parentprecedence, file) ); 936 } 937 else 938 { 939 /* default: <hdlrname>(<child1>, <child2>, ...) */ 940 switch( stage ) 941 { 942 case SCIP_EXPRITER_ENTEREXPR : 943 { 944 SCIPmessageFPrintInfo(messagehdlr, file, "%s", SCIPexprhdlrGetName(expr->exprhdlr)); 945 if( expr->nchildren > 0 ) 946 { 947 SCIPmessageFPrintInfo(messagehdlr, file, "("); 948 } 949 break; 950 } 951 952 case SCIP_EXPRITER_VISITEDCHILD : 953 { 954 assert(currentchild >= 0); 955 assert(currentchild < expr->nchildren); 956 if( currentchild < expr->nchildren-1 ) 957 { 958 SCIPmessageFPrintInfo(messagehdlr, file, ", "); 959 } 960 else 961 { 962 SCIPmessageFPrintInfo(messagehdlr, file, ")"); 963 } 964 965 break; 966 } 967 968 case SCIP_EXPRITER_VISITINGCHILD : 969 case SCIP_EXPRITER_LEAVEEXPR : 970 default: 971 break; 972 } 973 } 974 975 return SCIP_OKAY; 976 } 977 978 /** calls the parse callback of an expression handler 979 * 980 * The method parses an expression. 981 * It should be called when parsing an expression and an operator with the expr handler name is found. 982 * 983 * @see SCIP_DECL_EXPRPARSE 984 */ 985 SCIP_RETCODE SCIPexprhdlrParseExpr( 986 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 987 SCIP_SET* set, /**< global SCIP settings */ 988 const char* string, /**< string containing expression to be parse */ 989 const char** endstring, /**< buffer to store the position of string after parsing */ 990 SCIP_EXPR** expr, /**< buffer to store the parsed expression */ 991 SCIP_Bool* success, /**< buffer to store whether the parsing was successful or not */ 992 SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call on expression copy to create ownerdata */ 993 void* ownercreatedata /**< data to pass to ownercreate */ 994 ) 995 { 996 assert(exprhdlr != NULL); 997 assert(set != NULL); 998 assert(expr != NULL); 999 1000 *expr = NULL; 1001 1002 if( exprhdlr->parse == NULL ) 1003 { 1004 /* TODO we could just look for a comma separated list of operands and try to initialize the expr with this one? 1005 * That would be sufficient for sin, cos, exp, log, abs, for example. 1006 */ 1007 SCIPdebugMessage("Expression handler <%s> has no parsing method.\n", SCIPexprhdlrGetName(exprhdlr)); 1008 *success = FALSE; 1009 return SCIP_OKAY; 1010 } 1011 1012 /* give control to exprhdlr's parser */ 1013 SCIP_CALL( exprhdlr->parse(set->scip, exprhdlr, string, endstring, expr, success, ownercreate, ownercreatedata) ); 1014 1015 assert(*success || (*expr == NULL)); 1016 1017 return SCIP_OKAY; 1018 } 1019 1020 /** calls the curvature check callback of an expression handler 1021 * 1022 * @see SCIP_DECL_EXPRCURVATURE 1023 */ 1024 SCIP_RETCODE SCIPexprhdlrCurvatureExpr( 1025 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1026 SCIP_SET* set, /**< global SCIP settings */ 1027 SCIP_EXPR* expr, /**< expression to check the curvature for */ 1028 SCIP_EXPRCURV exprcurvature, /**< desired curvature of this expression */ 1029 SCIP_Bool* success, /**< buffer to store whether the desired curvature be obtained */ 1030 SCIP_EXPRCURV* childcurv /**< array to store required curvature for each child */ 1031 ) 1032 { 1033 assert(exprhdlr != NULL); 1034 assert(set != NULL); 1035 assert(expr != NULL); 1036 assert(expr->exprhdlr == exprhdlr); 1037 assert(success != NULL); 1038 1039 *success = FALSE; 1040 1041 if( exprhdlr->curvature != NULL ) 1042 { 1043 SCIP_CALL( exprhdlr->curvature(set->scip, expr, exprcurvature, success, childcurv) ); 1044 } 1045 1046 return SCIP_OKAY; 1047 } 1048 1049 /** calls the monotonicity check callback of an expression handler 1050 * 1051 * @see SCIP_DECL_EXPRMONOTONICITY 1052 */ 1053 SCIP_RETCODE SCIPexprhdlrMonotonicityExpr( 1054 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1055 SCIP_SET* set, /**< global SCIP settings */ 1056 SCIP_EXPR* expr, /**< expression to check the monotonicity for */ 1057 int childidx, /**< index of the considered child expression */ 1058 SCIP_MONOTONE* result /**< buffer to store the monotonicity */ 1059 ) 1060 { 1061 assert(exprhdlr != NULL); 1062 assert(set != NULL); 1063 assert(expr != NULL); 1064 assert(expr->exprhdlr == exprhdlr); 1065 assert(result != NULL); 1066 1067 *result = SCIP_MONOTONE_UNKNOWN; 1068 1069 /* check whether the expression handler implements the monotonicity callback */ 1070 if( exprhdlr->monotonicity != NULL ) 1071 { 1072 SCIP_CALL( exprhdlr->monotonicity(set->scip, expr, childidx, result) ); 1073 } 1074 1075 return SCIP_OKAY; 1076 } 1077 1078 /** calls the integrality check callback of an expression handler 1079 * 1080 * @see SCIP_DECL_EXPRINTEGRALITY 1081 */ 1082 SCIP_RETCODE SCIPexprhdlrIntegralityExpr( 1083 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1084 SCIP_SET* set, /**< global SCIP settings */ 1085 SCIP_EXPR* expr, /**< expression to check integrality for */ 1086 SCIP_Bool* isintegral /**< buffer to store whether expression is integral */ 1087 ) 1088 { 1089 assert(exprhdlr != NULL); 1090 assert(set != NULL); 1091 assert(expr != NULL); 1092 assert(expr->exprhdlr == exprhdlr); 1093 assert(isintegral != NULL); 1094 1095 *isintegral = FALSE; 1096 1097 /* check whether the expression handler implements the monotonicity callback */ 1098 if( exprhdlr->integrality != NULL ) 1099 { 1100 SCIP_CALL( exprhdlr->integrality(set->scip, expr, isintegral) ); 1101 } 1102 1103 return SCIP_OKAY; 1104 } 1105 1106 /** calls the hash callback of an expression handler 1107 * 1108 * The method hashes an expression by taking the hashes of its children into account. 1109 * 1110 * @see SCIP_DECL_EXPRHASH 1111 */ 1112 SCIP_RETCODE SCIPexprhdlrHashExpr( 1113 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1114 SCIP_SET* set, /**< global SCIP settings */ 1115 SCIP_EXPR* expr, /**< expression to be hashed */ 1116 unsigned int* hashkey, /**< buffer to store the hash value */ 1117 unsigned int* childrenhashes /**< array with hash values of children */ 1118 ) 1119 { 1120 assert(exprhdlr != NULL); 1121 assert(set != NULL); 1122 assert(expr != NULL); 1123 assert(expr->exprhdlr == exprhdlr); 1124 assert(hashkey != NULL); 1125 assert(childrenhashes != NULL || expr->nchildren == 0); 1126 1127 if( expr->exprhdlr->hash != NULL ) 1128 { 1129 SCIP_CALL( expr->exprhdlr->hash(set->scip, expr, hashkey, childrenhashes) ); 1130 } 1131 else 1132 { 1133 int i; 1134 1135 /* compute initial hash from expression handler name if callback is not implemented 1136 * this can lead to more collisions and thus a larger number of expensive expression compare calls 1137 */ 1138 *hashkey = 0; 1139 for( i = 0; expr->exprhdlr->name[i] != '\0'; i++ ) 1140 *hashkey += (unsigned int) expr->exprhdlr->name[i]; /*lint !e571*/ 1141 1142 *hashkey = SCIPcalcFibHash((SCIP_Real)*hashkey); 1143 1144 /* now make use of the hashkeys of the children */ 1145 for( i = 0; i < expr->nchildren; ++i ) 1146 *hashkey ^= childrenhashes[i]; 1147 } 1148 1149 return SCIP_OKAY; 1150 } 1151 1152 /** calls the compare callback of an expression handler 1153 * 1154 * The method receives two expressions, expr1 and expr2, and returns 1155 * - -1 if expr1 < expr2, 1156 * - 0 if expr1 = expr2, 1157 * - 1 if expr1 > expr2. 1158 * 1159 * @see SCIP_DECL_EXPRCOMPARE 1160 */ 1161 int SCIPexprhdlrCompareExpr( 1162 SCIP_SET* set, /**< global SCIP settings */ 1163 SCIP_EXPR* expr1, /**< first expression in comparison */ 1164 SCIP_EXPR* expr2 /**< second expression in comparison */ 1165 ) 1166 { 1167 int i; 1168 1169 assert(expr1 != NULL); 1170 assert(expr2 != NULL); 1171 assert(expr1->exprhdlr == expr2->exprhdlr); 1172 1173 if( expr1->exprhdlr->compare != NULL ) 1174 { 1175 /* enforces OR1-OR4 */ 1176 return expr1->exprhdlr->compare(set->scip, expr1, expr2); 1177 } 1178 1179 /* enforces OR5: default comparison method of expressions of the same type: 1180 * expr1 < expr2 if and only if expr1_i = expr2_i for all i < k and expr1_k < expr2_k. 1181 * if there is no such k, use number of children to decide 1182 * if number of children is equal, both expressions are equal 1183 * @note: Warning, this method doesn't know about expression data. So if your expressions have special data, 1184 * you must implement the compare callback: SCIP_DECL_EXPRCOMPARE 1185 */ 1186 for( i = 0; i < expr1->nchildren && i < expr2->nchildren; ++i ) 1187 { 1188 int compareresult = SCIPexprCompare(set, expr1->children[i], expr2->children[i]); 1189 if( compareresult != 0 ) 1190 return compareresult; 1191 } 1192 1193 return expr1->nchildren == expr2->nchildren ? 0 : expr1->nchildren < expr2->nchildren ? -1 : 1; 1194 } 1195 1196 /** calls the evaluation callback of an expression handler 1197 * 1198 * The method evaluates an expression by taking the values of its children into account. 1199 * 1200 * Further, allows to evaluate w.r.t. given expression and children values instead of those stored in children expressions. 1201 * 1202 * @see SCIP_DECL_EXPREVAL 1203 */ 1204 SCIP_RETCODE SCIPexprhdlrEvalExpr( 1205 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1206 SCIP_SET* set, /**< global SCIP settings */ 1207 BMS_BUFMEM* bufmem, /**< buffer memory, can be NULL if childrenvals is NULL */ 1208 SCIP_EXPR* expr, /**< expression to be evaluated */ 1209 SCIP_Real* val, /**< buffer to store value of expression */ 1210 SCIP_Real* childrenvals, /**< values for children, or NULL if values stored in children should be used */ 1211 SCIP_SOL* sol /**< solution that is evaluated (can be NULL) */ 1212 ) 1213 { 1214 SCIP_Real* origvals = NULL; 1215 1216 assert(exprhdlr != NULL); 1217 assert(set != NULL); 1218 assert(expr != NULL); 1219 assert(expr->exprhdlr == exprhdlr); 1220 assert(exprhdlr->eval != NULL); 1221 assert(val != NULL); 1222 1223 /* temporarily overwrite the evalvalue in all children with values from childrenvals */ 1224 if( childrenvals != NULL && expr->nchildren > 0 ) 1225 { 1226 int c; 1227 1228 assert(bufmem != NULL); 1229 1230 SCIP_ALLOC( BMSallocBufferMemoryArray(bufmem, &origvals, expr->nchildren) ); 1231 1232 for( c = 0; c < expr->nchildren; ++c ) 1233 { 1234 origvals[c] = expr->children[c]->evalvalue; 1235 expr->children[c]->evalvalue = childrenvals[c]; 1236 } 1237 } 1238 1239 /* call expression eval callback */ 1240 SCIP_CALL( exprhdlr->eval(set->scip, expr, val, sol) ); 1241 1242 /* if there was some evaluation error (e.g., overflow) that hasn't been caught yet, then do so now */ 1243 if( !SCIPisFinite(*val) ) 1244 *val = SCIP_INVALID; 1245 1246 /* restore original evalvalues in children */ 1247 if( origvals != NULL ) 1248 { 1249 int c; 1250 for( c = 0; c < expr->nchildren; ++c ) 1251 expr->children[c]->evalvalue = origvals[c]; 1252 1253 BMSfreeBufferMemoryArray(bufmem, &origvals); 1254 } 1255 1256 return SCIP_OKAY; 1257 } 1258 1259 /** calls the backward derivative evaluation callback of an expression handler 1260 * 1261 * The method should compute the partial derivative of expr w.r.t its child at childidx. 1262 * That is, it returns 1263 * \f[ 1264 * \frac{\partial \text{expr}}{\partial \text{child}_{\text{childidx}}} 1265 * \f] 1266 * 1267 * Further, allows to differentiate w.r.t. given expression and children values instead of those stored in children expressions. 1268 * 1269 * @see SCIP_DECL_EXPRBWDIFF 1270 */ 1271 SCIP_RETCODE SCIPexprhdlrBwDiffExpr( 1272 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1273 SCIP_SET* set, /**< global SCIP settings */ 1274 BMS_BUFMEM* bufmem, /**< buffer memory, can be NULL if childrenvals is NULL */ 1275 SCIP_EXPR* expr, /**< expression to be differentiated */ 1276 int childidx, /**< index of the child */ 1277 SCIP_Real* derivative, /**< buffer to store the partial derivative w.r.t. the i-th children */ 1278 SCIP_Real* childrenvals, /**< values for children, or NULL if values stored in children should be used */ 1279 SCIP_Real exprval /**< value for expression, used only if childrenvals is not NULL */ 1280 ) 1281 { 1282 SCIP_Real* origchildrenvals; 1283 SCIP_Real origexprval = SCIP_INVALID; 1284 int c; 1285 1286 assert(exprhdlr != NULL); 1287 assert(set != NULL); 1288 assert(expr != NULL); 1289 assert(expr->exprhdlr == exprhdlr); 1290 assert(derivative != NULL); 1291 1292 if( exprhdlr->bwdiff == NULL ) 1293 { 1294 *derivative = SCIP_INVALID; 1295 return SCIP_OKAY; 1296 } 1297 1298 if( childrenvals != NULL ) 1299 { 1300 /* temporarily overwrite the evalvalue in all children and expr with values from childrenvals and exprval, resp. */ 1301 if( expr->nchildren > 0 ) 1302 { 1303 assert(bufmem != NULL); 1304 SCIP_ALLOC( BMSallocBufferMemoryArray(bufmem, &origchildrenvals, expr->nchildren) ); 1305 1306 for( c = 0; c < expr->nchildren; ++c ) 1307 { 1308 origchildrenvals[c] = expr->children[c]->evalvalue; 1309 expr->children[c]->evalvalue = childrenvals[c]; 1310 } 1311 } 1312 1313 origexprval = expr->evalvalue; 1314 expr->evalvalue = exprval; 1315 } 1316 1317 SCIP_CALL( expr->exprhdlr->bwdiff(set->scip, expr, childidx, derivative) ); 1318 1319 /* if there was some evaluation error (e.g., overflow) that hasn't been caught yet, then do so now */ 1320 if( !SCIPisFinite(*derivative) ) 1321 *derivative = SCIP_INVALID; 1322 1323 /* restore original evalvalues in children */ 1324 if( childrenvals != NULL ) 1325 { 1326 if( expr->nchildren > 0 ) 1327 { 1328 for( c = 0; c < expr->nchildren; ++c ) 1329 expr->children[c]->evalvalue = origchildrenvals[c]; /*lint !e644*/ 1330 1331 BMSfreeBufferMemoryArray(bufmem, &origchildrenvals); 1332 } 1333 1334 expr->evalvalue = origexprval; 1335 } 1336 1337 return SCIP_OKAY; 1338 } 1339 1340 /** calls the forward differentiation callback of an expression handler 1341 * 1342 * @see SCIP_DECL_EXPRFWDIFF 1343 */ 1344 SCIP_RETCODE SCIPexprhdlrFwDiffExpr( 1345 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1346 SCIP_SET* set, /**< global SCIP settings */ 1347 SCIP_EXPR* expr, /**< expression to be differentiated */ 1348 SCIP_Real* dot, /**< buffer to store derivative value */ 1349 SCIP_SOL* direction /**< direction of the derivative (useful only for var expressions) */ 1350 ) 1351 { 1352 assert(exprhdlr != NULL); 1353 assert(set != NULL); 1354 assert(expr != NULL); 1355 assert(expr->exprhdlr == exprhdlr); 1356 assert(dot != NULL); 1357 1358 if( exprhdlr->fwdiff == NULL ) 1359 { 1360 *dot = SCIP_INVALID; 1361 return SCIP_OKAY; 1362 } 1363 1364 SCIP_CALL( exprhdlr->fwdiff(set->scip, expr, dot, direction) ); 1365 1366 /* if there was some evaluation error (e.g., overflow) that hasn't been caught yet, then do so now */ 1367 if( !SCIPisFinite(*dot) ) 1368 *dot = SCIP_INVALID; 1369 1370 return SCIP_OKAY; 1371 } 1372 1373 /** calls the evaluation and forward-differentiation callback of an expression handler 1374 * 1375 * The method evaluates an expression by taking the values of its children into account. 1376 * The method differentiates an expression by taking the values and directional derivatives of its children into account. 1377 * 1378 * Further, allows to evaluate and differentiate w.r.t. given values for children instead of those stored in children expressions. 1379 * 1380 * It probably doesn't make sense to call this function for a variable-expression if sol and/or direction are not given. 1381 * 1382 * @see SCIP_DECL_EXPREVAL 1383 * @see SCIP_DECL_EXPRFWDIFF 1384 */ 1385 SCIP_RETCODE SCIPexprhdlrEvalFwDiffExpr( 1386 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1387 SCIP_SET* set, /**< global SCIP settings */ 1388 BMS_BUFMEM* bufmem, /**< buffer memory, can be NULL if childrenvals is NULL */ 1389 SCIP_EXPR* expr, /**< expression to be evaluated */ 1390 SCIP_Real* val, /**< buffer to store value of expression */ 1391 SCIP_Real* dot, /**< buffer to store derivative value */ 1392 SCIP_Real* childrenvals, /**< values for children, or NULL if values stored in children should be used */ 1393 SCIP_SOL* sol, /**< solution that is evaluated (can be NULL) */ 1394 SCIP_Real* childrendirs, /**< directional derivatives for children, or NULL if dot-values stored in children should be used */ 1395 SCIP_SOL* direction /**< direction of the derivative (useful only for var expressions, can be NULL if childrendirs is given) */ 1396 ) 1397 { 1398 SCIP_Real origval; 1399 SCIP_Real* origvals = NULL; 1400 SCIP_Real* origdots = NULL; 1401 1402 assert(exprhdlr != NULL); 1403 assert(set != NULL); 1404 assert(expr != NULL); 1405 assert(expr->exprhdlr == exprhdlr); 1406 assert(exprhdlr->eval != NULL); 1407 assert(val != NULL); 1408 assert(dot != NULL); 1409 1410 /* temporarily overwrite the evalvalue in all children with values from childrenvals */ 1411 if( childrenvals != NULL && expr->nchildren > 0 ) 1412 { 1413 int c; 1414 1415 assert(bufmem != NULL); 1416 1417 SCIP_ALLOC( BMSallocBufferMemoryArray(bufmem, &origvals, expr->nchildren) ); 1418 1419 for( c = 0; c < expr->nchildren; ++c ) 1420 { 1421 origvals[c] = expr->children[c]->evalvalue; 1422 expr->children[c]->evalvalue = childrenvals[c]; 1423 } 1424 } 1425 1426 /* temporarily overwrite the dot in all children with values from childrendirs */ 1427 if( childrendirs != NULL && expr->nchildren > 0 ) 1428 { 1429 int c; 1430 1431 assert(bufmem != NULL); 1432 1433 SCIP_ALLOC( BMSallocBufferMemoryArray(bufmem, &origdots, expr->nchildren) ); 1434 1435 for( c = 0; c < expr->nchildren; ++c ) 1436 { 1437 origdots[c] = expr->children[c]->dot; 1438 expr->children[c]->dot = childrendirs[c]; 1439 } 1440 } 1441 1442 /* remember original value */ 1443 origval = expr->evalvalue; 1444 1445 /* call expression eval callback */ 1446 SCIP_CALL( exprhdlr->eval(set->scip, expr, val, sol) ); 1447 1448 /* if there was some evaluation error (e.g., overflow) that hasn't been caught yet, then do so now */ 1449 if( !SCIPisFinite(*val) ) 1450 *val = SCIP_INVALID; 1451 1452 /* temporarily overwrite evalvalue of expr, since some exprhdlr (e.g., product) access this value in fwdiff */ 1453 expr->evalvalue = *val; 1454 1455 /* call forward-differentiation callback (if available) */ 1456 SCIP_CALL( SCIPexprhdlrFwDiffExpr(exprhdlr, set, expr, dot, direction) ); 1457 1458 /* restore original value */ 1459 expr->evalvalue = origval; 1460 1461 /* restore original dots in children */ 1462 if( origdots != NULL ) 1463 { 1464 int c; 1465 for( c = 0; c < expr->nchildren; ++c ) 1466 expr->children[c]->dot = origdots[c]; 1467 1468 BMSfreeBufferMemoryArray(bufmem, &origdots); 1469 } 1470 1471 /* restore original evalvalues in children */ 1472 if( origvals != NULL ) 1473 { 1474 int c; 1475 for( c = 0; c < expr->nchildren; ++c ) 1476 expr->children[c]->evalvalue = origvals[c]; 1477 1478 BMSfreeBufferMemoryArray(bufmem, &origvals); 1479 } 1480 1481 return SCIP_OKAY; 1482 } 1483 1484 /** calls the evaluation callback for Hessian directions (backward over forward) of an expression handler 1485 * 1486 * @see SCIP_DECL_EXPRBWFWDIFF 1487 */ 1488 SCIP_RETCODE SCIPexprhdlrBwFwDiffExpr( 1489 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1490 SCIP_SET* set, /**< global SCIP settings */ 1491 SCIP_EXPR* expr, /**< expression to be differentiated */ 1492 int childidx, /**< index of the child */ 1493 SCIP_Real* bardot, /**< buffer to store derivative value */ 1494 SCIP_SOL* direction /**< direction of the derivative (useful only for var expressions) */ 1495 ) 1496 { 1497 assert(exprhdlr != NULL); 1498 assert(set != NULL); 1499 assert(expr != NULL); 1500 assert(expr->exprhdlr == exprhdlr); 1501 assert(childidx >= 0); 1502 assert(childidx < expr->nchildren); 1503 assert(bardot != NULL); 1504 1505 if( exprhdlr->bwfwdiff == NULL ) 1506 { 1507 *bardot = SCIP_INVALID; 1508 return SCIP_OKAY; 1509 } 1510 1511 SCIP_CALL( expr->exprhdlr->bwfwdiff(set->scip, expr, childidx, bardot, direction) ); 1512 1513 /* if there was some evaluation error (e.g., overflow) that hasn't been caught yet, then do so now */ 1514 if( !SCIPisFinite(*bardot) ) 1515 *bardot = SCIP_INVALID; 1516 1517 return SCIP_OKAY; 1518 } 1519 1520 /** calls the interval evaluation callback of an expression handler 1521 * 1522 * @see SCIP_DECL_EXPRINTEVAL 1523 */ 1524 SCIP_RETCODE SCIPexprhdlrIntEvalExpr( 1525 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1526 SCIP_SET* set, /**< global SCIP settings */ 1527 SCIP_EXPR* expr, /**< expression to be evaluated */ 1528 SCIP_INTERVAL* interval, /**< buffer where to store interval */ 1529 SCIP_DECL_EXPR_INTEVALVAR((*intevalvar)), /**< callback to be called when interval-evaluating a variable */ 1530 void* intevalvardata /**< data to be passed to intevalvar callback */ 1531 ) 1532 { 1533 assert(exprhdlr != NULL); 1534 assert(set != NULL); 1535 assert(expr != NULL); 1536 assert(expr->exprhdlr == exprhdlr); 1537 assert(interval != NULL); 1538 1539 if( exprhdlr->inteval != NULL ) 1540 { 1541 SCIPclockStart(exprhdlr->intevaltime, set); 1542 SCIP_CALL( exprhdlr->inteval(set->scip, expr, interval, intevalvar, intevalvardata) ); 1543 SCIPclockStop(exprhdlr->intevaltime, set); 1544 1545 ++exprhdlr->nintevalcalls; 1546 } 1547 1548 return SCIP_OKAY; 1549 } 1550 1551 /** calls the estimator callback of an expression handler 1552 * 1553 * @see SCIP_DECL_EXPRESTIMATE 1554 */ 1555 SCIP_RETCODE SCIPexprhdlrEstimateExpr( 1556 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1557 SCIP_SET* set, /**< global SCIP settings */ 1558 SCIP_EXPR* expr, /**< expression to be estimated */ 1559 SCIP_INTERVAL* localbounds, /**< current bounds for children */ 1560 SCIP_INTERVAL* globalbounds, /**< global bounds for children */ 1561 SCIP_Real* refpoint, /**< children values for the reference point where to estimate */ 1562 SCIP_Bool overestimate, /**< whether the expression needs to be over- or underestimated */ 1563 SCIP_Real targetvalue, /**< a value that the estimator shall exceed, can be +/-infinity */ 1564 SCIP_Real* coefs, /**< array to store coefficients of estimator */ 1565 SCIP_Real* constant, /**< buffer to store constant part of estimator */ 1566 SCIP_Bool* islocal, /**< buffer to store whether estimator is valid locally only */ 1567 SCIP_Bool* success, /**< buffer to indicate whether an estimator could be computed */ 1568 SCIP_Bool* branchcand /**< array to indicate which children (not) to consider for branching */ 1569 ) 1570 { 1571 assert(exprhdlr != NULL); 1572 assert(set != NULL); 1573 assert(expr != NULL); 1574 assert(expr->exprhdlr == exprhdlr); 1575 assert(coefs != NULL); 1576 assert(islocal != NULL); 1577 assert(success != NULL); 1578 1579 *success = FALSE; 1580 1581 if( exprhdlr->estimate != NULL ) 1582 { 1583 SCIPclockStart(exprhdlr->estimatetime, set); 1584 SCIP_CALL( exprhdlr->estimate(set->scip, expr, localbounds, globalbounds, refpoint, overestimate, targetvalue, 1585 coefs, constant, islocal, success, branchcand) ); 1586 SCIPclockStop(exprhdlr->estimatetime, set); 1587 1588 /* update statistics */ 1589 ++exprhdlr->nestimatecalls; 1590 } 1591 1592 return SCIP_OKAY; 1593 } 1594 1595 /** calls the intitial estimators callback of an expression handler 1596 * 1597 * @see SCIP_DECL_EXPRINITESTIMATES 1598 */ 1599 SCIP_RETCODE SCIPexprhdlrInitEstimatesExpr( 1600 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1601 SCIP_SET* set, /**< global SCIP settings */ 1602 SCIP_EXPR* expr, /**< expression to be estimated */ 1603 SCIP_INTERVAL* bounds, /**< bounds for children */ 1604 SCIP_Bool overestimate, /**< whether the expression shall be overestimated or underestimated */ 1605 SCIP_Real* coefs[SCIP_EXPR_MAXINITESTIMATES], /**< buffer to store coefficients of computed estimators */ 1606 SCIP_Real constant[SCIP_EXPR_MAXINITESTIMATES], /**< buffer to store constant of computed estimators */ 1607 int* nreturned /**< buffer to store number of estimators that have been computed */ 1608 ) 1609 { 1610 assert(exprhdlr != NULL); 1611 assert(set != NULL); 1612 assert(expr != NULL); 1613 assert(expr->exprhdlr == exprhdlr); 1614 assert(nreturned != NULL); 1615 1616 *nreturned = 0; 1617 1618 if( exprhdlr->initestimates ) 1619 { 1620 SCIPclockStart(expr->exprhdlr->estimatetime, set); 1621 SCIP_CALL( exprhdlr->initestimates(set->scip, expr, bounds, overestimate, coefs, constant, nreturned) ); 1622 SCIPclockStop(expr->exprhdlr->estimatetime, set); 1623 1624 ++exprhdlr->nestimatecalls; 1625 } 1626 1627 return SCIP_OKAY; 1628 } 1629 1630 /** calls the simplification callback of an expression handler 1631 * 1632 * @see SCIP_DECL_EXPRSIMPLIFY 1633 */ 1634 SCIP_RETCODE SCIPexprhdlrSimplifyExpr( 1635 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1636 SCIP_SET* set, /**< global SCIP settings */ 1637 SCIP_EXPR* expr, /**< expression to simplify */ 1638 SCIP_EXPR** simplifiedexpr, /**< buffer to store the simplified expression */ 1639 SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call on expression copy to create ownerdata */ 1640 void* ownercreatedata /**< data to pass to ownercreate */ 1641 ) 1642 { 1643 assert(exprhdlr != NULL); 1644 assert(set != NULL); 1645 assert(expr != NULL); 1646 assert(expr->exprhdlr == exprhdlr); 1647 assert(simplifiedexpr != NULL); 1648 1649 if( exprhdlr->simplify != NULL ) 1650 { 1651 SCIPclockStart(expr->exprhdlr->simplifytime, set); 1652 SCIP_CALL( exprhdlr->simplify(set->scip, expr, simplifiedexpr, ownercreate, ownercreatedata) ); 1653 SCIPclockStop(expr->exprhdlr->simplifytime, set); 1654 1655 /* update statistics */ 1656 ++exprhdlr->nsimplifycalls; 1657 if( expr != *simplifiedexpr ) 1658 ++exprhdlr->nsimplified; 1659 } 1660 else 1661 { 1662 *simplifiedexpr = expr; 1663 1664 /* if an expression handler doesn't implement simplify, we assume that it is already simplified 1665 * we have to capture it, since it must simulate a "normal" simplified call in which a new expression is created 1666 */ 1667 SCIPexprCapture(expr); 1668 } 1669 1670 return SCIP_OKAY; 1671 } 1672 1673 /** calls the reverse propagation callback of an expression handler 1674 * 1675 * The method propagates given bounds over the children of an expression. 1676 * 1677 * @see SCIP_DECL_EXPRREVERSEPROP 1678 */ 1679 SCIP_RETCODE SCIPexprhdlrReversePropExpr( 1680 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1681 SCIP_SET* set, /**< global SCIP settings */ 1682 SCIP_EXPR* expr, /**< expression to propagate */ 1683 SCIP_INTERVAL bounds, /**< the bounds on the expression that should be propagated */ 1684 SCIP_INTERVAL* childrenbounds, /**< array to store computed bounds for children, initialized with current activity */ 1685 SCIP_Bool* infeasible /**< buffer to store whether a children bounds were propagated to an empty interval */ 1686 ) 1687 { 1688 assert(exprhdlr != NULL); 1689 assert(set != NULL); 1690 assert(expr != NULL); 1691 assert(expr->exprhdlr == exprhdlr); 1692 assert(childrenbounds != NULL || expr->nchildren == 0); 1693 assert(infeasible != NULL); 1694 1695 *infeasible = FALSE; 1696 1697 if( exprhdlr->reverseprop != NULL ) 1698 { 1699 SCIPclockStart(exprhdlr->proptime, set); 1700 SCIP_CALL( exprhdlr->reverseprop(set->scip, expr, bounds, childrenbounds, infeasible) ); 1701 SCIPclockStop(exprhdlr->proptime, set); 1702 1703 /* update statistics */ 1704 if( *infeasible ) 1705 ++expr->exprhdlr->ncutoffs; 1706 ++expr->exprhdlr->npropcalls; 1707 } 1708 1709 return SCIP_OKAY; 1710 } 1711 1712 /**@name Expression Methods */ 1713 /**@{ */ 1714 1715 /* from expr.h */ 1716 1717 #ifdef NDEBUG 1718 #undef SCIPexprCapture 1719 #undef SCIPexprIsVar 1720 #undef SCIPexprIsValue 1721 #undef SCIPexprIsSum 1722 #undef SCIPexprIsProduct 1723 #undef SCIPexprIsPower 1724 #endif 1725 1726 /** creates and captures an expression with given expression data and children */ 1727 SCIP_RETCODE SCIPexprCreate( 1728 SCIP_SET* set, /**< global SCIP settings */ 1729 BMS_BLKMEM* blkmem, /**< block memory */ 1730 SCIP_EXPR** expr, /**< pointer where to store expression */ 1731 SCIP_EXPRHDLR* exprhdlr, /**< expression handler */ 1732 SCIP_EXPRDATA* exprdata, /**< expression data (expression assumes ownership) */ 1733 int nchildren, /**< number of children */ 1734 SCIP_EXPR** children, /**< children (can be NULL if nchildren is 0) */ 1735 SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call to create ownerdata */ 1736 void* ownercreatedata /**< data to pass to ownercreate */ 1737 ) 1738 { 1739 int c; 1740 1741 assert(set != NULL); 1742 assert(blkmem != NULL); 1743 assert(expr != NULL); 1744 assert(exprhdlr != NULL); 1745 assert(children != NULL || nchildren == 0); 1746 assert(exprdata == NULL || exprhdlr->copydata != NULL); /* copydata must be available if there is expression data */ 1747 assert(exprdata == NULL || exprhdlr->freedata != NULL); /* freedata must be available if there is expression data */ 1748 1749 SCIP_ALLOC( BMSallocClearBlockMemory(blkmem, expr) ); 1750 1751 (*expr)->exprhdlr = exprhdlr; 1752 (*expr)->exprdata = exprdata; 1753 (*expr)->activitytag = -1; /* to be less than initial domchgcount */ 1754 (*expr)->curvature = SCIP_EXPRCURV_UNKNOWN; 1755 1756 /* initialize activity to entire interval */ 1757 SCIPintervalSetEntire(SCIP_INTERVAL_INFINITY, &(*expr)->activity); 1758 1759 if( nchildren > 0 ) 1760 { 1761 SCIP_ALLOC( BMSduplicateBlockMemoryArray(blkmem, &(*expr)->children, children, nchildren) ); 1762 (*expr)->nchildren = nchildren; 1763 (*expr)->childrensize = nchildren; 1764 1765 for( c = 0; c < nchildren; ++c ) 1766 SCIPexprCapture((*expr)->children[c]); 1767 } 1768 1769 SCIPexprCapture(*expr); 1770 1771 ++exprhdlr->ncreated; 1772 1773 /* initializes the ownerdata */ 1774 if( ownercreate != NULL ) 1775 { 1776 SCIP_CALL( ownercreate(set->scip, *expr, &(*expr)->ownerdata, &(*expr)->ownerfree, &(*expr)->ownerprint, 1777 &(*expr)->ownerevalactivity, ownercreatedata) ); 1778 } 1779 1780 return SCIP_OKAY; 1781 } 1782 1783 /** appends child to the children list of expr */ 1784 SCIP_RETCODE SCIPexprAppendChild( 1785 SCIP_SET* set, /**< global SCIP settings */ 1786 BMS_BLKMEM* blkmem, /**< block memory */ 1787 SCIP_EXPR* expr, /**< expression */ 1788 SCIP_EXPR* child /**< expression to be appended */ 1789 ) 1790 { 1791 assert(set != NULL); 1792 assert(blkmem != NULL); 1793 assert(child != NULL); 1794 assert(expr->nchildren <= expr->childrensize); 1795 1796 if( expr->nchildren == expr->childrensize ) 1797 { 1798 expr->childrensize = SCIPsetCalcMemGrowSize(set, expr->nchildren+1); 1799 SCIP_ALLOC( BMSreallocBlockMemoryArray(blkmem, &expr->children, expr->nchildren, expr->childrensize) ); 1800 } 1801 1802 expr->children[expr->nchildren] = child; 1803 ++expr->nchildren; 1804 1805 /* capture child */ 1806 SCIPexprCapture(child); 1807 1808 return SCIP_OKAY; 1809 } 1810 1811 /** overwrites/replaces a child of an expressions 1812 * 1813 * @note the old child is released and the newchild is captured, unless they are the same (=same pointer) 1814 */ 1815 SCIP_RETCODE SCIPexprReplaceChild( 1816 SCIP_SET* set, /**< global SCIP settings */ 1817 SCIP_STAT* stat, /**< dynamic problem statistics */ 1818 BMS_BLKMEM* blkmem, /**< block memory */ 1819 SCIP_EXPR* expr, /**< expression where a child is going to be replaced */ 1820 int childidx, /**< index of child being replaced */ 1821 SCIP_EXPR* newchild /**< the new child */ 1822 ) 1823 { 1824 assert(set != NULL); 1825 assert(blkmem != NULL); 1826 assert(expr != NULL); 1827 assert(newchild != NULL); 1828 assert(childidx >= 0); 1829 assert(childidx < expr->nchildren); 1830 1831 /* do nothing if child is not changing */ 1832 if( newchild == expr->children[childidx] ) 1833 return SCIP_OKAY; 1834 1835 /* capture new child (do this before releasing the old child in case there are equal */ 1836 SCIPexprCapture(newchild); 1837 1838 SCIP_CALL( SCIPexprRelease(set, stat, blkmem, &(expr->children[childidx])) ); 1839 expr->children[childidx] = newchild; 1840 1841 return SCIP_OKAY; 1842 } 1843 1844 /** remove all children of expr */ 1845 SCIP_RETCODE SCIPexprRemoveChildren( 1846 SCIP_SET* set, /**< global SCIP settings */ 1847 SCIP_STAT* stat, /**< dynamic problem statistics */ 1848 BMS_BLKMEM* blkmem, /**< block memory */ 1849 SCIP_EXPR* expr /**< expression */ 1850 ) 1851 { 1852 int c; 1853 1854 assert(set != NULL); 1855 assert(blkmem != NULL); 1856 assert(expr != NULL); 1857 1858 for( c = 0; c < expr->nchildren; ++c ) 1859 { 1860 assert(expr->children[c] != NULL); 1861 SCIP_CALL( SCIPexprRelease(set, stat, blkmem, &(expr->children[c])) ); 1862 } 1863 1864 expr->nchildren = 0; 1865 1866 return SCIP_OKAY; 1867 } 1868 1869 /** copies an expression including subexpressions 1870 * 1871 * @note If copying fails due to an expression handler not being available in the targetscip, then *targetexpr will be set to NULL. 1872 * 1873 * For all or some expressions, a mapping to an existing expression can be specified via the mapexpr callback. 1874 * The mapped expression (including its children) will not be copied in this case and its ownerdata will not be touched. 1875 * If, however, the mapexpr callback returns NULL for the targetexpr, then the expr will be copied in the usual way. 1876 */ 1877 SCIP_RETCODE SCIPexprCopy( 1878 SCIP_SET* set, /**< global SCIP settings */ 1879 SCIP_STAT* stat, /**< dynamic problem statistics */ 1880 BMS_BLKMEM* blkmem, /**< block memory */ 1881 SCIP_SET* targetset, /**< global SCIP settings data structure where target expression will live */ 1882 SCIP_STAT* targetstat, /**< dynamic problem statistics in target SCIP */ 1883 BMS_BLKMEM* targetblkmem, /**< block memory in target SCIP */ 1884 SCIP_EXPR* sourceexpr, /**< expression to be copied */ 1885 SCIP_EXPR** targetexpr, /**< buffer to store pointer to copy of source expression */ 1886 SCIP_DECL_EXPR_MAPEXPR((*mapexpr)), /**< expression mapping function, or NULL for creating new expressions */ 1887 void* mapexprdata, /**< data of expression mapping function */ 1888 SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call on expression copy to create ownerdata */ 1889 void* ownercreatedata /**< data to pass to ownercreate */ 1890 ) 1891 { 1892 SCIP_EXPRITER* it; 1893 SCIP_EXPRITER_USERDATA expriteruserdata; 1894 SCIP_EXPR* expr; 1895 SCIP* sourcescip = set->scip; /* SCIP data structure corresponding to source expression */ 1896 SCIP* targetscip = targetset->scip; /* SCIP data structure where target expression will live */ 1897 1898 assert(set != NULL); 1899 assert(stat != NULL); 1900 assert(blkmem != NULL); 1901 assert(targetset != NULL); 1902 assert(sourceexpr != NULL); 1903 assert(targetexpr != NULL); 1904 assert(sourcescip != NULL); 1905 assert(targetscip != NULL); 1906 1907 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 1908 SCIP_CALL( SCIPexpriterInit(it, sourceexpr, SCIP_EXPRITER_DFS, TRUE) ); /*TODO use FALSE, i.e., don't duplicate common subexpr? */ 1909 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_ENTEREXPR | SCIP_EXPRITER_VISITEDCHILD); 1910 1911 expr = sourceexpr; 1912 while( !SCIPexpriterIsEnd(it) ) 1913 { 1914 switch( SCIPexpriterGetStageDFS(it) ) 1915 { 1916 case SCIP_EXPRITER_ENTEREXPR : 1917 { 1918 /* create expr that will hold the copy */ 1919 SCIP_EXPRHDLR* targetexprhdlr; 1920 SCIP_EXPRDATA* targetexprdata; 1921 SCIP_EXPR* exprcopy = NULL; 1922 1923 if( mapexpr != NULL ) 1924 { 1925 SCIP_CALL( mapexpr(targetscip, &exprcopy, sourcescip, expr, ownercreate, ownercreatedata, mapexprdata) ); 1926 if( exprcopy != NULL ) 1927 { 1928 /* map callback gave us an expression to use for the copy */ 1929 /* store targetexpr */ 1930 expriteruserdata.ptrval = exprcopy; 1931 SCIPexpriterSetCurrentUserData(it, expriteruserdata); 1932 1933 /* skip subexpression (assume that exprcopy is a complete copy) and continue */ 1934 expr = SCIPexpriterSkipDFS(it); 1935 continue; 1936 } 1937 } 1938 1939 /* get the exprhdlr of the target scip */ 1940 if( targetscip != sourcescip ) 1941 { 1942 targetexprhdlr = SCIPsetFindExprhdlr(targetset, expr->exprhdlr->name); 1943 1944 if( targetexprhdlr == NULL ) 1945 { 1946 /* expression handler not in target scip (probably did not have a copy callback) -> abort */ 1947 expriteruserdata.ptrval = NULL; 1948 SCIPexpriterSetCurrentUserData(it, expriteruserdata); 1949 1950 expr = SCIPexpriterSkipDFS(it); 1951 continue; 1952 } 1953 } 1954 else 1955 { 1956 targetexprhdlr = expr->exprhdlr; 1957 } 1958 assert(targetexprhdlr != NULL); 1959 1960 /* copy expression data */ 1961 if( expr->exprdata != NULL ) 1962 { 1963 assert(expr->exprhdlr->copydata != NULL); 1964 SCIP_CALL( expr->exprhdlr->copydata(targetscip, targetexprhdlr, &targetexprdata, sourcescip, expr) ); 1965 } 1966 else 1967 { 1968 targetexprdata = NULL; 1969 } 1970 1971 /* create in targetexpr an expression of the same type as expr, but without children for now */ 1972 SCIP_CALL( SCIPexprCreate(targetset, targetblkmem, &exprcopy, targetexprhdlr, targetexprdata, 0, NULL, 1973 ownercreate, ownercreatedata) ); 1974 1975 /* store targetexpr */ 1976 expriteruserdata.ptrval = exprcopy; 1977 SCIPexpriterSetCurrentUserData(it, expriteruserdata); 1978 1979 break; 1980 } 1981 1982 case SCIP_EXPRITER_VISITEDCHILD : 1983 { 1984 /* just visited child so a copy of himself should be available; append it */ 1985 SCIP_EXPR* exprcopy; 1986 SCIP_EXPR* childcopy; 1987 1988 exprcopy = (SCIP_EXPR*)SCIPexpriterGetCurrentUserData(it).ptrval; 1989 1990 /* get copy of child */ 1991 childcopy = (SCIP_EXPR*)SCIPexpriterGetChildUserDataDFS(it).ptrval; 1992 if( childcopy == NULL ) 1993 { 1994 /* abort */ 1995 /* release exprcopy (should free also the already copied children) */ 1996 SCIP_CALL( SCIPexprRelease(targetset, targetstat, targetblkmem, (SCIP_EXPR**)&exprcopy) ); 1997 1998 expriteruserdata.ptrval = NULL; 1999 SCIPexpriterSetCurrentUserData(it, expriteruserdata); 2000 2001 expr = SCIPexpriterSkipDFS(it); 2002 continue; 2003 } 2004 2005 /* append child to exprcopy */ 2006 SCIP_CALL( SCIPexprAppendChild(targetset, targetblkmem, exprcopy, childcopy) ); 2007 2008 /* release childcopy (still captured by exprcopy) */ 2009 SCIP_CALL( SCIPexprRelease(targetset, targetstat, targetblkmem, &childcopy) ); 2010 2011 break; 2012 } 2013 2014 default: 2015 /* we should never be called in this stage */ 2016 SCIPABORT(); 2017 break; 2018 } 2019 2020 expr = SCIPexpriterGetNext(it); 2021 } 2022 2023 /* the target expression should be stored in the userdata of the sourceexpr (can be NULL if aborted) */ 2024 *targetexpr = (SCIP_EXPR*)SCIPexpriterGetExprUserData(it, sourceexpr).ptrval; 2025 2026 SCIPexpriterFree(&it); 2027 2028 return SCIP_OKAY; 2029 } 2030 2031 /** duplicates the given expression without its children */ 2032 SCIP_RETCODE SCIPexprDuplicateShallow( 2033 SCIP_SET* set, /**< global SCIP settings */ 2034 BMS_BLKMEM* blkmem, /**< block memory */ 2035 SCIP_EXPR* expr, /**< original expression */ 2036 SCIP_EXPR** copyexpr, /**< buffer to store (shallow) duplicate of expr */ 2037 SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call on expression copy to create ownerdata */ 2038 void* ownercreatedata /**< data to pass to ownercreate */ 2039 ) 2040 { 2041 SCIP_EXPRDATA* exprdatacopy = NULL; 2042 2043 assert(set != NULL); 2044 assert(blkmem != NULL); 2045 assert(expr != NULL); 2046 assert(copyexpr != NULL); 2047 2048 /* copy expression data */ 2049 if( expr->exprdata != NULL ) 2050 { 2051 assert(expr->exprhdlr->copydata != NULL); 2052 SCIP_CALL( expr->exprhdlr->copydata(set->scip, expr->exprhdlr, &exprdatacopy, set->scip, expr) ); 2053 } 2054 2055 /* create expression with same handler and copied data, but without children */ 2056 SCIP_CALL( SCIPexprCreate(set, blkmem, copyexpr, expr->exprhdlr, exprdatacopy, 0, NULL, ownercreate, 2057 ownercreatedata) ); 2058 2059 return SCIP_OKAY; 2060 } 2061 2062 /** captures an expression (increments usage count) */ 2063 void SCIPexprCapture( 2064 SCIP_EXPR* expr /**< expression */ 2065 ) 2066 { 2067 assert(expr != NULL); 2068 2069 ++expr->nuses; 2070 } 2071 2072 /** releases an expression (decrements usage count and possibly frees expression) */ 2073 SCIP_RETCODE SCIPexprRelease( 2074 SCIP_SET* set, /**< global SCIP settings */ 2075 SCIP_STAT* stat, /**< dynamic problem statistics */ 2076 BMS_BLKMEM* blkmem, /**< block memory */ 2077 SCIP_EXPR** rootexpr /**< pointer to expression */ 2078 ) 2079 { 2080 SCIP_EXPRITER* it; 2081 SCIP_EXPR* expr; 2082 2083 assert(rootexpr != NULL); 2084 assert(*rootexpr != NULL); 2085 assert((*rootexpr)->nuses > 0); 2086 2087 if( (*rootexpr)->nuses > 1 ) 2088 { 2089 --(*rootexpr)->nuses; 2090 *rootexpr = NULL; 2091 2092 return SCIP_OKAY; 2093 } 2094 2095 /* handle the root expr separately: free ownerdata, quaddata, and exprdata first */ 2096 2097 /* call ownerfree callback, if given 2098 * we intentially call this also if ownerdata is NULL, so owner can be notified without storing data 2099 */ 2100 if( (*rootexpr)->ownerfree != NULL ) 2101 { 2102 SCIP_CALL( (*rootexpr)->ownerfree(set->scip, *rootexpr, &(*rootexpr)->ownerdata) ); 2103 assert((*rootexpr)->ownerdata == NULL); 2104 } 2105 2106 /* free quadratic info */ 2107 SCIPexprFreeQuadratic(blkmem, *rootexpr); 2108 2109 /* free expression data */ 2110 if( (*rootexpr)->exprdata != NULL ) 2111 { 2112 assert((*rootexpr)->exprhdlr->freedata != NULL); 2113 SCIP_CALL( (*rootexpr)->exprhdlr->freedata(set->scip, *rootexpr) ); 2114 } 2115 2116 /* now release and free children, where no longer in use */ 2117 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2118 SCIP_CALL( SCIPexpriterInit(it, *rootexpr, SCIP_EXPRITER_DFS, TRUE) ); 2119 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD | SCIP_EXPRITER_VISITEDCHILD); 2120 for( expr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it) ; ) 2121 { 2122 /* expression should be used by its parent and maybe by the iterator (only the root!) 2123 * in VISITEDCHILD we assert that expression is only used by its parent 2124 */ 2125 assert(expr != NULL); 2126 assert(0 <= expr->nuses && expr->nuses <= 2); 2127 2128 switch( SCIPexpriterGetStageDFS(it) ) 2129 { 2130 case SCIP_EXPRITER_VISITINGCHILD : 2131 { 2132 /* check whether a child needs to be visited (nuses == 1) 2133 * if not, then we still have to release it 2134 */ 2135 SCIP_EXPR* child; 2136 2137 child = SCIPexpriterGetChildExprDFS(it); 2138 if( child->nuses > 1 ) 2139 { 2140 /* child is not going to be freed: just release it */ 2141 SCIP_CALL( SCIPexprRelease(set, stat, blkmem, &child) ); 2142 expr = SCIPexpriterSkipDFS(it); 2143 continue; 2144 } 2145 2146 assert(child->nuses == 1); 2147 2148 /* free child's quaddata, ownerdata, and exprdata when entering child */ 2149 if( child->ownerfree != NULL ) 2150 { 2151 SCIP_CALL( child->ownerfree(set->scip, child, &child->ownerdata) ); 2152 assert(child->ownerdata == NULL); 2153 } 2154 2155 /* free quadratic info */ 2156 SCIPexprFreeQuadratic(blkmem, child); 2157 2158 /* free expression data */ 2159 if( child->exprdata != NULL ) 2160 { 2161 assert(child->exprhdlr->freedata != NULL); 2162 SCIP_CALL( child->exprhdlr->freedata(set->scip, child) ); 2163 assert(child->exprdata == NULL); 2164 } 2165 2166 break; 2167 } 2168 2169 case SCIP_EXPRITER_VISITEDCHILD : 2170 { 2171 /* free child after visiting it */ 2172 SCIP_EXPR* child; 2173 2174 child = SCIPexpriterGetChildExprDFS(it); 2175 /* child should only be used by its parent */ 2176 assert(child->nuses == 1); 2177 2178 /* child should have no data associated */ 2179 assert(child->exprdata == NULL); 2180 2181 /* free child expression */ 2182 SCIP_CALL( freeExpr(blkmem, &child) ); 2183 expr->children[SCIPexpriterGetChildIdxDFS(it)] = NULL; 2184 2185 break; 2186 } 2187 2188 default: 2189 SCIPABORT(); /* we should never be called in this stage */ 2190 break; 2191 } 2192 2193 expr = SCIPexpriterGetNext(it); 2194 } 2195 2196 SCIPexpriterFree(&it); 2197 2198 /* handle the root expr separately: free its children and itself here */ 2199 SCIP_CALL( freeExpr(blkmem, rootexpr) ); 2200 2201 return SCIP_OKAY; 2202 } 2203 2204 /** returns whether an expression is a variable expression */ 2205 SCIP_Bool SCIPexprIsVar( 2206 SCIP_SET* set, /**< global SCIP settings */ 2207 SCIP_EXPR* expr /**< expression */ 2208 ) 2209 { 2210 assert(set != NULL); 2211 assert(expr != NULL); 2212 2213 return expr->exprhdlr == set->exprhdlrvar; 2214 } 2215 2216 /** returns whether an expression is a value expression */ 2217 SCIP_Bool SCIPexprIsValue( 2218 SCIP_SET* set, /**< global SCIP settings */ 2219 SCIP_EXPR* expr /**< expression */ 2220 ) 2221 { 2222 assert(set != NULL); 2223 assert(expr != NULL); 2224 2225 return expr->exprhdlr == set->exprhdlrval; 2226 } 2227 2228 /** returns whether an expression is a sum expression */ 2229 SCIP_Bool SCIPexprIsSum( 2230 SCIP_SET* set, /**< global SCIP settings */ 2231 SCIP_EXPR* expr /**< expression */ 2232 ) 2233 { 2234 assert(set != NULL); 2235 assert(expr != NULL); 2236 2237 return expr->exprhdlr == set->exprhdlrsum; 2238 } 2239 2240 /** returns whether an expression is a product expression */ 2241 SCIP_Bool SCIPexprIsProduct( 2242 SCIP_SET* set, /**< global SCIP settings */ 2243 SCIP_EXPR* expr /**< expression */ 2244 ) 2245 { 2246 assert(set != NULL); 2247 assert(expr != NULL); 2248 2249 return expr->exprhdlr == set->exprhdlrproduct; 2250 } 2251 2252 /** returns whether an expression is a power expression */ 2253 SCIP_Bool SCIPexprIsPower( 2254 SCIP_SET* set, /**< global SCIP settings */ 2255 SCIP_EXPR* expr /**< expression */ 2256 ) 2257 { 2258 assert(set != NULL); 2259 assert(expr != NULL); 2260 2261 return expr->exprhdlr == set->exprhdlrpow; 2262 } 2263 2264 /** print an expression as info-message */ 2265 SCIP_RETCODE SCIPexprPrint( 2266 SCIP_SET* set, /**< global SCIP settings */ 2267 SCIP_STAT* stat, /**< dynamic problem statistics */ 2268 BMS_BLKMEM* blkmem, /**< block memory */ 2269 SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */ 2270 FILE* file, /**< file to print to, or NULL for stdout */ 2271 SCIP_EXPR* expr /**< expression to be printed */ 2272 ) 2273 { 2274 SCIP_EXPRITER* it; 2275 SCIP_EXPRITER_STAGE stage; 2276 int currentchild; 2277 unsigned int parentprecedence; 2278 2279 assert(set != NULL); 2280 assert(blkmem != NULL); 2281 assert(expr != NULL); 2282 2283 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2284 SCIP_CALL( SCIPexpriterInit(it, expr, SCIP_EXPRITER_DFS, TRUE) ); 2285 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_ALLSTAGES); 2286 2287 while( !SCIPexpriterIsEnd(it) ) 2288 { 2289 assert(expr->exprhdlr != NULL); 2290 stage = SCIPexpriterGetStageDFS(it); 2291 2292 if( stage == SCIP_EXPRITER_VISITEDCHILD || stage == SCIP_EXPRITER_VISITINGCHILD ) 2293 currentchild = SCIPexpriterGetChildIdxDFS(it); 2294 else 2295 currentchild = -1; 2296 2297 if( SCIPexpriterGetParentDFS(it) != NULL ) 2298 parentprecedence = SCIPexprhdlrGetPrecedence(SCIPexprGetHdlr(SCIPexpriterGetParentDFS(it))); 2299 else 2300 parentprecedence = 0; 2301 2302 SCIP_CALL( SCIPexprhdlrPrintExpr(expr->exprhdlr, set, messagehdlr, expr, stage, currentchild, 2303 parentprecedence, file) ); 2304 2305 expr = SCIPexpriterGetNext(it); 2306 } 2307 2308 SCIPexpriterFree(&it); 2309 2310 return SCIP_OKAY; 2311 } 2312 2313 /** initializes printing of expressions in dot format to a give FILE* pointer */ 2314 SCIP_RETCODE SCIPexprPrintDotInit( 2315 SCIP_SET* set, /**< global SCIP settings */ 2316 SCIP_STAT* stat, /**< dynamic problem statistics */ 2317 BMS_BLKMEM* blkmem, /**< block memory */ 2318 SCIP_EXPRPRINTDATA** printdata, /**< buffer to store dot printing data */ 2319 FILE* file, /**< file to print to, or NULL for stdout */ 2320 SCIP_EXPRPRINT_WHAT whattoprint /**< info on what to print for each expression */ 2321 ) 2322 { 2323 assert(set != NULL); 2324 assert(stat != NULL); 2325 assert(blkmem != NULL); 2326 assert(printdata != NULL); 2327 2328 if( file == NULL ) 2329 file = stdout; 2330 2331 SCIP_ALLOC( BMSallocBlockMemory(blkmem, printdata) ); 2332 2333 (*printdata)->file = file; 2334 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &(*printdata)->iterator) ); 2335 (*printdata)->closefile = FALSE; 2336 (*printdata)->whattoprint = whattoprint; 2337 SCIP_CALL( SCIPhashmapCreate(&(*printdata)->leaveexprs, blkmem, 100) ); 2338 2339 fputs("strict digraph exprgraph {\n", file); 2340 fputs("node [fontcolor=white, style=filled, rankdir=LR]\n", file); 2341 2342 return SCIP_OKAY; 2343 } 2344 2345 /** initializes printing of expressions in dot format to a file with given filename */ 2346 SCIP_RETCODE SCIPexprPrintDotInit2( 2347 SCIP_SET* set, /**< global SCIP settings */ 2348 SCIP_STAT* stat, /**< dynamic problem statistics */ 2349 BMS_BLKMEM* blkmem, /**< block memory */ 2350 SCIP_EXPRPRINTDATA** printdata, /**< buffer to store dot printing data */ 2351 const char* filename, /**< name of file to print to */ 2352 SCIP_EXPRPRINT_WHAT whattoprint /**< info on what to print for each expression */ 2353 ) 2354 { 2355 FILE* f; 2356 2357 assert(set != NULL); 2358 assert(stat != NULL); 2359 assert(blkmem != NULL); 2360 assert(printdata != NULL); 2361 assert(filename != NULL); 2362 2363 f = fopen(filename, "w"); 2364 if( f == NULL ) 2365 { 2366 SCIPerrorMessage("could not open file <%s> for writing\n", filename); /* error code would be in errno */ 2367 return SCIP_FILECREATEERROR; 2368 } 2369 2370 SCIP_CALL_FINALLY( SCIPexprPrintDotInit(set, stat, blkmem, printdata, f, whattoprint), 2371 fclose(f) ); 2372 (*printdata)->closefile = TRUE; 2373 2374 return SCIP_OKAY; 2375 } /*lint !e429*/ 2376 2377 /** main part of printing an expression in dot format */ 2378 SCIP_RETCODE SCIPexprPrintDot( 2379 SCIP_SET* set, /**< global SCIP settings */ 2380 SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */ 2381 SCIP_EXPRPRINTDATA* printdata, /**< data as initialized by \ref SCIPprintExprDotInit() */ 2382 SCIP_EXPR* expr /**< expression to be printed */ 2383 ) 2384 { 2385 SCIP_Real color; 2386 int c; 2387 2388 assert(set != NULL); 2389 assert(printdata != NULL); 2390 assert(expr != NULL); 2391 assert(expr->exprhdlr != NULL); 2392 2393 SCIP_CALL( SCIPexpriterInit(printdata->iterator, expr, SCIP_EXPRITER_DFS, FALSE) ); 2394 2395 while( !SCIPexpriterIsEnd(printdata->iterator) ) 2396 { 2397 /* print expression as dot node */ 2398 2399 if( expr->nchildren == 0 ) 2400 { 2401 SCIP_CALL( SCIPhashmapInsert(printdata->leaveexprs, (void*)expr, NULL) ); 2402 } 2403 2404 /* make up some color from the expression type (it's name) */ 2405 color = 0.0; 2406 for( c = 0; expr->exprhdlr->name[c] != '\0'; ++c ) 2407 color += (tolower(expr->exprhdlr->name[c]) - 'a') / 26.0; 2408 color = SCIPsetFrac(set, color); 2409 fprintf(printdata->file, "n%p [fillcolor=\"%g,%g,%g\", label=\"", (void*)expr, color, color, color); 2410 2411 if( printdata->whattoprint & SCIP_EXPRPRINT_EXPRHDLR ) 2412 { 2413 fprintf(printdata->file, "%s\\n", expr->exprhdlr->name); 2414 } 2415 2416 if( printdata->whattoprint & SCIP_EXPRPRINT_EXPRSTRING ) 2417 { 2418 SCIP_CALL( SCIPexprhdlrPrintExpr(expr->exprhdlr, set, messagehdlr, expr, SCIP_EXPRITER_ENTEREXPR, -1, 0, 2419 printdata->file) ); 2420 for( c = 0; c < expr->nchildren; ++c ) 2421 { 2422 SCIP_CALL( SCIPexprhdlrPrintExpr(expr->exprhdlr, set, messagehdlr, expr, SCIP_EXPRITER_VISITINGCHILD, 2423 c, 0, printdata->file) ); 2424 fprintf(printdata->file, "c%d", c); 2425 SCIP_CALL( SCIPexprhdlrPrintExpr(expr->exprhdlr, set, messagehdlr, expr, SCIP_EXPRITER_VISITEDCHILD, 2426 c, 0, printdata->file) ); 2427 } 2428 SCIP_CALL( SCIPexprhdlrPrintExpr(expr->exprhdlr, set, messagehdlr, expr, SCIP_EXPRITER_LEAVEEXPR, -1, 0, 2429 printdata->file) ); 2430 2431 fputs("\\n", printdata->file); 2432 } 2433 2434 if( printdata->whattoprint & SCIP_EXPRPRINT_NUSES ) 2435 { 2436 /* print number of uses */ 2437 fprintf(printdata->file, "%d uses\\n", expr->nuses); 2438 } 2439 2440 if( printdata->whattoprint & SCIP_EXPRPRINT_OWNER ) 2441 { 2442 /* print ownerdata */ 2443 if( expr->ownerprint != NULL ) 2444 { 2445 SCIP_CALL( expr->ownerprint(set->scip, printdata->file, expr, expr->ownerdata) ); 2446 } 2447 else if( expr->ownerdata != NULL ) 2448 { 2449 fprintf(printdata->file, "owner=%p\\n", (void*)expr->ownerdata); 2450 } 2451 } 2452 2453 if( printdata->whattoprint & SCIP_EXPRPRINT_EVALVALUE ) 2454 { 2455 /* print eval value */ 2456 fprintf(printdata->file, "val=%g", expr->evalvalue); 2457 2458 if( (printdata->whattoprint & SCIP_EXPRPRINT_EVALTAG) == SCIP_EXPRPRINT_EVALTAG ) 2459 { 2460 /* print also eval tag */ 2461 fprintf(printdata->file, " (%" SCIP_LONGINT_FORMAT ")", expr->evaltag); 2462 } 2463 fputs("\\n", printdata->file); 2464 } 2465 2466 if( printdata->whattoprint & SCIP_EXPRPRINT_ACTIVITY ) 2467 { 2468 /* print activity */ 2469 fprintf(printdata->file, "[%g,%g]", expr->activity.inf, expr->activity.sup); 2470 2471 if( (printdata->whattoprint & SCIP_EXPRPRINT_ACTIVITYTAG) == SCIP_EXPRPRINT_ACTIVITYTAG ) 2472 { 2473 /* print also activity eval tag */ 2474 fprintf(printdata->file, " (%" SCIP_LONGINT_FORMAT ")", expr->activitytag); 2475 } 2476 fputs("\\n", printdata->file); 2477 } 2478 2479 fputs("\"]\n", printdata->file); /* end of label and end of node */ 2480 2481 /* add edges from expr to its children */ 2482 for( c = 0; c < expr->nchildren; ++c ) 2483 fprintf(printdata->file, "n%p -> n%p [label=\"c%d\"]\n", (void*)expr, (void*)expr->children[c], c); 2484 2485 expr = SCIPexpriterGetNext(printdata->iterator); 2486 } 2487 2488 return SCIP_OKAY; 2489 } 2490 2491 /** finishes printing of expressions in dot format */ 2492 SCIP_RETCODE SCIPexprPrintDotFinal( 2493 SCIP_SET* set, /**< global SCIP settings */ 2494 SCIP_STAT* stat, /**< dynamic problem statistics */ 2495 BMS_BLKMEM* blkmem, /**< block memory */ 2496 SCIP_EXPRPRINTDATA** printdata /**< buffer where dot printing data has been stored */ 2497 ) 2498 { 2499 SCIP_EXPR* expr; 2500 SCIP_HASHMAPENTRY* entry; 2501 FILE* file; 2502 int i; 2503 2504 assert(set != NULL); 2505 assert(stat != NULL); 2506 assert(blkmem != NULL); 2507 assert(printdata != NULL); 2508 assert(*printdata != NULL); 2509 2510 file = (*printdata)->file; 2511 assert(file != NULL); 2512 2513 /* iterate through all entries of the map */ 2514 fputs("{rank=same;", file); 2515 for( i = 0; i < SCIPhashmapGetNEntries((*printdata)->leaveexprs); ++i ) 2516 { 2517 entry = SCIPhashmapGetEntry((*printdata)->leaveexprs, i); 2518 2519 if( entry != NULL ) 2520 { 2521 expr = (SCIP_EXPR*) SCIPhashmapEntryGetOrigin(entry); 2522 assert(expr != NULL); 2523 assert(expr->nchildren == 0); 2524 2525 fprintf(file, " n%p", (void*)expr); 2526 } 2527 } 2528 fprintf(file, "}\n"); 2529 2530 fprintf(file, "}\n"); 2531 2532 SCIPhashmapFree(&(*printdata)->leaveexprs); 2533 2534 SCIPexpriterFree(&(*printdata)->iterator); 2535 2536 if( (*printdata)->closefile ) 2537 fclose((*printdata)->file); 2538 2539 BMSfreeBlockMemory(blkmem, printdata); 2540 2541 return SCIP_OKAY; 2542 } 2543 2544 /** prints structure of an expression a la Maple's dismantle */ 2545 SCIP_RETCODE SCIPexprDismantle( 2546 SCIP_SET* set, /**< global SCIP settings */ 2547 SCIP_STAT* stat, /**< dynamic problem statistics */ 2548 BMS_BLKMEM* blkmem, /**< block memory */ 2549 SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */ 2550 FILE* file, /**< file to print to, or NULL for stdout */ 2551 SCIP_EXPR* expr /**< expression to dismantle */ 2552 ) 2553 { 2554 SCIP_EXPRITER* it; 2555 int depth = -1; 2556 2557 assert(set != NULL); 2558 assert(stat != NULL); 2559 assert(blkmem != NULL); 2560 assert(messagehdlr != NULL); 2561 assert(expr != NULL); 2562 2563 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2564 SCIP_CALL( SCIPexpriterInit(it, expr, SCIP_EXPRITER_DFS, TRUE) ); 2565 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_ENTEREXPR | SCIP_EXPRITER_VISITINGCHILD | SCIP_EXPRITER_LEAVEEXPR); 2566 2567 for( ; !SCIPexpriterIsEnd(it); expr = SCIPexpriterGetNext(it) ) 2568 { 2569 switch( SCIPexpriterGetStageDFS(it) ) 2570 { 2571 case SCIP_EXPRITER_ENTEREXPR: 2572 { 2573 int nspaces; 2574 2575 ++depth; 2576 nspaces = 3 * depth; 2577 2578 /* use depth of expression to align output */ 2579 SCIPmessageFPrintInfo(messagehdlr, file, "%*s[%s]: ", nspaces, "", expr->exprhdlr->name); 2580 2581 if( SCIPexprIsVar(set, expr) ) 2582 { 2583 SCIP_VAR* var; 2584 2585 var = SCIPgetVarExprVar(expr); 2586 SCIPmessageFPrintInfo(messagehdlr, file, "%s in [%g, %g]", SCIPvarGetName(var), SCIPvarGetLbLocal(var), 2587 SCIPvarGetUbLocal(var)); 2588 } 2589 else if( SCIPexprIsSum(set, expr) ) 2590 SCIPmessageFPrintInfo(messagehdlr, file, "%g", SCIPgetConstantExprSum(expr)); 2591 else if( SCIPexprIsProduct(set, expr) ) 2592 SCIPmessageFPrintInfo(messagehdlr, file, "%g", SCIPgetCoefExprProduct(expr)); 2593 else if( SCIPexprIsValue(set, expr) ) 2594 SCIPmessageFPrintInfo(messagehdlr, file, "%g", SCIPgetValueExprValue(expr)); 2595 else if( SCIPexprIsPower(set, expr) || strcmp(expr->exprhdlr->name, "signpower") == 0) 2596 SCIPmessageFPrintInfo(messagehdlr, file, "%g", SCIPgetExponentExprPow(expr)); 2597 2598 SCIPmessageFPrintInfo(messagehdlr, file, "\n"); 2599 2600 if( expr->ownerprint != NULL ) 2601 { 2602 SCIPmessageFPrintInfo(messagehdlr, file, "%*s ", nspaces, ""); 2603 SCIP_CALL( expr->ownerprint(set->scip, file, expr, expr->ownerdata) ); 2604 } 2605 2606 break; 2607 } 2608 2609 case SCIP_EXPRITER_VISITINGCHILD: 2610 { 2611 int nspaces = 3 * depth; 2612 2613 if( SCIPexprIsSum(set, expr) ) 2614 { 2615 SCIPmessageFPrintInfo(messagehdlr, file, "%*s ", nspaces, ""); 2616 SCIPmessageFPrintInfo(messagehdlr, file, "[coef]: %g\n", SCIPgetCoefsExprSum(expr)[SCIPexpriterGetChildIdxDFS(it)]); 2617 } 2618 2619 break; 2620 } 2621 2622 case SCIP_EXPRITER_LEAVEEXPR: 2623 { 2624 --depth; 2625 break; 2626 } 2627 2628 default: 2629 /* shouldn't be here */ 2630 SCIPABORT(); 2631 break; 2632 } 2633 } 2634 2635 SCIPexpriterFree(&it); 2636 2637 return SCIP_OKAY; 2638 } 2639 2640 /** evaluate an expression in a point 2641 * 2642 * Iterates over expressions to also evaluate children, if necessary. 2643 * Value can be received via SCIPexprGetEvalValue(). 2644 * If an evaluation error (division by zero, ...) occurs, this value will 2645 * be set to SCIP_INVALID. 2646 * 2647 * If a nonzero \p soltag is passed, then only (sub)expressions are 2648 * reevaluated that have a different solution tag. If a soltag of 0 2649 * is passed, then subexpressions are always reevaluated. 2650 * The tag is stored together with the value and can be received via 2651 * SCIPexprGetEvalTag(). 2652 */ 2653 SCIP_RETCODE SCIPexprEval( 2654 SCIP_SET* set, /**< global SCIP settings */ 2655 SCIP_STAT* stat, /**< dynamic problem statistics */ 2656 BMS_BLKMEM* blkmem, /**< block memory */ 2657 SCIP_EXPR* expr, /**< expression to be evaluated */ 2658 SCIP_SOL* sol, /**< solution to be evaluated */ 2659 SCIP_Longint soltag /**< tag that uniquely identifies the solution (with its values), or 0. */ 2660 ) 2661 { 2662 SCIP_EXPRITER* it; 2663 2664 assert(set != NULL); 2665 assert(stat != NULL); 2666 assert(blkmem != NULL); 2667 assert(expr != NULL); 2668 2669 /* if value is up-to-date, then nothing to do */ 2670 if( soltag != 0 && expr->evaltag == soltag ) 2671 return SCIP_OKAY; 2672 2673 /* assume we'll get a domain error, so we don't have to get this expr back if we abort the iteration 2674 * if there is no domain error, then we will overwrite the evalvalue in the last leaveexpr stage 2675 */ 2676 expr->evalvalue = SCIP_INVALID; 2677 expr->evaltag = soltag; 2678 2679 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2680 SCIP_CALL( SCIPexpriterInit(it, expr, SCIP_EXPRITER_DFS, TRUE) ); 2681 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD | SCIP_EXPRITER_LEAVEEXPR); 2682 2683 while( !SCIPexpriterIsEnd(it) ) 2684 { 2685 switch( SCIPexpriterGetStageDFS(it) ) 2686 { 2687 case SCIP_EXPRITER_VISITINGCHILD : 2688 { 2689 SCIP_EXPR* child; 2690 2691 if( soltag == 0 ) 2692 break; 2693 2694 /* check whether child has been evaluated for that solution already */ 2695 child = SCIPexpriterGetChildExprDFS(it); 2696 if( soltag == child->evaltag ) 2697 { 2698 if( child->evalvalue == SCIP_INVALID ) 2699 goto TERMINATE; 2700 2701 /* skip this child 2702 * this already returns the next one, so continue with loop 2703 */ 2704 expr = SCIPexpriterSkipDFS(it); 2705 continue; 2706 } 2707 2708 break; 2709 } 2710 2711 case SCIP_EXPRITER_LEAVEEXPR : 2712 { 2713 SCIP_CALL( SCIPexprhdlrEvalExpr(expr->exprhdlr, set, NULL , expr, &expr->evalvalue, NULL, sol) ); 2714 expr->evaltag = soltag; 2715 2716 if( expr->evalvalue == SCIP_INVALID ) 2717 goto TERMINATE; 2718 2719 break; 2720 } 2721 2722 default : 2723 /* we should never be here */ 2724 SCIPABORT(); 2725 break; 2726 } 2727 2728 expr = SCIPexpriterGetNext(it); 2729 } 2730 2731 TERMINATE: 2732 SCIPexpriterFree(&it); 2733 2734 return SCIP_OKAY; 2735 } 2736 2737 /** evaluates gradient of an expression for a given point 2738 * 2739 * Initiates an expression walk to also evaluate children, if necessary. 2740 * Value can be received via SCIPgetExprPartialDiffNonlinear(). 2741 * If an error (division by zero, ...) occurs, this value will 2742 * be set to SCIP_INVALID. 2743 */ 2744 SCIP_RETCODE SCIPexprEvalGradient( 2745 SCIP_SET* set, /**< global SCIP settings */ 2746 SCIP_STAT* stat, /**< dynamic problem statistics */ 2747 BMS_BLKMEM* blkmem, /**< block memory */ 2748 SCIP_EXPR* rootexpr, /**< expression to be evaluated */ 2749 SCIP_SOL* sol, /**< solution to be evaluated (NULL for the current LP solution) */ 2750 SCIP_Longint soltag /**< tag that uniquely identifies the solution (with its values), or 0. */ 2751 ) 2752 { 2753 SCIP_EXPRITER* it; 2754 SCIP_EXPR* expr; 2755 SCIP_EXPR* child; 2756 SCIP_Real derivative; 2757 SCIP_Longint difftag; 2758 2759 assert(set != NULL); 2760 assert(stat != NULL); 2761 assert(blkmem != NULL); 2762 assert(rootexpr != NULL); 2763 2764 /* ensure expression is evaluated */ 2765 SCIP_CALL( SCIPexprEval(set, stat, blkmem, rootexpr, sol, soltag) ); 2766 2767 /* check if expression could not be evaluated */ 2768 if( SCIPexprGetEvalValue(rootexpr) == SCIP_INVALID ) 2769 { 2770 rootexpr->derivative = SCIP_INVALID; 2771 return SCIP_OKAY; 2772 } 2773 2774 if( SCIPexprIsValue(set, rootexpr) ) 2775 { 2776 rootexpr->derivative = 0.0; 2777 return SCIP_OKAY; 2778 } 2779 2780 difftag = ++(stat->exprlastdifftag); 2781 2782 rootexpr->derivative = 1.0; 2783 rootexpr->difftag = difftag; 2784 2785 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2786 SCIP_CALL( SCIPexpriterInit(it, rootexpr, SCIP_EXPRITER_DFS, TRUE) ); 2787 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD); 2788 2789 for( expr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it); expr = SCIPexpriterGetNext(it) ) 2790 { 2791 assert(expr->evalvalue != SCIP_INVALID); 2792 2793 child = SCIPexpriterGetChildExprDFS(it); 2794 assert(child != NULL); 2795 2796 /* reset the value of the partial derivative w.r.t. a variable expression if we see it for the first time */ 2797 if( child->difftag != difftag && SCIPexprIsVar(set, child) ) 2798 child->derivative = 0.0; 2799 2800 /* update differentiation tag of the child */ 2801 child->difftag = difftag; 2802 2803 /* call backward differentiation callback */ 2804 if( SCIPexprIsValue(set, child) ) 2805 { 2806 derivative = 0.0; 2807 } 2808 else 2809 { 2810 derivative = SCIP_INVALID; 2811 SCIP_CALL( SCIPexprhdlrBwDiffExpr(expr->exprhdlr, set, NULL, expr, SCIPexpriterGetChildIdxDFS(it), 2812 &derivative, NULL, 0.0) ); 2813 2814 if( derivative == SCIP_INVALID ) 2815 { 2816 rootexpr->derivative = SCIP_INVALID; 2817 break; 2818 } 2819 } 2820 2821 /* update partial derivative stored in the child expression 2822 * for a variable, we have to sum up the partial derivatives of the root w.r.t. this variable over all parents 2823 * for other intermediate expressions, we only store the partial derivative of the root w.r.t. this expression 2824 */ 2825 if( !SCIPexprIsVar(set, child) ) 2826 child->derivative = expr->derivative * derivative; 2827 else 2828 child->derivative += expr->derivative * derivative; 2829 } 2830 2831 SCIPexpriterFree(&it); 2832 2833 return SCIP_OKAY; 2834 } 2835 2836 /** evaluates Hessian-vector product of an expression for a given point and direction 2837 * 2838 * Evaluates children, if necessary. 2839 * Value can be received via SCIPgetExprPartialDiffGradientDirNonlinear() 2840 * If an error (division by zero, ...) occurs, this value will 2841 * be set to SCIP_INVALID. 2842 */ 2843 SCIP_RETCODE SCIPexprEvalHessianDir( 2844 SCIP_SET* set, /**< global SCIP settings */ 2845 SCIP_STAT* stat, /**< dynamic problem statistics */ 2846 BMS_BLKMEM* blkmem, /**< block memory */ 2847 SCIP_EXPR* rootexpr, /**< expression to be evaluated */ 2848 SCIP_SOL* sol, /**< solution to be evaluated (NULL for the current LP solution) */ 2849 SCIP_Longint soltag, /**< tag that uniquely identifies the solution (with its values), or 0. */ 2850 SCIP_SOL* direction /**< direction */ 2851 ) 2852 { 2853 SCIP_EXPRITER* it; 2854 SCIP_EXPR* expr; 2855 SCIP_EXPR* child; 2856 SCIP_Real derivative; 2857 SCIP_Real hessiandir; 2858 2859 assert(set != NULL); 2860 assert(stat != NULL); 2861 assert(blkmem != NULL); 2862 assert(rootexpr != NULL); 2863 2864 if( SCIPexprIsValue(set, rootexpr) ) 2865 { 2866 rootexpr->dot = 0.0; 2867 rootexpr->bardot = 0.0; 2868 return SCIP_OKAY; 2869 } 2870 2871 /* evaluate expression and directional derivative */ 2872 SCIP_CALL( evalAndDiff(set, stat, blkmem, rootexpr, sol, soltag, direction) ); 2873 2874 if( rootexpr->evalvalue == SCIP_INVALID ) 2875 { 2876 rootexpr->derivative = SCIP_INVALID; 2877 rootexpr->bardot = SCIP_INVALID; 2878 return SCIP_OKAY; 2879 } 2880 2881 rootexpr->derivative = 1.0; 2882 2883 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2884 SCIP_CALL( SCIPexpriterInit(it, rootexpr, SCIP_EXPRITER_DFS, TRUE) ); 2885 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD); 2886 2887 /* compute reverse diff and bardots: i.e. hessian times direction */ 2888 for( expr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it); expr = SCIPexpriterGetNext(it) ) 2889 { 2890 assert(expr->evalvalue != SCIP_INVALID); 2891 2892 child = SCIPexpriterGetChildExprDFS(it); 2893 assert(child != NULL); 2894 2895 /* call backward and forward-backward differentiation callback */ 2896 if( SCIPexprIsValue(set, child) ) 2897 { 2898 derivative = 0.0; 2899 hessiandir = 0.0; 2900 } 2901 else 2902 { 2903 derivative = SCIP_INVALID; 2904 hessiandir = SCIP_INVALID; 2905 SCIP_CALL( SCIPexprhdlrBwDiffExpr(expr->exprhdlr, set, NULL, expr, SCIPexpriterGetChildIdxDFS(it), 2906 &derivative, NULL, SCIP_INVALID) ); 2907 SCIP_CALL( SCIPexprhdlrBwFwDiffExpr(expr->exprhdlr, set, expr, SCIPexpriterGetChildIdxDFS(it), 2908 &hessiandir, NULL) ); 2909 2910 if( derivative == SCIP_INVALID || hessiandir == SCIP_INVALID ) 2911 { 2912 rootexpr->derivative = SCIP_INVALID; 2913 rootexpr->bardot = SCIP_INVALID; 2914 break; 2915 } 2916 } 2917 2918 /* update partial derivative and hessian stored in the child expression 2919 * for a variable, we have to sum up the partial derivatives of the root w.r.t. this variable over all parents 2920 * for other intermediate expressions, we only store the partial derivative of the root w.r.t. this expression 2921 */ 2922 if( !SCIPexprIsVar(set, child) ) 2923 { 2924 child->derivative = expr->derivative * derivative; 2925 child->bardot = expr->bardot * derivative + expr->derivative * hessiandir; 2926 } 2927 else 2928 { 2929 child->derivative += expr->derivative * derivative; 2930 child->bardot += expr->bardot * derivative + expr->derivative * hessiandir; 2931 } 2932 } 2933 2934 SCIPexpriterFree(&it); 2935 2936 return SCIP_OKAY; 2937 } 2938 2939 /** possibly reevaluates and then returns the activity of the expression 2940 * 2941 * Reevaluate activity if currently stored is no longer uptodate. 2942 * If the expr owner provided a evalactivity-callback, then call this. 2943 * Otherwise, loop over descendants and compare activitytag with stat's domchgcount, i.e., 2944 * whether some bound was changed since last evaluation, to check whether exprhdlrs INTEVAL should be called. 2945 * 2946 * @note If expression is set to be integral, then activities are tightened to integral values. 2947 * Thus, ensure that the integrality information is valid (if set to TRUE; the default (FALSE) is always ok). 2948 */ 2949 SCIP_RETCODE SCIPexprEvalActivity( 2950 SCIP_SET* set, /**< global SCIP settings */ 2951 SCIP_STAT* stat, /**< dynamic problem statistics */ 2952 BMS_BLKMEM* blkmem, /**< block memory */ 2953 SCIP_EXPR* rootexpr /**< expression */ 2954 ) 2955 { 2956 SCIP_EXPRITER* it; 2957 SCIP_EXPR* expr; 2958 2959 assert(set != NULL); 2960 assert(stat != NULL); 2961 assert(blkmem != NULL); 2962 assert(rootexpr != NULL); 2963 2964 if( rootexpr->ownerevalactivity != NULL ) 2965 { 2966 /* call owner callback for activity-eval */ 2967 SCIP_CALL( rootexpr->ownerevalactivity(set->scip, rootexpr, rootexpr->ownerdata) ); 2968 2969 return SCIP_OKAY; 2970 } 2971 2972 /* fallback if no callback is given */ 2973 2974 assert(rootexpr->activitytag <= stat->domchgcount); 2975 2976 /* if value is up-to-date, then nothing to do */ 2977 if( rootexpr->activitytag == stat->domchgcount ) 2978 { 2979 #ifdef DEBUG_PROP 2980 SCIPsetDebugMsg(set, "activitytag of root expr equals domchgcount (%u), skip evalactivity\n", stat->domchgcount); 2981 #endif 2982 2983 return SCIP_OKAY; 2984 } 2985 2986 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 2987 SCIP_CALL( SCIPexpriterInit(it, rootexpr, SCIP_EXPRITER_DFS, TRUE) ); 2988 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD | SCIP_EXPRITER_LEAVEEXPR); 2989 2990 for( expr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it); ) 2991 { 2992 switch( SCIPexpriterGetStageDFS(it) ) 2993 { 2994 case SCIP_EXPRITER_VISITINGCHILD : 2995 { 2996 /* skip child if it has been evaluated already */ 2997 SCIP_EXPR* child; 2998 2999 child = SCIPexpriterGetChildExprDFS(it); 3000 if( child->activitytag == stat->domchgcount ) 3001 { 3002 expr = SCIPexpriterSkipDFS(it); 3003 continue; 3004 } 3005 3006 break; 3007 } 3008 3009 case SCIP_EXPRITER_LEAVEEXPR : 3010 { 3011 /* we should not have entered this expression if its activity was already uptodate */ 3012 assert(expr->activitytag < stat->domchgcount); 3013 3014 /* reset activity to entire if invalid, so we can use it as starting point below */ 3015 SCIPintervalSetEntire(SCIP_INTERVAL_INFINITY, &expr->activity); 3016 3017 #ifdef DEBUG_PROP 3018 SCIPsetDebugMsg(set, "interval evaluation of expr %p ", (void*)expr); 3019 SCIP_CALL( SCIPprintExpr(set->scip, expr, NULL) ); 3020 SCIPsetDebugMsgPrint(set, "\n"); 3021 #endif 3022 3023 /* call the inteval callback of the exprhdlr */ 3024 SCIP_CALL( SCIPexprhdlrIntEvalExpr(expr->exprhdlr, set, expr, &expr->activity, NULL, NULL) ); 3025 #ifdef DEBUG_PROP 3026 SCIPsetDebugMsg(set, " exprhdlr <%s>::inteval = [%.20g, %.20g]", expr->exprhdlr->name, expr->activity.inf, 3027 expr->activity.sup); 3028 #endif 3029 3030 /* if expression is integral, then we try to tighten the interval bounds a bit 3031 * this should undo the addition of some unnecessary safety added by use of nextafter() in interval 3032 * arithmetics, e.g., when doing pow() it would be ok to use ceil() and floor(), but for safety we 3033 * use SCIPceil and SCIPfloor for now the default intevalVar does not relax variables, so can omit 3034 * expressions without children (constants should be ok, too) 3035 */ 3036 if( expr->isintegral && expr->nchildren > 0 ) 3037 { 3038 if( expr->activity.inf > -SCIP_INTERVAL_INFINITY ) 3039 expr->activity.inf = SCIPsetCeil(set, expr->activity.inf); 3040 if( expr->activity.sup < SCIP_INTERVAL_INFINITY ) 3041 expr->activity.sup = SCIPsetFloor(set, expr->activity.sup); 3042 #ifdef DEBUG_PROP 3043 SCIPsetDebugMsg(set, " applying integrality: [%.20g, %.20g]\n", expr->activity.inf, expr->activity.sup); 3044 #endif 3045 } 3046 3047 /* mark activity as empty if either the lower/upper bound is above/below +/- SCIPinfinity() 3048 * TODO this is a problem if dual-presolve fixed a variable to +/- infinity 3049 */ 3050 if( SCIPsetIsInfinity(set, expr->activity.inf) || SCIPsetIsInfinity(set, -expr->activity.sup) ) 3051 { 3052 SCIPsetDebugMsg(set, "treat activity [%g,%g] as empty as beyond infinity\n", expr->activity.inf, expr->activity.sup); 3053 SCIPintervalSetEmpty(&expr->activity); 3054 } 3055 3056 /* remember that activity is uptodate now */ 3057 expr->activitytag = stat->domchgcount; 3058 3059 break; 3060 } 3061 3062 default: 3063 /* you should never be here */ 3064 SCIPABORT(); 3065 break; 3066 } 3067 3068 expr = SCIPexpriterGetNext(it); 3069 } 3070 3071 SCIPexpriterFree(&it); 3072 3073 return SCIP_OKAY; 3074 } 3075 3076 /** compare expressions 3077 * 3078 * @return -1, 0 or 1 if expr1 <, =, > expr2, respectively 3079 * @note The given expressions are assumed to be simplified. 3080 */ 3081 int SCIPexprCompare( 3082 SCIP_SET* set, /**< global SCIP settings */ 3083 SCIP_EXPR* expr1, /**< first expression */ 3084 SCIP_EXPR* expr2 /**< second expression */ 3085 ) 3086 { 3087 SCIP_EXPRHDLR* exprhdlr1; 3088 SCIP_EXPRHDLR* exprhdlr2; 3089 int retval; 3090 3091 exprhdlr1 = expr1->exprhdlr; 3092 exprhdlr2 = expr2->exprhdlr; 3093 3094 /* expressions are of the same kind/type; use compare callback or default method */ 3095 if( exprhdlr1 == exprhdlr2 ) 3096 return SCIPexprhdlrCompareExpr(set, expr1, expr2); 3097 3098 /* expressions are of different kind/type */ 3099 /* enforces OR6 */ 3100 if( SCIPexprIsValue(set, expr1) ) 3101 { 3102 return -1; 3103 } 3104 /* enforces OR12 */ 3105 if( SCIPexprIsValue(set, expr2) ) 3106 return -SCIPexprCompare(set, expr2, expr1); 3107 3108 /* enforces OR7 */ 3109 if( SCIPexprIsSum(set, expr1) ) 3110 { 3111 int compareresult; 3112 int nchildren; 3113 3114 nchildren = expr1->nchildren; 3115 compareresult = SCIPexprCompare(set, expr1->children[nchildren-1], expr2); 3116 3117 if( compareresult != 0 ) 3118 return compareresult; 3119 3120 /* "base" of the largest expression of the sum is equal to expr2, coefficient might tell us that 3121 * expr2 is larger */ 3122 if( SCIPgetCoefsExprSum(expr1)[nchildren-1] < 1.0 ) 3123 return -1; 3124 3125 /* largest expression of sum is larger or equal than expr2 => expr1 > expr2 */ 3126 return 1; 3127 } 3128 /* enforces OR12 */ 3129 if( SCIPexprIsSum(set, expr2) ) 3130 return -SCIPexprCompare(set, expr2, expr1); 3131 3132 /* enforces OR8 */ 3133 if( SCIPexprIsProduct(set, expr1) ) 3134 { 3135 int compareresult; 3136 int nchildren; 3137 3138 nchildren = expr1->nchildren; 3139 compareresult = SCIPexprCompare(set, expr1->children[nchildren-1], expr2); 3140 3141 if( compareresult != 0 ) 3142 return compareresult; 3143 3144 /* largest expression of product is larger or equal than expr2 => expr1 > expr2 */ 3145 return 1; 3146 } 3147 /* enforces OR12 */ 3148 if( SCIPexprIsProduct(set, expr2) ) 3149 return -SCIPexprCompare(set, expr2, expr1); 3150 3151 /* enforces OR9 */ 3152 if( SCIPexprIsPower(set, expr1) ) 3153 { 3154 int compareresult; 3155 3156 compareresult = SCIPexprCompare(set, expr1->children[0], expr2); 3157 3158 if( compareresult != 0 ) 3159 return compareresult; 3160 3161 /* base equal to expr2, exponent might tell us that expr2 is larger */ 3162 if( SCIPgetExponentExprPow(expr1) < 1.0 ) 3163 return -1; 3164 3165 /* power expression is larger => expr1 > expr2 */ 3166 return 1; 3167 } 3168 /* enforces OR12 */ 3169 if( SCIPexprIsPower(set, expr2) ) 3170 return -SCIPexprCompare(set, expr2, expr1); 3171 3172 /* enforces OR10 */ 3173 if( SCIPexprIsVar(set, expr1) ) 3174 return -1; 3175 /* enforces OR12 */ 3176 if( SCIPexprIsVar(set, expr2) ) 3177 return -SCIPexprCompare(set, expr2, expr1); 3178 3179 /* enforces OR11 */ 3180 retval = strcmp(SCIPexprhdlrGetName(exprhdlr1), SCIPexprhdlrGetName(exprhdlr2)); 3181 return retval == 0 ? 0 : retval < 0 ? -1 : 1; 3182 } 3183 3184 /** simplifies an expression 3185 * 3186 * @see SCIPsimplifyExpr 3187 */ 3188 SCIP_RETCODE SCIPexprSimplify( 3189 SCIP_SET* set, /**< global SCIP settings */ 3190 SCIP_STAT* stat, /**< dynamic problem statistics */ 3191 BMS_BLKMEM* blkmem, /**< block memory */ 3192 SCIP_EXPR* rootexpr, /**< expression to be simplified */ 3193 SCIP_EXPR** simplified, /**< buffer to store simplified expression */ 3194 SCIP_Bool* changed, /**< buffer to store if rootexpr actually changed */ 3195 SCIP_Bool* infeasible, /**< buffer to store whether infeasibility has been detected */ 3196 SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call to create ownerdata */ 3197 void* ownercreatedata /**< data to pass to ownercreate */ 3198 ) 3199 { 3200 SCIP_EXPR* expr; 3201 SCIP_EXPRITER* it; 3202 3203 assert(rootexpr != NULL); 3204 assert(simplified != NULL); 3205 assert(changed != NULL); 3206 assert(infeasible != NULL); 3207 3208 /* simplify bottom up 3209 * when leaving an expression it simplifies it and stores the simplified expr in its iterators expression data 3210 * after the child was visited, it is replaced with the simplified expr 3211 */ 3212 SCIP_CALL( SCIPexpriterCreate(stat, blkmem, &it) ); 3213 SCIP_CALL( SCIPexpriterInit(it, rootexpr, SCIP_EXPRITER_DFS, TRUE) ); /* TODO can we set allowrevisited to FALSE?*/ 3214 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITEDCHILD | SCIP_EXPRITER_LEAVEEXPR); 3215 3216 *changed = FALSE; 3217 *infeasible = FALSE; 3218 for( expr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it); expr = SCIPexpriterGetNext(it) ) 3219 { 3220 switch( SCIPexpriterGetStageDFS(it) ) 3221 { 3222 case SCIP_EXPRITER_VISITEDCHILD: 3223 { 3224 SCIP_EXPR* newchild; 3225 SCIP_EXPR* child; 3226 3227 newchild = (SCIP_EXPR*)SCIPexpriterGetChildUserDataDFS(it).ptrval; 3228 child = SCIPexpriterGetChildExprDFS(it); 3229 assert(newchild != NULL); 3230 3231 /* if child got simplified, replace it with the new child */ 3232 if( newchild != child ) 3233 { 3234 SCIP_CALL( SCIPexprReplaceChild(set, stat, blkmem, expr, SCIPexpriterGetChildIdxDFS(it), newchild) ); 3235 } 3236 3237 /* we do not need to hold newchild anymore */ 3238 SCIP_CALL( SCIPexprRelease(set, stat, blkmem, &newchild) ); 3239 3240 break; 3241 } 3242 3243 case SCIP_EXPRITER_LEAVEEXPR: 3244 { 3245 SCIP_EXPR* refexpr = NULL; 3246 SCIP_EXPRITER_USERDATA iterdata; 3247 3248 /* TODO we should do constant folding (handle that all children are value-expressions) here in a generic way 3249 * instead of reimplementing it in every handler 3250 */ 3251 3252 /* use simplification of expression handlers */ 3253 SCIP_CALL( SCIPexprhdlrSimplifyExpr(expr->exprhdlr, set, expr, &refexpr, ownercreate, ownercreatedata) ); 3254 assert(refexpr != NULL); 3255 if( expr != refexpr ) 3256 *changed = TRUE; 3257 3258 iterdata.ptrval = (void*) refexpr; 3259 SCIPexpriterSetCurrentUserData(it, iterdata); 3260 3261 break; 3262 } 3263 3264 default: 3265 SCIPABORT(); /* we should never be called in this stage */ 3266 break; 3267 } 3268 } 3269 3270 *simplified = (SCIP_EXPR*)SCIPexpriterGetExprUserData(it, rootexpr).ptrval; 3271 assert(*simplified != NULL); 3272 3273 SCIPexpriterFree(&it); 3274 3275 return SCIP_OKAY; 3276 } 3277 3278 /** method to retrieve symmetry information from an expression 3279 * 3280 * @see SCIPgetSymDataExpr 3281 */ 3282 SCIP_RETCODE SCIPexprGetSymData( 3283 SCIP_SET* set, /**< global SCIP settings */ 3284 SCIP_EXPR* expr, /**< expression from which information is retrieved */ 3285 SYM_EXPRDATA** symdata /**< buffer to store symmetry information */ 3286 ) 3287 { 3288 SCIP_EXPRHDLR* exprhdlr; 3289 3290 assert(set != NULL); 3291 assert(expr != NULL); 3292 assert(symdata != NULL); 3293 3294 exprhdlr = SCIPexprGetHdlr(expr); 3295 assert(exprhdlr != NULL); 3296 3297 if( exprhdlr->getsymdata != NULL ) 3298 { 3299 SCIP_CALL( exprhdlr->getsymdata(set->scip, expr, symdata) ); 3300 } 3301 else 3302 *symdata = NULL; 3303 3304 return SCIP_OKAY; 3305 } 3306 3307 /** checks whether an expression is quadratic 3308 * 3309 * An expression is quadratic if it is either a power expression with exponent 2.0, a product of two expressions, 3310 * or a sum of terms where at least one is a square or a product of two. 3311 * 3312 * Use \ref SCIPexprGetQuadraticData to get data about the representation as quadratic. 3313 */ 3314 SCIP_RETCODE SCIPexprCheckQuadratic( 3315 SCIP_SET* set, /**< global SCIP settings */ 3316 BMS_BLKMEM* blkmem, /**< block memory */ 3317 SCIP_EXPR* expr, /**< expression */ 3318 SCIP_Bool* isquadratic /**< buffer to store result */ 3319 ) 3320 { 3321 SCIP_HASHMAP* expr2idx; 3322 SCIP_HASHMAP* seenexpr = NULL; 3323 int nquadterms = 0; 3324 int nlinterms = 0; 3325 int nbilinterms = 0; 3326 int c; 3327 3328 assert(set != NULL); 3329 assert(blkmem != NULL); 3330 assert(expr != NULL); 3331 assert(isquadratic != NULL); 3332 3333 if( expr->quadchecked ) 3334 { 3335 *isquadratic = expr->quaddata != NULL; 3336 return SCIP_OKAY; 3337 } 3338 assert(expr->quaddata == NULL); 3339 3340 expr->quadchecked = TRUE; 3341 *isquadratic = FALSE; 3342 3343 /* check if expression is a quadratic expression */ 3344 SCIPsetDebugMsg(set, "checking if expr %p is quadratic\n", (void*)expr); 3345 3346 /* handle single square term */ 3347 if( SCIPexprIsPower(set, expr) && SCIPgetExponentExprPow(expr) == 2.0 ) 3348 { 3349 SCIPsetDebugMsg(set, "expr %p looks like square: fill data structures\n", (void*)expr); 3350 SCIP_ALLOC( BMSallocClearBlockMemory(blkmem, &expr->quaddata) ); 3351 3352 expr->quaddata->nquadexprs = 1; 3353 SCIP_ALLOC( BMSallocClearBlockMemoryArray(blkmem, &expr->quaddata->quadexprterms, 1) ); 3354 expr->quaddata->quadexprterms[0].expr = expr->children[0]; 3355 expr->quaddata->quadexprterms[0].sqrcoef = 1.0; 3356 3357 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, expr->quaddata->quadexprterms[0].expr); 3358 3359 *isquadratic = TRUE; 3360 return SCIP_OKAY; 3361 } 3362 3363 /* handle single bilinear term */ 3364 if( SCIPexprIsProduct(set, expr) && SCIPexprGetNChildren(expr) == 2 ) 3365 { 3366 SCIPsetDebugMsg(set, "expr %p looks like bilinear product: fill data structures\n", (void*)expr); 3367 SCIP_ALLOC( BMSallocClearBlockMemory(blkmem, &expr->quaddata) ); 3368 expr->quaddata->nquadexprs = 2; 3369 3370 SCIP_ALLOC( BMSallocClearBlockMemoryArray(blkmem, &expr->quaddata->quadexprterms, 2) ); 3371 expr->quaddata->quadexprterms[0].expr = SCIPexprGetChildren(expr)[0]; 3372 expr->quaddata->quadexprterms[0].nadjbilin = 1; 3373 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &expr->quaddata->quadexprterms[0].adjbilin, 1) ); 3374 expr->quaddata->quadexprterms[0].adjbilin[0] = 0; 3375 3376 expr->quaddata->quadexprterms[1].expr = SCIPexprGetChildren(expr)[1]; 3377 expr->quaddata->quadexprterms[1].nadjbilin = 1; 3378 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &expr->quaddata->quadexprterms[1].adjbilin, 1) ); 3379 expr->quaddata->quadexprterms[1].adjbilin[0] = 0; 3380 3381 expr->quaddata->nbilinexprterms = 1; 3382 SCIP_ALLOC( BMSallocClearBlockMemoryArray(blkmem, &expr->quaddata->bilinexprterms, 1) ); 3383 expr->quaddata->bilinexprterms[0].expr1 = SCIPexprGetChildren(expr)[1]; 3384 expr->quaddata->bilinexprterms[0].expr2 = SCIPexprGetChildren(expr)[0]; 3385 expr->quaddata->bilinexprterms[0].coef = 1.0; 3386 3387 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, expr->quaddata->quadexprterms[0].expr) 3388 && SCIPexprIsVar(set, expr->quaddata->quadexprterms[1].expr); 3389 3390 *isquadratic = TRUE; 3391 return SCIP_OKAY; 3392 } 3393 3394 /* neither a sum, nor a square, nor a bilinear term */ 3395 if( !SCIPexprIsSum(set, expr) ) 3396 return SCIP_OKAY; 3397 3398 SCIP_CALL( SCIPhashmapCreate(&seenexpr, blkmem, 2*SCIPexprGetNChildren(expr)) ); 3399 for( c = 0; c < SCIPexprGetNChildren(expr); ++c ) 3400 { 3401 SCIP_EXPR* child; 3402 3403 child = SCIPexprGetChildren(expr)[c]; 3404 assert(child != NULL); 3405 3406 if( SCIPexprIsPower(set, child) && SCIPgetExponentExprPow(child) == 2.0 ) /* quadratic term */ 3407 { 3408 SCIP_CALL( quadDetectProcessExpr(SCIPexprGetChildren(child)[0], seenexpr, &nquadterms, &nlinterms) ); 3409 } 3410 else if( SCIPexprIsProduct(set, child) && SCIPexprGetNChildren(child) == 2 ) /* bilinear term */ 3411 { 3412 ++nbilinterms; 3413 SCIP_CALL( quadDetectProcessExpr(SCIPexprGetChildren(child)[0], seenexpr, &nquadterms, &nlinterms) ); 3414 SCIP_CALL( quadDetectProcessExpr(SCIPexprGetChildren(child)[1], seenexpr, &nquadterms, &nlinterms) ); 3415 } 3416 else 3417 { 3418 /* first time seen linearly --> assign -1; ++nlinterms 3419 * not first time --> assign +=1; 3420 */ 3421 if( SCIPhashmapExists(seenexpr, (void*)child) ) 3422 { 3423 assert(SCIPhashmapGetImageInt(seenexpr, (void*)child) > 0); 3424 3425 SCIP_CALL( SCIPhashmapSetImageInt(seenexpr, (void*)child, SCIPhashmapGetImageInt(seenexpr, (void*)child) + 1) ); 3426 } 3427 else 3428 { 3429 ++nlinterms; 3430 SCIP_CALL( SCIPhashmapInsertInt(seenexpr, (void*)child, -1) ); 3431 } 3432 } 3433 } 3434 3435 if( nquadterms == 0 ) 3436 { 3437 /* only linear sum */ 3438 SCIPhashmapFree(&seenexpr); 3439 return SCIP_OKAY; 3440 } 3441 3442 SCIPsetDebugMsg(set, "expr %p looks quadratic: fill data structures\n", (void*)expr); 3443 3444 /* expr2idx maps expressions to indices; if index > 0, it is its index in the linexprs array, otherwise -index-1 is 3445 * its index in the quadexprterms array 3446 */ 3447 SCIP_CALL( SCIPhashmapCreate(&expr2idx, blkmem, nquadterms + nlinterms) ); 3448 3449 /* allocate memory, etc */ 3450 SCIP_ALLOC( BMSallocClearBlockMemory(blkmem, &expr->quaddata) ); 3451 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &expr->quaddata->quadexprterms, nquadterms) ); 3452 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &expr->quaddata->linexprs, nlinterms) ); 3453 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &expr->quaddata->lincoefs, nlinterms) ); 3454 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &expr->quaddata->bilinexprterms, nbilinterms) ); 3455 3456 expr->quaddata->constant = SCIPgetConstantExprSum(expr); 3457 3458 expr->quaddata->allexprsarevars = TRUE; 3459 /* for every term of the sum-expr */ 3460 for( c = 0; c < SCIPexprGetNChildren(expr); ++c ) 3461 { 3462 SCIP_EXPR* child; 3463 SCIP_Real coef; 3464 3465 child = SCIPexprGetChildren(expr)[c]; 3466 coef = SCIPgetCoefsExprSum(expr)[c]; 3467 3468 assert(child != NULL); 3469 assert(coef != 0.0); 3470 3471 if( SCIPexprIsPower(set, child) && SCIPgetExponentExprPow(child) == 2.0 ) /* quadratic term */ 3472 { 3473 SCIP_QUADEXPR_QUADTERM* quadexprterm; 3474 assert(SCIPexprGetNChildren(child) == 1); 3475 3476 child = SCIPexprGetChildren(child)[0]; 3477 assert(SCIPhashmapGetImageInt(seenexpr, (void *)child) > 0); 3478 3479 SCIP_CALL( quadDetectGetQuadexprterm(blkmem, child, expr2idx, seenexpr, expr->quaddata, &quadexprterm) ); 3480 assert(quadexprterm->expr == child); 3481 quadexprterm->sqrcoef = coef; 3482 quadexprterm->sqrexpr = SCIPexprGetChildren(expr)[c]; 3483 3484 if( expr->quaddata->allexprsarevars ) 3485 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, quadexprterm->expr); 3486 } 3487 else if( SCIPexprIsProduct(set, child) && SCIPexprGetNChildren(child) == 2 ) /* bilinear term */ 3488 { 3489 SCIP_QUADEXPR_BILINTERM* bilinexprterm; 3490 SCIP_QUADEXPR_QUADTERM* quadexprterm; 3491 SCIP_EXPR* expr1; 3492 SCIP_EXPR* expr2; 3493 3494 assert(SCIPgetCoefExprProduct(child) == 1.0); 3495 3496 expr1 = SCIPexprGetChildren(child)[0]; 3497 expr2 = SCIPexprGetChildren(child)[1]; 3498 assert(expr1 != NULL && expr2 != NULL); 3499 3500 bilinexprterm = &expr->quaddata->bilinexprterms[expr->quaddata->nbilinexprterms]; 3501 3502 bilinexprterm->coef = coef; 3503 if( SCIPhashmapGetImageInt(seenexpr, (void*)expr1) >= SCIPhashmapGetImageInt(seenexpr, (void*)expr2) ) 3504 { 3505 bilinexprterm->expr1 = expr1; 3506 bilinexprterm->expr2 = expr2; 3507 } 3508 else 3509 { 3510 bilinexprterm->expr1 = expr2; 3511 bilinexprterm->expr2 = expr1; 3512 } 3513 bilinexprterm->prodexpr = child; 3514 3515 SCIP_CALL( quadDetectGetQuadexprterm(blkmem, expr1, expr2idx, seenexpr, expr->quaddata, &quadexprterm) ); 3516 assert(quadexprterm->expr == expr1); 3517 quadexprterm->adjbilin[quadexprterm->nadjbilin] = expr->quaddata->nbilinexprterms; 3518 ++quadexprterm->nadjbilin; 3519 3520 if( expr->quaddata->allexprsarevars ) 3521 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, quadexprterm->expr); 3522 3523 SCIP_CALL( quadDetectGetQuadexprterm(blkmem, expr2, expr2idx, seenexpr, expr->quaddata, &quadexprterm) ); 3524 assert(quadexprterm->expr == expr2); 3525 quadexprterm->adjbilin[quadexprterm->nadjbilin] = expr->quaddata->nbilinexprterms; 3526 ++quadexprterm->nadjbilin; 3527 3528 if( expr->quaddata->allexprsarevars ) 3529 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, quadexprterm->expr); 3530 3531 ++expr->quaddata->nbilinexprterms; 3532 3533 /* store position of second factor in quadexprterms */ 3534 bilinexprterm->pos2 = SCIPhashmapGetImageInt(expr2idx, (void*)bilinexprterm->expr2); 3535 } 3536 else /* linear term */ 3537 { 3538 if( SCIPhashmapGetImageInt(seenexpr, (void*)child) < 0 ) 3539 { 3540 assert(SCIPhashmapGetImageInt(seenexpr, (void*)child) == -1); 3541 3542 /* expression only appears linearly */ 3543 expr->quaddata->linexprs[expr->quaddata->nlinexprs] = child; 3544 expr->quaddata->lincoefs[expr->quaddata->nlinexprs] = coef; 3545 expr->quaddata->nlinexprs++; 3546 3547 if( expr->quaddata->allexprsarevars ) 3548 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, child); 3549 } 3550 else 3551 { 3552 /* expression appears non-linearly: set lin coef */ 3553 SCIP_QUADEXPR_QUADTERM* quadexprterm; 3554 assert(SCIPhashmapGetImageInt(seenexpr, (void*)child) > 0); 3555 3556 SCIP_CALL( quadDetectGetQuadexprterm(blkmem, child, expr2idx, seenexpr, expr->quaddata, &quadexprterm) ); 3557 assert(quadexprterm->expr == child); 3558 quadexprterm->lincoef = coef; 3559 3560 if( expr->quaddata->allexprsarevars ) 3561 expr->quaddata->allexprsarevars = SCIPexprIsVar(set, quadexprterm->expr); 3562 } 3563 } 3564 } 3565 assert(expr->quaddata->nquadexprs == nquadterms); 3566 assert(expr->quaddata->nlinexprs == nlinterms); 3567 assert(expr->quaddata->nbilinexprterms == nbilinterms); 3568 3569 SCIPhashmapFree(&seenexpr); 3570 SCIPhashmapFree(&expr2idx); 3571 3572 *isquadratic = TRUE; 3573 3574 return SCIP_OKAY; 3575 } 3576 3577 /** frees information on quadratic representation of an expression 3578 * 3579 * Reverts SCIPexprCheckQuadratic(). 3580 * Before doing changes to an expression, it can be useful to call this function. 3581 */ 3582 void SCIPexprFreeQuadratic( 3583 BMS_BLKMEM* blkmem, /**< block memory */ 3584 SCIP_EXPR* expr /**< expression */ 3585 ) 3586 { 3587 int i; 3588 int n; 3589 3590 assert(blkmem != NULL); 3591 assert(expr != NULL); 3592 3593 expr->quadchecked = FALSE; 3594 3595 if( expr->quaddata == NULL ) 3596 return; 3597 3598 n = expr->quaddata->nquadexprs; 3599 3600 BMSfreeBlockMemoryArrayNull(blkmem, &expr->quaddata->linexprs, expr->quaddata->nlinexprs); 3601 BMSfreeBlockMemoryArrayNull(blkmem, &expr->quaddata->lincoefs, expr->quaddata->nlinexprs); 3602 BMSfreeBlockMemoryArrayNull(blkmem, &expr->quaddata->bilinexprterms, expr->quaddata->nbilinexprterms); 3603 BMSfreeBlockMemoryArrayNull(blkmem, &expr->quaddata->eigenvalues, n); 3604 if( expr->quaddata->eigenvectors != NULL ) /* check for NULL here before calculating n*n to avoid (harmless) overflow in case of large n */ 3605 BMSfreeBlockMemoryArray(blkmem, &expr->quaddata->eigenvectors, n * n); /*lint !e647*/ 3606 3607 for( i = 0; i < n; ++i ) 3608 { 3609 BMSfreeBlockMemoryArrayNull(blkmem, &expr->quaddata->quadexprterms[i].adjbilin, 3610 expr->quaddata->quadexprterms[i].adjbilinsize); 3611 } 3612 BMSfreeBlockMemoryArrayNull(blkmem, &expr->quaddata->quadexprterms, n); 3613 3614 BMSfreeBlockMemory(blkmem, &expr->quaddata); 3615 } 3616 3617 /** Checks the curvature of the quadratic function stored in quaddata 3618 * 3619 * For this, it builds the matrix Q of quadratic coefficients and computes its eigenvalues using LAPACK. 3620 * If Q is 3621 * - semidefinite positive -> curv is set to convex, 3622 * - semidefinite negative -> curv is set to concave, 3623 * - otherwise -> curv is set to unknown. 3624 * 3625 * If `assumevarfixed` is given and some expressions in quadratic terms correspond to variables present in 3626 * this hashmap, then the corresponding rows and columns are ignored in the matrix Q. 3627 */ 3628 SCIP_RETCODE SCIPexprComputeQuadraticCurvature( 3629 SCIP_SET* set, /**< global SCIP settings */ 3630 BMS_BLKMEM* blkmem, /**< block memory */ 3631 BMS_BUFMEM* bufmem, /**< buffer memory */ 3632 SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */ 3633 SCIP_EXPR* expr, /**< quadratic expression */ 3634 SCIP_EXPRCURV* curv, /**< pointer to store the curvature of quadratics */ 3635 SCIP_HASHMAP* assumevarfixed, /**< hashmap containing variables that should be assumed to be fixed, or NULL */ 3636 SCIP_Bool storeeigeninfo /**< whether the eigenvalues and eigenvectors should be stored */ 3637 ) 3638 { 3639 SCIP_QUADEXPR* quaddata; 3640 SCIP_HASHMAP* expr2matrix; 3641 double* matrix; 3642 double* alleigval; 3643 int nvars; 3644 int nn; 3645 int n; 3646 int i; 3647 3648 assert(set != NULL); 3649 assert(blkmem != NULL); 3650 assert(bufmem != NULL); 3651 assert(messagehdlr != NULL); 3652 assert(expr != NULL); 3653 assert(curv != NULL); 3654 3655 quaddata = expr->quaddata; 3656 assert(quaddata != NULL); 3657 3658 /* do not store eigen information if we are not considering full matrix */ 3659 if( assumevarfixed != NULL ) 3660 storeeigeninfo = FALSE; 3661 3662 if( quaddata->eigeninfostored || (quaddata->curvaturechecked && !storeeigeninfo) ) 3663 { 3664 *curv = quaddata->curvature; 3665 /* if we are convex or concave on the full set of variables, then we will also be so on a subset */ 3666 if( assumevarfixed == NULL || quaddata->curvature != SCIP_EXPRCURV_UNKNOWN ) 3667 return SCIP_OKAY; 3668 } 3669 assert(quaddata->curvature == SCIP_EXPRCURV_UNKNOWN || assumevarfixed != NULL 3670 || (storeeigeninfo && !quaddata->eigeninfostored)); 3671 3672 *curv = SCIP_EXPRCURV_UNKNOWN; 3673 3674 n = quaddata->nquadexprs; 3675 3676 /* do not check curvature if nn will be too large 3677 * we want nn * sizeof(real) to fit into an unsigned int, so n must be <= sqrt(unit_max/sizeof(real)) 3678 * sqrt(2*214748364/8) = 7327.1475350234 3679 */ 3680 if( n > 7000 ) 3681 { 3682 SCIPmessageFPrintVerbInfo(messagehdlr, set->disp_verblevel, SCIP_VERBLEVEL_FULL, NULL, 3683 "number of quadratic variables is too large (%d) to check the curvature\n", n); 3684 return SCIP_OKAY; 3685 } 3686 3687 /* TODO do some simple tests first; like diagonal entries don't change sign, etc */ 3688 3689 if( ! SCIPlapackIsAvailable() ) 3690 return SCIP_OKAY; 3691 3692 nn = n * n; 3693 assert(nn > 0); 3694 assert((unsigned)nn < UINT_MAX / sizeof(SCIP_Real)); 3695 3696 if( storeeigeninfo ) 3697 { 3698 SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, &quaddata->eigenvalues, n)); 3699 SCIP_ALLOC( BMSallocClearBlockMemoryArray(blkmem, &quaddata->eigenvectors, nn)); 3700 3701 alleigval = quaddata->eigenvalues; 3702 matrix = quaddata->eigenvectors; 3703 } 3704 else 3705 { 3706 SCIP_ALLOC( BMSallocBufferMemoryArray(bufmem, &alleigval, n) ); 3707 SCIP_ALLOC( BMSallocClearBufferMemoryArray(bufmem, &matrix, nn) ); 3708 } 3709 3710 SCIP_CALL( SCIPhashmapCreate(&expr2matrix, blkmem, n) ); 3711 3712 /* fill matrix's diagonal */ 3713 nvars = 0; 3714 for( i = 0; i < n; ++i ) 3715 { 3716 SCIP_QUADEXPR_QUADTERM quadexprterm; 3717 3718 quadexprterm = quaddata->quadexprterms[i]; 3719 3720 assert(!SCIPhashmapExists(expr2matrix, (void*)quadexprterm.expr)); 3721 3722 /* skip expr if it is a variable mentioned in assumevarfixed */ 3723 if( assumevarfixed != NULL && SCIPexprIsVar(set, quadexprterm.expr) 3724 && SCIPhashmapExists(assumevarfixed, (void*)SCIPgetVarExprVar(quadexprterm.expr)) ) 3725 continue; 3726 3727 if( quadexprterm.sqrcoef == 0.0 && ! storeeigeninfo ) 3728 { 3729 assert(quadexprterm.nadjbilin > 0); 3730 /* SCIPdebugMsg(scip, "var <%s> appears in bilinear term but is not squared 3731 * --> indefinite quadratic\n", SCIPvarGetName(quadexprterm.var)); */ 3732 goto CLEANUP; 3733 } 3734 3735 matrix[nvars * n + nvars] = quadexprterm.sqrcoef; 3736 3737 /* remember row of variable in matrix */ 3738 SCIP_CALL( SCIPhashmapInsert(expr2matrix, (void *)quadexprterm.expr, (void *)(size_t)nvars) ); 3739 nvars++; 3740 } 3741 3742 /* fill matrix's upper-diagonal */ 3743 for( i = 0; i < quaddata->nbilinexprterms; ++i ) 3744 { 3745 SCIP_QUADEXPR_BILINTERM bilinexprterm; 3746 int col; 3747 int row; 3748 3749 bilinexprterm = quaddata->bilinexprterms[i]; 3750 3751 /* each factor should have been added to expr2matrix unless it corresponds to a variable mentioned in assumevarfixed */ 3752 assert(SCIPhashmapExists(expr2matrix, (void*)bilinexprterm.expr1) 3753 || (assumevarfixed != NULL && SCIPexprIsVar(set, bilinexprterm.expr1) 3754 && SCIPhashmapExists(assumevarfixed, (void*)SCIPgetVarExprVar(bilinexprterm.expr1)))); 3755 assert(SCIPhashmapExists(expr2matrix, (void*)bilinexprterm.expr2) 3756 || (assumevarfixed != NULL && SCIPexprIsVar(set, bilinexprterm.expr2) 3757 && SCIPhashmapExists(assumevarfixed, (void*)SCIPgetVarExprVar(bilinexprterm.expr2)))); 3758 3759 /* skip bilinear terms where at least one of the factors should be assumed to be fixed 3760 * (i.e., not present in expr2matrix map) */ 3761 if( !SCIPhashmapExists(expr2matrix, (void*)bilinexprterm.expr1) 3762 || !SCIPhashmapExists(expr2matrix, (void*)bilinexprterm.expr2) ) 3763 continue; 3764 3765 row = (int)(size_t)SCIPhashmapGetImage(expr2matrix, bilinexprterm.expr1); 3766 col = (int)(size_t)SCIPhashmapGetImage(expr2matrix, bilinexprterm.expr2); 3767 3768 assert(row != col); 3769 3770 if( row < col ) 3771 matrix[row * n + col] = bilinexprterm.coef / 2.0; 3772 else 3773 matrix[col * n + row] = bilinexprterm.coef / 2.0; 3774 } 3775 3776 /* compute eigenvalues */ 3777 if( SCIPlapackComputeEigenvalues(bufmem, storeeigeninfo, n, matrix, alleigval) != SCIP_OKAY ) 3778 { 3779 SCIPmessagePrintWarning(messagehdlr, "Failed to compute eigenvalues of quadratic coefficient " 3780 "matrix --> don't know curvature\n"); 3781 goto CLEANUP; 3782 } 3783 3784 /* check convexity */ 3785 if( !SCIPsetIsNegative(set, alleigval[0]) ) 3786 *curv = SCIP_EXPRCURV_CONVEX; 3787 else if( !SCIPsetIsPositive(set, alleigval[n-1]) ) 3788 *curv = SCIP_EXPRCURV_CONCAVE; 3789 3790 CLEANUP: 3791 SCIPhashmapFree(&expr2matrix); 3792 3793 if( !storeeigeninfo ) 3794 { 3795 BMSfreeBufferMemoryArray(bufmem, &matrix); 3796 BMSfreeBufferMemoryArray(bufmem, &alleigval); 3797 } 3798 else 3799 { 3800 assert(!quaddata->eigeninfostored); 3801 quaddata->eigeninfostored = TRUE; 3802 } 3803 3804 /* if checked convexity on full Q matrix, then remember it 3805 * if indefinite on submatrix, then it will also be indefinite on full matrix, so can remember that, too */ 3806 if( assumevarfixed == NULL || *curv == SCIP_EXPRCURV_UNKNOWN ) 3807 { 3808 quaddata->curvature = *curv; 3809 quaddata->curvaturechecked = TRUE; 3810 } 3811 3812 return SCIP_OKAY; 3813 } 3814 3815 3816 /* from pub_expr.h */ 3817 3818 #ifdef NDEBUG 3819 #undef SCIPexprGetNUses 3820 #undef SCIPexprGetNChildren 3821 #undef SCIPexprGetChildren 3822 #undef SCIPexprGetHdlr 3823 #undef SCIPexprGetData 3824 #undef SCIPexprSetData 3825 #undef SCIPexprGetOwnerData 3826 #undef SCIPexprGetEvalValue 3827 #undef SCIPexprGetEvalTag 3828 #undef SCIPexprGetDerivative 3829 #undef SCIPexprGetDot 3830 #undef SCIPexprGetBardot 3831 #undef SCIPexprGetDiffTag 3832 #undef SCIPexprGetActivity 3833 #undef SCIPexprGetActivityTag 3834 #undef SCIPexprSetActivity 3835 #undef SCIPexprGetCurvature 3836 #undef SCIPexprSetCurvature 3837 #undef SCIPexprIsIntegral 3838 #undef SCIPexprSetIntegrality 3839 #undef SCIPexprAreQuadraticExprsVariables 3840 #endif 3841 3842 /** gets the number of times the expression is currently captured */ 3843 int SCIPexprGetNUses( 3844 SCIP_EXPR* expr /**< expression */ 3845 ) 3846 { 3847 assert(expr != NULL); 3848 3849 return expr->nuses; 3850 } 3851 3852 /** gives the number of children of an expression */ 3853 int SCIPexprGetNChildren( 3854 SCIP_EXPR* expr /**< expression */ 3855 ) 3856 { 3857 assert(expr != NULL); 3858 3859 return expr->nchildren; 3860 } 3861 3862 /** gives the children of an expression (can be NULL if no children) */ 3863 SCIP_EXPR** SCIPexprGetChildren( 3864 SCIP_EXPR* expr /**< expression */ 3865 ) 3866 { 3867 assert(expr != NULL); 3868 3869 return expr->children; 3870 } 3871 3872 /** gets the expression handler of an expression 3873 * 3874 * This identifies the type of the expression (sum, variable, ...). 3875 */ 3876 SCIP_EXPRHDLR* SCIPexprGetHdlr( 3877 SCIP_EXPR* expr /**< expression */ 3878 ) 3879 { 3880 assert(expr != NULL); 3881 3882 return expr->exprhdlr; 3883 } 3884 3885 /** gets the expression data of an expression */ 3886 SCIP_EXPRDATA* SCIPexprGetData( 3887 SCIP_EXPR* expr /**< expression */ 3888 ) 3889 { 3890 assert(expr != NULL); 3891 3892 return expr->exprdata; 3893 } 3894 3895 /** sets the expression data of an expression 3896 * 3897 * The pointer to possible old data is overwritten and the 3898 * freedata-callback is not called before. 3899 * This function is intended to be used by expression handler only. 3900 */ 3901 void SCIPexprSetData( 3902 SCIP_EXPR* expr, /**< expression */ 3903 SCIP_EXPRDATA* exprdata /**< expression data to be set (can be NULL) */ 3904 ) 3905 { 3906 assert(expr != NULL); 3907 assert(exprdata == NULL || expr->exprhdlr->copydata != NULL); /* copydata must be available if there is expression data */ 3908 assert(exprdata == NULL || expr->exprhdlr->freedata != NULL); /* freedata must be available if there is expression data */ 3909 3910 expr->exprdata = exprdata; 3911 } 3912 3913 /** gets the data that the owner of an expression has stored in an expression */ 3914 SCIP_EXPR_OWNERDATA* SCIPexprGetOwnerData( 3915 SCIP_EXPR* expr /**< expression */ 3916 ) 3917 { 3918 assert(expr != NULL); 3919 3920 return expr->ownerdata; 3921 } 3922 3923 /** gives the value from the last evaluation of an expression (or SCIP_INVALID if there was an eval error) 3924 * 3925 * @see SCIPevalExpr to evaluate the expression at a given solution. 3926 */ 3927 SCIP_Real SCIPexprGetEvalValue( 3928 SCIP_EXPR* expr /**< expression */ 3929 ) 3930 { 3931 assert(expr != NULL); 3932 3933 return expr->evalvalue; 3934 } 3935 3936 /** gives the evaluation tag from the last evaluation, or 0 3937 * 3938 * @see SCIPevalExpr 3939 */ 3940 SCIP_Longint SCIPexprGetEvalTag( 3941 SCIP_EXPR* expr /**< expression */ 3942 ) 3943 { 3944 assert(expr != NULL); 3945 3946 return expr->evaltag; 3947 } 3948 3949 /** returns the derivative stored in an expression (or SCIP_INVALID if there was an evaluation error) 3950 * 3951 * @see SCIPevalExprGradient 3952 */ 3953 SCIP_Real SCIPexprGetDerivative( 3954 SCIP_EXPR* expr /**< expression */ 3955 ) 3956 { 3957 assert(expr != NULL); 3958 3959 return expr->derivative; 3960 } 3961 3962 /** gives the value of directional derivative from the last evaluation of a directional derivative of 3963 * expression (or SCIP_INVALID if there was an error) 3964 * 3965 * @see SCIPevalExprHessianDir 3966 */ 3967 SCIP_Real SCIPexprGetDot( 3968 SCIP_EXPR* expr /**< expression */ 3969 ) 3970 { 3971 assert(expr != NULL); 3972 3973 return expr->dot; 3974 } 3975 3976 /** gives the value of directional derivative from the last evaluation of a directional derivative of 3977 * derivative of root (or SCIP_INVALID if there was an error) 3978 * 3979 * @see SCIPevalExprHessianDir 3980 */ 3981 SCIP_Real SCIPexprGetBardot( 3982 SCIP_EXPR* expr /**< expression */ 3983 ) 3984 { 3985 assert(expr != NULL); 3986 3987 return expr->bardot; 3988 } 3989 3990 /** returns the difftag stored in an expression 3991 * 3992 * can be used to check whether partial derivative value is valid 3993 * 3994 * @see SCIPevalExprGradient 3995 */ 3996 SCIP_Longint SCIPexprGetDiffTag( 3997 SCIP_EXPR* expr /**< expression */ 3998 ) 3999 { 4000 assert(expr != NULL); 4001 4002 return expr->difftag; 4003 } 4004 4005 /** returns the activity that is currently stored for an expression 4006 * 4007 * @see SCIPevalExprActivity 4008 */ 4009 SCIP_INTERVAL SCIPexprGetActivity( 4010 SCIP_EXPR* expr /**< expression */ 4011 ) 4012 { 4013 assert(expr != NULL); 4014 4015 return expr->activity; 4016 } 4017 4018 /** returns the tag associated with the activity of the expression 4019 * 4020 * It can depend on the owner of the expression how to interpret this tag. 4021 * SCIPevalExprActivity() compares with `stat->domchgcount`. 4022 * 4023 * @see SCIPevalExprActivity 4024 */ 4025 SCIP_Longint SCIPexprGetActivityTag( 4026 SCIP_EXPR* expr /**< expression */ 4027 ) 4028 { 4029 assert(expr != NULL); 4030 4031 return expr->activitytag; 4032 } 4033 4034 /** set the activity with tag for an expression */ 4035 void SCIPexprSetActivity( 4036 SCIP_EXPR* expr, /**< expression */ 4037 SCIP_INTERVAL activity, /**< new activity */ 4038 SCIP_Longint activitytag /**< tag associated with activity */ 4039 ) 4040 { 4041 assert(expr != NULL); 4042 4043 expr->activity = activity; 4044 expr->activitytag = activitytag; 4045 } 4046 4047 /** returns the curvature of an expression 4048 * 4049 * @note Call SCIPcomputeExprCurvature() before calling this function. 4050 */ 4051 SCIP_EXPRCURV SCIPexprGetCurvature( 4052 SCIP_EXPR* expr /**< expression */ 4053 ) 4054 { 4055 assert(expr != NULL); 4056 4057 return expr->curvature; 4058 } 4059 4060 /** sets the curvature of an expression */ 4061 void SCIPexprSetCurvature( 4062 SCIP_EXPR* expr, /**< expression */ 4063 SCIP_EXPRCURV curvature /**< curvature of the expression */ 4064 ) 4065 { 4066 assert(expr != NULL); 4067 4068 expr->curvature = curvature; 4069 } 4070 4071 /** returns whether an expression is integral */ 4072 SCIP_Bool SCIPexprIsIntegral( 4073 SCIP_EXPR* expr /**< expression */ 4074 ) 4075 { 4076 assert(expr != NULL); 4077 4078 return expr->isintegral; 4079 } 4080 4081 /** sets the integrality flag of an expression */ 4082 void SCIPexprSetIntegrality( 4083 SCIP_EXPR* expr, /**< expression */ 4084 SCIP_Bool isintegral /**< integrality of the expression */ 4085 ) 4086 { 4087 assert(expr != NULL); 4088 4089 expr->isintegral = isintegral; 4090 } 4091 4092 /** gives the coefficients and expressions that define a quadratic expression 4093 * 4094 * It can return the constant part, the number, arguments, and coefficients of the purely linear part 4095 * and the number of quadratic terms and bilinear terms. 4096 * Note that for arguments that appear in the quadratic part, a linear coefficient is 4097 * stored with the quadratic term. 4098 * Use SCIPexprGetQuadraticQuadTerm() and SCIPexprGetQuadraticBilinTerm() 4099 * to access the data for a quadratic or bilinear term. 4100 * 4101 * It can also return the eigenvalues and the eigenvectors of the matrix \f$Q\f$ when the quadratic is written 4102 * as \f$x^T Q x + b^T x + c^T y + d\f$, where \f$c^T y\f$ defines the purely linear part. 4103 * Note, however, that to have access to them one needs to call SCIPcomputeExprQuadraticCurvature() 4104 * with `storeeigeninfo=TRUE`. If the eigen information was not stored or it failed to be computed, 4105 * `eigenvalues` and `eigenvectors` will be set to NULL. 4106 * 4107 * This function returns pointers to internal data in linexprs and lincoefs. 4108 * The user must not change this data. 4109 * 4110 * @attention SCIPcheckExprQuadratic() needs to be called first to check whether expression is quadratic and initialize the data of the quadratic representation. 4111 */ 4112 void SCIPexprGetQuadraticData( 4113 SCIP_EXPR* expr, /**< quadratic expression */ 4114 SCIP_Real* constant, /**< buffer to store constant term, or NULL */ 4115 int* nlinexprs, /**< buffer to store number of expressions that appear linearly, or NULL */ 4116 SCIP_EXPR*** linexprs, /**< buffer to store pointer to array of expressions that appear linearly, or NULL */ 4117 SCIP_Real** lincoefs, /**< buffer to store pointer to array of coefficients of expressions that appear linearly, or NULL */ 4118 int* nquadexprs, /**< buffer to store number of expressions in quadratic terms, or NULL */ 4119 int* nbilinexprs, /**< buffer to store number of bilinear expressions terms, or NULL */ 4120 SCIP_Real** eigenvalues, /**< buffer to store pointer to array of eigenvalues of Q, or NULL */ 4121 SCIP_Real** eigenvectors /**< buffer to store pointer to array of eigenvectors of Q, or NULL */ 4122 ) 4123 { 4124 SCIP_QUADEXPR* quaddata; 4125 4126 assert(expr != NULL); 4127 4128 quaddata = expr->quaddata; 4129 assert(quaddata != NULL); 4130 4131 if( constant != NULL ) 4132 *constant = quaddata->constant; 4133 if( nlinexprs != NULL ) 4134 *nlinexprs = quaddata->nlinexprs; 4135 if( linexprs != NULL ) 4136 *linexprs = quaddata->linexprs; 4137 if( lincoefs != NULL ) 4138 *lincoefs = quaddata->lincoefs; 4139 if( nquadexprs != NULL ) 4140 *nquadexprs = quaddata->nquadexprs; 4141 if( nbilinexprs != NULL ) 4142 *nbilinexprs = quaddata->nbilinexprterms; 4143 if( eigenvalues != NULL ) 4144 *eigenvalues = quaddata->eigenvalues; 4145 if( eigenvectors != NULL ) 4146 *eigenvectors = quaddata->eigenvectors; 4147 } 4148 4149 /** gives the data of a quadratic expression term 4150 * 4151 * For a term \f$a \cdot \text{expr}^2 + b \cdot \text{expr} + \sum_i (c_i \cdot \text{expr} \cdot \text{otherexpr}_i)\f$, returns 4152 * `expr`, \f$a\f$, \f$b\f$, the number of summands, and indices of bilinear terms in the quadratic expressions `bilinexprterms`. 4153 * 4154 * This function returns pointers to internal data in adjbilin. 4155 * The user must not change this data. 4156 */ 4157 void SCIPexprGetQuadraticQuadTerm( 4158 SCIP_EXPR* quadexpr, /**< quadratic expression */ 4159 int termidx, /**< index of quadratic term */ 4160 SCIP_EXPR** expr, /**< buffer to store pointer to argument expression (the 'x') of this term, or NULL */ 4161 SCIP_Real* lincoef, /**< buffer to store linear coefficient of variable, or NULL */ 4162 SCIP_Real* sqrcoef, /**< buffer to store square coefficient of variable, or NULL */ 4163 int* nadjbilin, /**< buffer to store number of bilinear terms this variable is involved in, or NULL */ 4164 int** adjbilin, /**< buffer to store pointer to indices of associated bilinear terms, or NULL */ 4165 SCIP_EXPR** sqrexpr /**< buffer to store pointer to square expression (the 'x^2') of this term or NULL if no square expression, or NULL */ 4166 ) 4167 { 4168 SCIP_QUADEXPR_QUADTERM* quadexprterm; 4169 4170 assert(quadexpr != NULL); 4171 assert(quadexpr->quaddata != NULL); 4172 assert(quadexpr->quaddata->quadexprterms != NULL); 4173 assert(termidx >= 0); 4174 assert(termidx < quadexpr->quaddata->nquadexprs); 4175 4176 quadexprterm = &quadexpr->quaddata->quadexprterms[termidx]; 4177 4178 if( expr != NULL ) 4179 *expr = quadexprterm->expr; 4180 if( lincoef != NULL ) 4181 *lincoef = quadexprterm->lincoef; 4182 if( sqrcoef != NULL ) 4183 *sqrcoef = quadexprterm->sqrcoef; 4184 if( nadjbilin != NULL ) 4185 *nadjbilin = quadexprterm->nadjbilin; 4186 if( adjbilin != NULL ) 4187 *adjbilin = quadexprterm->adjbilin; 4188 if( sqrexpr != NULL ) 4189 *sqrexpr = quadexprterm->sqrexpr; 4190 } 4191 4192 /** gives the data of a bilinear expression term 4193 * 4194 * For a term a*expr1*expr2, returns expr1, expr2, a, and 4195 * the position of the quadratic expression term that uses expr2 in the quadratic expressions `quadexprterms`. 4196 */ 4197 void SCIPexprGetQuadraticBilinTerm( 4198 SCIP_EXPR* expr, /**< quadratic expression */ 4199 int termidx, /**< index of bilinear term */ 4200 SCIP_EXPR** expr1, /**< buffer to store first factor, or NULL */ 4201 SCIP_EXPR** expr2, /**< buffer to store second factor, or NULL */ 4202 SCIP_Real* coef, /**< buffer to coefficient, or NULL */ 4203 int* pos2, /**< buffer to position of expr2 in quadexprterms array of quadratic expression, or NULL */ 4204 SCIP_EXPR** prodexpr /**< buffer to store pointer to expression that is product if first and second factor, or NULL */ 4205 ) 4206 { 4207 SCIP_QUADEXPR_BILINTERM* bilinexprterm; 4208 4209 assert(expr != NULL); 4210 assert(expr->quaddata != NULL); 4211 assert(expr->quaddata->bilinexprterms != NULL); 4212 assert(termidx >= 0); 4213 assert(termidx < expr->quaddata->nbilinexprterms); 4214 4215 bilinexprterm = &expr->quaddata->bilinexprterms[termidx]; 4216 4217 if( expr1 != NULL ) 4218 *expr1 = bilinexprterm->expr1; 4219 if( expr2 != NULL ) 4220 *expr2 = bilinexprterm->expr2; 4221 if( coef != NULL ) 4222 *coef = bilinexprterm->coef; 4223 if( pos2 != NULL ) 4224 *pos2 = bilinexprterm->pos2; 4225 if( prodexpr != NULL ) 4226 *prodexpr = bilinexprterm->prodexpr; 4227 } 4228 4229 /** returns whether all expressions that are used in a quadratic expression are variable expressions 4230 * 4231 * @return TRUE iff all `linexprs` and `quadexprterms[.].expr` are variable expressions 4232 */ 4233 SCIP_Bool SCIPexprAreQuadraticExprsVariables( 4234 SCIP_EXPR* expr /**< quadratic expression */ 4235 ) 4236 { 4237 assert(expr != NULL); 4238 assert(expr->quaddata != NULL); 4239 4240 return expr->quaddata->allexprsarevars; 4241 } 4242 4243 /** returns a monomial representation of a product expression 4244 * 4245 * The array to store all factor expressions needs to be of size the number of 4246 * children in the expression which is given by SCIPexprGetNChildren(). 4247 * 4248 * Given a non-trivial monomial expression, the function finds its representation as \f$cx^\alpha\f$, where 4249 * \f$c\f$ is a real coefficient, \f$x\f$ is a vector of auxiliary or original variables (where some entries can 4250 * be NULL is the auxiliary variable has not been created yet), and \f$\alpha\f$ is a real vector of exponents. 4251 * 4252 * A non-trivial monomial is a product of a least two expressions. 4253 */ 4254 SCIP_RETCODE SCIPexprGetMonomialData( 4255 SCIP_SET* set, /**< global SCIP settings */ 4256 BMS_BLKMEM* blkmem, /**< block memory */ 4257 SCIP_EXPR* expr, /**< expression */ 4258 SCIP_Real* coef, /**< coefficient \f$c\f$ */ 4259 SCIP_Real* exponents, /**< exponents \f$\alpha\f$ */ 4260 SCIP_EXPR** factors /**< factor expressions \f$x\f$ */ 4261 ) 4262 { 4263 SCIP_EXPR* child; 4264 int c; 4265 int nexprs; 4266 4267 assert(set != NULL); 4268 assert(blkmem != NULL); 4269 assert(expr != NULL); 4270 assert(coef != NULL); 4271 assert(exponents != NULL); 4272 assert(factors != NULL); 4273 4274 assert(SCIPexprIsProduct(set, expr)); 4275 4276 *coef = SCIPgetCoefExprProduct(expr); 4277 nexprs = SCIPexprGetNChildren(expr); 4278 4279 for( c = 0; c < nexprs; ++c ) 4280 { 4281 child = SCIPexprGetChildren(expr)[c]; 4282 4283 if( SCIPexprIsPower(set, child) ) 4284 { 4285 exponents[c] = SCIPgetExponentExprPow(child); 4286 factors[c] = SCIPexprGetChildren(child)[0]; 4287 } 4288 else 4289 { 4290 exponents[c] = 1.0; 4291 factors[c] = child; 4292 } 4293 } 4294 4295 return SCIP_OKAY; 4296 } 4297 4298 /**@} */ 4299