Changeset 1591


Ignore:
Timestamp:
Feb 7, 2011 7:48:46 AM (9 years ago)
Author:
forrest
Message:

modifications to heuristics and allow missing out some printout

Location:
trunk/Cbc
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/Cbc/examples/Makefile.in

    r1565 r1591  
    1919# CHANGEME: Here is the name of all object files corresponding to the source
    2020#           code that you wrote in order to define the problem statement
    21 OBJS =  $(DRIVER).@OBJEXT@ \
    22         CbcCompareUser.@OBJEXT@ \
    23         CbcSolver2.@OBJEXT@ \
    24         CbcSolver3.@OBJEXT@ \
    25         ClpQuadInterface.@OBJEXT@ \
    26         CbcBranchFollow2.@OBJEXT@ \
    27         CbcBranchUser.@OBJEXT@ \
    28         CbcBranchLink.@OBJEXT@
     21OBJS =  $(DRIVER).@OBJEXT@
     22
     23#       CbcCompareUser.@OBJEXT@ \
     24#       CbcSolver2.@OBJEXT@ \
     25#       CbcSolver3.@OBJEXT@ \
     26#       ClpQuadInterface.@OBJEXT@ \
     27#       CbcBranchFollow2.@OBJEXT@ \
     28#       CbcBranchUser.@OBJEXT@ \
     29#       CbcBranchLink.@OBJEXT@
    2930
    3031# CHANGEME: Additional libraries

  • trunk/Cbc/src/CbcHeuristicGreedy.cpp

    r1589 r1591  
    990990    const int * columnLength = matrix_.getVectorLengths();
    991991    int * sosRow = new int [numberColumns];
     992    int nonSOS=0;
    992993    // If bit set then use current
    993994    if ((algorithm_&1)!=0) {
    … …  
    10051006          // SOS
    10061007          rhs[iRow]=-1.0;
    1007         } else if (rowLower[iRow] > 0.0) {
     1008        } else if (rowLower[iRow] > 0.0 && rowUpper[iRow] < 1.0e10) {
    10081009          good = false;
    10091010        } else if (rowUpper[iRow] < 0.0) {
    10101011          good = false;
     1012        } else if (rowUpper[iRow] < 1.0e10) {
     1013          rhs[iRow]=rowUpper[iRow];
    10111014        } else {
    1012           rhs[iRow]=rowUpper[iRow];
     1015          rhs[iRow]=rowLower[iRow];
    10131016        }
    10141017      }
    10151018      for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
     1019        if (!columnLength[iColumn])
     1020          continue;
    10161021        if (columnLower[iColumn] < 0.0 || columnUpper[iColumn] > 1.0)
    10171022          good = false;
    … …  
    10191024        int nSOS=0;
    10201025        int iSOS=-1;
     1026        if (!solver->isInteger(iColumn))
     1027          good = false;
    10211028        for (j = columnStart[iColumn];
    10221029             j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    … …  
    10311038          }
    10321039        }
    1033         if (nSOS>1||!solver->isBinary(iColumn))
     1040        if (nSOS>1)
    10341041          good = false;
     1042        else if (!nSOS)
     1043          nonSOS++;
    10351044        sosRow[iColumn] = iSOS;
    10361045      }
    10371046      if (!good) {
     1047        delete [] sosRow;
    10381048        delete [] rhs;
    1039         delete [] sosRow;
    10401049        setWhen(0); // switch off
    10411050        return 0;
    … …  
    10461055    const double * solution = solver->getColSolution();
    10471056    const double * objective = solver->getObjCoefficients();
     1057    const double * rowLower = solver->getRowLower();
     1058    const double * rowUpper = solver->getRowUpper();
    10481059    double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);
    10491060    double primalTolerance;
    … …  
    10521063    numRuns_++;
    10531064    assert (numberRows == matrix_.getNumRows());
     1065    // set up linked list for sets
     1066    int * firstGub = new int [numberRows];
     1067    int * nextGub = new int [numberColumns];
    10541068    int iRow, iColumn;
    10551069    double direction = solver->getObjSense();
    10561070    double * slackCost = new double [numberRows];
    10571071    double * modifiedCost = CoinCopyOfArray(objective,numberColumns);
    1058     for (int iRow = 0;iRow < numberRows; iRow++)
     1072    for (int iRow = 0;iRow < numberRows; iRow++) {
    10591073      slackCost[iRow]=1.0e30;
     1074      firstGub[iRow]=-1;
     1075    }
    10601076    // Take off cost of gub slack
    10611077    for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
     1078      nextGub[iColumn]=-1;
    10621079      int iRow = sosRow[iColumn];
    10631080      if (columnLength[iColumn] == 1&&iRow>=0) {
    … …  
    10751092        sos[iRow]=1;
    10761093        rhs[iRow]=1.0;
     1094      } else if (rhs[iRow] != rowUpper[iRow]) {
     1095        // G row
     1096        sos[iRow]=-1;
    10771097      }
    10781098      if( slackCost[iRow] == 1.0e30) {
    … …  
    10891109           j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    10901110        int iRow = row[j];
    1091         if (sos[iRow]) {
     1111        if (sos[iRow]>0) {
    10921112          cost -= slackCost[iRow];
     1113          if (firstGub[iRow]<0) {
     1114            firstGub[iRow]=iColumn;
     1115          } else {
     1116            int jColumn = firstGub[iRow];
     1117            while (nextGub[jColumn]>=0)
     1118              jColumn=nextGub[jColumn];
     1119            nextGub[jColumn]=iColumn;
     1120          }
    10931121        }
    10941122      }
    … …  
    11051133    double * newSolution = new double [numberColumns];
    11061134    double * rowActivity = new double[numberRows];
    1107     memset(rowActivity, 0, numberRows*sizeof(double));
     1135    double * contribution = new double [numberColumns];
     1136    int * which = new int [numberColumns];
     1137    double * newSolution0 = new double [numberColumns];
    11081138    if ((algorithm_&(2|4))==0) {
    11091139      // get solution as small as possible
    11101140      for (iColumn = 0; iColumn < numberColumns; iColumn++)
    1111         newSolution[iColumn] = columnLower[iColumn];
     1141        newSolution0[iColumn] = columnLower[iColumn];
    11121142    } else {
    11131143      // Get rounded down solution
    11141144      for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1115         double value = solution[iColumn];
     1145        double value = solution[iColumn];
    11161146        // Round down integer
    11171147        if (fabs(floor(value + 0.5) - value) < integerTolerance) {
    … …  
    11201150          value = CoinMax(floor(value), columnLower[iColumn]);
    11211151        }
    1122         // make sure clean
    1123         value = CoinMin(value, columnUpper[iColumn]);
    1124         value = CoinMax(value, columnLower[iColumn]);
    1125         newSolution[iColumn] = value;
     1152        // make sure clean
     1153        value = CoinMin(value, columnUpper[iColumn]);
     1154        value = CoinMax(value, columnLower[iColumn]);
     1155        newSolution0[iColumn] = value;
    11261156      }
    11271157    }
    1128     // get row activity
    1129     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1130       CoinBigIndex j;
    1131       double value = newSolution[iColumn];
    1132       for (j = columnStart[iColumn];
    1133            j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    1134         int iRow = row[j];
    1135         rowActivity[iRow] += value * element[j];
     1158    double * rowWeight = new double [numberRows];
     1159    for (int i=0;i<numberRows;i++)
     1160      rowWeight[i]=1.0;
     1161    double costBias = 0.0;
     1162    int nPass = ((algorithm_&4)!=0) ? 1 : 10;
     1163    for (int iPass=0;iPass<nPass;iPass++) {
     1164      newSolutionValue = -offset+offset2;
     1165      memcpy(newSolution,newSolution0,numberColumns*sizeof(double));
     1166      // get row activity
     1167      memset(rowActivity, 0, numberRows*sizeof(double));
     1168      for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1169        CoinBigIndex j;
     1170        double value = newSolution[iColumn];
     1171        for (j = columnStart[iColumn];
     1172             j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1173          int iRow = row[j];
     1174          rowActivity[iRow] += value * element[j];
     1175        }
    11361176      }
    1137     }
    1138     double * contribution = new double [numberColumns];
    1139     int * which = new int [numberColumns];
    1140     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1141       CoinBigIndex j;
    1142       double value = newSolution[iColumn];
    1143       double cost =  modifiedCost[iColumn];
    1144       double forSort = 0.0;
    1145       bool hasSlack=false;
    1146       bool willFit=true;
    1147       newSolutionValue += value * cost;
    1148       for (j = columnStart[iColumn];
    1149            j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    1150         int iRow = row[j];
    1151         if (sos[iRow] == 2)
    1152           hasSlack = true;
    1153         forSort += element[j];
    1154         double gap = rhs[iRow] - rowActivity[iRow]+1.0e-8;
    1155         if (gap<element[j])
    1156           willFit = false;
     1177      if (!rowWeight) {
     1178        rowWeight = CoinCopyOfArray(rowActivity,numberRows);
    11571179      }
    1158       bool isSlack = hasSlack && (columnLength[iColumn]==1);
    1159       if ((algorithm_&4)!=0)
    1160         forSort=1.0;
    1161       // Use smallest cost if will fit
    1162       if (willFit /*&& hasSlack*/ &&
    1163           value == 0.0 && columnUpper[iColumn]) {
    1164         if (hasSlack) {
    1165           if (cost>0.0) {
    1166             forSort = 2.0e30;
    1167           } else if (cost==0.0) {
    1168             if (!isSlack)
    1169               forSort = 1.0e29;
     1180      for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1181        CoinBigIndex j;
     1182        double value = newSolution[iColumn];
     1183        double cost =  modifiedCost[iColumn];
     1184        double forSort = 1.0e-24;
     1185        bool hasSlack=false;
     1186        bool willFit=true;
     1187        bool gRow=false;
     1188        newSolutionValue += value * cost;
     1189        cost += 1.0e-12;
     1190        for (j = columnStart[iColumn];
     1191             j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1192          int iRow = row[j];
     1193          int type = sos[iRow];
     1194          double gap = rhs[iRow] - rowActivity[iRow]+1.0e-8;
     1195          switch (type) {
     1196          case -1:
     1197            // G row
     1198            gRow = true;
     1199#if 0
     1200            if (rhs[iRow]>rowWeight[iRow]||(algorithm_&(2|4))!=0)
     1201              forSort += element[j];
    11701202            else
    1171               forSort = 1.0e28;
     1203              forSort += 0.1*element[j];
     1204#else
     1205            forSort += rowWeight[iRow]*element[j];
     1206#endif
     1207            break;
     1208          case 0:
     1209            // L row
     1210            if (gap<element[j]) {
     1211              willFit = false;
     1212            } else {
     1213              forSort += element[j];
     1214            }
     1215            break;
     1216          case 1:
     1217            // SOS without slack
     1218            if (gap<element[j]) {
     1219              willFit = false;
     1220            }
     1221            break;
     1222          case 2:
     1223            // SOS with slack
     1224            hasSlack = true;
     1225            if (gap<element[j]) {
     1226              willFit = false;
     1227            }
     1228            break;
     1229          }
     1230        }
     1231        bool isSlack = hasSlack && (columnLength[iColumn]==1);
     1232        if (forSort<1.0e-24)
     1233          forSort = 1.0e-12;
     1234        if ((algorithm_&4)!=0 && forSort > 1.0e-24)
     1235          forSort=1.0;
     1236        // Use smallest cost if will fit
     1237        if (willFit && (hasSlack||gRow) &&
     1238            value == 0.0 && columnUpper[iColumn]) {
     1239          if (hasSlack && !gRow) {
     1240            if (cost>1.0e-12) {
     1241              forSort = 2.0e30;
     1242            } else if (cost==1.0e-12) {
     1243              if (!isSlack)
     1244                forSort = 1.0e29;
     1245              else
     1246                forSort = 1.0e28;
     1247            } else {
     1248              forSort = cost/forSort;
     1249            }
    11721250          } else {
    1173             forSort = cost/forSort;
     1251            if (!gRow||true)
     1252              forSort = (cost+costBias)/forSort;
     1253            else
     1254              forSort = 1.0e-12/forSort;
    11741255          }
    11751256        } else {
    1176           forSort = cost/forSort;
     1257          // put at end
     1258          forSort = 1.0e30;
    11771259        }
    1178       } else {
    1179         // put at end
    1180         forSort = 1.0e30;
     1260        which[iColumn]=iColumn;
     1261        contribution[iColumn]= forSort;
    11811262      }
    1182       which[iColumn]=iColumn;
    1183       contribution[iColumn]= forSort;
    1184     }
    1185     CoinSort_2(contribution,contribution+numberColumns,which);
    1186     // Go through columns
    1187     for (int jColumn = 0; jColumn < numberColumns; jColumn++) {
    1188       int iColumn = which[jColumn];
    1189       double value = newSolution[iColumn];
    1190       if (value)
    1191         continue;
    1192       bool possible = true;
    1193       CoinBigIndex j;
    1194       for (j = columnStart[iColumn];
    1195            j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    1196         int iRow = row[j];
    1197         if (sos[iRow]&&rowActivity[iRow]) {
    1198           possible = false;
    1199         } else {
    1200           double gap = rhs[iRow] - rowActivity[iRow]+1.0e-8;
    1201           if (gap<element[j])
     1263      CoinSort_2(contribution,contribution+numberColumns,which);
     1264      // Go through columns
     1265      int nAdded=0;
     1266      int nSlacks=0;
     1267      for (int jColumn = 0; jColumn < numberColumns; jColumn++) {
     1268        if (contribution[jColumn]>=1.0e30)
     1269          break;
     1270        int iColumn = which[jColumn];
     1271        double value = newSolution[iColumn];
     1272        if (value)
     1273          continue;
     1274        bool possible = true;
     1275        CoinBigIndex j;
     1276        for (j = columnStart[iColumn];
     1277             j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1278          int iRow = row[j];
     1279          if (sos[iRow]>0&&rowActivity[iRow]) {
    12021280            possible = false;
     1281          } else {
     1282            double gap = rhs[iRow] - rowActivity[iRow]+1.0e-8;
     1283            if (gap<element[j]&&sos[iRow]>=0)
     1284              possible = false;
     1285          }
     1286        }
     1287        if (possible) {
     1288          //#define REPORT 1
     1289#ifdef REPORT
     1290          if ((nAdded%1000)==0) {
     1291            double gap=0.0;
     1292            for (int i=0;i<numberRows;i++) {
     1293              if (rowUpper[i]>1.0e20)
     1294                gap += CoinMax(rowLower[i]-rowActivity[i],0.0);
     1295            }
     1296            if (gap)
     1297              printf("after %d added gap %g - %d slacks\n",
     1298                     nAdded,gap,nSlacks);
     1299          }
     1300#endif
     1301          nAdded++;
     1302          if (columnLength[iColumn]==1)
     1303            nSlacks++;
     1304          // Increase chosen column
     1305          newSolution[iColumn] = 1.0;
     1306          double cost =  modifiedCost[iColumn];
     1307          newSolutionValue += cost;
     1308          for (CoinBigIndex j = columnStart[iColumn];
     1309               j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1310            int iRow = row[j];
     1311            rowActivity[iRow] += element[j];
     1312          }
    12031313        }
    12041314      }
    1205       if (possible) {
    1206         // Increase chosen column
    1207         newSolution[iColumn] = 1.0;
    1208         double cost =  modifiedCost[iColumn];
    1209         newSolutionValue += cost;
    1210         for (CoinBigIndex j = columnStart[iColumn];
    1211              j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    1212           int iRow = row[j];
    1213           rowActivity[iRow] += element[j];
    1214         }
     1315#ifdef REPORT
     1316      {
     1317        double under=0.0;
     1318        double over=0.0;
     1319        double gap = 0.0;
     1320        int nUnder=0;
     1321        int nOver=0;
     1322        int nGap=0;
     1323        for (iRow = 0; iRow < numberRows; iRow++) {
     1324          if (rowActivity[iRow] < rowLower[iRow] - 10.0*primalTolerance) {
     1325            double value = rowLower[iRow]-rowActivity[iRow];
     1326#if REPORT>1
     1327            printf("below on %d is %g - activity %g lower %g\n",
     1328                   iRow,value,rowActivity[iRow],rowLower[iRow]);
     1329#endif
     1330            under += value;
     1331            nUnder++;
     1332          } else if (rowActivity[iRow] > rowUpper[iRow] + 10.0*primalTolerance) {
     1333            double value = rowActivity[iRow]-rowUpper[iRow];
     1334#if REPORT>1
     1335            printf("above on %d is %g - activity %g upper %g\n",
     1336                   iRow,value,rowActivity[iRow],rowUpper[iRow]);
     1337#endif
     1338            over += value;
     1339            nOver++;
     1340          } else {
     1341            double value = rowActivity[iRow]-rowLower[iRow];
     1342            if (value && value < 1.0e20) {
     1343#if REPORT>1
     1344              printf("gap on %d is %g - activity %g lower %g\n",
     1345                     iRow,value,rowActivity[iRow],rowLower[iRow]);
     1346#endif
     1347              gap += value;
     1348              nGap++;
     1349            }
     1350          }
     1351        }
     1352        printf("final under %g (%d) - over %g (%d) - free %g (%d) - %d added - solvalue %g\n",
     1353               under,nUnder,over,nOver,gap,nGap,nAdded,newSolutionValue);
    12151354      }
    1216     }
     1355#endif
     1356      double gap = 0.0;
     1357      double over = 0.0;
     1358      int nL=0;
     1359      int nG=0;
     1360      int iUnder=-1;
     1361      int nUnder=0;
     1362      for (iRow = 0; iRow < numberRows; iRow++) {
     1363        if (rowLower[iRow]<-1.0e20)
     1364          nL++;
     1365        if (rowUpper[iRow]>1.0e20)
     1366          nG++;
     1367        if (rowActivity[iRow] < rowLower[iRow] - 10.0*primalTolerance) {
     1368          gap += rowLower[iRow]-rowActivity[iRow];
     1369          nUnder++;
     1370          iUnder=iRow;
     1371          rowWeight[iRow] *= 1.1;
     1372        } else if (rowActivity[iRow] > rowUpper[iRow] + 10.0*primalTolerance) {
     1373          gap += rowActivity[iRow]-rowUpper[iRow];
     1374        } else {
     1375          over += rowActivity[iRow]-rowLower[iRow];
     1376          //rowWeight[iRow] *= 0.9;
     1377        }
     1378      }
     1379      if (nG&&!nL) {
     1380        // can we fix
     1381        // get list of columns which can go down without making
     1382        // things much worse
     1383        int nPossible=0;
     1384        int nEasyDown=0;
     1385        int nSlackDown=0;
     1386        for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     1387          if (newSolution[iColumn]&&
     1388              columnUpper[iColumn]>columnLower[iColumn]) {
     1389            bool canGoDown=true;
     1390            bool under = false;
     1391            int iSos=-1;
     1392            for (CoinBigIndex j = columnStart[iColumn];
     1393                 j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1394              int iRow = row[j];
     1395              if (sos[iRow]<0) {
     1396                double over = rowActivity[iRow]-rowLower[iRow];
     1397                if (over>=0.0&&element[j]>over+1.0e-12) {
     1398                  canGoDown=false;
     1399                  break;
     1400                } else if (over<0.0) {
     1401                  under = true;
     1402                }
     1403              } else {
     1404                iSos=iRow;
     1405              }
     1406            }
     1407            if (canGoDown) {
     1408              if (!under) {
     1409                if (iSos>=0) {
     1410                  // find cheapest
     1411                  double cheapest=modifiedCost[iColumn];
     1412                  int iCheapest = -1;
     1413                  int jColumn = firstGub[iSos];
     1414                  assert (jColumn>=0);
     1415                  while (jColumn>=0) {
     1416                    if (modifiedCost[jColumn]<cheapest) {
     1417                      cheapest=modifiedCost[jColumn];
     1418                      iCheapest=jColumn;
     1419                    }
     1420                    jColumn = nextGub[jColumn];
     1421                  }
     1422                  if (iCheapest>=0) {
     1423                    // Decrease column
     1424                    newSolution[iColumn] = 0.0;
     1425                    newSolutionValue -= modifiedCost[iColumn];
     1426                    for (CoinBigIndex j = columnStart[iColumn];
     1427                         j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1428                      int iRow = row[j];
     1429                      rowActivity[iRow] -= element[j];
     1430                    }
     1431                    // Increase chosen column
     1432                    newSolution[iCheapest] = 1.0;
     1433                    newSolutionValue += modifiedCost[iCheapest];
     1434                    for (CoinBigIndex j = columnStart[iCheapest];
     1435                         j < columnStart[iCheapest] + columnLength[iCheapest]; j++) {
     1436                      int iRow = row[j];
     1437                      rowActivity[iRow] += element[j];
     1438                    }
     1439                    nEasyDown++;
     1440                    if (columnLength[iColumn]>1) {
     1441                      //printf("%d is easy down\n",iColumn);
     1442                    } else {
     1443                      nSlackDown++;
     1444                    }
     1445                  }
     1446                } else if (modifiedCost[iColumn]>0.0) {
     1447                  // easy down
     1448                  // Decrease column
     1449                  newSolution[iColumn] = 0.0;
     1450                  newSolutionValue -= modifiedCost[iColumn];
     1451                  for (CoinBigIndex j = columnStart[iColumn];
     1452                       j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1453                    int iRow = row[j];
     1454                    rowActivity[iRow] -= element[j];
     1455                  }
     1456                  nEasyDown++;
     1457                }
     1458              } else {
     1459                which[nPossible++]=iColumn;
     1460              }
     1461            }
     1462          }
     1463        }
     1464#ifdef REPORT
     1465        printf("%d possible down, %d easy down of which %d are slacks\n",
     1466               nPossible,nEasyDown,nSlackDown);
     1467#endif
     1468        double * needed = new double [numberRows];
     1469        for (int i=0;i<numberRows;i++) {
     1470          double value = rowLower[i] - rowActivity[i];
     1471          if (value<1.0e-8)
     1472            value=0.0;
     1473          needed[i]=value;
     1474        }
     1475        if (gap && /*nUnder==1 &&*/ nonSOS) {
     1476          double * weight = new double [numberColumns];
     1477          int * sort = new int [numberColumns];
     1478          // look at ones not in set
     1479          int nPossible=0;
     1480          for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     1481            if (!newSolution[iColumn]&&
     1482                columnUpper[iColumn]>columnLower[iColumn]) {
     1483              int iSos=-1;
     1484              double value=0.0;
     1485              for (CoinBigIndex j = columnStart[iColumn];
     1486                   j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1487                int iRow = row[j];
     1488                if (sos[iRow]<0) {
     1489                  if (needed[iRow])
     1490                    value += CoinMin(element[j]/needed[iRow],1.0);
     1491                } else {
     1492                  iSos=iRow;
     1493                }
     1494              }
     1495              if (value && iSos<0) {
     1496                weight[nPossible]=-value;
     1497                sort[nPossible++]=iColumn;
     1498              }
     1499            }
     1500          }
     1501          CoinSort_2(weight,weight+nPossible,sort);
     1502          for (int i=0;i<nPossible;i++) {
     1503            int iColumn = sort[i];
     1504            double helps=0.0;
     1505            for (CoinBigIndex j = columnStart[iColumn];
     1506                 j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1507              int iRow = row[j];
     1508              if (needed[iRow])
     1509                helps += CoinMin(needed[iRow],element[j]);
     1510            }
     1511            if (helps) {
     1512              newSolution[iColumn] = 1.0;
     1513              newSolutionValue += modifiedCost[iColumn];
     1514              for (CoinBigIndex j = columnStart[iColumn];
     1515                   j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1516                int iRow = row[j];
     1517                rowActivity[iRow] += element[j];
     1518                if (needed[iRow]) {
     1519                  needed[iRow] -= element[j];
     1520                  if (needed[iRow]<1.0e-8)
     1521                    needed[iRow]=0.0;
     1522                }
     1523              }
     1524              gap -= helps;
     1525#ifdef REPORT
     1526              {
     1527                double gap2 = 0.0;
     1528                for (iRow = 0; iRow < numberRows; iRow++) {
     1529                  if (rowActivity[iRow] < rowLower[iRow] - 10.0*primalTolerance) {
     1530                    gap2 += rowLower[iRow]-rowActivity[iRow];
     1531                  }
     1532                }
     1533                printf("estimated gap (nonsos) %g - computed %g\n",
     1534                       gap,gap2);
     1535              }
     1536#endif
     1537              if (gap<1.0e-12)
     1538                break;
     1539            }
     1540          }
     1541          delete [] weight;
     1542          delete [] sort;
     1543        }
     1544        if (gap&&nPossible/*&&nUnder==1*/&&true&&model_->bestSolution()) {
     1545          double * weight = new double [numberColumns];
     1546          int * sort = new int [numberColumns];
     1547          // look at ones in sets
     1548          const double * goodSolution = model_->bestSolution();
     1549          int nPossible=0;
     1550          double largestWeight=0.0;
     1551          for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     1552            if (!newSolution[iColumn]&&goodSolution[iColumn]&&
     1553                columnUpper[iColumn]>columnLower[iColumn]) {
     1554              int iSos=-1;
     1555              double value=0.0;
     1556              for (CoinBigIndex j = columnStart[iColumn];
     1557                   j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1558                int iRow = row[j];
     1559                if (sos[iRow]<0) {
     1560                  if (needed[iRow])
     1561                    value += CoinMin(element[j]/needed[iRow],1.0);
     1562                } else {
     1563                  iSos=iRow;
     1564                }
     1565              }
     1566              if (value&&iSos>=0) {
     1567                // see if value bigger than current
     1568                int jColumn = firstGub[iSos];
     1569                assert (jColumn>=0);
     1570                while (jColumn>=0) {
     1571                  if (newSolution[jColumn])
     1572                    break;
     1573                  jColumn = nextGub[jColumn];
     1574                }
     1575                assert (jColumn>=0);
     1576                double value2=0.0;
     1577                for (CoinBigIndex j = columnStart[jColumn];
     1578                     j < columnStart[jColumn] + columnLength[jColumn]; j++) {
     1579                  int iRow = row[j];
     1580                  if (needed[iRow])
     1581                    value2 += CoinMin(element[j]/needed[iRow],1.0);
     1582                }
     1583                if (value>value2) {
     1584                  weight[nPossible]=-(value-value2);
     1585                  largestWeight = CoinMax(largestWeight,(value-value2));
     1586                  sort[nPossible++]=iColumn;
     1587                }
     1588              }
     1589            }
     1590          }
     1591          if (nPossible) {
     1592            double * temp = new double [numberRows];
     1593            int * which2 = new int [numberRows];
     1594            memset(temp,0,numberRows*sizeof(double));
     1595            // modify so ones just more than gap best
     1596            if (largestWeight>gap&&nUnder==1) {
     1597              double offset = 4*largestWeight;
     1598              for (int i=0;i<nPossible;i++) {
     1599                double value = -weight[i];
     1600                if (value>gap-1.0e-12)
     1601                  weight[i] = -(offset-(value-gap));
     1602              }
     1603            }
     1604            CoinSort_2(weight,weight+nPossible,sort);
     1605            for (int i=0;i<nPossible;i++) {
     1606              int iColumn = sort[i];
     1607              int n=0;
     1608              // find jColumn
     1609              int iSos=-1;
     1610              for (CoinBigIndex j = columnStart[iColumn];
     1611                   j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     1612                int iRow = row[j];
     1613                temp[iRow]=element[j];
     1614                which2[n++]=iRow;
     1615                if (sos[iRow]>=0) {
     1616                  iSos=iRow;
     1617                }
     1618              }
     1619              int jColumn = firstGub[iSos];
     1620              assert (jColumn>=0);
     1621              while (jColumn>=0) {
     1622                if (newSolution[jColumn])
     1623                  break;
     1624                jColumn = nextGub[jColumn];
     1625              }
     1626              assert (jColumn>=0);
     1627              for (CoinBigIndex j = columnStart[jColumn];
     1628                   j < columnStart[jColumn] + columnLength[jColumn]; j++) {
     1629                int iRow = row[j];
     1630                if (!temp[iRow])
     1631                  which2[n++]=iRow;
     1632                temp[iRow] -= element[j];
     1633              }
     1634              double helps = 0.0;
     1635              for (int i=0;i<n;i++) {
     1636                int iRow = which2[i];
     1637                double newValue = rowActivity[iRow]+temp[iRow];
     1638                if (temp[iRow]>1.0e-8) {
     1639                  if (rowActivity[iRow]<rowLower[iRow]-1.0e-8) {
     1640                    helps += CoinMin(temp[iRow],
     1641                                     rowLower[iRow]-rowActivity[iRow]);
     1642                  }
     1643                } else if (temp[iRow]<-1.0e-8) {
     1644                  if (newValue<rowLower[iRow]-1.0e-12) {
     1645                    helps -= CoinMin(-temp[iRow],
     1646                                     1.0*(rowLower[iRow]-newValue));
     1647                  }
     1648                }
     1649              }
     1650              if (helps>0.0) {
     1651                newSolution[iColumn]=1.0;
     1652                newSolution[jColumn]=0.0;
     1653                newSolutionValue += modifiedCost[iColumn]-modifiedCost[jColumn];
     1654                for (int i=0;i<n;i++) {
     1655                  int iRow = which2[i];
     1656                  double newValue = rowActivity[iRow]+temp[iRow];
     1657                  rowActivity[iRow] = newValue;
     1658                  temp[iRow]=0.0;
     1659                }
     1660                gap -= helps;
     1661#ifdef REPORT
     1662                {
     1663                  double gap2 = 0.0;
     1664                  for (iRow = 0; iRow < numberRows; iRow++) {
     1665                    if (rowActivity[iRow] < rowLower[iRow] - 10.0*primalTolerance) {
     1666                      gap2 += rowLower[iRow]-rowActivity[iRow];
     1667                    }
     1668                  }
     1669                  printf("estimated gap %g - computed %g\n",
     1670                         gap,gap2);
     1671                }
     1672#endif
     1673                if (gap<1.0e-8)
     1674                  break;
     1675              } else {
     1676                for (int i=0;i<n;i++)
     1677                  temp[which2[i]]=0.0;
     1678              }
     1679            }
     1680            delete [] which2;
     1681            delete [] temp;
     1682          }
     1683          delete [] weight;
     1684          delete [] sort;
     1685        }
     1686        delete [] needed;
     1687      }
     1688#ifdef REPORT
     1689      {
     1690        double gap=0.0;
     1691        double over = 0.0;
     1692        for (iRow = 0; iRow < numberRows; iRow++) {
     1693          if (rowActivity[iRow] < rowLower[iRow] - 10.0*primalTolerance) {
     1694            double value = rowLower[iRow]-rowActivity[iRow];
     1695#if REPORT>1
     1696            printf("below on %d is %g - activity %g lower %g\n",
     1697                   iRow,value,rowActivity[iRow],rowLower[iRow]);
     1698#endif
     1699            gap += value;
     1700          } else if (rowActivity[iRow] > rowUpper[iRow] + 10.0*primalTolerance) {
     1701            double value = rowActivity[iRow]-rowUpper[iRow];
     1702#if REPORT>1
     1703            printf("above on %d is %g - activity %g upper %g\n",
     1704                   iRow,value,rowActivity[iRow],rowUpper[iRow]);
     1705#endif
     1706            gap += value;
     1707          } else {
     1708            double value = rowActivity[iRow]-rowLower[iRow];
     1709            if (value) {
     1710#if REPORT>1
     1711              printf("over on %d is %g - activity %g lower %g\n",
     1712                     iRow,value,rowActivity[iRow],rowLower[iRow]);
     1713#endif
     1714              over += value;
     1715            }
     1716          }
     1717        }
     1718        printf("modified final gap %g - over %g - %d added - solvalue %g\n",
     1719               gap,over,nAdded,newSolutionValue);
     1720      }
     1721#endif
     1722      if (!gap) {
     1723        break;
     1724      } else {
     1725        if (iPass==0) {
     1726          costBias = 10.0*newSolutionValue/static_cast<double>(nAdded);
     1727        } else {
     1728          costBias *= 10.0;
     1729        }
     1730      }
     1731    }
     1732    delete [] newSolution0;
     1733    delete [] rowWeight;
    12171734    delete [] sos;
     1735    delete [] firstGub;
     1736    delete [] nextGub;
    12181737    if (newSolutionValue < solutionValue) {
    12191738        // check feasible
    1220       const double * rowLower = solver->getRowLower();
    1221       const double * rowUpper = solver->getRowUpper();
    12221739      memset(rowActivity, 0, numberRows*sizeof(double));
    12231740        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    … …  
    12931810                setWhen(0);
    12941811        }
    1295         // Only works if coefficients positive and all rows L or SOS
     1812        // Only works if coefficients positive and all rows L/G or SOS
    12961813        OsiSolverInterface * solver = model_->solver();
    12971814        const double * columnUpper = solver->getColUpper();
    … …  
    13121829            // SOS
    13131830            originalRhs_[iRow]=-1.0;
    1314           } else if (rowLower[iRow] > 0.0) {
     1831          } else if (rowLower[iRow] > 0.0 && rowUpper[iRow] < 1.0e10) {
    13151832                good = false;
    13161833          } else if (rowUpper[iRow] < 0.0) {
    13171834                good = false;
     1835          } else if (rowUpper[iRow] < 1.0e10) {
     1836            originalRhs_[iRow]=rowUpper[iRow];
    13181837          } else {
    1319             originalRhs_[iRow]=rowUpper[iRow];
     1838            originalRhs_[iRow]=rowLower[iRow];
    13201839          }
    13211840        }
    13221841        int numberColumns = solver->getNumCols();
    13231842        for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
     1843          if (!columnLength[iColumn])
     1844            continue;
    13241845            if (columnLower[iColumn] < 0.0 || columnUpper[iColumn] > 1.0)
    13251846                good = false;
    13261847            CoinBigIndex j;
    13271848            int nSOS=0;
     1849            if (!solver->isInteger(iColumn))
     1850              good = false;
    13281851            for (j = columnStart[iColumn];
    13291852                    j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    … …  
    13371860                }
    13381861            }
    1339             if (nSOS>1||!solver->isBinary(iColumn))
     1862            if (nSOS > 1)
    13401863              good = false;
    13411864        }

  • trunk/Cbc/src/CbcHeuristicRENS.cpp

    r1587 r1591  
    570570            }
    571571            char line[100];
    572             sprintf(line,"Heuristic %s fixed %d to one",
     572            sprintf(line,"Heuristic %s fixed %d to one (%d sets)",
    573573                    heuristicName(),
    574                     nFixed);
     574                    nFixed,nSOS);
    575575            model_->messageHandler()->message(CBC_FPUMP1, model_->messages())
    576576              << line

  • trunk/Cbc/src/CbcSolver.cpp

    r1585 r1591  
    16211621        }
    16221622        std::string field;
     1623#ifndef CBC_QUIET
    16231624        if (!noPrinting_) {
    16241625           sprintf(generalPrint,
    … …  
    16581659            }
    16591660        }
     1661#endif
    16601662        while (1) {
    16611663            // next command
    … …  
    26682670                                           "Enumerated nodes: 0\n");
    26692671                                   sprintf(generalPrint + strlen(generalPrint),
    2670                                            "Total iterations: 0\n");
     2672                                           "Total iterations: 0\n");
     2673#ifndef CBC_QUIET
    26712674                                   sprintf(generalPrint + strlen(generalPrint),
    26722675                                           "Time (seconds):   %.2f\n",
    26732676                                           CoinCpuTime() - time0);
     2677#endif
    26742678                                   generalMessageHandler->message(CLP_GENERAL, generalMessages)
    26752679                                      << generalPrint
    … …  
    60746078                                            "Total iterations: %d\n",
    60756079                                            babModel_->getIterationCount());
     6080#ifndef CBC_QUIET
    60766081                                    sprintf(generalPrint + strlen(generalPrint),
    60776082                                            "Time (seconds):   %.2f\n",
    … …  
    60806085                                            "Total time:       %.2f\n",
    60816086                                            time2 - time0);
     6087#endif
    60826088                                    generalMessageHandler->message(CLP_GENERAL, generalMessages)
    60836089                                       << generalPrint
Note: See TracChangeset for help on using the changeset viewer.