- Timestamp:
- Aug 28, 2007 9:22:47 AM (14 years ago)
- Location:
- trunk/Cbc/src
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/Cbc/src/CbcLinked.cpp
r719 r771 276 276 void OsiSolverLink::resolve() 277 277 { 278 if (false) { 279 bool takeHint; 280 OsiHintStrength strength; 281 // Switch off printing if asked to 282 bool gotHint = (getHintParam(OsiDoReducePrint,takeHint,strength)); 283 assert (gotHint); 284 if (strength!=OsiHintIgnore&&takeHint) { 285 printf("no printing\n"); 286 } else { 287 printf("printing\n"); 288 } 289 } 278 290 specialOptions_ =0; 279 291 modelPtr_->setWhatsChanged(0); … … 680 692 upper[i]=upper2[i]; 681 693 } 694 qpTemp.setLogLevel(modelPtr_->logLevel()); 682 695 qpTemp.primal(); 683 696 assert (!qpTemp.problemStatus()); … … 1214 1227 element[n++]=-1.0; 1215 1228 double offset = - coinModel.objectiveOffset(); 1216 assert (!offset); // get sign right if happens 1229 //assert (!offset); // get sign right if happens 1230 printf("***** offset %g\n",offset); 1217 1231 coinModel.setObjectiveOffset(0.0); 1218 1232 double lowerBound = -COIN_DBL_MAX; -
trunk/Cbc/src/CbcSolver.cpp
r768 r771 134 134 #include "CbcCompareActual.hpp" 135 135 #include "CbcBranchActual.hpp" 136 #include "CbcBranchLotsize.hpp" 136 137 #include "CbcOrClpParam.hpp" 137 138 #include "CbcCutGenerator.hpp" … … 496 497 4 do BAB and just set best solution 497 498 -2 cleanup afterwards if using 2 498 On ou put - number fixed499 On output - number fixed 499 500 */ 500 501 static OsiClpSolverInterface * fixVubs(CbcModel & model, int numberLevels, 501 502 int & doAction, CoinMessageHandler * generalMessageHandler) 502 503 { 504 double time1 = CoinCpuTime(); 503 505 OsiSolverInterface * solver = model.solver()->clone(); 504 506 OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver); … … 507 509 lpSolver->tightenPrimalBounds(0.0,11,true); 508 510 char generalPrint[200]; 509 //const double *objective = solver->getObjCoefficients() ;511 const double *objective = lpSolver->getObjCoefficients() ; 510 512 double *columnLower = lpSolver->columnLower() ; 511 513 double *columnUpper = lpSolver->columnUpper() ; … … 548 550 double tolerance = lpSolver->primalTolerance(); 549 551 int * check = new int[numberRows]; 550 for (iRow=0;iRow<numberRows;iRow++) 551 check[iRow]=-1; 552 for (iRow=0;iRow<numberRows;iRow++) { 553 check[iRow]=-2; // don't check 554 if (rowLower[iRow]<-1.0e6&&rowUpper[iRow]>1.0e6) 555 continue;// unlikely 556 // possible row 557 int numberPositive=0; 558 int iPositive=-1; 559 int numberNegative=0; 560 int iNegative=-1; 561 CoinBigIndex rStart = rowStart[iRow]; 562 CoinBigIndex rEnd = rowStart[iRow]+rowLength[iRow]; 563 CoinBigIndex j; 564 int kColumn; 565 for (j = rStart; j < rEnd; ++j) { 566 double value=elementByRow[j]; 567 kColumn = column[j]; 568 if (columnUpper[kColumn]>columnLower[kColumn]) { 569 if (value > 0.0) { 570 numberPositive++; 571 iPositive=kColumn; 572 } else { 573 numberNegative++; 574 iNegative=kColumn; 575 } 576 } 577 } 578 if (numberPositive==1&&numberNegative==1) 579 check[iRow]=-1; // try both 580 if (numberPositive==1&&rowLower[iRow]>-1.0e20) 581 check[iRow]=iPositive; 582 else if (numberNegative==1&&rowUpper[iRow]<1.0e20) 583 check[iRow]=iNegative; 584 } 552 585 for (iColumn=0;iColumn<numberColumns;iColumn++) { 553 586 fix[iColumn]=-1; … … 564 597 i<columnStart[iColumn]+columnLength[iColumn];i++) { 565 598 iRow = row[i]; 566 if ( rowLower[iRow]<-1.0e6&&rowUpper[iRow]>1.0e6)599 if (check[iRow]!=-1&&check[iRow]!=iColumn) 567 600 continue; // unlikely 568 //==569 601 // possible row 570 602 int infiniteUpper = 0; … … 801 833 fixColumn2[iColumn+1]=numberOther; 802 834 } 835 // Create other way and convert fixColumn to counts 803 836 for ( iColumn=0;iColumn<numberColumns;iColumn++) { 804 837 for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { … … 813 846 while (true) { 814 847 int numberLayered=0; 815 int numberInteger=0;816 848 for ( iColumn=0;iColumn<numberColumns;iColumn++) { 817 849 if (fix[iColumn]==kLayer) { … … 819 851 int jColumn=otherColumn2[i]; 820 852 if (fix[jColumn]==kLayer) { 821 fix[iColumn]=kLayer+1 ;853 fix[iColumn]=kLayer+100; 822 854 } 823 855 } … … 825 857 if (fix[iColumn]==kLayer) { 826 858 numberLayered++; 827 if (solver->isInteger(iColumn))828 numberInteger++;829 859 } 830 860 } 831 861 if (numberLayered) { 832 printf("%d (%d integer) at priority %d\n",numberLayered,numberInteger,kLayer); 833 kLayer++; 862 kLayer+=100; 834 863 } else { 835 864 break; 836 865 } 837 866 } 867 for (int iPass=0;iPass<2;iPass++) { 868 for (int jLayer=0;jLayer<kLayer;jLayer++) { 869 int check[]={-1,0,1,2,3,4,5,10,50,100,500,1000,5000,10000,INT_MAX}; 870 int nCheck = (int) (sizeof(check)/sizeof(int)); 871 int countsI[20]; 872 int countsC[20]; 873 assert (nCheck<=20); 874 memset(countsI,0,nCheck*sizeof(int)); 875 memset(countsC,0,nCheck*sizeof(int)); 876 check[nCheck-1]=numberColumns; 877 int numberLayered=0; 878 int numberInteger=0; 879 for ( iColumn=0;iColumn<numberColumns;iColumn++) { 880 if (fix[iColumn]==jLayer) { 881 numberLayered++; 882 int nFix = fixColumn[iColumn+1]-fixColumn[iColumn]; 883 if (iPass) { 884 // just integers 885 nFix=0; 886 for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { 887 int jColumn=otherColumn[i]; 888 if (solver->isInteger(jColumn)) 889 nFix++; 890 } 891 } 892 int iFix; 893 for (iFix=0;iFix<nCheck;iFix++) { 894 if (nFix<=check[iFix]) 895 break; 896 } 897 assert (iFix<nCheck); 898 if (solver->isInteger(iColumn)) { 899 numberInteger++; 900 countsI[iFix]++; 901 } else { 902 countsC[iFix]++; 903 } 904 } 905 } 906 if (numberLayered) { 907 printf("%d (%d integer) at priority %d\n",numberLayered,numberInteger,1+(jLayer/100)); 908 char buffer[50]; 909 for (int i=1;i<nCheck;i++) { 910 if (countsI[i]||countsC[i]) { 911 if (i==1) 912 sprintf(buffer," == zero "); 913 else if (i<nCheck-1) 914 sprintf(buffer,"> %6d and <= %6d ",check[i-1],check[i]); 915 else 916 sprintf(buffer,"> %6d ",check[i-1]); 917 printf("%s %8d integers and %8d continuous\n",buffer,countsI[i],countsC[i]); 918 } 919 } 920 } 921 } 922 } 838 923 delete [] counts; 839 delete [] fixColumn;840 delete [] otherColumn;841 924 delete [] otherColumn2; 842 925 delete [] fixColumn2; 843 926 // Now do fixing 927 { 928 // switch off presolve and up weight 929 ClpSolve solveOptions; 930 //solveOptions.setPresolveType(ClpSolve::presolveOff,0); 931 solveOptions.setSolveType(ClpSolve::usePrimalorSprint); 932 //solveOptions.setSpecialOption(1,3,30); // sprint 933 int numberColumns = lpSolver->numberColumns(); 934 int iColumn; 935 bool allSlack=true; 936 for (iColumn=0;iColumn<numberColumns;iColumn++) { 937 if (lpSolver->getColumnStatus(iColumn)==ClpSimplex::basic) { 938 allSlack=false; 939 break; 940 } 941 } 942 if (allSlack) 943 solveOptions.setSpecialOption(1,2,50); // idiot 944 lpSolver->setInfeasibilityCost(1.0e11); 945 lpSolver->defaultFactorizationFrequency(); 946 lpSolver->initialSolve(solveOptions); 947 double * columnLower = lpSolver->columnLower(); 948 double * columnUpper = lpSolver->columnUpper(); 949 double * fullSolution = lpSolver->primalColumnSolution(); 950 const double * dj = lpSolver->dualColumnSolution(); 951 // First squash small and fix big 952 //double small=1.0e-7; 953 //double djTol=20.0; 954 int iPass=0; 955 #define MAXPROB 3 956 ClpSimplex models[MAXPROB]; 957 int pass[MAXPROB]; 958 int kPass=-1; 959 int kLayer=0; 960 while (true) { 961 double largest=-0.1; 962 double smallest=1.1; 963 int iLargest=-1; 964 int iSmallest=-1; 965 int atZero=0; 966 int atOne=0; 967 int toZero=0; 968 int toOne=0; 969 int numberFree=0; 970 int numberGreater=0; 971 columnLower = lpSolver->columnLower(); 972 columnUpper = lpSolver->columnUpper(); 973 fullSolution = lpSolver->primalColumnSolution(); 974 for ( iColumn=0;iColumn<numberColumns;iColumn++) { 975 if (!solver->isInteger(iColumn)||fix[iColumn]>kLayer) 976 continue; 977 double value = fullSolution[iColumn]; 978 if (value>1.00001) { 979 numberGreater++; 980 continue; 981 } 982 double lower = columnLower[iColumn]; 983 double upper = columnUpper[iColumn]; 984 if (lower==upper) { 985 if (lower) 986 atOne++; 987 else 988 atZero++; 989 continue; 990 } 991 if (value<1.0e-7) { 992 toZero++; 993 columnUpper[iColumn]=0.0; 994 continue; 995 } 996 if (value>1.0-1.0e-7) { 997 toOne++; 998 columnLower[iColumn]=1.0; 999 continue; 1000 } 1001 numberFree++; 1002 if (value<smallest) { 1003 smallest=value; 1004 iSmallest=iColumn; 1005 } 1006 if (value>largest) { 1007 largest=value; 1008 iLargest=iColumn; 1009 } 1010 } 1011 if (toZero||toOne) 1012 printf("%d at 0 fixed and %d at one fixed\n",toZero,toOne); 1013 printf("%d variables free, %d fixed to 0, %d to 1 - smallest %g, largest %g\n", 1014 numberFree,atZero,atOne,smallest,largest); 1015 if (numberGreater&&!iPass) 1016 printf("%d variables have value > 1.0\n",numberGreater); 1017 int jLayer=0; 1018 int nFixed=-1; 1019 int nTotalFixed=0; 1020 while (nFixed) { 1021 nFixed=0; 1022 for ( iColumn=0;iColumn<numberColumns;iColumn++) { 1023 if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) { 1024 for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { 1025 int jColumn=otherColumn[i]; 1026 columnUpper[jColumn]=0.0; 1027 nFixed++; 1028 } 1029 } 1030 } 1031 nTotalFixed += nFixed; 1032 jLayer++; 1033 } 1034 printf("This fixes %d variables in lower priorities\n",nTotalFixed); 1035 if (iLargest<0) 1036 break; 1037 double movement; 1038 int way; 1039 if (smallest<=1.0-largest&&smallest<0.2) { 1040 columnUpper[iSmallest]=0.0; 1041 movement=smallest; 1042 way=-1; 1043 } else { 1044 columnLower[iLargest]=1.0; 1045 movement=1.0-largest; 1046 way=1; 1047 } 1048 double saveObj = lpSolver->objectiveValue(); 1049 iPass++; 1050 kPass = iPass%MAXPROB; 1051 models[kPass]=*lpSolver; 1052 pass[kPass]=iPass; 1053 double maxCostUp = COIN_DBL_MAX; 1054 if (way==-1) 1055 maxCostUp= (1.0-movement)*objective[iSmallest]; 1056 lpSolver->setDualObjectiveLimit(saveObj+maxCostUp); 1057 lpSolver->dual(); 1058 double moveObj = lpSolver->objectiveValue()-saveObj; 1059 printf("movement %s was %g costing %g\n", 1060 (smallest<=1.0-largest) ? "down" : "up",movement,moveObj); 1061 if (way==-1&&(moveObj>=(1.0-movement)*objective[iSmallest]||lpSolver->status())) { 1062 // go up 1063 columnLower = models[kPass].columnLower(); 1064 columnUpper = models[kPass].columnUpper(); 1065 columnLower[iSmallest]=1.0; 1066 columnUpper[iSmallest]=lpSolver->columnUpper()[iSmallest]; 1067 *lpSolver=models[kPass]; 1068 columnLower = lpSolver->columnLower(); 1069 columnUpper = lpSolver->columnUpper(); 1070 fullSolution = lpSolver->primalColumnSolution(); 1071 dj = lpSolver->dualColumnSolution(); 1072 columnLower[iSmallest]=1.0; 1073 columnUpper[iSmallest]=lpSolver->columnUpper()[iSmallest]; 1074 lpSolver->dual(); 1075 } 1076 models[kPass]=*lpSolver; 1077 } 1078 } 1079 printf("Fixing took %g seconds\n",CoinCpuTime()-time1); 844 1080 delete [] fix; 1081 delete [] fixColumn; 1082 delete [] otherColumn; 845 1083 delete solver; 846 1084 return NULL; … … 2185 2423 // need some relative granularity 2186 2424 si->setDefaultBound(100.0); 2187 si->setDefaultMeshSize(0.001); 2425 double dextra3 = parameters[whichParam(DEXTRA3,numberParameters,parameters)].doubleValue(); 2426 if (dextra3) 2427 si->setDefaultMeshSize(dextra3); 2188 2428 si->setDefaultBound(100000.0); 2189 2429 si->setIntegerPriority(1000); … … 3350 3590 #endif 3351 3591 } else { 3592 #ifndef DISALLOW_PRINTING 3352 3593 std::cout<<"** Current model not valid"<<std::endl; 3594 #endif 3353 3595 } 3354 3596 break; … … 3386 3628 delete model2; 3387 3629 } else { 3630 #ifndef DISALLOW_PRINTING 3388 3631 std::cout<<"** Current model not valid"<<std::endl; 3632 #endif 3389 3633 } 3390 3634 break; … … 3395 3639 std::cout<<"** Analysis indicates model infeasible"<<std::endl; 3396 3640 } else { 3641 #ifndef DISALLOW_PRINTING 3397 3642 std::cout<<"** Current model not valid"<<std::endl; 3643 #endif 3398 3644 } 3399 3645 break; … … 3416 3662 } 3417 3663 } else { 3664 #ifndef DISALLOW_PRINTING 3418 3665 std::cout<<"** Current model not valid"<<std::endl; 3666 #endif 3419 3667 } 3420 3668 break; … … 3432 3680 <<CoinMessageEol; 3433 3681 } else { 3682 #ifndef DISALLOW_PRINTING 3434 3683 std::cout<<"** Current model not valid"<<std::endl; 3684 #endif 3435 3685 } 3436 3686 break; … … 3453 3703 } 3454 3704 } else { 3705 #ifndef DISALLOW_PRINTING 3455 3706 std::cout<<"** Current model not valid"<<std::endl; 3707 #endif 3456 3708 } 3457 3709 break; … … 3482 3734 goodModel=true; 3483 3735 /* 3484 Run branch-and-cut. First set a few options -- node comparison, scaling. If 3485 the solver is Clp, consider running some presolve code (not yet converted 3486 this to generic OSI) with branch-and-cut. If presolve is disabled, or the 3487 solver is not Clp, simply run branch-and-cut. Print elapsed time at the end. 3736 Run branch-and-cut. First set a few options -- node comparison, scaling. 3737 Print elapsed time at the end. 3488 3738 */ 3489 3739 case BAB: // branchAndBound … … 3522 3772 // need some relative granularity 3523 3773 si->setDefaultBound(100.0); 3524 si->setDefaultMeshSize(0.001); 3774 double dextra3 = parameters[whichParam(DEXTRA3,numberParameters,parameters)].doubleValue(); 3775 if (dextra3) 3776 si->setDefaultMeshSize(dextra3); 3525 3777 si->setDefaultBound(1000.0); 3526 3778 si->setIntegerPriority(1000); … … 3650 3902 cbcModel->initialSolve(); 3651 3903 if (clpModel->tightenPrimalBounds()!=0) { 3904 #ifndef DISALLOW_PRINTING 3652 3905 std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl; 3653 exit(1); 3906 #endif 3907 break; 3654 3908 } 3655 3909 clpModel->dual(); // clean up … … 3754 4008 } 3755 4009 if (!complicatedInteger&&clpSolver->tightenPrimalBounds()!=0) { 4010 #ifndef DISALLOW_PRINTING 3756 4011 std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl; 3757 exit(1); 4012 #endif 4013 model.setProblemStatus(0); 4014 model.setSecondaryStatus(1); 4015 // and in babModel if exists 4016 if (babModel) { 4017 babModel->setProblemStatus(0); 4018 babModel->setSecondaryStatus(1); 4019 } 4020 break; 3758 4021 } 3759 4022 if (clpSolver->dualBound()==1.0e10) { … … 3825 4088 if (tightenFactor&&!complicatedInteger) { 3826 4089 if (modelC->tightenPrimalBounds(tightenFactor)!=0) { 4090 #ifndef DISALLOW_PRINTING 3827 4091 std::cout<<"Problem is infeasible!"<<std::endl; 4092 #endif 4093 model.setProblemStatus(0); 4094 model.setSecondaryStatus(1); 4095 // and in babModel if exists 4096 if (babModel) { 4097 babModel->setProblemStatus(0); 4098 babModel->setSecondaryStatus(1); 4099 } 3828 4100 break; 3829 4101 } … … 3973 4245 prohibited[iColumn]=1; 3974 4246 } 4247 } 4248 CbcLotsize * obj2 = 4249 dynamic_cast <CbcLotsize *>(oldObjects[iObj]) ; 4250 if (obj2) { 4251 int iColumn = obj2->columnNumber(); 4252 prohibited[iColumn]=1; 3975 4253 } 3976 4254 } … … 4089 4367 //modelC->setLogLevel(0); 4090 4368 if (!complicatedInteger&&modelC->tightenPrimalBounds()!=0) { 4369 #ifndef DISALLOW_PRINTING 4091 4370 std::cout<<"Problem is infeasible!"<<std::endl; 4371 #endif 4372 model.setProblemStatus(0); 4373 model.setSecondaryStatus(1); 4374 // and in babModel if exists 4375 if (babModel) { 4376 babModel->setProblemStatus(0); 4377 babModel->setSecondaryStatus(1); 4378 } 4092 4379 break; 4093 4380 } … … 4982 5269 serendipity.setHeuristicName("linked"); 4983 5270 babModel->addHeuristic(&serendipity); 5271 double dextra3 = parameters[whichParam(DEXTRA3,numberParameters,parameters)].doubleValue(); 5272 if (dextra3) 5273 solver3->setMeshSizes(dextra3); 4984 5274 int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000; 4985 5275 CglStored stored; … … 5603 5893 //babModel=NULL; 5604 5894 } else { 5895 #ifndef DISALLOW_PRINTING 5605 5896 std::cout << "** Current model not valid" << std::endl ; 5897 #endif 5606 5898 } 5607 5899 break ; … … 5843 6135 // need some relative granularity 5844 6136 si->setDefaultBound(100.0); 5845 si->setDefaultMeshSize(0.001); 6137 double dextra3 = parameters[whichParam(DEXTRA3,numberParameters,parameters)].doubleValue(); 6138 if (dextra3) 6139 si->setDefaultMeshSize(dextra3); 5846 6140 si->setDefaultBound(100.0); 5847 6141 si->setIntegerPriority(1000); … … 5998 6292 } 5999 6293 } else { 6294 #ifndef DISALLOW_PRINTING 6000 6295 std::cout<<"** Current model not valid"<<std::endl; 6296 #endif 6001 6297 } 6002 6298 break; … … 6053 6349 } 6054 6350 } else { 6351 #ifndef DISALLOW_PRINTING 6055 6352 std::cout<<"** Current model not valid"<<std::endl; 6353 #endif 6056 6354 } 6057 6355 break; … … 6324 6622 } 6325 6623 } else { 6624 #ifndef DISALLOW_PRINTING 6326 6625 std::cout<<"** Current model not valid"<<std::endl; 6626 #endif 6327 6627 } 6328 6628 break; … … 6380 6680 } 6381 6681 } else { 6682 #ifndef DISALLOW_PRINTING 6382 6683 std::cout<<"** Current model not valid"<<std::endl; 6684 #endif 6383 6685 } 6384 6686 break; … … 6439 6741 } 6440 6742 } else { 6743 #ifndef DISALLOW_PRINTING 6441 6744 std::cout<<"** Current model not valid"<<std::endl; 6745 #endif 6442 6746 } 6443 6747 break; … … 7112 7416 } 7113 7417 } else { 7418 #ifndef DISALLOW_PRINTING 7114 7419 std::cout<<"** Current model not valid"<<std::endl; 7115 7420 #endif 7116 7421 } 7117 7422 break; … … 7145 7450 saveSolution(lpSolver,fileName); 7146 7451 } else { 7452 #ifndef DISALLOW_PRINTING 7147 7453 std::cout<<"** Current model not valid"<<std::endl; 7148 7454 #endif 7149 7455 } 7150 7456 break;
Note: See TracChangeset
for help on using the changeset viewer.