Changeset 968


Ignore:
Timestamp:
Feb 17, 2011 2:27:01 PM (11 years ago)
Author:
lou
Message:

(snapshot) Final commit prior to excising clique code. The sole purpose of
this commit is so that there's a single revision that can be identified when
someone comes back to try and recover clique functionality.

Location:
branches/CglWorking101215/src/CglProbing
Files:
2 added
6 edited

Legend:

Unmodified
Added
Removed
  • branches/CglWorking101215/src/CglProbing/CglProbing.cpp

    r964 r968  
    4040
    4141namespace {
     42
     43/*
     44  Compare two sets of upper and lower bound arrays and generate column cuts.
     45*/
     46
     47void makeColCuts (int nCols, OsiCuts &cs,
     48                  const char *const intVar,
     49                  const double *const colsol,
     50                  const double *const origlbs,
     51                  const double *const origubs,
     52                  const double *const soln,
     53                  double *const newlbs,
     54                  double *const newubs
     55                 )
     56{
     57  int numberChanged = 0 ;
     58  int ifCut = 0 ;
     59  CoinPackedVector lbs ;
     60  CoinPackedVector ubs ;
     61/*
     62  Scan the bound arrays, recording changed bounds and counting changes and
     63  actual cuts.
     64*/
     65  for (int i = 0 ; i < nCols ; ++i) {
     66    if (intVar[i]) {
     67      newubs[i] = CoinMin(origubs[i],floor(newubs[i]+1.0e-4)) ;
     68      if (newubs[i] < origubs[i]-1.0e-8) {
     69        if (newubs[i] < colsol[i]-1.0e-8) ifCut++ ;
     70        ubs.insert(i,newubs[i]) ;
     71        numberChanged++ ;
     72      }
     73      newlbs[i] = CoinMax(origlbs[i],ceil(newlbs[i]-1.0e-4)) ;
     74      if (newlbs[i] > origlbs[i]+1.0e-8) {
     75        if (newlbs[i] > colsol[i]+1.0e-8) ifCut++ ;
     76        lbs.insert(i,newlbs[i]) ;
     77        numberChanged++ ;
     78      }
     79    }
     80  }
     81/*
     82  Stash the changes in a column cut. If the changes actually cut off some part
     83  of the current solution, boost the effectiveness.
     84*/
     85  if (numberChanged) {
     86    OsiColCut cc ;
     87    cc.setUbs(ubs) ;
     88    cc.setLbs(lbs) ;
     89    if (ifCut) {
     90      cc.setEffectiveness(100.0) ;
     91    } else {
     92      cc.setEffectiveness(1.0e-5) ;
     93    }
     94    cs.insert(cc) ;
     95  }
     96}
     97
     98
     99
    42100/*
    43101  This method will attempt to identify variables that will naturally take on
     
    51109    * Why not execute repeatedly, as the problem is simplified due to fixed
    52110      variables?
     111  This would involve calling analyze as part of integer preprocessing.
    53112
    54113  The expression abs(x - floor(x+.5)) evaluates to zero when x is integer and
     
    57116  coefficients aij = 1/k, k integer. This guarantees that b/aij is integer
    58117  as long as b is integer.
     118
     119  TODO: Should this routine return a count, instead of a boolean? No
     120        particular reason other than more information.  -- lh, 110211 --
    59121*/
    60122
     
    844906// Generate disaggregation cuts
    845907//-------------------------------------------------------------------
    846 void CglProbing::generateCuts(const OsiSolverInterface & si, OsiCuts & cs,
     908
     909/*
     910  The traditional const method. Row and column cuts will be returned in cs.
     911*/
     912void CglProbing::generateCuts(const OsiSolverInterface &si, OsiCuts &cs,
    847913                              const CglTreeInfo info2) const
    848914{
     
    9441010/*
    9451011  Delete the working arrays and restore the row cut mode.
     1012
     1013  TODO: Why are we fiddling with colLower_, colUpper_? Track this down.
     1014        -- lh, 110212 --
    9461015*/
    9471016  delete [] rowLower ;
     
    9641033
    9651034/*
    966   TODO: Note that this is not a const method (contrast with generateCuts).
    967         -- lh, 100924
     1035  Note that this is not a const method (contrast with generateCuts).
    9681036 
    969   TODO: Note that the CglTreeInfo object is not const here. In generateCuts,
     1037  TODO: Note that the CglTreeInfo object is not const either. In generateCuts,
    9701038        it's passed in as a const and a nonconst clone is created. I need to
    9711039        look into that.  -- lh, 100924 --
     1040
     1041        Possible reason is use of CglTreeInfo to pass back strengthened rows
     1042        for row replacement during preprocessing.  -- lh, 110211 --
    9721043
    9731044  Generally, check generateCuts for additional comment. I went through it first
     
    10251096/*
    10261097  Create working arrays for row and column bounds.
     1098
     1099  We need nRows+1 because we may use the objective function as a constraint
     1100  in the working system within gutsOfGenerateCuts.
    10271101*/
    10281102  int nRows = si.getNumRows() ;
    1029   double * rowLower = new double[nRows+1] ;
    1030   double * rowUpper = new double[nRows+1] ;
     1103  double *rowLower = new double[nRows+1] ;
     1104  double *rowUpper = new double[nRows+1] ;
    10311105
    10321106  int nCols = si.getNumCols() ;
    1033   double * colLower = new double[nCols] ;
    1034   double * colUpper = new double[nCols] ;
     1107  double *colLower = new double[nCols] ;
     1108  double *colUpper = new double[nCols] ;
    10351109/*
    10361110  Do the work.
     
    10701144
    10711145  If we're debugging, check whether this is a plausible outcome.
    1072 
    1073   TODO: Remind me --- Why can't I use the debugger created at the top of the
    1074         method. And don't I need to dispose of both of these?
    10751146*/
    10761147  if (ninfeas) {
     
    10951166  Mode 3 is just like mode 2, except that the tightened row bounds will be
    10961167  returned to the client.
    1097 */
    1098   if (mode_ == 3) {
     1168
     1169  Don't replace bounds supplied by the client if we've gone infeasible.
     1170  There's no guarantee of consistency in the values reported back from
     1171  gutsOfGenerateCuts.
     1172*/
     1173  if (mode_ == 3 && !ninfeas) {
    10991174    delete [] rowLower_ ;
    11001175    delete [] rowUpper_ ;
     
    11051180    delete [] rowUpper ;
    11061181  }
    1107 /*
    1108   But given that this is generateCutsAndModify, we certainly want to return the
    1109   tightened column bounds.
    1110 */
    1111   delete [] colLower_ ;
    1112   delete [] colUpper_ ;
    1113   colLower_ = colLower ;
    1114   colUpper_ = colUpper ;
     1182  if (!ninfeas) {
     1183    delete [] colLower_ ;
     1184    delete [] colUpper_ ;
     1185    colLower_ = colLower ;
     1186    colUpper_ = colUpper ;
     1187  } else {
     1188    delete [] colUpper ;
     1189    delete [] colLower ;
     1190  }
    11151191/*
    11161192  jjf: Setup information
     
    11281204# endif
    11291205
    1130   return ninfeas ;
     1206  return (ninfeas) ;
    11311207}
    11321208
     
    11381214
    11391215  rowLower, rowUpper, colLower, and colUpper should be allocated by the caller
    1140   and will be filled in by gutsOfGenerateCuts.
     1216  and will be filled in by gutsOfGenerateCuts. If the system is found to be
     1217  infeasible, there's no particular guarantee of consistency for these arrays
     1218  and the content should not be used.
    11411219*/
    11421220int CglProbing::gutsOfGenerateCuts (const OsiSolverInterface &si,
     
    12011279  const double *colsol = si.getColSolution() ;
    12021280/*
    1203   Stage 1: Preliminary processing.
     1281  Stage 1: Preliminary Processing (Common setup?)
    12041282
    12051283  During Stage 1 we'll set `reasonable' bounds on variables and constraints,
     
    12341312  If we lose feasibility here, cut our losses and bail out while it's still
    12351313  uncomplicated.
    1236 
    1237   TODO: Really, analyze should return a count instead of a boolean.
    1238         -- lh, 110204 --
    12391314*/
    12401315  if (!info->inTree && !info->pass) {
     
    12601335
    12611336  TODO: Down in probeClique, there's a clause that forces the cutoff to
    1262         infty for mode 0.
    1263 */
    1264   double offset ;
     1337        infty for mode 0. Why not here?  -- lh, 110211 --
     1338*/
     1339  double offset = 0.0 ;
    12651340  si.getDblParam(OsiObjOffset,offset) ;
    1266   double cutoff ;
     1341  double cutoff = 0.0 ;
    12671342  bool cutoff_available = si.getDblParam(OsiDualObjectiveLimit,cutoff) ;
    12681343  if (!cutoff_available) {
     
    12961371/*
    12971372  Create the working row copy from the si's model or from a snapshot (only if
    1298   the user specified mode 0 and a snapshot exists).
     1373  the user specified mode 0 and a snapshot exists). rowLower and rowUpper are
     1374  valid at completion, but if we initialise the local row copy from a snapshot
     1375  the number of valid entries may not match the number of rows for the system
     1376  held in the si.
    12991377*/
    13001378  CoinPackedMatrix *rowCopy = setupRowCopy(mode,useObj,cutoff,
     
    13131391  realRows is the cross-reference between rows in rowCopy and rows of the
    13141392  original model. rowStartPos is the first positive coefficient in a row.
     1393
     1394  Note that the number of rows in the local copy is almost certainly not the
     1395  same as the number of rows in the system held in the si. You can't just copy
     1396  rowLower and rowUpper; they must be translated if realRows exists.
    13151397*/
    13161398  int *realRows = 0 ;
     
    13281410/*
    13291411  Make a column-major copy, after we've winnowed the rows to the set we want
    1330   to work with.
     1412  to work with. The break out the internal arrays for subsequent use.
    13311413
    13321414  TODO: I'm asking myself, ``Why does the snapshot keep a colum-major
    1333         copy? -- lh, 101021 --
     1415        copy? It doesn't seem to be used. -- lh, 101021 --
    13341416*/
    13351417  CoinPackedMatrix * columnCopy = new CoinPackedMatrix(*rowCopy,0,0,true) ;
    1336 
    1337 /*
    1338   End of Stage 1.
    1339 */
    1340 
    1341 # if CGL_DEBUG > 0
    1342   const OsiRowCutDebugger *debugger = si.getRowCutDebugger() ;
    1343   if (debugger && !debugger->onOptimalPath(si))
    1344     debugger = NULL ;
    1345 # else
    1346   const OsiRowCutDebugger * debugger = NULL ;
    1347 # endif
    1348 
    1349 /*
    1350   Begin Stage 2
    1351 */
    1352 
    13531418  const int *column = rowCopy->getIndices() ;
    13541419  const CoinBigIndex *rowStart = rowCopy->getVectorStarts() ;
     
    13561421  const double *rowElements = rowCopy->getElements() ;
    13571422/*
    1358   jjf: Arrays so user can find out what happened
    1359 
    1360   TODO: Misleading? lookedAt_ seems to be used to hold the set of variables
    1361         selected for probing. The question is whether the selection activity
    1362         is duplicated in probe() and/or probeCliques().
     1423  Allocate an array to hold indices of the columns we will probe.
    13631424*/
    13641425  if (!lookedAt_) {
     
    13671428  numberThisTime_ = 0 ;
    13681429/*
    1369   jjf: Let us never add more than twice the number of rows worth of row cuts
    1370   jjf: Keep cuts out of cs until end so we can find duplicates quickly
    1371 
    1372   TODO: Well, ok, but none of the values for nRowsFake agree with that comment
    1373         (assuming nRowsFake is the limit on the number of row cuts). And the
    1374         row_cut constructor also adjusts the value.
    1375         -- lh, 100923 --
     1430  Allocate a container to hold the cuts we'll generate locally. We'll only use
     1431  this in mode 0.
    13761432
    13771433  TODO: There's a similar calculation in probe() to size the local container
    1378         used there. It allocates much more space.
     1434        used there. It allocates much more space. But note that the container
     1435        passed to probe() is the OsiCuts container cs, given as a parameter to
     1436        this method. rowCut is going to be used down in mode 0 for something.
    13791437*/
    13801438  int nRowsSafe = CoinMin(nRows,si.getNumRows()) ;
    1381   int nRowsFake = info->inTree ? nRowsSafe/3 : nRowsSafe ;
    1382   if (!info->inTree && !info->pass)
    1383     nRowsFake *= 5 ;
    1384   CglProbingRowCut rowCut(nRowsFake,!info->inTree) ;
    1385 
     1439  int cutCapacity = info->inTree ? nRowsSafe/3 : nRowsSafe ;
     1440  if (!info->inTree && !info->pass) cutCapacity *= 5 ;
     1441  CglProbingRowCut rowCut(cutCapacity,!info->inTree) ;
     1442/*
     1443  Allocate the row processing record (markR) and the arrays that will hold
     1444  row lhs bounds (minR, maxR).
     1445*/
    13861446  int * markR = new int [nRows] ;
    13871447  double * minR = new double [nRows] ;
    13881448  double * maxR = new double [nRows] ;
    13891449/*
     1450  End of Stage 1: Setup.
     1451
     1452  Begin Stage 2: Mode-dependent Setup and Probing
     1453*/
     1454/*
    13901455  In modes 1 and 2, we're allowed to tighten bounds. Call tighten() to do bound
    13911456  propagation (two pass limit), then create `column cuts' to record the
    13921457  tightened bounds.
    1393 
    1394   TODO: Surely we could skip this if we've already lost feasibility! Or does
    1395         tighten() bail immediately? Not immediately obvious. This method
    1396         really needs a cleanupAndGoHome() method that can be called at the
    1397         point where we discover infeasibility.
    1398         -- lh, 100923 --
    13991458
    14001459  NOTE: Split tighten() into tighten() and tightenClique() -- lh, 101124 --
     
    14121471                        intVar,2,primalTolerance_);
    14131472    }
    1414     if (!feasible)
    1415       ninfeas = 1 ;
    14161473    if (!ninfeas) {   // block extends to end of if (mode)
    14171474/*
    1418   jjf: create column cuts where integer bounds have changed
    1419 
    1420   TODO: Looks like another good candidate for a method. Already a block.
    1421         -- lh, 100923 --
    1422  
    1423   TODO: ifCut could just as well be a boolean.  -- lh, 100923 --
    1424 */
    1425       {
    1426         const double *lower = si.getColLower() ;
    1427         const double *upper = si.getColUpper() ;
    1428         const double *colsol = si.getColSolution() ;
    1429         int numberChanged = 0 ;
    1430         int ifCut = 0 ;
    1431         CoinPackedVector lbs ;
    1432         CoinPackedVector ubs ;
    1433         for (int i = 0 ; i < nCols ; ++i) {
    1434           if (intVar[i]) {
    1435             colUpper[i] = CoinMin(upper[i],floor(colUpper[i]+1.0e-4)) ;
    1436             if (colUpper[i] < upper[i]-1.0e-8) {
    1437               if (colUpper[i] < colsol[i]-1.0e-8)
    1438                 ifCut = 1 ;
    1439               ubs.insert(i,colUpper[i]) ;
    1440               numberChanged++ ;
    1441             }
    1442             colLower[i] = CoinMax(lower[i],ceil(colLower[i]-1.0e-4)) ;
    1443             if (colLower[i] > lower[i]+1.0e-8) {
    1444               if (colLower[i] > colsol[i]+1.0e-8)
    1445                 ifCut = 1 ;
    1446               lbs.insert(i,colLower[i]) ;
    1447               numberChanged++ ;
    1448             }
    1449           }
    1450         }
    1451         if (numberChanged) {
    1452           OsiColCut cc ;
    1453           cc.setUbs(ubs) ;
    1454           cc.setLbs(lbs) ;
    1455           if (ifCut) {
    1456             cc.setEffectiveness(100.0) ;
    1457           } else {
    1458             cc.setEffectiveness(1.0e-5) ;
    1459           }
    1460 #         if CGL_DEBUG > 0
    1461           CglProbingDebug::checkBounds(si,cc) ;
    1462 #         endif
    1463           cs.insert(cc) ;
    1464         }
    1465       }  // end create column cuts block
    1466 
    1467       // This next block extends to end of mode 1,2 block
    1468       if (maxProbe > 0) {
     1475  Create column cuts where integer bounds have changed.
     1476*/
     1477      makeColCuts(nCols,cs,intVar,colsol,si.getColLower(),si.getColUpper(),
     1478                  si.getColSolution(),colLower,colUpper) ;
     1479#     if CGL_DEBUG > 0
     1480      CglProbingDebug::checkBounds(si,cs.getColCut(cs.sizeColCuts())) ;
     1481#     endif
     1482/*
     1483  We've calculated changes in column bounds due to propagation. To go further,
     1484  we need to be able to probe.
     1485*/
     1486      if (maxProbe > 0) {    // block extends to end of if (mode)
    14691487        numberThisTime_ = 0 ;
    14701488/*
     
    14771495  TODO: Why do we need a second run for this? Why not do it as part of
    14781496        tighten()? My old notes indicate that tighten calculates the row
    1479         bounds (and indeed it must to propagate column bounds). On a purely
    1480         mechanical level, tighten does not fill in minR, maxR, markR.
    1481         -- lh, 100923 --
     1497        bounds (and indeed it must to propagate column bounds). It's not
     1498        like we're saving any work, other than transcribing the bounds into
     1499        minR and maxR. On a purely mechanical level, tighten does not fill
     1500        in minR, maxR, markR.  -- lh, 100923 --
    14821501*/
    14831502        calcRowBounds(colLower,colUpper,column,rowElements,
     
    15581577        }
    15591578/*
    1560   lookedAt_ now contains a list of indices of variables to probe. Rows worth
    1561   processing are tagged with -1 in markR
    1562 */
    1563 /*
    1564   TODO: Eh? Sort by index for repeatability? Seems pointless, and it is
    1565         commented out. -- lh, 100924 --
    1566 */
    1567         // sort to be clean
    1568         // std::sort(lookedAt_,lookedAt_+numberThisTime_) ;
    1569 /*
    1570   And do the probing, with or without cliques.
     1579  lookedAt_ now contains a list of indices of variables to probe. Rows
     1580  worth processing are tagged with -1 in markR. Do the probing, with or
     1581  without cliques.
    15711582*/
    15721583        if (!numberCliques_) {
     
    15761587                          useObj,useCutoff,cutoff) ;
    15771588        } else {
    1578           ninfeas = probeCliques(si,debugger,cs,
    1579                                  colLower,colUpper,rowCopy,columnCopy,
     1589          ninfeas = probeCliques(si,cs,colLower,colUpper,rowCopy,columnCopy,
    15801590                                 realRows,rowLower,rowUpper,
    15811591                                 intVar,minR,maxR,markR,info,
     
    15881598  End of block for mode = 1 or 2; begin mode = 0.
    15891599 
    1590   Mode 0 works with the snapshot and does not do bound propagation.
    1591 
    1592   TODO: In fact, it's unclear to me just what we're doing here. The activity
    1593         is completely unsymmetric with the previous case. -- lh, 100924 --
     1600  Mode 0 works with the snapshot and does not do bound propagation. But note
     1601  that much of the preparatory work that's done above for modes 1 and 2 in
     1602  Stage 1 is done in snapshot() when the snapshot is created, including
     1603  bound propagation. It really should be possible to unify a lot of this.
     1604
     1605  TODO: Note that this code has been effectively disabled since 080428 (r629)
     1606        when the only use (in CglPreProcess::modified) was guarded with
     1607        #if 0 / #endif. Expect oddities.  -- lh, 110212 --
     1608
     1609  TODO: And now that I've worked through it, I have real questions. The key is
     1610        the table down at ACTIONS.
     1611       
     1612        We start out by calling probeClique and get back a private container
     1613        of cuts. The code is clear that it's expecting only disaggregation
     1614        cuts. We go through a lot of work to classify these cuts and process
     1615        them into an alternate form which is stored in cutVector_. Then
     1616        we work through cutVector_, looking for the situations identified
     1617        in ACTIONS. Which are exactly the things that should be picked up
     1618        by monotone actions, disaggregation cuts, and implication cuts. So
     1619        if this code ever finds anything, it indicates an algorithm error
     1620        someplace else. (Just possibly, we might pick up something across
     1621        calls to this CglProbing object. cutVector_ is preserved for as
     1622        long as the object is in existence.) Then we go through a lot
     1623        more work to go back over cutVector_ one more time, checking for
     1624        violations and reconstituting the cuts we processed away at the
     1625        start. I'm still shaking my head.
     1626
     1627        What's worth keeping? Well, the basic concept of the snapshot: We
     1628        don't have to recreate the row and column copies with each call. And
     1629        there's a simple placeholder where John's comment says `could augment
     1630        cliques'. Which is true, but I wonder if that, too, isn't already done
     1631        somewhere else, or could be done somewhere else more efficiently.
     1632
     1633        So my conclusion, at this point, is that the only distinction worth
     1634        keeping between mode 0 and modes 1 -- 3 is that mode 0 works from the
     1635        snapshot. I'll see if I notice more along the way.
     1636
     1637        -- lh, 110212 --
     1638
     1639        Essentially confirmed in correspondence with JJF. Clique code
     1640        imported from OSL but never fully functional.  -- lh, 110216 --
    15941641*/
    15951642  else
     
    16001647       make up list of new variables to look at
    16011648
    1602   Initial hypothesis: We're working from a snapshot, so we've done at least
    1603   one round of probing and cut generation. In particular, we should have
    1604   disaggregation cuts, stashed in cutVector_. So we're only going to probe
    1605   variables where we have no information. We're looking at binary
    1606   variables because that's the natural index for a collection of
    1607   disaggregation cuts generated by splitting a constraint where a single
    1608   binary variable controls whether a bunch of other variables can be
    1609   nonzero.
    1610 
    1611   The comment in CglProbing for mode 0: `if no information exists for
    1612   that variable then probing will be done' makes more sense now. And if you
    1613   take that together with the comment `only with snapshot; otherwise as mode
    1614   1' then the tweaking of mode at the head of generateCutsAndModify almost
    1615   makes sense.
    1616 
    1617   -- lh, 101021 --
    1618 
    1619   TODO: Check out how we can get away with reordering the major operations
    1620         here. In the mode != 0 block, it's calcRowBounds; select & sort; probe.
    1621         Here, it's select and sort; calcRowBounds; probe.  -- lh, 101021 --
    1622 
    1623         Now I think about it, there's no problem. calcRowBounds flags
    1624         rows worth processing, which is independent of columns worth
    1625         processing. Order makes no difference.  -- lh, 101125 --
    1626 
     1649  We're working from a snapshot.  In particular, we should have information
     1650  about probing stashed in cutVector_ (for each x<p> = sequence, a vector of
     1651  affected variables x<t> = index). So we're only going to probe variables
     1652  where we have no information.
     1653
     1654  When the snapshot is created, cutVector_ is initialised with all binary
     1655  variables. So the snapshot documentation is a little deceptive. On first use
     1656  after creation, we'll probe all binary variables, in order of distance from
     1657  integrality, subject to whatever limits have been set.
    16271658*/
    16281659    assert(numberIntegers == numberIntegers_) ;
    16291660    numberThisTime_ = 0 ;
    16301661    const double *colsol = si.getColSolution() ;
    1631     double_int_pair *array = new double_int_pair [nCols] ;
     1662    double_int_pair *array = new double_int_pair [number01Integers_] ;
    16321663#   ifdef ZEROFAULT
    1633     std::memset(array,0,sizeof(double_int_pair)*nCols) ;
     1664    std::memset(array,0,sizeof(double_int_pair)*number01Integers_) ;
    16341665#   endif
    16351666    for (int i = 0 ; i < number01Integers_ ; i++) {
     
    16461677      lookedAt_[i] = array[i].sequence ;
    16471678    }
    1648     // sort to be clean
    1649     // std::sort(lookedAt_,lookedAt_+numberThisTime_) ;
    16501679    delete [] array ;
    1651     // get min max etc for rows
    1652     calcRowBounds(colLower, colUpper, column, rowElements,
    1653              rowStart, rowLength, rowLower, rowUpper,
    1654              minR , maxR , markR, nRows) ;
     1680/*
     1681  lookedAt_ now contains the indices of the binary variables for which we have
     1682  no probing information. Calculate row lhs min and max, and call probeClique
     1683  to do the probing (note that we will have clique information!)
     1684*/
     1685    calcRowBounds(colLower,colUpper,column,rowElements,rowStart,rowLength,
     1686                  rowLower,rowUpper,minR,maxR,markR,nRows) ;
    16551687/*
    16561688  jjf: don't do cuts at all if 0 (i.e. we are just checking bounds)
    1657 
    1658   TODO: The comment is wrong, at least for generateCuts. What happens there is
    1659         that any negative value of rowCuts_ is converted to 0x4, which is
    1660         supposedly column cuts only. But it's definitely not 0, so this block
    1661         will still execute.  -- lh, 100924 --
    16621689
    16631690  TODO: Why don't we get a choice between probe and probeCliques? And why did
     
    16671694    OsiCuts csNew ;
    16681695    if (rowCuts_) {
    1669       ninfeas = probeCliques(si,debugger,csNew,colLower,colUpper,
     1696      ninfeas = probeCliques(si,csNew,colLower,colUpper,
    16701697                             rowCopy,columnCopy,
    16711698                             realRows,rowLower,rowUpper,
     
    16741701    }
    16751702/*
    1676   TODO: And I have no idea what we're getting into here --- completely
    1677         asymmetric to the mode 1 and 2 case. I'm going to bail from here for a
    1678         bit and look over the startup code and probing methods, then come back
    1679         to this.  -- lh, 100924 --
    1680 
    1681         Made it back this far with an understanding of cliques, but have been
    1682         spending far too much time at home with Windows. As a working theory,
    1683         it's this next bit that does disaggregation cuts, or strengthens
    1684         disaggregation cuts with implication cuts from probing.
    1685 
    1686         If I interpret the slides from the 3/97 presentation correctly,
    1687         probing should do implication cuts, and disaggregation would be a
    1688         separate activity.  So I'm off to look at the probing code and
    1689         understand what it does. There's 500 lines of uncommented code coming
    1690         up, followed by another 200 lines with the comment `see if any
    1691         disaggregation cuts are violated.' I'll be back. -- lh, 101021 --
     1703  TODO: I've made it to this point twice now (100924, 101021) and each time
     1704        bailed out to sort out other parts. Now it's time.
     1705
     1706        The first thing that hits me is `If we have no snapshot to work from
     1707        (rowCuts_ = 0), csNew is surely empty. Why are we even executing
     1708        this next bit of code? For that matter, why are we here?! Mode 0
     1709        should have been converted to mode 1.
     1710
     1711        After reading ahead about 350 lines (look for END CLASSIFICATION) it's
     1712        clear what we're doing in the next bit of code. We take the
     1713        disaggregation cuts returned by probeCliques (and gods help us if
     1714        there's anything else) and do a complicated case analysis to try and
     1715        classify the type of cut we're looking at. We have about 150 lines of
     1716        classification code whose sole purpose is to do a count, followed by
     1717        the same 150 lines of classification code, this time filling in the
     1718        entries and labelling them with the classification. This whole code
     1719        block is highly suspect.
     1720
     1721        -- lh, 110211 --
     1722
     1723  AFTER PROBECLIQUE
    16921724*/
    16931725    if (!ninfeas) {
    1694       // go through row cuts
    16951726      int nCuts = csNew.sizeRowCuts() ;
    16961727      int iCut ;
    1697       // need space for backward lookup
    1698       // just for ones being looked at
     1728/*
     1729  backward will allow us to go from a probe index p to the corresponding in
     1730  cutVector_.
     1731  onList is an indicator that x<p> is in cutVector_
     1732
     1733  Presumably onList will become more? It's not clear below why we need it.
     1734*/
    16991735      int *backward = new int [2*nCols] ;
    17001736      int *onList = backward+nCols ;
     
    17081744        onList[j] = 1 ;
    17091745      }
    1710       // first do counts
    1711       // we know initialized to zero
    1712       for (iCut=0;iCut<nCuts;iCut++) {
     1746/*
     1747  jjf: first do counts; we know initialized to zero
     1748
     1749  We know this because we only probed variables that didn't have information.
     1750  (Of course, if we guess wrong below while trying to determine the type of
     1751  constraint, that goes out the window.)
     1752
     1753  TODO: I've fixed up a few obvious errors (redundant tests, incorrect
     1754        semantics) in the next two loops to classify cuts, but there are
     1755        lots more. I'm going to wait 'til I have a grasp of the whole
     1756        activity before trying to do more detail fixes. The possibility
     1757        of classification error remains. Not to mention the `can do
     1758        something?' comments.  -- lh, 110212 --
     1759
     1760  TODO: There were a number of instances in the next two loops where there's a
     1761        test for backward[k] < 0 in one place and, close by, a test for
     1762        onList[k]. Given initialisation immediately above, these tests should
     1763        be equivalent. It'd be nice if the code made that apparent, and I've
     1764        settled for backward[.] < 0.
     1765        -- lh, 110212 --
     1766*/
     1767      for (iCut = 0 ; iCut < nCuts ; iCut++) {
    17131768        OsiRowCut rcut ;
    17141769        CoinPackedVector rpv ;
    17151770        rcut = csNew.rowCut(iCut) ;
    17161771        rpv = rcut.row() ;
    1717         assert(rpv.getNumElements()==2) ;
    1718         const int * indices = rpv.getIndices() ;
    1719         double* elements = rpv.getElements() ;
    1720         double lb=rcut.lb() ;
    1721         // find out which integer
    1722         int which=0 ;
    1723         int i=backward[indices[0]] ;
    1724         if (i<0||!onList[indices[0]]) {
    1725           which=1 ;
    1726           i=backward[indices[1]] ;
    1727           // Just possible variable was general integer but now 0-1
    1728           if (!onList[indices[which]])
    1729             continue ;
     1772/*
     1773  Hmmm ... probeClique looks to be capable of generating coefficient
     1774  strengthening cuts. The code is there, at least. But this assert is only
     1775  guaranteed for disaggregation cuts.  -- lh, 110211 --
     1776
     1777  Examination of the Cgl timeline says that mode 0 has been effectively
     1778  disabled since 080428 r629, so this assert may well predate the introduction
     1779  of coefficient strengthening.  -- lh, 110212 --
     1780*/
     1781        assert(rpv.getNumElements() == 2) ;
     1782/*
     1783  Get the index of the probing variable x<p> from the cut and use it to locate
     1784  the entry in cutVector_. The other index will be x<t>.
     1785
     1786  TODO: The probe variable x<p> could be either index, apparently. But
     1787        it's entirely possible that x<p> in one cut could be x<t> in
     1788        another. What happens then?  -- lh, 110211 --
     1789
     1790  TODO: John allows that neither variable might be in cutVector_, because it's
     1791        possible for a variable to be converted to binary if its bounds are
     1792        reduced to 0-1. So we have the freak accident of a disaggregation
     1793        constraint where x<t> is general integer and x<p> was converted from
     1794        general integer to binary.  -- lh, 110212 --
     1795*/
     1796        const int *indices = rpv.getIndices() ;
     1797        int which = 0 ;
     1798        int i = backward[indices[0]] ;
     1799        if (i < 0) {
     1800          which = 1 ;
     1801          i = backward[indices[1]] ;
     1802          if (i < 0) continue ;
    17301803        }
    17311804        int other = indices[1-which] ;
    1732         if (lb==-DBL_MAX) {
     1805        double *elements = rpv.getElements() ;
     1806        double lb = rcut.lb() ;
     1807/*
     1808  Sort out what kind of cut we're looking at. It's a bit hard to say
     1809  whether this is correct for John's disaggregation cut generation code
     1810  (I think not), and it's a cinch it's wrong for my disaggCuts method. But
     1811  the details can be worked out.
     1812
     1813  TODO: If this is going to be done easily and efficiently, it should be done
     1814        when the disaggregation cut is generated. The critical questions
     1815        (identity of x<p>, probe direction, identity of x<t>, bound change)
     1816        are all known at that point. It's far more difficult to recover that
     1817        information here. Having said that, the current implementation always
     1818        constructs disaggregation cuts in a canonical form, with x<p> in the
     1819        second position.   -- lh, 110212 --
     1820
     1821  In-code comments for this next bit are JJF.
     1822*/
     1823        if (lb == -DBL_MAX) {
    17331824          if (!rcut.ub()) {
    17341825            // UB
    1735             if (elements[which]<0.0) {
    1736               //assert (elements[1-which]>0.0) ;
     1826            double elWhich = elements[which] ;
     1827            double elOther = elements[1-which] ;
     1828            if (elWhich < 0.0) {
     1829              // assert (elOther > 0.0) ;
    17371830              // delta to 0 => x to 0.0
    17381831              cutVector_[i].length++ ;
    17391832            } else {
    1740               if (elements[1-which]<0.0 && fabs(elements[which]/elements[1-which]-
    1741                                               colUpper[other])<1.0e-5) {
     1833              if (elOther < 0.0 &&
     1834                  fabs(elWhich/elOther-colUpper[other]) < 1.0e-5) {
    17421835                // delta to 1 => x to upper bound
    17431836                cutVector_[i].length++ ;
     
    17461839                  double value0 = elements[0] ;
    17471840                  double value1 = elements[1] ;
    1748                   if (value0*value1==-1.0) {
     1841                  if (value0*value1 == -1.0) {
    17491842                    // can do something ?
    1750                     int j=backward[other] ;
     1843                    // (lh) is x<t> binary? if not, j < 0 here.
     1844                    int j = backward[other] ;
    17511845                    cutVector_[i].length++ ;
    17521846                    cutVector_[j].length++ ;
     
    17591853          } else {
    17601854            if (onList[other]) {
    1761               if (elements[0]==1.0 && elements[1]==1.0&&rcut.ub()==1.0) {
     1855              if (elements[0] == 1.0 &&
     1856                  elements[1] == 1.0 && rcut.ub() == 1.0) {
    17621857                // can do something ?
    1763                 int j=backward[other] ;
     1858                // (lh) is x<t> binary? if not, j < 0 here.
     1859                int j = backward[other] ;
    17641860                cutVector_[i].length++ ;
    17651861                cutVector_[j].length++ ;
     
    17701866          }
    17711867        } else {
    1772           assert(rcut.ub()==DBL_MAX) ;
     1868/*
     1869  End (stuff) <= b, begin b <= (stuff). Note that this case is not symmetric
     1870  with the above case --- there's nothing for lb != 0.
     1871*/
     1872          assert(rcut.ub() == DBL_MAX) ;
    17731873          if (!lb) {
    17741874            // LB
    1775             if (elements[which]>0.0) {
    1776               //assert (elements[1-which]<0.0) ;
     1875            double elWhich = elements[which] ;
     1876            double elOther = elements[1-which] ;
     1877            if (elWhich > 0.0) {
     1878              // assert (elOther < 0.0) ;
    17771879              // delta to 0 => x to 0.0
    17781880              // flip so same as UB
    17791881              cutVector_[i].length++;
    17801882            } else {
    1781               if (elements[1-which]<0.0 && fabs(elements[which]/elements[1-which]-
    1782                                               colUpper[other])<1.0e-5) {
     1883              if (elOther < 0.0 &&
     1884                  fabs(elWhich/elOther-colUpper[other]) < 1.0e-5) {
    17831885                // delta to 1 => x to upper bound
    17841886                cutVector_[i].length++ ;
     
    17871889                  double value0 = elements[0] ;
    17881890                  double value1 = elements[1] ;
    1789                   if (value0*value1==-1.0) {
     1891                  if (value0*value1 == -1.0) {
    17901892                    // can do something ?
    1791                     int j=backward[other] ;
     1893                    int j = backward[other] ;
    17921894                    cutVector_[i].length++ ;
    17931895                    cutVector_[j].length++ ;
     
    17981900              }
    17991901            }
    1800           }
    1801         }
    1802       }
    1803       // allocate space
    1804       for (int i=0;i<number01Integers_;i++) {
    1805         int j=cutVector_[i].sequence ;
     1902          }        // end of b == 0
     1903        }       // end of b <= (stuff)
     1904      }      // end of counting loop
     1905/*
     1906  We have a set of counts attached to entries in cutVector_ that tell us how
     1907  long the index arrays will need to be. Allocate the arrays.
     1908*/
     1909      for (int i = 0 ; i < number01Integers_ ; i++) {
     1910        int j = cutVector_[i].sequence ;
    18061911        if (onList[j] && !cutVector_[i].index) {
    1807           disaggregation thisOne=cutVector_[i] ;
    1808           cutVector_[i].index=new disaggregationAction [thisOne.length] ;
    1809           cutVector_[i].length=0 ;
     1912          disaggregation thisOne = cutVector_[i] ;
     1913          cutVector_[i].index = new disaggregationAction [thisOne.length] ;
     1914          cutVector_[i].length = 0 ;
    18101915        }
    18111916      }
    1812       // now put in
    1813       for (iCut=0;iCut<nCuts;iCut++) {
     1917/*
     1918  And load up the arrays.
     1919
     1920  The classification structure here is cut-and-paste from above; all the
     1921  comments about the case analysis apply again. The Agree/Disagree forms
     1922  assume binary probe and target.   -- lh, 110212 --
     1923
     1924  TODO: Note that the index stored in disaggregation.index[] can be either an
     1925        index for cutVector_ (zeroOneInDisagg true) or the original
     1926        column index (false). Presence in cutVector_ implies binary; it's
     1927        still unclear to me if the reverse is always true. cutVector_
     1928        may not be loaded with every binary variable.  -- lh, 110212 --
     1929
     1930  TODO: What you can see, after a bit of thought, is that for a disaggregation
     1931        cut on two binary variables, there's a symmetry. For example,
     1932        x<p> -> 1 ==> x<t> -> 0 gives x<p> <= (1-x<t>), which also arises from
     1933        x<t> -> 1 ==> x<p> -> 0 (with the roles of x<p> and x<t> exchanged).
     1934        It looks like the classification scheme is trying to exploit this, but
     1935        not getting the decision tree right. We need to check for two binary
     1936        variables right at the top.  -- lh, 110212 --
     1937
     1938*/
     1939      for (iCut = 0 ; iCut < nCuts ; iCut++) {
    18141940        OsiRowCut rcut ;
    18151941        CoinPackedVector rpv ;
     
    18171943        rcut = csNew.rowCut(iCut) ;
    18181944        rpv = rcut.row() ;
    1819         assert(rpv.getNumElements()==2) ;
    1820         const int * indices = rpv.getIndices() ;
    1821         double* elements = rpv.getElements() ;
    1822         double lb=rcut.lb() ;
    1823         // find out which integer
    1824         // find out which integer
    1825         int which=0 ;
    1826         int i=backward[indices[0]] ;
    1827         if (i<0||!onList[indices[0]]) {
    1828           which=1 ;
    1829           i=backward[indices[1]] ;
    1830           // Just possible variable was general integer but now 0-1
    1831           if (!onList[indices[which]])
    1832             continue ;
     1945        assert(rpv.getNumElements() == 2) ;
     1946        const int *indices = rpv.getIndices() ;
     1947        double *elements = rpv.getElements() ;
     1948        double lb = rcut.lb() ;
     1949        // (lh) assume x<p> is the first binary variable we find?
     1950        int which = 0 ;
     1951        int i = backward[indices[0]] ;
     1952        if (i < 0) {
     1953          which = 1 ;
     1954          i = backward[indices[1]] ;
     1955          if (i < 0) continue ;
    18331956        }
     1957        disaggregation *cVi = &cutVector_[i] ;
    18341958        int other = indices[1-which] ;
    1835         int j = other ? backward[other] : -1 ;
    1836         if (lb==-DBL_MAX) {
    1837           if (!rcut.ub()) {
    1838             // UB
    1839             if (elements[which]<0.0) {
    1840               iput=cutVector_[i].length ;
    1841               if (j>=0)
    1842                 setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
     1959        int j = backward[other] ;
     1960        disaggregation *cVj = 0 ;
     1961        if (j >= 0) cVj = &cutVector_[j] ;
     1962
     1963        if (lb == -DBL_MAX) {   // <= constraint
     1964          if (!rcut.ub()) {     // b = 0
     1965            double elWhich = elements[which] ;
     1966            double elOther = elements[1-which] ;
     1967            if (elWhich < 0.0) {
     1968              // delta to 0 => x to 0.0
     1969              // Agree: -x<p>+x<t> <= 0
     1970              iput = cVi->length ;
     1971              if (j >= 0)
     1972                setAffectedInDisagg(cVi->index[iput],j) ;
    18431973              else
    1844                 setAffectedInDisaggregation(cutVector_[i].index[iput],other) ;
    1845               setWhenAtUBInDisaggregation(cutVector_[i].index[iput],false) ;
    1846               setAffectedToUBInDisaggregation(cutVector_[i].index[iput],false) ;
    1847               setZeroOneInDisaggregation(cutVector_[i].index[iput],onList[other]!=0) ;
    1848               cutVector_[i].length++ ;
     1974                setAffectedInDisagg(cVi->index[iput],other) ;
     1975              setWhenAtUBInDisagg(cVi->index[iput],false) ;
     1976              setAffectedToUBInDisagg(cVi->index[iput],false) ;
     1977              setZeroOneInDisagg(cVi->index[iput],(j >= 0)) ;
     1978              cVi->length++ ;
    18491979            } else {
    1850               if (elements[1-which]<0.0 && fabs(elements[which]/elements[1-which]-
    1851                                               colUpper[other])<1.0e-5) {
     1980              /*
     1981                Wrong classification check. Assume x<p> binary, x<t> general
     1982                integer, l<p> = l<t> = 0. l'<p> = 1. Canonical form is
     1983                (l<t>-l'<t>)x<p> + (l'<p>-l<p>)x<t> >= l<t>l'<p>-l'<t>l<p>
     1984                reduces to
     1985                l'<t>x<p> - x<t> <= 0
     1986                l'<t> is tightened colLower, equivalent to colUpper only for
     1987                binary variables or fixed general integer (shades of
     1988                goingToTrueBound). Also, we've carefully tested to ensure
     1989                that a<p> > 0 and a<t> < 0.  -- lh, 110212 --
     1990              */
     1991              if (elOther < 0.0 &&
     1992                  fabs(elWhich/elOther-colUpper[other]) < 1.0e-5) {
    18521993                // delta to 1 => x to upper bound
    1853                 iput=cutVector_[i].length ;
    1854                 if (j>=0)
    1855                   setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
     1994                // Agree: x<p>-x<t> <= 0
     1995                iput = cVi->length ;
     1996                if (j >= 0)
     1997                  setAffectedInDisagg(cVi->index[iput],j) ;
    18561998                else
    1857                   setAffectedInDisaggregation(cutVector_[i].index[iput],other) ;
    1858                 setWhenAtUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1859                 setAffectedToUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1860                 setZeroOneInDisaggregation(cutVector_[i].index[iput],onList[other]!=0) ;
    1861                 cutVector_[i].length++ ;
     1999                  setAffectedInDisagg(cVi->index[iput],other) ;
     2000                setWhenAtUBInDisagg(cVi->index[iput],true) ;
     2001                setAffectedToUBInDisagg(cVi->index[iput],true) ;
     2002                setZeroOneInDisagg(cVi->index[iput],(j >= 0)) ;
     2003                cVi->length++ ;
    18622004              } else {
    1863                 if (onList[other]) {
     2005                /*
     2006                  a<p> > 0, a<t> < 0, b == 0, x<p>  binary. Now add x<t>
     2007                  binary (j >= 0) and a<p>*a<t> = -1.0.  Effectively, we're
     2008                  looking for x<p> - x<t> <= 0. Sure, we can do something. We
     2009                  just did it, immediately above. Symmetry here says that
     2010                  x<p> -> 1 ==> x<t> -> 1 and x<t> -> 0 ==> x<p> -> 0
     2011                  -- lh, 110212 --
     2012                */
     2013                if (j >= 0) {
    18642014                  double value0 = elements[0] ;
    18652015                  double value1 = elements[1] ;
    1866                   if (value0*value1==-1.0) {
     2016                  if (value0*value1 == -1.0) {
    18672017                    // can do something ?
    1868                     int j=backward[other] ;
    1869                     assert (j>=0) ;
    18702018                    // flip so value0 1.0
    1871                     if (value1==1.0) {
    1872                       j=i ;
    1873                       i=backward[other] ;
    1874                       value1=value0 ;
    1875                       value0=1.0 ;
     2019                    if (value1 == 1.0) {
     2020                      j = i ;
     2021                      i = backward[other] ;
     2022                      value1 = value0 ;
     2023                      value0 = 1.0 ;
    18762024                    }
    1877                     assert (value0==1.0) ;
    1878                     assert (value1==-1.0) ;
    1879                     iput=cutVector_[i].length ;
    1880                     setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
    1881                     setWhenAtUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1882                     setAffectedToUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1883                     setZeroOneInDisaggregation(cutVector_[i].index[iput],true) ;
    1884                     cutVector_[i].length++ ;
    1885                     iput=cutVector_[j].length ;
    1886                     setAffectedInDisaggregation(cutVector_[j].index[iput],i) ;
    1887                     setWhenAtUBInDisaggregation(cutVector_[j].index[iput],false) ;
    1888                     setAffectedToUBInDisaggregation(cutVector_[j].index[iput],false) ;
    1889                     setZeroOneInDisaggregation(cutVector_[j].index[iput],true) ;
    1890                     cutVector_[j].length++ ;
     2025                    assert (value0 == 1.0) ;
     2026                    assert (value1 == -1.0) ;
     2027                    iput = cVi->length ;
     2028                    setAffectedInDisagg(cVi->index[iput],j) ;
     2029                    setWhenAtUBInDisagg(cVi->index[iput],true) ;
     2030                    setAffectedToUBInDisagg(cVi->index[iput],true) ;
     2031                    setZeroOneInDisagg(cVi->index[iput],true) ;
     2032                    cVi->length++ ;
     2033                    iput = cVj->length ;
     2034                    setAffectedInDisagg(cVj->index[iput],i) ;
     2035                    setWhenAtUBInDisagg(cVj->index[iput],false) ;
     2036                    setAffectedToUBInDisagg(cVj->index[iput],false) ;
     2037                    setZeroOneInDisagg(cVj->index[iput],true) ;
     2038                    cVj->length++ ;
    18912039                  }
    18922040                }
     
    18942042            }
    18952043          } else {
     2044            /*
     2045              stuff <= b, b != 0. Additionally require x<t> binary,
     2046              a<p> = a<t> = 1.0, b = 1. So we're looking for x<p> + x<t> <= 1.
     2047              This is the binary disaggregation for x<p> to 1 ==> x<t> to 0,
     2048              and x<t> -> 1 ==> x<p> -> 0.
     2049            */
    18962050            if (onList[other]) {
    1897               if (elements[0]==1.0 && elements[1]==1.0&&rcut.ub()==1.0) {
     2051              if (elements[0] == 1.0 &&
     2052                  elements[1] == 1.0 && rcut.ub() == 1.0) {
    18982053                // can do something ?
    1899                 int j=backward[other] ;
    1900                 assert (j>=0) ;
    1901                 iput=cutVector_[i].length ;
    1902                 setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
    1903                 setWhenAtUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1904                 setAffectedToUBInDisaggregation(cutVector_[i].index[iput],false) ;
    1905                 setZeroOneInDisaggregation(cutVector_[i].index[iput],true) ;
    1906                 cutVector_[i].length++ ;
    1907                 iput=cutVector_[j].length ;
    1908                 setAffectedInDisaggregation(cutVector_[j].index[iput],i) ;
    1909                 setWhenAtUBInDisaggregation(cutVector_[j].index[iput],true) ;
    1910                 setAffectedToUBInDisaggregation(cutVector_[j].index[iput],false) ;
    1911                 setZeroOneInDisaggregation(cutVector_[j].index[iput],true) ;
    1912                 cutVector_[j].length++ ;
     2054                int j = backward[other] ;
     2055                assert ( j >= 0) ;
     2056                iput = cVi->length ;
     2057                setAffectedInDisagg(cVi->index[iput],j) ;
     2058                setWhenAtUBInDisagg(cVi->index[iput],true) ;
     2059                setAffectedToUBInDisagg(cVi->index[iput],false) ;
     2060                setZeroOneInDisagg(cVi->index[iput],true) ;
     2061                cVi->length++ ;
     2062                iput = cVj->length ;
     2063                setAffectedInDisagg(cVj->index[iput],i) ;
     2064                setWhenAtUBInDisagg(cVj->index[iput],true) ;
     2065                setAffectedToUBInDisagg(cVj->index[iput],false) ;
     2066                setZeroOneInDisagg(cVj->index[iput],true) ;
     2067                cVj->length++ ;
    19132068              } else {
    19142069#ifdef COIN_DEVELOP
     
    19202075          }
    19212076        } else {
    1922           assert(rcut.ub()==DBL_MAX) ;
     2077/*
     2078  More of the same, this time for constraints of the form b <= (stuff).
     2079*/
     2080          assert(rcut.ub() == DBL_MAX) ;
    19232081          if (!lb) {
    19242082            // LB
    1925             if (elements[which]>0.0) {
    1926               iput=cutVector_[i].length ;
    1927               if (j>=0)
    1928                 setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
     2083            double elWhich = elements[which] ;
     2084            double elOther = elements[1-which] ;
     2085            if (elWhich > 0.0) {
     2086              iput = cVi->length ;
     2087              if (j >= 0)
     2088                setAffectedInDisagg(cVi->index[iput],j) ;
    19292089              else
    1930                 setAffectedInDisaggregation(cutVector_[i].index[iput],other) ;
    1931               setWhenAtUBInDisaggregation(cutVector_[i].index[iput],false) ;
    1932               setAffectedToUBInDisaggregation(cutVector_[i].index[iput],false) ;
    1933               setZeroOneInDisaggregation(cutVector_[i].index[iput],onList[other]!=0) ;
    1934               cutVector_[i].length++ ;
     2090                setAffectedInDisagg(cVi->index[iput],other) ;
     2091              setWhenAtUBInDisagg(cVi->index[iput],false) ;
     2092              setAffectedToUBInDisagg(cVi->index[iput],false) ;
     2093              setZeroOneInDisagg(cVi->index[iput],onList[other]!=0) ;
     2094              cVi->length++ ;
    19352095            } else {
    1936               if (elements[1-which]<0.0 && fabs(elements[which]/elements[1-which]-
    1937                                               colUpper[other])<1.0e-5) {
    1938                 iput=cutVector_[i].length ;
    1939                 if (j>=0)
    1940                   setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
     2096              if (elOther < 0.0 &&
     2097                  fabs(elWhich/elOther-colUpper[other]) < 1.0e-5) {
     2098                iput = cVi->length ;
     2099                if (j >= 0)
     2100                  setAffectedInDisagg(cVi->index[iput],j) ;
    19412101                else
    1942                   setAffectedInDisaggregation(cutVector_[i].index[iput],other) ;
    1943                 setWhenAtUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1944                 setAffectedToUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1945                 setZeroOneInDisaggregation(cutVector_[i].index[iput],onList[other]!=0) ;
    1946                 cutVector_[i].length++ ;
     2102                  setAffectedInDisagg(cVi->index[iput],other) ;
     2103                setWhenAtUBInDisagg(cVi->index[iput],true) ;
     2104                setAffectedToUBInDisagg(cVi->index[iput],true) ;
     2105                setZeroOneInDisagg(cVi->index[iput],onList[other]!=0) ;
     2106                cVi->length++ ;
    19472107              } else {
    19482108                if (onList[other]) {
    19492109                  double value0 = elements[0] ;
    19502110                  double value1 = elements[1] ;
    1951                   if (value0*value1==-1.0) {
     2111                  if (value0*value1 == -1.0) {
    19522112                    // can do something ?
    1953                     int j=backward[other] ;
     2113                    int j = backward[other] ;
    19542114                    assert (j>=0) ;
    19552115                    // flip so value0 -1.0
    1956                     if (value1==-1.0) {
    1957                       j=i ;
    1958                       i=backward[other] ;
    1959                       value1=value0 ;
    1960                       value0=-1.0 ;
     2116                    if (value1 == -1.0) {
     2117                      j = i ;
     2118                      i = backward[other] ;
     2119                      value1 = value0 ;
     2120                      value0 = -1.0 ;
    19612121                    }
    1962                     assert (value0==-1.0) ;
    1963                     assert (value1==1.0) ;
    1964                     iput=cutVector_[i].length ;
    1965                     setAffectedInDisaggregation(cutVector_[i].index[iput],j) ;
    1966                     setWhenAtUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1967                     setAffectedToUBInDisaggregation(cutVector_[i].index[iput],true) ;
    1968                     setZeroOneInDisaggregation(cutVector_[i].index[iput],true) ;
    1969                     cutVector_[i].length++ ;
    1970                     iput=cutVector_[j].length ;
    1971                     setAffectedInDisaggregation(cutVector_[j].index[iput],i) ;
    1972                     setWhenAtUBInDisaggregation(cutVector_[j].index[iput],false) ;
    1973                     setAffectedToUBInDisaggregation(cutVector_[j].index[iput],false) ;
    1974                     setZeroOneInDisaggregation(cutVector_[j].index[iput],true) ;
    1975                     cutVector_[j].length++ ;
     2122                    assert (value0 == -1.0) ;
     2123                    assert (value1 == 1.0) ;
     2124                    iput = cVi->length ;
     2125                    setAffectedInDisagg(cVi->index[iput],j) ;
     2126                    setWhenAtUBInDisagg(cVi->index[iput],true) ;
     2127                    setAffectedToUBInDisagg(cVi->index[iput],true) ;
     2128                    setZeroOneInDisagg(cVi->index[iput],true) ;
     2129                    cVi->length++ ;
     2130                    iput = cVj->length ;
     2131                    setAffectedInDisagg(cVj->index[iput],i) ;
     2132                    setWhenAtUBInDisagg(cVj->index[iput],false) ;
     2133                    setAffectedToUBInDisagg(cVj->index[iput],false) ;
     2134                    setZeroOneInDisagg(cVj->index[iput],true) ;
     2135                    cVj->length++ ;
    19762136                  }
    19772137                }
     
    19822142      }
    19832143      delete [] backward ;
    1984       // Now sort and get rid of duplicates
    1985       // could also see if any are cliques
    1986       int longest=0 ;
    1987       for (int i=0;i<number01Integers_;i++)
    1988         longest = CoinMax(longest, cutVector_[i].length) ;
    1989       unsigned int * sortit = new unsigned int[longest] ;
    1990       for (int i=0;i<number01Integers_;i++) {
    1991         disaggregation & thisOne=cutVector_[i] ;
     2144/*
     2145  Note that onList is just a pointer into the block of space headed by
     2146  backward.
     2147
     2148  END CLASSIFICATION
     2149
     2150  We know have a vector of indices attached to each entry in cutVector_,
     2151  with bits set to indicate just what sort of disaggregation cut we're
     2152  dealing with.
     2153
     2154  jjf: Now sort and get rid of duplicates;  could also see if any are cliques
     2155
     2156  To do the sort, compose a unique signature for each disaggregation entry, as
     2157  (index)@(status bits), and sort the resulting vector. Then compare adjacent
     2158  entries, squeezing out duplicates.
     2159*/
     2160      int longest = 0 ;
     2161      for (int i = 0 ; i < number01Integers_ ; i++)
     2162        longest = CoinMax(longest,cutVector_[i].length) ;
     2163      unsigned int *sortit = new unsigned int[longest] ;
     2164
     2165      for (int i = 0 ; i < number01Integers_ ; i++) {
     2166
     2167        disaggregation &thisOne = cutVector_[i] ;
    19922168        int k ;
    19932169        int number = thisOne.length ;
    1994         for (k=0;k<number;k++) {
    1995           int affected = affectedInDisaggregation(thisOne.index[k]) ;
    1996           int zeroOne = zeroOneInDisaggregation(thisOne.index[k]) ? 1 : 0 ;
    1997           int whenAtUB = whenAtUBInDisaggregation(thisOne.index[k]) ? 1 : 0 ;
    1998           int affectedToUB = affectedToUBInDisaggregation(thisOne.index[k]) ? 1: 0 ;
    1999           sortit[k]=(affected<<3)|(zeroOne<<2)|(whenAtUB<<1)|affectedToUB ;
     2170        for (k = 0 ; k < number ; k++) {
     2171          int affected = affectedInDisagg(thisOne.index[k]) ;
     2172          int zeroOne = zeroOneInDisagg(thisOne.index[k]) ? 1 : 0 ;
     2173          int whenAtUB = whenAtUBInDisagg(thisOne.index[k]) ? 1 : 0 ;
     2174          int affectedToUB = affectedToUBInDisagg(thisOne.index[k]) ? 1 : 0 ;
     2175          sortit[k] = (affected<<3)|(zeroOne<<2)|(whenAtUB<<1)|affectedToUB ;
    20002176        }
    20012177        std::sort(sortit,sortit+number) ;
     2178
    20022179        int affectedLast = 0xffffffff ;
    20032180        int zeroOneLast = 0 ;
    20042181        int whenAtUBLast = 0 ;
    2005         int affectedToUBLast = 0;
    2006         int put=0 ;
    2007         for (k=0;k<number;k++) {
     2182        int affectedToUBLast = 0 ;
     2183        int put = 0 ;
     2184        for (k = 0 ; k < number ; k++) {
    20082185          int affected = sortit[k]>>3 ;
    20092186          int zeroOne = (sortit[k]&4)>>2 ;
     
    20112188          int affectedToUB = sortit[k]&1 ;
    20122189          disaggregationAction action ;
    2013           action.affected=0 ;
    2014           setAffectedInDisaggregation(action,affected) ;
    2015           setZeroOneInDisaggregation(action,zeroOne!=0) ;
    2016           setWhenAtUBInDisaggregation(action,whenAtUB!=0) ;
    2017           setAffectedToUBInDisaggregation(action,affectedToUB!=0) ;
    2018           if (affected!=affectedLast||zeroOne!=zeroOneLast) {
     2190          action.affected = 0 ;
     2191          setAffectedInDisagg(action,affected) ;
     2192          setZeroOneInDisagg(action,(zeroOne != 0)) ;
     2193          setWhenAtUBInDisagg(action,(whenAtUB != 0)) ;
     2194          setAffectedToUBInDisagg(action,(affectedToUB != 0)) ;
     2195          if (affected != affectedLast || zeroOne != zeroOneLast) {
    20192196            // new variable
    2020             thisOne.index[put++]=action ;
    2021           } else if (whenAtUB!=whenAtUBLast||affectedToUB!=affectedToUBLast) {
     2197            thisOne.index[put++] = action ;
     2198          } else if (whenAtUB != whenAtUBLast ||
     2199                     affectedToUB != affectedToUBLast) {
    20222200            // new action - what can we discover
    2023             thisOne.index[put++]=action ;
    2024             int j=cutVector_[i].sequence ;
    2025             int k=affected ;
     2201            thisOne.index[put++] = action ;
     2202            int j = cutVector_[i].sequence ;
     2203            int k = affected ;
    20262204            if (zeroOne) {
    2027               k=cutVector_[k].sequence ;
    2028               if (logLevel_>1)
    2029                 printf("For %d %d 0-1 pair",j,k) ;
     2205              k = cutVector_[k].sequence ;
     2206              if (logLevel_ > 1)
     2207                std::cout
     2208                  << "  for 0-1 pair " << j << " and "
     2209                  << k << ":" << std::endl ;
    20302210            } else {
    2031               if (logLevel_>1)
    2032                 printf("For %d %d pair",j,k) ;
     2211              if (logLevel_ > 1)
     2212                std::cout
     2213                  << "  for pair " << j << " and " << k << ":" << std::endl ;
    20332214            }
    2034             if (logLevel_>1)
    2035               printf(" old whenAtUB, affectedToUB %d %d, new whenAtUB, affectedToUB %d %d\n",
    2036                      whenAtUBLast, affectedToUBLast,whenAtUB, affectedToUB) ;
     2215            if (logLevel_ > 1)
     2216              std::cout
     2217                << "    previous: " << ((whenAtUBLast)?"up":"down") << " probe "
     2218                << "tightened " << ((affectedToUBLast)?"lower":"upper")
     2219                << " bound" << std::endl ;
     2220              std::cout
     2221                << "    current: " << ((whenAtUB)?"up":"down") << " probe "
     2222                << "tightened " << ((affectedToUB)?"lower":"upper")
     2223                << " bound" << std::endl ;
    20372224          }
    2038           affectedLast=affected ;
    2039           zeroOneLast=zeroOne ;
    2040           whenAtUBLast=whenAtUB ;
    2041           affectedToUBLast=affectedToUB ;
     2225          affectedLast = affected ;
     2226          zeroOneLast = zeroOne ;
     2227          whenAtUBLast = whenAtUB ;
     2228          affectedToUBLast = affectedToUB ;
    20422229        }
    2043         if (put<number) {
    2044           //printf("%d reduced from %d to %d\n",i,number,put) ;
    2045           thisOne.length=put ;
     2230        if (put < number) {
     2231          if (logLevel_ > 1)
     2232            std::cout
     2233              << "  actions for " << i << "reduced from " << number
     2234              << " to " << put << "." << std::endl ;
     2235          thisOne.length = put ;
    20462236        }
    20472237      }
    2048       // And look at all where two 0-1 variables involved
    2049       for (int i=0;i<number01Integers_;i++) {
    2050         disaggregation & thisOne=cutVector_[i] ;
    2051         int k ;
     2238/*
     2239  End sort. Now look over the lists of disaggregation cuts and see what we can
     2240  learn from the ones where both x<p> and x<t> are binary.
     2241
     2242  TODO: Unless I totally miss my guess here, we're about to derive implication
     2243        cuts. And indeed, if I'd read ahead 20 lines, that's exactly what
     2244        we're up to. Not to mention looking again for monotone actions and
     2245        disaggregation cuts. So the real question is `Why aren't we doing all
     2246        this by calling the appropriate methods in probeClique?
     2247
     2248        What could be done here is to try and extend cliques. But there's only
     2249        a placeholder case statement. Which is consistent with the comment up
     2250        in CglPreProcess at the call to snapshot: `out for now - think about
     2251        cliques'.
     2252       
     2253        -- lh, 110212 --
     2254*/
     2255      for (int i = 0 ; i < number01Integers_ ; i++) {
     2256        disaggregation &thisOne = cutVector_[i] ;
    20522257        int number = thisOne.length ;
    2053         for (k=0;k<number;k++) {
    2054           int affected = affectedInDisaggregation(thisOne.index[k]) ;
    2055           bool zeroOne = zeroOneInDisaggregation(thisOne.index[k]) ;
    2056           if (zeroOne && static_cast<int>(affected)>i) {
    2057             bool whenAtUB = whenAtUBInDisaggregation(thisOne.index[k]) ;
    2058             bool affectedToUB = affectedToUBInDisaggregation(thisOne.index[k]) ;
    2059             disaggregation otherOne=cutVector_[affected] ;
     2258        for (int k = 0 ; k < number ; k++) {
     2259          int affected = affectedInDisagg(thisOne.index[k]) ;
     2260          bool zeroOne = zeroOneInDisagg(thisOne.index[k]) ;
     2261          /*
     2262            affected < i means we've already processed the equivalent cut
     2263            when we processed the list associated with affected, with the
     2264            roles of affected and i exchanged.
     2265          */
     2266          if (zeroOne && affected > i) {
     2267            bool whenAtUB = whenAtUBInDisagg(thisOne.index[k]) ;
     2268            bool affectedToUB = affectedToUBInDisagg(thisOne.index[k]) ;
     2269            disaggregation otherOne = cutVector_[affected] ;
    20602270            int numberOther = otherOne.length ;
    20612271            // Could do binary search if a lot
    2062             int lastAction=-1 ;
    2063             for (int j=0;j<numberOther;j++) {
    2064               if (affectedInDisaggregation(otherOne.index[j])==i) {
    2065                 bool whenAtUBOther = whenAtUBInDisaggregation(otherOne.index[j]) ;
    2066                 bool affectedToUBOther = affectedToUBInDisaggregation(otherOne.index[j]) ;
    2067                 /* action -
     2272            int lastAction = -1 ;
     2273            for (int j = 0 ; j < numberOther ; j++) {
     2274              if (affectedInDisagg(otherOne.index[j])==i) {
     2275                bool whenAtUBOther = whenAtUBInDisagg(otherOne.index[j]) ;
     2276                bool affectedToUBOther = affectedToUBInDisagg(otherOne.index[j]) ;
     2277                /* ACTIONS
     2278               
    20682279                   0 -> x + y <=1 (1,1 impossible)
    20692280                   1 -> x - y <=0 (1,0 impossible)
     
    20762287                  22 -> y == 0
    20772288                  23 -> y == 1
     2289
     2290                  (lh) 0 -- 3 are standard disaggregation forms. 10 -- 11 are
     2291                  implication cuts. 20 -- 23 are caught as monotone actions or
     2292                  implication cuts. We're just treading old ground. If this
     2293                  code discovers anything, it indicates an algorithm error
     2294                  elsewhere.  -- lh, 110212 --
    20782295                */
    20792296                int action=-1 ;
     
    22232440      }
    22242441      delete [] sortit ;
    2225     }
     2442    }           // end AFTER PROBECLIQUE
     2443/*
     2444  End of the block that starts with if (!ninfeas), right after the call to
     2445  probeClique.
     2446
     2447  Now have a look at the cuts in cutVector_ and see if any are violated. Only
     2448  look at the entries for unfixed variables, and then only if the current
     2449  solution value is not integral.
     2450
     2451  TODO: The immediate question is ``Why are we doing this if we're
     2452        infeasible?' I suppose you can argue that cutVector_ will be mostly
     2453        empty (we haven't loaded in the disaggregation cuts), but it will be
     2454        present. And didn't we check these cuts for violation when we created
     2455        them? Apparently not -- a quick look ahead says we're finally going to
     2456        add the cut to cs if we find it's violated.   -- lh, 110212 --
     2457
     2458  TODO: Good grief! We're about to reconstitute the cuts that we processed
     2459        down into the entries in cutVector_. And we're not even going to do
     2460        a complete job of it. Surely there must be a better way.
     2461        -- lh, 110212 --
     2462*/
    22262463    if (cutVector_) {
    2227       // now see if any disaggregation cuts are violated
    2228       for (int i=0;i<number01Integers_;i++) {
    2229         int j=cutVector_[i].sequence ;
    2230         double solInt=colsol[j] ;
    2231         double  upper, solValue ;
    2232         int icol ;
     2464      for (int i = 0 ; i < number01Integers_ ; i++) {
     2465        int j = cutVector_[i].sequence ;
     2466        double solInt = colsol[j] ;
     2467        double  upper ;
     2468        double solValue ;
     2469        int icol = -1 ;
    22332470        int index[2] ;
    22342471        double element[2] ;
    2235         if (colUpper[j]-colLower[j]>1.0e-8) {
     2472        if (colUpper[j]-colLower[j] > 1.0e-8) {
    22362473          double away = fabs(0.5-(solInt-floor(solInt))) ;
    2237           if (away<0.4999999) {
    2238             disaggregation thisOne=cutVector_[i] ;
     2474          if (away < 0.4999999) {
     2475            disaggregation thisOne = cutVector_[i] ;
    22392476            int k ;
    22402477            OsiRowCut rc ;
    2241             for (k=0;k<thisOne.length;k++) {
    2242               icol = affectedInDisaggregation(thisOne.index[k]) ;
    2243               if (zeroOneInDisaggregation(thisOne.index[k]))
     2478            for (k = 0 ; k < thisOne.length ; k++) {
     2479              icol = affectedInDisagg(thisOne.index[k]) ;
     2480              if (zeroOneInDisagg(thisOne.index[k]))
    22442481                icol = cutVector_[icol].sequence ;
    2245               solValue=colsol[icol] ;
    2246               upper=colUpper_[icol] ;
    2247               double infeasibility=0.0 ;
    2248               if (!whenAtUBInDisaggregation(thisOne.index[k])) {
    2249                 if (!affectedToUBInDisaggregation(thisOne.index[k])) {
     2482              solValue = colsol[icol] ;
     2483              upper = colUpper_[icol] ;
     2484              double infeasibility = 0.0 ;
     2485              if (!whenAtUBInDisagg(thisOne.index[k])) {
     2486                if (!affectedToUBInDisagg(thisOne.index[k])) {
    22502487                  // delta -> 0 => x to lb (at present just 0)
    22512488                  infeasibility = solValue - upper * solInt ;
     
    22532490                    rc.setLb(-DBL_MAX) ;
    22542491                    rc.setUb(0.0) ;
    2255                     index[0]=icol ;
    2256                     element[0]=1.0 ;
    2257                     index[1]=j ;
    2258                     element[1]= -upper ;
     2492                    index[0] = icol ;
     2493                    element[0] = 1.0 ;
     2494                    index[1] = j ;
     2495                    element[1] = -upper ;
    22592496                  } else {
    2260                     infeasibility=0.0 ;
     2497                    infeasibility = 0.0 ;
    22612498                  }
    22622499                } else {
     
    22652502                }
    22662503              } else {
    2267                 if (affectedToUBInDisaggregation(thisOne.index[k])) {
     2504                if (affectedToUBInDisagg(thisOne.index[k])) {
    22682505                  // delta -> 1 => x to ub (?)
    2269                   icol = affectedInDisaggregation(thisOne.index[k]) ;
    2270                   if (zeroOneInDisaggregation(thisOne.index[k]))
     2506                  icol = affectedInDisagg(thisOne.index[k]) ;
     2507                  if (zeroOneInDisagg(thisOne.index[k]))
    22712508                    icol = cutVector_[icol].sequence ;
    2272                   solValue=colsol[icol] ;
    2273                   upper=colUpper_[icol] ;
     2509                  solValue = colsol[icol] ;
     2510                  upper = colUpper_[icol] ;
    22742511                  if (!colLower[icol]) {
    22752512                    infeasibility = upper * solInt - solValue ;
     
    22772514                      rc.setLb(-DBL_MAX) ;
    22782515                      rc.setUb(0.0) ;
    2279                       index[0]=icol ;
    2280                       element[0]=-1.0 ;
    2281                       index[1]=j ;
    2282                       element[1]= upper ;
     2516                      index[0] = icol ;
     2517                      element[0] = -1.0 ;
     2518                      index[1] = j ;
     2519                      element[1] = upper ;
    22832520                    } else {
    2284                       infeasibility=0.0 ;
     2521                      infeasibility = 0.0 ;
    22852522                    }
    22862523                  } else {
    2287                     assert (upper==colLower[icol]) ;
    2288                     infeasibility=0.0 ;
     2524                    assert (upper == colLower[icol]) ;
     2525                    infeasibility = 0.0 ;
    22892526                  }
    22902527                } else {
    22912528                  // delta + delta2 <= 1
    2292                   assert (zeroOneInDisaggregation(thisOne.index[k])) ;
     2529                  assert (zeroOneInDisagg(thisOne.index[k])) ;
    22932530                  // delta -> 1 => delta2 -> 0
    2294                   icol = affectedInDisaggregation(thisOne.index[k]) ;
     2531                  icol = affectedInDisagg(thisOne.index[k]) ;
    22952532                  icol = cutVector_[icol].sequence ;
    22962533                  // only do if icol > j
    2297                   if (icol >j && colUpper[icol] ) {
    2298                     solValue=colsol[icol] ;
     2534                  if (icol > j && colUpper[icol] ) {
     2535                    solValue = colsol[icol] ;
    22992536                    if (!colLower[icol]) {
    23002537                      infeasibility = solInt + solValue - 1.0 ;
     
    23022539                        rc.setLb(-DBL_MAX) ;
    23032540                        rc.setUb(1.0) ;
    2304                         index[0]=icol ;
    2305                         element[0]=1.0 ;
    2306                         index[1]=j ;
    2307                         element[1]= 1.0 ;
     2541                        index[0] = icol ;
     2542                        element[0] = 1.0 ;
     2543                        index[1] = j ;
     2544                        element[1] = 1.0 ;
    23082545                      } else {
    2309                         infeasibility=0.0 ;
     2546                        infeasibility = 0.0 ;
    23102547                      }
    23112548                    } else {
    2312                       assert (upper==colLower[icol]) ;
    2313                       infeasibility=0.0 ;
     2549                      assert (upper == colLower[icol]) ;
     2550                      infeasibility = 0.0 ;
    23142551                    }
    23152552                  }
     
    23212558                if (logLevel_>1)
    23222559                  printf("%g <= %g * x%d + %g * x%d <= %g\n",
    2323                          rc.lb(),element[0],index[0],element[1],index[1],rc.ub()) ;
     2560                         rc.lb(),element[0],index[0],element[1],
     2561                         index[1],rc.ub()) ;
    23242562#if CGL_DEBUG > 0
    23252563                if (debugger) assert(!debugger->invalidCut(rc));
     
    23972635
    23982636/*
    2399   Parameters useObj, useCutoff, and cutoff are unused as of 110208. Added to
    2400   probe(), added here for symmetry.
    2401 */
    2402 
    2403 // Does probing and adding cuts
    2404 int CglProbing::probeCliques( const OsiSolverInterface & si,
    2405                               const OsiRowCutDebugger *
    2406 #if CGL_DEBUG > 0
    2407                          debugger
    2408 #endif
    2409                               ,OsiCuts & cs,
    2410                               double * colLower, double * colUpper,
    2411                        CoinPackedMatrix *rowCopy,
    2412                               CoinPackedMatrix *columnCopy, const int * realRows,
    2413                        double * rowLower, double * rowUpper,
    2414                        char * intVar, double * minR, double * maxR,
    2415                        int * markR,
    2416                        CglTreeInfo * info,
    2417                        bool useObj, bool useCutoff, double cutoff) const
    2418 {
    2419   // Set up maxes
    2420   int maxStack = info->inTree ? maxStack_ : maxStackRoot_;
    2421   int nRows=rowCopy->getNumRows();
    2422   int nCols=rowCopy->getNumCols();
    2423   double * colsol = new double[nCols];
    2424   double * djs = new double[nCols];
    2425   const double * currentColLower = si.getColLower();
    2426   const double * currentColUpper = si.getColUpper();
    2427   double * tempL = new double [nCols];
    2428   double * tempU = new double [nCols];
    2429   int * markC = new int [nCols];
    2430   int * stackC = new int [2*nCols];
    2431   int * stackR = new int [nRows];
    2432   double * saveL = new double [2*nCols];
    2433   double * saveU = new double [2*nCols];
    2434   double * saveMin = new double [nRows];
    2435   double * saveMax = new double [nRows];
    2436   double * element = new double[nCols];
    2437   int * index = new int[nCols];
    2438   // For trying to extend cliques
    2439   int * cliqueStack=NULL;
    2440   int * cliqueCount=NULL;
    2441   int * to_01=NULL;
    2442   if (!mode_) {
    2443     to_01 = new int[nCols];
    2444     cliqueStack = new int[numberCliques_];
    2445     cliqueCount = new int[numberCliques_];
    2446     int i;
    2447     for (i=0;i<numberCliques_;i++) {
    2448       cliqueCount[i]=cliqueStart_[i+1]-cliqueStart_[i];
    2449     }
    2450     for (i=0;i<nCols;i++)
    2451       to_01[i]=-1;
    2452     for (i=0;i<number01Integers_;i++) {
    2453       int j=cutVector_[i].sequence;
    2454       to_01[j]=i;
    2455     }
    2456   }
    2457   // Let us never add more than twice the number of rows worth of row cuts
    2458   // Keep cuts out of cs until end so we can find duplicates quickly
    2459   int nRowsFake = info->inTree ? nRows/3 : nRows;
    2460   CglProbingRowCut rowCut(nRowsFake, !info->inTree);
    2461   const int * column = rowCopy->getIndices();
    2462   const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    2463   const int * rowLength = rowCopy->getVectorLengths();
    2464   const double * rowElements = rowCopy->getElements();
    2465   const int * row = columnCopy->getIndices();
    2466   const CoinBigIndex * columnStart = columnCopy->getVectorStarts();
    2467   const int * columnLength = columnCopy->getVectorLengths();
    2468   const double * columnElements = columnCopy->getElements();
    2469   double movement;
    2470   int i, j, k,kk,jj;
    2471   int kcol,krow;
    2472   bool anyColumnCuts=false;
    2473   double dbound, value, value2;
    2474   int ninfeas=0;
    2475   int rowCuts;
    2476   double disaggEffectiveness;
    2477   if (mode_) {
    2478     /* clean up djs and solution */
    2479     CoinMemcpyN(si.getReducedCost(),nCols,djs);
    2480     CoinMemcpyN( si.getColSolution(),nCols,colsol);
    2481     disaggEffectiveness=1.0e-3;
    2482     rowCuts=rowCuts_;
    2483   } else {
    2484     // need to go from a neutral place
    2485     memset(djs,0,nCols*sizeof(double));
    2486     CoinMemcpyN( si.getColSolution(),nCols,colsol);
    2487     disaggEffectiveness=-1.0e10;
    2488     if (rowCuts_!=4)
    2489       rowCuts=1;
    2490     else
    2491       rowCuts=4;
    2492   }
    2493   for (i = 0; i < nCols; ++i) {
    2494     /* was if (intVar[i]) */
    2495     if (1) {
    2496       if (colUpper[i]-colLower[i]>1.0e-8) {
    2497         if (colsol[i]<colLower[i]+primalTolerance_) {
    2498           colsol[i]=colLower[i];
    2499           djs[i] = CoinMax(0.0,djs[i]);
    2500         } else if (colsol[i]>colUpper[i]-primalTolerance_) {
    2501           colsol[i]=colUpper[i];
    2502           djs[i] = CoinMin(0.0,djs[i]);
    2503         } else {
    2504           djs[i]=0.0;
    2505         }
    2506         /*if (fabs(djs[i])<1.0e-5)
    2507           djs[i]=0.0;*/
    2508       }
    2509     }
    2510   }
    2511 
    2512   int ipass=0,nfixed=-1;
    2513 /*
    2514   Chop now that useCutoff and cutoff come in as parameters.
    2515   double cutoff;
    2516   bool cutoff_available = si.getDblParam(OsiDualObjectiveLimit,cutoff);
    2517   if (!cutoff_available||usingObjective_<0) { // cut off isn't set or isn't valid
    2518     cutoff = si.getInfinity();
    2519   }
    2520   cutoff *= si.getObjSense();
    2521   if (fabs(cutoff)>1.0e30)
    2522     assert (cutoff>1.0e30);
    2523 */
    2524   double current = si.getObjValue();
    2525   // make irrelevant if mode is 0
    2526   if (!mode_)
    2527     cutoff=DBL_MAX;
    2528   /* for both way coding */
    2529   int nstackC0=-1;
    2530   int * stackC0 = new int[maxStack];
    2531   double * lo0 = new double[maxStack];
    2532   double * up0 = new double[maxStack];
    2533   int nstackR,nstackC;
    2534   for (i=0;i<nCols;i++) {
    2535     if (colUpper[i]-colLower[i]<1.0e-8) {
    2536       markC[i]=3;
    2537     } else {
    2538       markC[i]=0;
    2539     }
    2540   }
    2541   double tolerance = 1.0e1*primalTolerance_;
    2542   int maxPass = info->inTree ? maxPass_ : maxPassRoot_;
    2543   // If we are going to replace coefficient then we don't need to be effective
    2544   double needEffectiveness = info->strengthenRow ? -1.0e10 : 1.0e-3;
    2545   while (ipass<maxPass&&nfixed) {
    2546     int iLook;
    2547     ipass++;
    2548     nfixed=0;
    2549     for (iLook=0;iLook<numberThisTime_;iLook++) {
    2550       double solval;
    2551       double down;
    2552       double up;
    2553       int awayFromBound=1;
    2554       j=lookedAt_[iLook];
    2555       solval=colsol[j];
    2556       down = floor(solval+tolerance);
    2557       up = ceil(solval-tolerance);
    2558       if(colUpper[j]-colLower[j]<1.0e-8) markC[j]=3;
    2559       if (markC[j]||!intVar[j]) continue;
    2560       double saveSolval = solval;
    2561       if (solval>=colUpper[j]-tolerance||solval<=colLower[j]+tolerance||up==down) {
    2562         awayFromBound=0;
    2563         if (solval<=colLower[j]+2.0*tolerance) {
    2564           solval = colLower[j]+1.0e-1;
    2565           down=colLower[j];
    2566           up=down+1.0;
    2567         } else if (solval>=colUpper[j]-2.0*tolerance) {
    2568           solval = colUpper[j]-1.0e-1;
    2569           up=colUpper[j];
    2570           down=up-1;
    2571         } else {
    2572           // odd
    2573           up=down+1.0;
    2574           solval = down+1.0e-1;
    2575         }
    2576       }
    2577       assert (up<=colUpper[j]);
    2578       assert (down>=colLower[j]);
    2579       assert (up>down);
    2580       if ((solval-down>1.0e-6&&up-solval>1.0e-6)||mode_!=1) {
    2581         int istackC,iway, istackR;
    2582         int way[]={1,2,1};
    2583         int feas[]={1,2,4};
    2584         int feasible=0;
    2585         int notFeasible;
    2586         for (iway=0;iway<3;iway ++) {
    2587           int fixThis=0;
    2588           double objVal=current;
    2589           int goingToTrueBound=0;
    2590           stackC[0]=j;
    2591           markC[j]=way[iway];
    2592           double solMovement;
    2593           if (way[iway]==1) {
    2594             movement=down-colUpper[j];
    2595             solMovement = down-colsol[j];
    2596             assert(movement<-0.99999);
    2597             if (fabs(down-colLower[j])<1.0e-7) {
    2598               goingToTrueBound=2;
    2599               down=colLower[j];
    2600             }
    2601           } else {
    2602             movement=up-colLower[j];
    2603             solMovement = up-colsol[j];
    2604             assert(movement>0.99999);
    2605             if (fabs(up-colUpper[j])<1.0e-7) {
    2606               goingToTrueBound=2;
    2607               up=colUpper[j];
    2608             }
    2609           }
    2610           if (goingToTrueBound&&(colUpper[j]-colLower[j]>1.5||colLower[j]))
    2611             goingToTrueBound=1;
    2612           // switch off disaggregation if not wanted
    2613           if ((rowCuts&1)==0)
    2614             goingToTrueBound=0;
    2615 #ifdef PRINT_DEBUG
    2616           if (fabs(movement)>1.01) {
    2617             printf("big %d %g %g %g\n",j,colLower[j],solval,colUpper[j]);
    2618           }
    2619 #endif
    2620           if (solMovement*djs[j]>0.0)
    2621             objVal += solMovement*djs[j];
    2622           nstackC=1;
    2623           nstackR=0;
    2624           saveL[0]=colLower[j];
    2625           saveU[0]=colUpper[j];
    2626           assert (saveU[0]>saveL[0]);
    2627           notFeasible=0;
    2628           if (movement<0.0) {
    2629             colUpper[j] += movement;
    2630             colUpper[j] = floor(colUpper[j]+0.5);
    2631 #ifdef PRINT_DEBUG
    2632             printf("** Trying %d down to 0\n",j);
    2633 #endif
    2634           } else {
    2635             colLower[j] += movement;
    2636             colLower[j] = floor(colLower[j]+0.5);
    2637 #ifdef PRINT_DEBUG
    2638             printf("** Trying %d up to 1\n",j);
    2639 #endif
    2640           }
    2641           if (fabs(colUpper[j]-colLower[j])<1.0e-6)
    2642             markC[j]=3; // say fixed
    2643           istackC=0;
    2644           /* update immediately */
    2645           for (k=columnStart[j];k<columnStart[j]+columnLength[j];k++) {
    2646             int irow = row[k];
    2647             value = columnElements[k];
    2648             assert (markR[irow]!=-2);
    2649             if (markR[irow]==-1) {
    2650               stackR[nstackR]=irow;
    2651               markR[irow]=nstackR;
    2652               saveMin[nstackR]=minR[irow];
    2653               saveMax[nstackR]=maxR[irow];
    2654               nstackR++;
    2655 #if 0
    2656             } else if (markR[irow]==-2) {
    2657               continue;
    2658 #endif
    2659             }
    2660             /* could check immediately if violation */
    2661             if (movement>0.0) {
    2662               /* up */
    2663               if (value>0.0) {
    2664                 /* up does not change - down does */
    2665                 if (minR[irow]>-1.0e10)
    2666                   minR[irow] += value;
    2667                 if (minR[irow]>rowUpper[irow]+1.0e-5) {
    2668                   notFeasible=1;
    2669                   istackC=1;
    2670                   break;
    2671                 }
    2672               } else {
    2673                 /* down does not change - up does */
    2674                 if (maxR[irow]<1.0e10)
    2675                   maxR[irow] += value;
    2676                 if (maxR[irow]<rowLower[irow]-1.0e-5) {
    2677                   notFeasible=1;
    2678                   istackC=1;
    2679                   break;
    2680                 }
    2681               }
    2682             } else {
    2683               /* down */
    2684               if (value<0.0) {
    2685                 /* up does not change - down does */
    2686                 if (minR[irow]>-1.0e10)
    2687                   minR[irow] -= value;
    2688                 if (minR[irow]>rowUpper[irow]+1.0e-5) {
    2689                   notFeasible=1;
    2690                   istackC=1;
    2691                   break;
    2692                 }
    2693               } else {
    2694                 /* down does not change - up does */
    2695                 if (maxR[irow]<1.0e10)
    2696                   maxR[irow] -= value;
    2697                 if (maxR[irow]<rowLower[irow]-1.0e-5) {
    2698                   notFeasible=1;
    2699                   istackC=1;
    2700                   break;
    2701                 }
    2702               }
    2703             }
    2704           }
    2705           while (istackC<nstackC&&nstackC<maxStack) {
    2706             int jway;
    2707             int jcol =stackC[istackC];
    2708             jway=markC[jcol];
    2709             // If not first and fixed then skip
    2710             if (jway==3&&istackC) {
    2711               //istackC++;
    2712               //continue;
    2713               //printf("fixed %d on stack\n",jcol);
    2714             }
    2715             // Do cliques
    2716             if (oneFixStart_&&oneFixStart_[jcol]>=0) {
    2717               int start;
    2718               int end;
    2719               if (colLower[jcol]>saveL[istackC]) {
    2720                 // going up
    2721                 start = oneFixStart_[jcol];
    2722                 end = zeroFixStart_[jcol];
    2723               } else {
    2724                 assert (colUpper[jcol]<saveU[istackC]);
    2725                 // going down
    2726                 start = zeroFixStart_[jcol];
    2727                 end = endFixStart_[jcol];
    2728               }
    2729               for (int i=start;i<end;i++) {
    2730                 int iClique = whichClique_[i];
    2731                 for (int k=cliqueStart_[iClique];k<cliqueStart_[iClique+1];k++) {
    2732                   int kcol = sequenceInCliqueEntry(cliqueEntry_[k]);
    2733                   if (jcol==kcol)
    2734                     continue;
    2735                   int kway = oneFixesInCliqueEntry(cliqueEntry_[k]);
    2736                   if (kcol!=jcol) {
    2737                     if (!markC[kcol]) {
    2738                       // not on list yet
    2739                       if (nstackC<2*maxStack) {
    2740                         markC[kcol] = 3; // say fixed
    2741                         fixThis++;
    2742                         stackC[nstackC]=kcol;
    2743                         saveL[nstackC]=colLower[kcol];
    2744                         saveU[nstackC]=colUpper[kcol];
    2745                         assert (saveU[nstackC]>saveL[nstackC]);
    2746                         nstackC++;
    2747                         if (!kway) {
    2748                           // going up
    2749                           double solMovement=1.0-colsol[kcol];
    2750                           if (solMovement>0.0001) {
    2751                             assert (djs[kcol]>=0.0);
    2752                             objVal += djs[kcol]*solMovement;
    2753                           }
    2754                           colLower[kcol]=1.0;
    2755                           /* update immediately */
    2756                           for (int jj =columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    2757                             krow = row[jj];
    2758                             value = columnElements[jj];
    2759                             assert (markR[krow]!=-2);
    2760                             if (markR[krow]==-1) {
    2761                               stackR[nstackR]=krow;
    2762                               markR[krow]=nstackR;
    2763                               saveMin[nstackR]=minR[krow];
    2764                               saveMax[nstackR]=maxR[krow];
    2765                               nstackR++;
    2766 #if 0
    2767                             } else if (markR[krow]==-2) {
    2768                               continue;
    2769 #endif
    2770                             }
    2771                             /* could check immediately if violation */
    2772                             /* up */
    2773                             if (value>0.0) {
    2774                               /* up does not change - down does */
    2775                               if (minR[krow]>-1.0e10)
    2776                                 minR[krow] += value;
    2777                               if (minR[krow]>rowUpper[krow]+1.0e-5) {
    2778                                 colUpper[kcol]=-1.0e50; /* force infeasible */
    2779                                 break;
    2780                               }
    2781                             } else {
    2782                               /* down does not change - up does */
    2783                               if (maxR[krow]<1.0e10)
    2784                                 maxR[krow] += value;
    2785                               if (maxR[krow]<rowLower[krow]-1.0e-5) {
    2786                                 notFeasible=1;
    2787                                 break;
    2788                               }
    2789                             }
    2790                           }
    2791                         } else {
    2792                           // going down
    2793                           double solMovement=0.0-colsol[kcol];
    2794                           if (solMovement<-0.0001) {
    2795                             assert (djs[kcol]<=0.0);
    2796                             objVal += djs[kcol]*solMovement;
    2797                           }
    2798                           colUpper[kcol]=0.0;
    2799                           /* update immediately */
    2800                           for (int jj =columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    2801                             krow = row[jj];
    2802                             value = columnElements[jj];
    2803                             assert (markR[krow]!=-2);
    2804                             if (markR[krow]==-1) {
    2805                               stackR[nstackR]=krow;
    2806                               markR[krow]=nstackR;
    2807                               saveMin[nstackR]=minR[krow];
    2808                               saveMax[nstackR]=maxR[krow];
    2809                               nstackR++;
    2810 #if 0
    2811                             } else if (markR[krow]==-2) {
    2812                               continue;
    2813 #endif
    2814                             }
    2815                             /* could check immediately if violation */
    2816                             /* down */
    2817                             if (value<0.0) {
    2818                               /* up does not change - down does */
    2819                               if (minR[krow]>-1.0e10)
    2820                                 minR[krow] -= value;
    2821                               if (minR[krow]>rowUpper[krow]+1.0e-5) {
    2822                                 notFeasible=1;
    2823                                 break;
    2824                               }
    2825                             } else {
    2826                               /* down does not change - up does */
    2827                               if (maxR[krow]<1.0e10)
    2828                                 maxR[krow] -= value;
    2829                               if (maxR[krow]<rowLower[krow]-1.0e-5) {
    2830                                 notFeasible=1;
    2831                                 break;
    2832                               }
    2833                             }
    2834                           }
    2835                         }
    2836                       }
    2837                     } else if (markC[kcol]==1) {
    2838                       // marked as going to 0
    2839                       assert (!colUpper[kcol]);
    2840                       if (!kway) {
    2841                         // contradiction
    2842                         notFeasible=1;
    2843                         break;
    2844                       }
    2845                     } else if (markC[kcol]==2) {
    2846                       // marked as going to 1
    2847                       assert (colLower[kcol]);
    2848                       if (kway) {
    2849                         // contradiction
    2850                         notFeasible=1;
    2851                         break;
    2852                       }
    2853                     } else {
    2854                       // marked as fixed
    2855                       assert (markC[kcol]==3);
    2856                       int jkway;
    2857                       if (colLower[kcol])
    2858                         jkway=1;
    2859                       else
    2860                         jkway=0;
    2861                       if (kway==jkway) {
    2862                         // contradiction
    2863                         notFeasible=1;
    2864                         break;
    2865                       }
    2866                     }
    2867                   }
    2868                 }
    2869                 if (notFeasible)
    2870                   break;
    2871               }
    2872               if (notFeasible)
    2873                 istackC=nstackC+1;
    2874             }
    2875             for (k=columnStart[jcol];k<columnStart[jcol]+columnLength[jcol];k++) {
    2876               // break if found not feasible
    2877               if (notFeasible)
    2878                 break;
    2879               int irow = row[k];
    2880               /*value = columnElements[k];*/
    2881               assert (markR[irow]!=-2);
    2882 #if 0
    2883               if (markR[irow]!=-2) {
    2884 #endif
    2885                 /* see if anything forced */
    2886                 for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];kk++) {
    2887                   double moveUp=0.0;
    2888                   double moveDown=0.0;
    2889                   double newUpper=-1.0,newLower=1.0;
    2890                   kcol=column[kk];
    2891                   bool onList = (markC[kcol]!=0);
    2892                   if (markC[kcol]!=3) {
    2893                     value2=rowElements[kk];
    2894                     int markIt=markC[kcol];
    2895                     if (value2 < 0.0) {
    2896                       if (colUpper[kcol] < 1e10 && (markIt&2)==0 &&
    2897                           rowUpper[irow]<1.0e10) {
    2898                         dbound = colUpper[kcol]+
    2899                           (rowUpper[irow]-minR[irow])/value2;
    2900                         if (dbound > colLower[kcol] + primalTolerance_) {
    2901                           if (intVar[kcol]) {
    2902                             markIt |= 2;
    2903                             newLower = ceil(dbound-primalTolerance_);
    2904                           } else {
    2905                             newLower=dbound;
    2906                             if (newLower+primalTolerance_>colUpper[kcol]&&
    2907                                 newLower-primalTolerance_<=colUpper[kcol]) {
    2908                               newLower=colUpper[kcol];
    2909                               markIt |= 2;
    2910                               //markIt=3;
    2911                             } else {
    2912                               // avoid problems - fix later ?
    2913                               markIt=3;
    2914                             }
    2915                           }
    2916                           moveUp = newLower-colLower[kcol];
    2917                         }
    2918                       }
    2919                       if (colLower[kcol] > -1e10 && (markIt&1)==0 &&
    2920                           rowLower[irow]>-1.0e10) {
    2921                         dbound = colLower[kcol] +
    2922                           (rowLower[irow]-maxR[irow])/value2;
    2923                         if (dbound < colUpper[kcol] - primalTolerance_) {
    2924                           if (intVar[kcol]) {
    2925                             markIt |= 1;
    2926                             newUpper = floor(dbound+primalTolerance_);
    2927                           } else {
    2928                             newUpper=dbound;
    2929                             if (newUpper-primalTolerance_<colLower[kcol]&&
    2930                                 newUpper+primalTolerance_>=colLower[kcol]) {
    2931                               newUpper=colLower[kcol];
    2932                               markIt |= 1;
    2933                               //markIt=3;
    2934                             } else {
    2935                               // avoid problems - fix later ?
    2936                               markIt=3;
    2937                             }
    2938                           }
    2939                           moveDown = newUpper-colUpper[kcol];
    2940                         }
    2941                       }
    2942                     } else {
    2943                       /* positive element */
    2944                       if (colUpper[kcol] < 1e10 && (markIt&2)==0 &&
    2945                           rowLower[irow]>-1.0e10) {
    2946                         dbound = colUpper[kcol] +
    2947                           (rowLower[irow]-maxR[irow])/value2;
    2948                         if (dbound > colLower[kcol] + primalTolerance_) {
    2949                           if (intVar[kcol]) {
    2950                             markIt |= 2;
    2951                             newLower = ceil(dbound-primalTolerance_);
    2952                           } else {
    2953                             newLower=dbound;
    2954                             if (newLower+primalTolerance_>colUpper[kcol]&&
    2955                                 newLower-primalTolerance_<=colUpper[kcol]) {
    2956                               newLower=colUpper[kcol];
    2957                               markIt |= 2;
    2958                               //markIt=3;
    2959                             } else {
    2960                               // avoid problems - fix later ?
    2961                               markIt=3;
    2962                             }
    2963                           }
    2964                           moveUp = newLower-colLower[kcol];
    2965                         }
    2966                       }
    2967                       if (colLower[kcol] > -1e10 && (markIt&1)==0 &&
    2968                           rowUpper[irow]<1.0e10) {
    2969                         dbound = colLower[kcol] +
    2970                           (rowUpper[irow]-minR[irow])/value2;
    2971                         if (dbound < colUpper[kcol] - primalTolerance_) {
    2972                           if (intVar[kcol]) {
    2973                             markIt |= 1;
    2974                             newUpper = floor(dbound+primalTolerance_);
    2975                           } else {
    2976                             newUpper=dbound;
    2977                             if (newUpper-primalTolerance_<colLower[kcol]&&
    2978                                 newUpper+primalTolerance_>=colLower[kcol]) {
    2979                               newUpper=colLower[kcol];
    2980                               markIt |= 1;
    2981                               //markIt=3;
    2982                             } else {
    2983                               // avoid problems - fix later ?
    2984                               markIt=3;
    2985                             }
    2986                           }
    2987                           moveDown = newUpper-colUpper[kcol];
    2988                         }
    2989                       }
    2990                     }
    2991                     if (nstackC<2*maxStack) {
    2992                       markC[kcol] = markIt;
    2993                     }
    2994                     if (moveUp&&nstackC<2*maxStack) {
    2995                       fixThis++;
    2996 #ifdef PRINT_DEBUG
    2997                       printf("lower bound on %d increased from %g to %g by row %d %g %g\n",kcol,colLower[kcol],newLower,irow,rowLower[irow],rowUpper[irow]);
    2998                       value=0.0;
    2999                       for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    3000                         int ii=column[jj];
    3001                         if (colUpper[ii]-colLower[ii]>primalTolerance_) {
    3002                           printf("(%d, %g) ",ii,rowElements[jj]);
    3003                         } else {
    3004                           value += rowElements[jj]*colLower[ii];
    3005                         }
    3006                       }
    3007                       printf(" - fixed %g\n",value);
    3008                       for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    3009                         int ii=column[jj];
    3010                         if (colUpper[ii]-colLower[ii]<primalTolerance_) {
    3011                           printf("(%d, %g, %g) ",ii,rowElements[jj],colLower[ii]);
    3012                         }
    3013                       }
    3014                       printf("\n");
    3015 #endif
    3016                       if (!onList) {
    3017                         stackC[nstackC]=kcol;
    3018                         saveL[nstackC]=colLower[kcol];
    3019                         saveU[nstackC]=colUpper[kcol];
    3020                         assert (saveU[nstackC]>saveL[nstackC]);
    3021                         nstackC++;
    3022                         onList=true;
    3023                       }
    3024                       if (newLower>colsol[kcol]) {
    3025                         if (djs[kcol]<0.0) {
    3026                           /* should be infeasible */
    3027                           assert (newLower>colUpper[kcol]+primalTolerance_);
    3028                         } else {
    3029                           objVal += moveUp*djs[kcol];
    3030                         }
    3031                       }
    3032                       if (intVar[kcol])
    3033                         newLower = CoinMax(colLower[kcol],ceil(newLower-1.0e-4));
    3034                       colLower[kcol]=newLower;
    3035                       if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
    3036                         markC[kcol]=3; // say fixed
    3037                       }
    3038                       /* update immediately */
    3039                       for (jj=columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    3040                         krow = row[jj];
    3041                         value = columnElements[jj];
    3042                         assert (markR[krow]!=-2);
    3043                         if (markR[krow]==-1) {
    3044                           stackR[nstackR]=krow;
    3045                           markR[krow]=nstackR;
    3046                           saveMin[nstackR]=minR[krow];
    3047                           saveMax[nstackR]=maxR[krow];
    3048                           nstackR++;
    3049 #if 0
    3050                         } else if (markR[krow]==-2) {
    3051                           continue;
    3052 #endif
    3053                         }
    3054                         /* could check immediately if violation */
    3055                         /* up */
    3056                         if (value>0.0) {
    3057                           /* up does not change - down does */
    3058                           if (minR[krow]>-1.0e10)
    3059                             minR[krow] += value*moveUp;
    3060                           if (minR[krow]>rowUpper[krow]+1.0e-5) {
    3061                             colUpper[kcol]=-1.0e50; /* force infeasible */
    3062                             break;
    3063                           }
    3064                         } else {
    3065                           /* down does not change - up does */
    3066                           if (maxR[krow]<1.0e10)
    3067                             maxR[krow] += value*moveUp;
    3068                           if (maxR[krow]<rowLower[krow]-1.0e-5) {
    3069                             colUpper[kcol]=-1.0e50; /* force infeasible */
    3070                             break;
    3071                           }
    3072                         }
    3073                       }
    3074                     }
    3075                     if (moveDown&&nstackC<2*maxStack) {
    3076                       fixThis++;
    3077 #ifdef PRINT_DEBUG
    3078                       printf("upper bound on %d decreased from %g to %g by row %d %g %g\n",kcol,colUpper[kcol],newUpper,irow,rowLower[irow],rowUpper[irow]);
    3079                       value=0.0;
    3080                       for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    3081                         int ii=column[jj];
    3082                         if (colUpper[ii]-colLower[ii]>primalTolerance_) {
    3083                           printf("(%d, %g) ",ii,rowElements[jj]);
    3084                         } else {
    3085                           value += rowElements[jj]*colLower[ii];
    3086                         }
    3087                       }
    3088                       printf(" - fixed %g\n",value);
    3089                       for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    3090                         int ii=column[jj];
    3091                         if (colUpper[ii]-colLower[ii]<primalTolerance_) {
    3092                           printf("(%d, %g, %g) ",ii,rowElements[jj],colLower[ii]);
    3093                         }
    3094                       }
    3095                       printf("\n");
    3096 #endif
    3097                       if (!onList) {
    3098                         stackC[nstackC]=kcol;
    3099                         saveL[nstackC]=colLower[kcol];
    3100                         saveU[nstackC]=colUpper[kcol];
    3101                         assert (saveU[nstackC]>saveL[nstackC]);
    3102                         nstackC++;
    3103                         onList=true;
    3104                       }
    3105                       if (newUpper<colsol[kcol]) {
    3106                         if (djs[kcol]>0.0) {
    3107                           /* should be infeasible */
    3108                           assert (colLower[kcol]>newUpper+primalTolerance_);
    3109                         } else {
    3110                           objVal += moveDown*djs[kcol];
    3111                         }
    3112                       }
    3113                       if (intVar[kcol])
    3114                         newUpper = CoinMin(colUpper[kcol],floor(newUpper+1.0e-4));
    3115                       colUpper[kcol]=newUpper;
    3116                       if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
    3117                         markC[kcol]=3; // say fixed
    3118                       }
    3119                       /* update immediately */
    3120                       for (jj=columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    3121                         krow = row[jj];
    3122                         value = columnElements[jj];
    3123                         assert (markR[krow]!=-2);
    3124                         if (markR[krow]==-1) {
    3125                           stackR[nstackR]=krow;
    3126                           markR[krow]=nstackR;
    3127                           saveMin[nstackR]=minR[krow];
    3128                           saveMax[nstackR]=maxR[krow];
    3129                           nstackR++;
    3130 #if 0
    3131                         } else if (markR[krow]==-2) {
    3132 #endif
    3133                           continue;
    3134                         }
    3135                         /* could check immediately if violation */
    3136                         /* down */
    3137                         if (value<0.0) {
    3138                           /* up does not change - down does */
    3139                           if (minR[krow]>-1.0e10)
    3140                             minR[krow] += value*moveDown;
    3141                           if (minR[krow]>rowUpper[krow]+1.0e-5) {
    3142                             colUpper[kcol]=-1.0e50; /* force infeasible */
    3143                             break;
    3144                           }
    3145                         } else {
    3146                           /* down does not change - up does */
    3147                           if (maxR[krow]<1.0e10)
    3148                             maxR[krow] += value*moveDown;
    3149                           if (maxR[krow]<rowLower[krow]-1.0e-5) {
    3150                             colUpper[kcol]=-1.0e50; /* force infeasible */
    3151                             break;
    3152                           }
    3153                         }
    3154                       }
    3155                     }
    3156                     if (colLower[kcol]>colUpper[kcol]+primalTolerance_) {
    3157                       notFeasible=1;;
    3158                       k=columnStart[jcol]+columnLength[jcol];
    3159                       istackC=nstackC+1;
    3160 #ifdef PRINT_DEBUG
    3161                       printf("** not feasible this way\n");
    3162 #endif
    3163                       break;
    3164                     }
    3165                   }
    3166                 }
    3167 #if 0
    3168               }
    3169 #endif
    3170             }
    3171             istackC++;
    3172           }
    3173           if (!notFeasible) {
    3174             if (objVal<=cutoff) {
    3175               feasible |= feas[iway];
    3176             } else {
    3177 #ifdef PRINT_DEBUG
    3178               printf("not feasible on dj\n");
    3179 #endif
    3180               notFeasible=1;
    3181               if (iway==1&&feasible==0) {
    3182                 /* not feasible at all */
    3183                 ninfeas=1;
    3184                 j=nCols-1;
    3185                 break;
    3186               }
    3187             }
    3188           } else if (iway==1&&feasible==0) {
    3189             /* not feasible at all */
    3190             ninfeas=1;
    3191             j=nCols-1;
    3192             iLook=numberThisTime_;
    3193             ipass=maxPass;
    3194             break;
    3195           }
    3196           if (notFeasible)
    3197             goingToTrueBound=0;
    3198           if (iway==2||(iway==1&&feasible==2)) {
    3199             /* keep */
    3200             iway=3;
    3201             nfixed++;
    3202             if (mode_) {
    3203               OsiColCut cc;
    3204               int nTot=0,nFix=0,nInt=0;
    3205               bool ifCut=false;
    3206               for (istackC=0;istackC<nstackC;istackC++) {
    3207                 int icol=stackC[istackC];
    3208                 if (intVar[icol]) {
    3209                   if (colUpper[icol]<currentColUpper[icol]-1.0e-4) {
    3210                     element[nFix]=colUpper[icol];
    3211                     index[nFix++]=icol;
    3212                     nInt++;
    3213                     if (colsol[icol]>colUpper[icol]+primalTolerance_) {
    3214                       ifCut=true;
    3215                       anyColumnCuts=true;
    3216                     }
    3217                   }
    3218                 }
    3219               }
    3220               if (nFix) {
    3221                 nTot=nFix;
    3222                 cc.setUbs(nFix,index,element);
    3223                 nFix=0;
    3224               }
    3225               for (istackC=0;istackC<nstackC;istackC++) {
    3226                 int icol=stackC[istackC];
    3227                 if (intVar[icol]) {
    3228                   if (colLower[icol]>currentColLower[icol]+1.0e-4) {
    3229                     element[nFix]=colLower[icol];
    3230                     index[nFix++]=icol;
    3231                     nInt++;
    3232                     if (colsol[icol]<colLower[icol]-primalTolerance_) {
    3233                       ifCut=true;
    3234                       anyColumnCuts=true;
    3235                     }
    3236                   }
    3237                 }
    3238               }
    3239               if (nFix) {
    3240                 nTot+=nFix;
    3241                 cc.setLbs(nFix,index,element);
    3242               }
    3243               // could tighten continuous as well
    3244               if (nInt) {
    3245                 if (ifCut) {
    3246                   cc.setEffectiveness(100.0);
    3247                 } else {
    3248                   cc.setEffectiveness(1.0e-5);
    3249                 }
    3250 #if CGL_DEBUG > 0
    3251                 CglProbingDebug::checkBounds(si,cc);
    3252 #endif
    3253                 cs.insert(cc);
    3254               }
    3255             }
    3256             for (istackC=0;istackC<nstackC;istackC++) {
    3257               int icol=stackC[istackC];
    3258               if (colUpper[icol]-colLower[icol]>primalTolerance_) {
    3259                 markC[icol]=0;
    3260               } else {
    3261                 markC[icol]=3;
    3262               }
    3263             }
    3264             for (istackR=0;istackR<nstackR;istackR++) {
    3265               int irow=stackR[istackR];
    3266               markR[irow]=-1;
    3267             }
    3268           } else {
    3269             /* is it worth seeing if can increase coefficients
    3270                or maybe better see if it is a cut */
    3271             if (iway==0) {
    3272               nstackC0=CoinMin(nstackC,maxStack);
    3273               double solMove = saveSolval-down;
    3274               double boundChange;
    3275               if (notFeasible) {
    3276                 nstackC0=0;
    3277               } else {
    3278                 for (istackC=0;istackC<nstackC0;istackC++) {
    3279                   int icol=stackC[istackC];
    3280                   stackC0[istackC]=icol;
    3281                   lo0[istackC]=colLower[icol];
    3282                   up0[istackC]=colUpper[icol];
    3283                 }
    3284               }
    3285               /* restore all */
    3286               int nCliquesAffected=0;
    3287               assert (iway==0);
    3288               for (istackC=nstackC-1;istackC>=0;istackC--) {
    3289                 int icol=stackC[istackC];
    3290                 double oldU=saveU[istackC];
    3291                 double oldL=saveL[istackC];
    3292                 if(goingToTrueBound==2&&istackC) {
    3293                   // Work for extending cliques
    3294                   if (!mode_&&numberCliques_) {
    3295                     int i_01 = to_01[icol];
    3296                     if (i_01>=0) {
    3297                       int start;
    3298                       int end;
    3299                       if (colLower[icol]) {
    3300                         // going up - but we want weak way
    3301                         start = zeroFixStart_[icol];
    3302                         end = endFixStart_[icol];
    3303                       } else {
    3304                         // going down - but we want weak way
    3305                         start = oneFixStart_[icol];
    3306                         end = zeroFixStart_[icol];
    3307                       }
    3308                       //if (end>start)
    3309                       //printf("j %d, other %d is in %d cliques\n",
    3310                       //     j,i_01,end-start);
    3311                       for (int i=start;i<end;i++) {
    3312                         int iClique = whichClique_[i];
    3313                         int size = cliqueStart_[iClique+1]-cliqueStart_[iClique];
    3314                         if (cliqueCount[iClique]==size) {
    3315                           // first time
    3316                           cliqueStack[nCliquesAffected++]=iClique;
    3317                         }
    3318                         // decrement counts
    3319                         cliqueCount[iClique]--;
    3320                       }
    3321                     }
    3322                   }
    3323                   // upper disaggregation cut would be
    3324                   // xval < upper + (old_upper-upper) (jval-down)
    3325                   boundChange = oldU-colUpper[icol];
    3326                   if (boundChange>0.0&&oldU<1.0e10&&
    3327                       (!mode_||colsol[icol]>colUpper[icol]
    3328                       + boundChange*solMove+primalTolerance_)) {
    3329                     // create cut
    3330                     OsiRowCut rc;
    3331                     rc.setLb(-DBL_MAX);
    3332                     rc.setUb(colUpper[icol]-down*boundChange);
    3333                     index[0]=icol;
    3334                     element[0]=1.0;
    3335                     index[1]=j;
    3336                     element[1]= - boundChange;
    3337                     // effectiveness is how far j moves
    3338                     double newSol = (colsol[icol]-colUpper[icol])/
    3339                       boundChange;
    3340                     if (mode_)
    3341                       assert(newSol>solMove);
    3342                     rc.setEffectiveness(newSol-solMove);
    3343                     if (rc.effectiveness()>disaggEffectiveness) {
    3344                       rc.setRow(2,index,element,false);
    3345 #if CGL_DEBUG > 0
    3346                       if (debugger) assert(!debugger->invalidCut(rc));
    3347 #endif
    3348                       rowCut.addCutIfNotDuplicate(rc);
    3349                     }
    3350                   }
    3351                   // lower disaggregation cut would be
    3352                   // xval > lower + (old_lower-lower) (jval-down)
    3353                   boundChange = oldL-colLower[icol];
    3354                   if (boundChange<0.0&&oldL>-1.0e10&&
    3355                       (!mode_||colsol[icol]<colLower[icol]
    3356                       + boundChange*solMove-primalTolerance_)) {
    3357                     // create cut
    3358                     OsiRowCut rc;
    3359                     rc.setLb(colLower[icol]-down*boundChange);
    3360                     rc.setUb(DBL_MAX);
    3361                     index[0]=icol;
    3362                     element[0]=1.0;
    3363                     index[1]=j;
    3364                     element[1]=- boundChange;
    3365                     // effectiveness is how far j moves
    3366                     double newSol = (colsol[icol]-colLower[icol])/
    3367                       boundChange;
    3368                     if (mode_)
    3369                       assert(newSol>solMove);
    3370                     rc.setEffectiveness(newSol-solMove);
    3371                     if (rc.effectiveness()>disaggEffectiveness) {
    3372                       rc.setRow(2,index,element,false);
    3373 #if CGL_DEBUG > 0
    3374                       if (debugger) assert(!debugger->invalidCut(rc));
    3375 #endif
    3376                       rowCut.addCutIfNotDuplicate(rc);
    3377 #if 0
    3378                       printf("%d original bounds %g, %g new Lo %g sol= %g int %d sol= %g\n",icol,oldL,oldU,colLower[icol],colsol[icol], j, colsol[j]);
    3379                       printf("-1.0 * x(%d) + %g * y(%d) <= %g\n",
    3380                              icol,boundChange,j,rc.ub());
    3381 #endif
    3382                     }
    3383                   }
    3384                 }
    3385                 colUpper[icol]=oldU;
    3386                 colLower[icol]=oldL;
    3387                 markC[icol]=0;
    3388               }
    3389               if (nCliquesAffected) {
    3390                 for (int i=0;i<nCliquesAffected;i++) {
    3391                   int iClique = cliqueStack[i];
    3392                   int size = cliqueCount[iClique];
    3393                   // restore
    3394                   cliqueCount[iClique]= cliqueStart_[iClique+1]-cliqueStart_[iClique];
    3395                   if (!size) {
    3396                     if (logLevel_>1)
    3397                       printf("** could extend clique by adding j!\n");
    3398                   }
    3399                 }
    3400               }
    3401               for (istackR=0;istackR<nstackR;istackR++) {
    3402                 int irow=stackR[istackR];
    3403                 // switch off strengthening if not wanted
    3404                 if ((rowCuts&2)!=0&&goingToTrueBound) {
    3405                   bool ifCut=anyColumnCuts;
    3406                   double gap = rowUpper[irow]-maxR[irow];
    3407                   double sum=0.0;
    3408                   if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    3409                     // see if the strengthened row is a cut
    3410                     for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3411                          kk++) {
    3412                       sum += rowElements[kk]*colsol[column[kk]];
    3413                     }
    3414                     if (sum-gap*colsol[j]>maxR[irow]+primalTolerance_||(info->strengthenRow&&rowLower[irow]<-1.0e20)) {
    3415                       // can be a cut
    3416                       // subtract gap from upper and integer coefficient
    3417                       // saveU and saveL spare
    3418                       int * index = reinterpret_cast<int *>(saveL);
    3419                       double * element = saveU;
    3420                       int n=0;
    3421                       bool coefficientExists=false;
    3422                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3423                            kk++) {
    3424                         if (column[kk]!=j) {
    3425                           index[n]=column[kk];
    3426                           element[n++]=rowElements[kk];
    3427                         } else {
    3428                           double value=rowElements[kk]-gap;
    3429                           if (fabs(value)>1.0e-12) {
    3430                             index[n]=column[kk];
    3431                             element[n++]=value;
    3432                           }
    3433                           coefficientExists=true;
    3434                         }
    3435                       }
    3436                       if (!coefficientExists) {
    3437                         index[n]=j;
    3438                         element[n++]=-gap;
    3439                       }
    3440                       OsiRowCut rc;
    3441                       rc.setLb(-DBL_MAX);
    3442                       rc.setUb(rowUpper[irow]-gap*(colLower[j]+1.0));
    3443                       // effectiveness is how far j moves
    3444                       rc.setEffectiveness((sum-gap*colsol[j]-maxR[irow])/gap);
    3445                       if (rc.effectiveness()>needEffectiveness) {
    3446                         rc.setRow(n,index,element,false);
    3447 #if CGL_DEBUG > 0
    3448                         if (debugger) assert(!debugger->invalidCut(rc));
    3449 #endif
    3450                         // If strengthenRow point to row
    3451                         //if(info->strengthenRow)
    3452                         //printf("a point to row %d\n",irow);
    3453 #ifdef STRENGTHEN_PRINT
    3454                       if (rowLower[irow]<-1.0e20) {
    3455                         printf("5Cut %g <= ",rc.lb());
    3456                         int k;
    3457                         for ( k=0;k<n;k++) {
    3458                           int iColumn = index[k];
    3459                           printf("%g*",element[k]);
    3460                           if (si.isInteger(iColumn))
    3461                             printf("i%d ",iColumn);
    3462                           else
    3463                             printf("x%d ",iColumn);
    3464                         }
    3465                         printf("<= %g\n",rc.ub());
    3466                         printf("Row %g <= ",rowLower[irow]);
    3467                         for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    3468                           int iColumn = column[k];
    3469                           printf("%g*",rowElements[k]);
    3470                           if (si.isInteger(iColumn))
    3471                             printf("i%d ",iColumn);
    3472                           else
    3473                             printf("x%d ",iColumn);
    3474                         }
    3475                         printf("<= %g\n",rowUpper[irow]);
    3476                       }
    3477 #endif
    3478                       int realRow = (rowLower[irow]<-1.0e20) ? irow : -1;
    3479                       if (realRows&&realRow>0)
    3480                         realRow=realRows[realRow];
    3481                       rc.setWhichRow(realRow) ;
    3482                       rowCut.addCutIfNotDuplicate(rc);
    3483                       }
    3484                     }
    3485                   }
    3486                   gap = minR[irow]-rowLower[irow];
    3487                   if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    3488                     // see if the strengthened row is a cut
    3489                     if (!sum) {
    3490                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3491                            kk++) {
    3492                         sum += rowElements[kk]*colsol[column[kk]];
    3493                       }
    3494                     }
    3495                     if (sum+gap*colsol[j]<minR[irow]+primalTolerance_||(info->strengthenRow&&rowUpper[irow]>1.0e20)) {
    3496                       // can be a cut
    3497                       // add gap to lower and integer coefficient
    3498                       // saveU and saveL spare
    3499                       int * index = reinterpret_cast<int *>(saveL);
    3500                       double * element = saveU;
    3501                       int n=0;
    3502                       bool coefficientExists=false;
    3503                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3504                            kk++) {
    3505                         if (column[kk]!=j) {
    3506                           index[n]=column[kk];
    3507                           element[n++]=rowElements[kk];
    3508                         } else {
    3509                           double value=rowElements[kk]+gap;
    3510                           if (fabs(value)>1.0e-12) {
    3511                             index[n]=column[kk];
    3512                             element[n++]=value;
    3513                           }
    3514                           coefficientExists=true;
    3515                         }
    3516                       }
    3517                       if (!coefficientExists) {
    3518                         index[n]=j;
    3519                         element[n++]=gap;
    3520                       }
    3521                       OsiRowCut rc;
    3522                       rc.setLb(rowLower[irow]+gap*(colLower[j]+1.0));
    3523                       rc.setUb(DBL_MAX);
    3524                       // effectiveness is how far j moves
    3525                       rc.setEffectiveness((minR[irow]-sum-gap*colsol[j])/gap);
    3526                       if (rc.effectiveness()>needEffectiveness) {
    3527                         rc.setRow(n,index,element,false);
    3528 #if CGL_DEBUG > 0
    3529                         if (debugger) assert(!debugger->invalidCut(rc));
    3530 #endif
    3531                         //if(info->strengthenRow)
    3532                         //printf("b point to row %d\n",irow);
    3533 #ifdef STRENGTHEN_PRINT
    3534                       if (rowUpper[irow]>1.0e20) {
    3535                         printf("6Cut %g <= ",rc.lb());
    3536                         int k;
    3537                         for ( k=0;k<n;k++) {
    3538                           int iColumn = index[k];
    3539                           printf("%g*",element[k]);
    3540                           if (si.isInteger(iColumn))
    3541                             printf("i%d ",iColumn);
    3542                           else
    3543                             printf("x%d ",iColumn);
    3544                         }
    3545                         printf("<= %g\n",rc.ub());
    3546                         printf("Row %g <= ",rowLower[irow]);
    3547                         for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    3548                           int iColumn = column[k];
    3549                           printf("%g*",rowElements[k]);
    3550                           if (si.isInteger(iColumn))
    3551                             printf("i%d ",iColumn);
    3552                           else
    3553                             printf("x%d ",iColumn);
    3554                         }
    3555                         printf("<= %g\n",rowUpper[irow]);
    3556                       }
    3557 #endif
    3558                       int realRow = (rowUpper[irow]>1.0e20) ? irow : -1;
    3559                       if (realRows&&realRow>0)
    3560                         realRow=realRows[realRow];
    3561                       rc.setWhichRow(realRow) ;
    3562                       rowCut.addCutIfNotDuplicate(rc);
    3563                       }
    3564                     }
    3565                   }
    3566                 }
    3567                 minR[irow]=saveMin[istackR];
    3568                 maxR[irow]=saveMax[istackR];
    3569                 markR[irow]=-1;
    3570               }
    3571             } else {
    3572               if (iway==1&&feasible==3) {
    3573                 iway=3;
    3574                 /* point back to stack */
    3575                 for (istackC=nstackC-1;istackC>=0;istackC--) {
    3576                   int icol=stackC[istackC];
    3577                   markC[icol]=istackC+1000;
    3578                 }
    3579                 if (mode_) {
    3580                   OsiColCut cc;
    3581                   int nTot=0,nFix=0,nInt=0;
    3582                   bool ifCut=false;
    3583                   for (istackC=0;istackC<nstackC0;istackC++) {
    3584                     int icol=stackC0[istackC];
    3585                     int istackC1=markC[icol]-1000;
    3586                     if (istackC1>=0) {
    3587                       if (CoinMin(lo0[istackC],colLower[icol])>saveL[istackC1]+1.0e-4) {
    3588                         saveL[istackC1]=CoinMin(lo0[istackC],colLower[icol]);
    3589                         if (intVar[icol]) {
    3590                           element[nFix]=saveL[istackC1];
    3591                           index[nFix++]=icol;
    3592                           nInt++;
    3593                           if (colsol[icol]<saveL[istackC1]-primalTolerance_)
    3594                             ifCut=true;
    3595                         }
    3596                         nfixed++;
    3597                       }
    3598                     }
    3599                   }
    3600                   if (nFix) {
    3601                     nTot=nFix;
    3602                     cc.setLbs(nFix,index,element);
    3603                     nFix=0;
    3604                   }
    3605                   for (istackC=0;istackC<nstackC0;istackC++) {
    3606                     int icol=stackC0[istackC];
    3607                     int istackC1=markC[icol]-1000;
    3608                     if (istackC1>=0) {
    3609                       if (CoinMax(up0[istackC],colUpper[icol])<saveU[istackC1]-1.0e-4) {
    3610                         saveU[istackC1]=CoinMax(up0[istackC],colUpper[icol]);
    3611                         if (intVar[icol]) {
    3612                           element[nFix]=saveU[istackC1];
    3613                           index[nFix++]=icol;
    3614                           nInt++;
    3615                           if (colsol[icol]>saveU[istackC1]+primalTolerance_)
    3616                             ifCut=true;
    3617                         }
    3618                         nfixed++;
    3619                       }
    3620                     }
    3621                   }
    3622                   if (nFix) {
    3623                     nTot+=nFix;
    3624                     cc.setUbs(nFix,index,element);
    3625                   }
    3626                   // could tighten continuous as well
    3627                   if (nInt) {
    3628                     if (ifCut) {
    3629                       cc.setEffectiveness(100.0);
    3630                     } else {
    3631                       cc.setEffectiveness(1.0e-5);
    3632                     }
    3633 #if CGL_DEBUG > 0
    3634                     CglProbingDebug::checkBounds(si,cc);
    3635 #endif
    3636                     cs.insert(cc);
    3637                   }
    3638                 }
    3639               } else {
    3640                 goingToTrueBound=0;
    3641               }
    3642               double solMove = up-saveSolval;
    3643               double boundChange;
    3644               /* restore all */
    3645               int nCliquesAffected=0;
    3646               for (istackC=nstackC-1;istackC>=0;istackC--) {
    3647                 int icol=stackC[istackC];
    3648                 double oldU=saveU[istackC];
    3649                 double oldL=saveL[istackC];
    3650                 if(goingToTrueBound==2&&istackC) {
    3651                   // Work for extending cliques
    3652                   if (!mode_&&numberCliques_&&iway==3) {
    3653                     int i_01 = to_01[icol];
    3654                     if (i_01>=0) {
    3655                       int start;
    3656                       int end;
    3657                       if (colLower[icol]) {
    3658                         // going up - but we want weak way
    3659                         start = zeroFixStart_[icol];
    3660                         end = endFixStart_[icol];
    3661                       } else {
    3662                         // going down - but we want weak way
    3663                         start = oneFixStart_[icol];
    3664                         end = zeroFixStart_[icol];
    3665                       }
    3666                       //if (end>start)
    3667                       //printf("up j %d, other %d is in %d cliques\n",
    3668                       //     j,i_01,end-start);
    3669                       for (int i=start;i<end;i++) {
    3670                         int iClique = whichClique_[i];
    3671                         int size = cliqueStart_[iClique+1]-cliqueStart_[iClique];
    3672                         if (cliqueCount[iClique]==size) {
    3673                           // first time
    3674                           cliqueStack[nCliquesAffected++]=iClique;
    3675                         }
    3676                         // decrement counts
    3677                         cliqueCount[iClique]--;
    3678                       }
    3679                     }
    3680                   }
    3681                   // upper disaggregation cut would be
    3682                   // xval < upper + (old_upper-upper) (up-jval)
    3683                   boundChange = oldU-colUpper[icol];
    3684                   if (boundChange>0.0&&oldU<1.0e10&&
    3685                       (!mode_||colsol[icol]>colUpper[icol]
    3686                       + boundChange*solMove+primalTolerance_)) {
    3687                     // create cut
    3688                     OsiRowCut rc;
    3689                     rc.setLb(-DBL_MAX);
    3690                     rc.setUb(colUpper[icol]+up*boundChange);
    3691                     index[0]=icol;
    3692                     element[0]=1.0;
    3693                     index[1]=j;
    3694                     element[1]= + boundChange;
    3695                     // effectiveness is how far j moves
    3696                     double newSol = (colsol[icol]-colUpper[icol])/
    3697                       boundChange;
    3698                     if (mode_)
    3699                       assert(newSol>solMove);
    3700                     rc.setEffectiveness(newSol-solMove);
    3701                     if (rc.effectiveness()>disaggEffectiveness) {
    3702                       rc.setRow(2,index,element,false);
    3703 #if CGL_DEBUG > 0
    3704                       if (debugger) assert(!debugger->invalidCut(rc));
    3705 #endif
    3706                       rowCut.addCutIfNotDuplicate(rc);
    3707                     }
    3708                   }
    3709                   // lower disaggregation cut would be
    3710                   // xval > lower + (old_lower-lower) (up-jval)
    3711                   boundChange = oldL-colLower[icol];
    3712                   if (boundChange<0.0&&oldL>-1.0e10&&
    3713                       (!mode_||colsol[icol]<colLower[icol]
    3714                       + boundChange*solMove-primalTolerance_)) {
    3715                     // create cut
    3716                     OsiRowCut rc;
    3717                     rc.setLb(colLower[icol]+up*boundChange);
    3718                     rc.setUb(DBL_MAX);
    3719                     index[0]=icol;
    3720                     element[0]=1.0;
    3721                     index[1]=j;
    3722                     element[1]= + boundChange;
    3723                     // effectiveness is how far j moves
    3724                     double newSol = (colsol[icol]-colLower[icol])/
    3725                       boundChange;
    3726                     if (mode_)
    3727                       assert(newSol>solMove);
    3728                     rc.setEffectiveness(newSol-solMove);
    3729                     if (rc.effectiveness()>disaggEffectiveness) {
    3730                       rc.setRow(2,index,element,false);
    3731 #if CGL_DEBUG > 0
    3732                       if (debugger) assert(!debugger->invalidCut(rc));
    3733 #endif
    3734                       rowCut.addCutIfNotDuplicate(rc);
    3735                     }
    3736                   }
    3737                 }
    3738                 colUpper[icol]=oldU;
    3739                 colLower[icol]=oldL;
    3740                 if (oldU>oldL+1.0e-4)
    3741                   markC[icol]=0;
    3742                 else
    3743                   markC[icol]=3;
    3744               }
    3745               if (nCliquesAffected) {
    3746                 for (int i=0;i<nCliquesAffected;i++) {
    3747                   int iClique = cliqueStack[i];
    3748                   int size = cliqueCount[iClique];
    3749                   // restore
    3750                   cliqueCount[iClique]= cliqueStart_[iClique+1]-cliqueStart_[iClique];
    3751                   if (!size) {
    3752                     if (logLevel_>1)
    3753                       printf("** could extend clique by adding j!\n");
    3754                   }
    3755                 }
    3756               }
    3757               for (istackR=0;istackR<nstackR;istackR++) {
    3758                 int irow=stackR[istackR];
    3759                 // switch off strengthening if not wanted
    3760                 if ((rowCuts&2)!=0&&goingToTrueBound) {
    3761                   bool ifCut=anyColumnCuts;
    3762                   double gap = rowUpper[irow]-maxR[irow];
    3763                   double sum=0.0;
    3764                   if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    3765                     // see if the strengthened row is a cut
    3766                     for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3767                          kk++) {
    3768                       sum += rowElements[kk]*colsol[column[kk]];
    3769                     }
    3770                     if (sum+gap*colsol[j]>rowUpper[irow]+primalTolerance_||(info->strengthenRow&&rowLower[irow]<-1.0e20)) {
    3771                       // can be a cut
    3772                       // add gap to integer coefficient
    3773                       // saveU and saveL spare
    3774                       int * index = reinterpret_cast<int *>(saveL);
    3775                       double * element = saveU;
    3776                       int n=0;
    3777                       bool coefficientExists=false;
    3778                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3779                            kk++) {
    3780                         if (column[kk]!=j) {
    3781                           index[n]=column[kk];
    3782                           element[n++]=rowElements[kk];
    3783                         } else {
    3784                           double value=rowElements[kk]+gap;
    3785                           if (fabs(value)>1.0e-12) {
    3786                             index[n]=column[kk];
    3787                             element[n++]=value;
    3788                           }
    3789                           coefficientExists=true;
    3790                         }
    3791                       }
    3792                       if (!coefficientExists) {
    3793                         index[n]=j;
    3794                         element[n++]=gap;
    3795                       }
    3796                       OsiRowCut rc;
    3797                       rc.setLb(-DBL_MAX);
    3798                       rc.setUb(rowUpper[irow]+gap*(colUpper[j]-1.0));
    3799                       // effectiveness is how far j moves
    3800                       rc.setEffectiveness((sum+gap*colsol[j]-rowUpper[irow])/gap);
    3801                       if (rc.effectiveness()>needEffectiveness) {
    3802                         rc.setRow(n,index,element,false);
    3803 #if CGL_DEBUG > 0
    3804                         if (debugger) assert(!debugger->invalidCut(rc));
    3805 #endif
    3806                         //if(info->strengthenRow)
    3807                         //printf("c point to row %d\n",irow);
    3808 #ifdef STRENGTHEN_PRINT
    3809                         if (rowLower[irow]<-1.0e20) {
    3810                           printf("7Cut %g <= ",rc.lb());
    3811                           int k;
    3812                           for ( k=0;k<n;k++) {
    3813                             int iColumn = index[k];
    3814                             printf("%g*",element[k]);
    3815                             if (si.isInteger(iColumn))
    3816                               printf("i%d ",iColumn);
    3817                             else
    3818                               printf("x%d ",iColumn);
    3819                           }
    3820                           printf("<= %g\n",rc.ub());
    3821                           printf("Row %g <= ",rowLower[irow]);
    3822                           for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    3823                             int iColumn = column[k];
    3824                             printf("%g*",rowElements[k]);
    3825                             if (si.isInteger(iColumn))
    3826                               printf("i%d ",iColumn);
    3827                             else
    3828                               printf("x%d ",iColumn);
    3829                           }
    3830                           printf("<= %g\n",rowUpper[irow]);
    3831                         }
    3832   #endif
    3833                         int realRow = (rowLower[irow]<-1.0e20) ? irow : -1;
    3834                         if (realRows&&realRow>0)
    3835                           realRow=realRows[realRow];
    3836                         rc.setWhichRow(realRow) ;
    3837                         rowCut.addCutIfNotDuplicate(rc);
    3838                       }
    3839                     }
    3840                   }
    3841                   gap = minR[irow]-rowLower[irow];
    3842                   if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    3843                     // see if the strengthened row is a cut
    3844                     if (!sum) {
    3845                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3846                            kk++) {
    3847                         sum += rowElements[kk]*colsol[column[kk]];
    3848                       }
    3849                     }
    3850                     if (sum-gap*colsol[j]<rowLower[irow]+primalTolerance_||(info->strengthenRow&&rowUpper[irow]>1.0e20)) {
    3851                       // can be a cut
    3852                       // subtract gap from integer coefficient
    3853                       // saveU and saveL spare
    3854                       int * index = reinterpret_cast<int *>(saveL);
    3855                       double * element = saveU;
    3856                       int n=0;
    3857                       bool coefficientExists=false;
    3858                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    3859                            kk++) {
    3860                         if (column[kk]!=j) {
    3861                           index[n]=column[kk];
    3862                           element[n++]=rowElements[kk];
    3863                         } else {
    3864                           double value=rowElements[kk]-gap;
    3865                           if (fabs(value)>1.0e-12) {
    3866                             index[n]=column[kk];
    3867                             element[n++]=value;
    3868                           }
    3869                           coefficientExists=true;
    3870                         }
    3871                       }
    3872                       if (!coefficientExists) {
    3873                         index[n]=j;
    3874                         element[n++]=-gap;
    3875                       }
    3876                       OsiRowCut rc;
    3877                       rc.setLb(rowLower[irow]-gap*(colUpper[j]-1));
    3878                       rc.setUb(DBL_MAX);
    3879                       // effectiveness is how far j moves
    3880                       rc.setEffectiveness((rowLower[irow]-sum+gap*colsol[j])/gap);
    3881                       if (rc.effectiveness()>needEffectiveness) {
    3882                         rc.setRow(n,index,element,false);
    3883 #if CGL_DEBUG > 0
    3884                         if (debugger) assert(!debugger->invalidCut(rc));
    3885 #endif
    3886                         //if(info->strengthenRow)
    3887                         //printf("d point to row %d\n",irow);
    3888 #ifdef STRENGTHEN_PRINT
    3889                         if (rowUpper[irow]>1.0e20) {
    3890                           printf("8Cut %g <= ",rc.lb());
    3891                           int k;
    3892                           for ( k=0;k<n;k++) {
    3893                             int iColumn = index[k];
    3894                             printf("%g*",element[k]);
    3895                             if (si.isInteger(iColumn))
    3896                               printf("i%d ",iColumn);
    3897                             else
    3898                               printf("x%d ",iColumn);
    3899                           }
    3900                           printf("<= %g\n",rc.ub());
    3901                           printf("Row %g <= ",rowLower[irow]);
    3902                           for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    3903                             int iColumn = column[k];
    3904                             printf("%g*",rowElements[k]);
    3905                             if (si.isInteger(iColumn))
    3906                               printf("i%d ",iColumn);
    3907                             else
    3908                               printf("x%d ",iColumn);
    3909                           }
    3910                           printf("<= %g\n",rowUpper[irow]);
    3911                         }
    3912   #endif
    3913                         int realRow = (rowUpper[irow]>1.0e20) ? irow : -1;
    3914                         if (realRows&&realRow>0)
    3915                           realRow=realRows[realRow];
    3916                         rc.setWhichRow(realRow) ;
    3917                         rowCut.addCutIfNotDuplicate(rc);
    3918                       }
    3919                     }
    3920                   }
    3921                 }
    3922                 minR[irow]=saveMin[istackR];
    3923                 maxR[irow]=saveMax[istackR];
    3924                 markR[irow]=-1;
    3925               }
    3926             }
    3927           }
    3928         }
    3929       }
    3930     }
    3931   }
    3932   delete [] cliqueStack;
    3933   delete [] cliqueCount;
    3934   delete [] to_01;
    3935   delete [] stackC0;
    3936   delete [] lo0;
    3937   delete [] up0;
    3938   delete [] tempL;
    3939   delete [] tempU;
    3940   delete [] markC;
    3941   delete [] stackC;
    3942   delete [] stackR;
    3943   delete [] saveL;
    3944   delete [] saveU;
    3945   delete [] saveMin;
    3946   delete [] saveMax;
    3947   delete [] index;
    3948   delete [] element;
    3949   delete [] djs;
    3950   delete [] colsol;
    3951   // Add in row cuts
    3952   if (!ninfeas) {
    3953     rowCut.addCuts(cs,info->strengthenRow,0);
    3954   }
    3955   return (ninfeas);
    3956 }
    3957 
    3958 
    3959 // Does probing and adding cuts for clique slacks
    3960 int
    3961 CglProbing::probeSlacks( const OsiSolverInterface & si,
    3962                           const OsiRowCutDebugger *
    3963 #if CGL_DEBUG > 0
    3964                          debugger
    3965 #endif
    3966                          ,OsiCuts & cs,
    3967                           double * colLower, double * colUpper, CoinPackedMatrix *rowCopy,
    3968                          CoinPackedMatrix *columnCopy,
    3969                           double * rowLower, double * rowUpper,
    3970                           char * intVar, double * minR, double * maxR,int * markR,
    3971                           CglTreeInfo * info) const
    3972 {
    3973   // Check if this is ever used.
    3974   assert(false) ;
    3975   if (!numberCliques_)
    3976     return 0;
    3977   // Set up maxes
    3978   int maxProbe = info->inTree ? maxProbe_ : maxProbeRoot_;
    3979   int maxStack = info->inTree ? maxStack_ : maxStackRoot_;
    3980   int nRows=rowCopy->getNumRows();
    3981   int nCols=rowCopy->getNumCols();
    3982   double * colsol = new double[nCols];
    3983   CoinMemcpyN( si.getColSolution(),nCols,colsol);
    3984   int rowCuts=rowCuts_;
    3985   double_int_pair * array = new double_int_pair [numberCliques_];
    3986   // look at <= cliques
    3987   int iClique;
    3988   int nLook=0;
    3989   for (iClique=0;iClique<numberCliques_;iClique++) {
    3990     if (!cliqueType_[iClique].equality) {
    3991       double sum=0.0;
    3992       for (int j=cliqueStart_[iClique];j<cliqueStart_[iClique+1];j++) {
    3993         int iColumn = sequenceInCliqueEntry(cliqueEntry_[j]);
    3994         double value = colsol[iColumn];
    3995         if (oneFixesInCliqueEntry(cliqueEntry_[j]))
    3996           sum += value;
    3997         else
    3998           sum -= value;
    3999       }
    4000       double away = fabs(0.5-(sum-floor(sum)));
    4001       if (away<0.49999) {
    4002         array[nLook].infeasibility=away;
    4003         array[nLook++].sequence=iClique;
    4004       }
    4005     }
    4006   }
    4007   std::sort(array,array+nLook,double_int_pair_compare());
    4008   nLook=CoinMin(nLook,maxProbe);
    4009   const double * currentColLower = si.getColLower();
    4010   const double * currentColUpper = si.getColUpper();
    4011   double * tempL = new double [nCols];
    4012   double * tempU = new double [nCols];
    4013   int * markC = new int [nCols];
    4014   int * stackC = new int [2*nCols];
    4015   int * stackR = new int [nRows];
    4016   double * saveL = new double [2*nCols];
    4017   double * saveU = new double [2*nCols];
    4018   double * saveMin = new double [nRows];
    4019   double * saveMax = new double [nRows];
    4020   double * element = new double[nCols];
    4021   int * index = new int[nCols];
    4022   // Let us never add more than twice the number of rows worth of row cuts
    4023   // Keep cuts out of cs until end so we can find duplicates quickly
    4024   int nRowsFake = info->inTree ? nRows/3 : nRows;
    4025   CglProbingRowCut rowCut(nRowsFake, !info->inTree);
    4026   const int * column = rowCopy->getIndices();
    4027   const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    4028   const int * rowLength = rowCopy->getVectorLengths();
    4029   const double * rowElements = rowCopy->getElements();
    4030   const int * row = columnCopy->getIndices();
    4031   const CoinBigIndex * columnStart = columnCopy->getVectorStarts();
    4032   const int * columnLength = columnCopy->getVectorLengths();
    4033   const double * columnElements = columnCopy->getElements();
    4034   double movement;
    4035   int i, j, k,kk,jj;
    4036   int kcol,irow,krow;
    4037   bool anyColumnCuts=false;
    4038   double dbound, value, value2;
    4039   int ninfeas=0;
    4040   for (i = 0; i < nCols; ++i) {
    4041     if (colUpper[i]-colLower[i]>1.0e-8) {
    4042       if (colsol[i]<colLower[i]+primalTolerance_) {
    4043         colsol[i]=colLower[i];
    4044       } else if (colsol[i]>colUpper[i]-primalTolerance_) {
    4045         colsol[i]=colUpper[i];
    4046       }
    4047     }
    4048   }
    4049 
    4050   int ipass=0,nfixed=-1;
    4051 
    4052   /* for both way coding */
    4053   int nstackC0=-1;
    4054   int * stackC0 = new int[maxStack];
    4055   double * lo0 = new double[maxStack];
    4056   double * up0 = new double[maxStack];
    4057   int nstackR,nstackC;
    4058   for (i=0;i<nCols;i++) {
    4059     if (colUpper[i]-colLower[i]<1.0e-8) {
    4060       markC[i]=3;
    4061     } else {
    4062       markC[i]=0;
    4063     }
    4064   }
    4065   double tolerance = 1.0e1*primalTolerance_;
    4066   // If we are going to replace coefficient then we don't need to be effective
    4067   int maxPass = info->inTree ? maxPass_ : maxPassRoot_;
    4068   double needEffectiveness = info->strengthenRow ? -1.0e10 : 1.0e-3;
    4069   while (ipass<maxPass&&nfixed) {
    4070     int iLook;
    4071     ipass++;
    4072     nfixed=0;
    4073     for (iLook=0;iLook<nLook;iLook++) {
    4074       double solval;
    4075       double down;
    4076       double up;
    4077       int iClique=array[iLook].sequence;
    4078       solval=0.0;
    4079       j=0;
    4080       for (j=cliqueStart_[iClique];j<cliqueStart_[iClique+1];j++) {
    4081         int iColumn = sequenceInCliqueEntry(cliqueEntry_[j]);
    4082         double value = colsol[iColumn];
    4083         if (oneFixesInCliqueEntry(cliqueEntry_[j]))
    4084           solval += value;
    4085         else
    4086           solval -= value;
    4087       }
    4088       down = floor(solval+tolerance);
    4089       up = ceil(solval-tolerance);
    4090       int istackC,iway, istackR;
    4091       int way[]={1,2,1};
    4092       int feas[]={1,2,4};
    4093       int feasible=0;
    4094       int notFeasible;
    4095       for (iway=0;iway<3;iway ++) {
    4096         int fixThis=0;
    4097         stackC[0]=j;
    4098         markC[j]=way[iway];
    4099         if (way[iway]==1) {
    4100           movement=down-colUpper[j];
    4101           assert(movement<-0.99999);
    4102           down=colLower[j];
    4103         } else {
    4104           movement=up-colLower[j];
    4105           assert(movement>0.99999);
    4106           up=colUpper[j];
    4107         }
    4108         nstackC=1;
    4109         nstackR=0;
    4110         saveL[0]=colLower[j];
    4111         saveU[0]=colUpper[j];
    4112         assert (saveU[0]>saveL[0]);
    4113         notFeasible=0;
    4114         if (movement<0.0) {
    4115           colUpper[j] += movement;
    4116           colUpper[j] = floor(colUpper[j]+0.5);
    4117 #ifdef PRINT_DEBUG
    4118           printf("** Trying %d down to 0\n",j);
    4119 #endif
    4120         } else {
    4121           colLower[j] += movement;
    4122           colLower[j] = floor(colLower[j]+0.5);
    4123 #ifdef PRINT_DEBUG
    4124           printf("** Trying %d up to 1\n",j);
    4125 #endif
    4126         }
    4127         if (fabs(colUpper[j]-colLower[j])<1.0e-6)
    4128           markC[j]=3; // say fixed
    4129         istackC=0;
    4130         /* update immediately */
    4131         for (k=columnStart[j];k<columnStart[j]+columnLength[j];k++) {
    4132           int irow = row[k];
    4133           value = columnElements[k];
    4134           if (markR[irow]==-1) {
    4135             stackR[nstackR]=irow;
    4136             markR[irow]=nstackR;
    4137             saveMin[nstackR]=minR[irow];
    4138             saveMax[nstackR]=maxR[irow];
    4139             nstackR++;
    4140           } else if (markR[irow]==-2) {
    4141             continue;
    4142           }
    4143           /* could check immediately if violation */
    4144           if (movement>0.0) {
    4145             /* up */
    4146             if (value>0.0) {
    4147               /* up does not change - down does */
    4148               if (minR[irow]>-1.0e10)
    4149                 minR[irow] += value;
    4150               if (minR[irow]>rowUpper[irow]+1.0e-5) {
    4151                 notFeasible=1;
    4152                 istackC=1;
    4153                 break;
    4154               }
    4155             } else {
    4156               /* down does not change - up does */
    4157               if (maxR[irow]<1.0e10)
    4158                 maxR[irow] += value;
    4159               if (maxR[irow]<rowLower[irow]-1.0e-5) {
    4160                 notFeasible=1;
    4161                 istackC=1;
    4162                 break;
    4163               }
    4164             }
    4165           } else {
    4166             /* down */
    4167             if (value<0.0) {
    4168               /* up does not change - down does */
    4169               if (minR[irow]>-1.0e10)
    4170                 minR[irow] -= value;
    4171               if (minR[irow]>rowUpper[irow]+1.0e-5) {
    4172                 notFeasible=1;
    4173                 istackC=1;
    4174                 break;
    4175               }
    4176             } else {
    4177               /* down does not change - up does */
    4178               if (maxR[irow]<1.0e10)
    4179                 maxR[irow] -= value;
    4180               if (maxR[irow]<rowLower[irow]-1.0e-5) {
    4181                 notFeasible=1;
    4182                 istackC=1;
    4183                 break;
    4184               }
    4185             }
    4186           }
    4187         }
    4188         while (istackC<nstackC&&nstackC<maxStack) {
    4189           int jway;
    4190           int jcol =stackC[istackC];
    4191           jway=markC[jcol];
    4192           // If not first and fixed then skip
    4193           if (jway==3&&istackC) {
    4194             //istackC++;
    4195             //continue;
    4196             //printf("fixed %d on stack\n",jcol);
    4197           }
    4198           // Do cliques
    4199           if (oneFixStart_&&oneFixStart_[jcol]>=0) {
    4200             int start;
    4201             int end;
    4202             if (colLower[jcol]>saveL[istackC]) {
    4203               // going up
    4204               start = oneFixStart_[jcol];
    4205               end = zeroFixStart_[jcol];
    4206             } else {
    4207               assert (colUpper[jcol]<saveU[istackC]);
    4208               // going down
    4209               start = zeroFixStart_[jcol];
    4210               end = endFixStart_[jcol];
    4211             }
    4212             for (int i=start;i<end;i++) {
    4213               int iClique = whichClique_[i];
    4214               for (int k=cliqueStart_[iClique];k<cliqueStart_[iClique+1];k++) {
    4215                 int kcol = sequenceInCliqueEntry(cliqueEntry_[k]);
    4216                 if (jcol==kcol)
    4217                   continue;
    4218                 int kway = oneFixesInCliqueEntry(cliqueEntry_[k]);
    4219                 if (kcol!=jcol) {
    4220                   if (!markC[kcol]) {
    4221                     // not on list yet
    4222                     if (nstackC<2*maxStack) {
    4223                       markC[kcol] = 3; // say fixed
    4224                       fixThis++;
    4225                       stackC[nstackC]=kcol;
    4226                       saveL[nstackC]=colLower[kcol];
    4227                       saveU[nstackC]=colUpper[kcol];
    4228                       assert (saveU[nstackC]>saveL[nstackC]);
    4229                       nstackC++;
    4230                       if (!kway) {
    4231                         // going up
    4232                         colLower[kcol]=1.0;
    4233                         /* update immediately */
    4234                         for (jj=columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    4235                           krow = row[jj];
    4236                           value = columnElements[jj];
    4237                           if (markR[krow]==-1) {
    4238                             stackR[nstackR]=krow;
    4239                             markR[krow]=nstackR;
    4240                             saveMin[nstackR]=minR[krow];
    4241                             saveMax[nstackR]=maxR[krow];
    4242                             nstackR++;
    4243                           } else if (markR[krow]==-2) {
    4244                             continue;
    4245                           }
    4246                           /* could check immediately if violation */
    4247                           /* up */
    4248                           if (value>0.0) {
    4249                             /* up does not change - down does */
    4250                             if (minR[krow]>-1.0e10)
    4251                               minR[krow] += value;
    4252                             if (minR[krow]>rowUpper[krow]+1.0e-5) {
    4253                               colUpper[kcol]=-1.0e50; /* force infeasible */
    4254                               break;
    4255                             }
    4256                           } else {
    4257                             /* down does not change - up does */
    4258                             if (maxR[krow]<1.0e10)
    4259                               maxR[krow] += value;
    4260                             if (maxR[krow]<rowLower[krow]-1.0e-5) {
    4261                               notFeasible=1;
    4262                               break;
    4263                             }
    4264                           }
    4265                         }
    4266                       } else {
    4267                         // going down
    4268                         colUpper[kcol]=0.0;
    4269                         /* update immediately */
    4270                         for (jj=columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    4271                           krow = row[jj];
    4272                           value = columnElements[jj];
    4273                           if (markR[krow]==-1) {
    4274                             stackR[nstackR]=krow;
    4275                             markR[krow]=nstackR;
    4276                             saveMin[nstackR]=minR[krow];
    4277                             saveMax[nstackR]=maxR[krow];
    4278                             nstackR++;
    4279                           } else if (markR[krow]==-2) {
    4280                             continue;
    4281                           }
    4282                           /* could check immediately if violation */
    4283                           /* down */
    4284                           if (value<0.0) {
    4285                             /* up does not change - down does */
    4286                             if (minR[krow]>-1.0e10)
    4287                               minR[krow] -= value;
    4288                             if (minR[krow]>rowUpper[krow]+1.0e-5) {
    4289                               notFeasible=1;
    4290                               break;
    4291                             }
    4292                           } else {
    4293                             /* down does not change - up does */
    4294                             if (maxR[krow]<1.0e10)
    4295                               maxR[krow] -= value;
    4296                             if (maxR[krow]<rowLower[krow]-1.0e-5) {
    4297                               notFeasible=1;
    4298                               break;
    4299                             }
    4300                           }
    4301                         }
    4302                       }
    4303                     }
    4304                   } else if (markC[kcol]==1) {
    4305                     // marked as going to 0
    4306                     assert (!colUpper[kcol]);
    4307                     if (!kway) {
    4308                       // contradiction
    4309                       notFeasible=1;
    4310                       break;
    4311                     }
    4312                   } else if (markC[kcol]==2) {
    4313                     // marked as going to 1
    4314                     assert (colLower[kcol]);
    4315                     if (kway) {
    4316                       // contradiction
    4317                       notFeasible=1;
    4318                       break;
    4319                     }
    4320                   } else {
    4321                     // marked as fixed
    4322                     assert (markC[kcol]==3);
    4323                     int jkway;
    4324                     if (colLower[kcol])
    4325                       jkway=1;
    4326                     else
    4327                       jkway=0;
    4328                     if (kway==jkway) {
    4329                       // contradiction
    4330                       notFeasible=1;
    4331                       break;
    4332                     }
    4333                   }
    4334                 }
    4335               }
    4336               if (notFeasible)
    4337                 break;
    4338             }
    4339             if (notFeasible)
    4340               istackC=nstackC+1;
    4341           }
    4342           for (k=columnStart[jcol];k<columnStart[jcol]+columnLength[jcol];k++) {
    4343             // break if found not feasible
    4344             if (notFeasible)
    4345               break;
    4346             irow = row[k];
    4347             /*value = columnElements[k];*/
    4348             if (markR[irow]!=-2) {
    4349               /* see if anything forced */
    4350               for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];kk++) {
    4351                 double moveUp=0.0;
    4352                 double moveDown=0.0;
    4353                 double newUpper=-1.0,newLower=1.0;
    4354                 kcol=column[kk];
    4355                 bool onList = (markC[kcol]!=0);
    4356                 if (markC[kcol]!=3) {
    4357                   value2=rowElements[kk];
    4358                   int markIt=markC[kcol];
    4359                   if (value2 < 0.0) {
    4360                     if (colUpper[kcol] < 1e10 && (markIt&2)==0 &&
    4361                         rowUpper[irow]<1.0e10) {
    4362                       dbound = colUpper[kcol]+
    4363                         (rowUpper[irow]-minR[irow])/value2;
    4364                       if (dbound > colLower[kcol] + primalTolerance_) {
    4365                         if (intVar[kcol]) {
    4366                           markIt |= 2;
    4367                           newLower = ceil(dbound-primalTolerance_);
    4368                         } else {
    4369                           newLower=dbound;
    4370                           if (newLower+primalTolerance_>colUpper[kcol]&&
    4371                               newLower-primalTolerance_<=colUpper[kcol]) {
    4372                             newLower=colUpper[kcol];
    4373                             markIt |= 2;
    4374                             markIt=3;
    4375                           } else {
    4376                             // avoid problems - fix later ?
    4377                             markIt=3;
    4378                           }
    4379                         }
    4380                         moveUp = newLower-colLower[kcol];
    4381                       }
    4382                     }
    4383                     if (colLower[kcol] > -1e10 && (markIt&1)==0 &&
    4384                         rowLower[irow]>-1.0e10) {
    4385                       dbound = colLower[kcol] +
    4386                         (rowLower[irow]-maxR[irow])/value2;
    4387                       if (dbound < colUpper[kcol] - primalTolerance_) {
    4388                         if (intVar[kcol]) {
    4389                           markIt |= 1;
    4390                           newUpper = floor(dbound+primalTolerance_);
    4391                         } else {
    4392                           newUpper=dbound;
    4393                           if (newUpper-primalTolerance_<colLower[kcol]&&
    4394                               newUpper+primalTolerance_>=colLower[kcol]) {
    4395                             newUpper=colLower[kcol];
    4396                             markIt |= 1;
    4397                             markIt=3;
    4398                           } else {
    4399                             // avoid problems - fix later ?
    4400                             markIt=3;
    4401                           }
    4402                         }
    4403                         moveDown = newUpper-colUpper[kcol];
    4404                       }
    4405                     }
    4406                   } else {
    4407                     /* positive element */
    4408                     if (colUpper[kcol] < 1e10 && (markIt&2)==0 &&
    4409                         rowLower[irow]>-1.0e10) {
    4410                       dbound = colUpper[kcol] +
    4411                         (rowLower[irow]-maxR[irow])/value2;
    4412                       if (dbound > colLower[kcol] + primalTolerance_) {
    4413                         if (intVar[kcol]) {
    4414                           markIt |= 2;
    4415                           newLower = ceil(dbound-primalTolerance_);
    4416                         } else {
    4417                           newLower=dbound;
    4418                           if (newLower+primalTolerance_>colUpper[kcol]&&
    4419                               newLower-primalTolerance_<=colUpper[kcol]) {
    4420                             newLower=colUpper[kcol];
    4421                             markIt |= 2;
    4422                             markIt=3;
    4423                           } else {
    4424                             // avoid problems - fix later ?
    4425                             markIt=3;
    4426                           }
    4427                         }
    4428                         moveUp = newLower-colLower[kcol];
    4429                       }
    4430                     }
    4431                     if (colLower[kcol] > -1e10 && (markIt&1)==0 &&
    4432                         rowUpper[irow]<1.0e10) {
    4433                       dbound = colLower[kcol] +
    4434                         (rowUpper[irow]-minR[irow])/value2;
    4435                       if (dbound < colUpper[kcol] - primalTolerance_) {
    4436                         if (intVar[kcol]) {
    4437                           markIt |= 1;
    4438                           newUpper = floor(dbound+primalTolerance_);
    4439                         } else {
    4440                           newUpper=dbound;
    4441                           if (newUpper-primalTolerance_<colLower[kcol]&&
    4442                               newUpper+primalTolerance_>=colLower[kcol]) {
    4443                             newUpper=colLower[kcol];
    4444                             markIt |= 1;
    4445                             markIt=3;
    4446                           } else {
    4447                             // avoid problems - fix later ?
    4448                             markIt=3;
    4449                           }
    4450                         }
    4451                         moveDown = newUpper-colUpper[kcol];
    4452                       }
    4453                     }
    4454                   }
    4455                   if (nstackC<2*maxStack) {
    4456                     markC[kcol] = markIt;
    4457                   }
    4458                   if (moveUp&&nstackC<2*maxStack) {
    4459                     fixThis++;
    4460 #ifdef PRINT_DEBUG
    4461                     printf("lower bound on %d increased from %g to %g by row %d %g %g\n",kcol,colLower[kcol],newLower,irow,rowLower[irow],rowUpper[irow]);
    4462                     value=0.0;
    4463                     for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    4464                       int ii=column[jj];
    4465                       if (colUpper[ii]-colLower[ii]>primalTolerance_) {
    4466                         printf("(%d, %g) ",ii,rowElements[jj]);
    4467                       } else {
    4468                         value += rowElements[jj]*colLower[ii];
    4469                       }
    4470                     }
    4471                     printf(" - fixed %g\n",value);
    4472                     for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    4473                       int ii=column[jj];
    4474                       if (colUpper[ii]-colLower[ii]<primalTolerance_) {
    4475                         printf("(%d, %g, %g) ",ii,rowElements[jj],colLower[ii]);
    4476                       }
    4477                     }
    4478                     printf("\n");
    4479 #endif
    4480                     if (!onList) {
    4481                       stackC[nstackC]=kcol;
    4482                       saveL[nstackC]=colLower[kcol];
    4483                       saveU[nstackC]=colUpper[kcol];
    4484                       assert (saveU[nstackC]>saveL[nstackC]);
    4485                       nstackC++;
    4486                       onList=true;
    4487                     }
    4488                     if (intVar[kcol])
    4489                       newLower = CoinMax(colLower[kcol],ceil(newLower-1.0e-4));
    4490                     colLower[kcol]=newLower;
    4491                     if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
    4492                       markC[kcol]=3; // say fixed
    4493                     }
    4494                     /* update immediately */
    4495                     for (jj=columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    4496                       krow = row[jj];
    4497                       value = columnElements[jj];
    4498                       if (markR[krow]==-1) {
    4499                         stackR[nstackR]=krow;
    4500                         markR[krow]=nstackR;
    4501                         saveMin[nstackR]=minR[krow];
    4502                         saveMax[nstackR]=maxR[krow];
    4503                         nstackR++;
    4504                       } else if (markR[krow]==-2) {
    4505                         continue;
    4506                       }
    4507                       /* could check immediately if violation */
    4508                       /* up */
    4509                       if (value>0.0) {
    4510                         /* up does not change - down does */
    4511                         if (minR[krow]>-1.0e10)
    4512                           minR[krow] += value*moveUp;
    4513                         if (minR[krow]>rowUpper[krow]+1.0e-5) {
    4514                           colUpper[kcol]=-1.0e50; /* force infeasible */
    4515                           break;
    4516                         }
    4517                       } else {
    4518                         /* down does not change - up does */
    4519                         if (maxR[krow]<1.0e10)
    4520                           maxR[krow] += value*moveUp;
    4521                         if (maxR[krow]<rowLower[krow]-1.0e-5) {
    4522                           colUpper[kcol]=-1.0e50; /* force infeasible */
    4523                           break;
    4524                         }
    4525                       }
    4526                     }
    4527                   }
    4528                   if (moveDown&&nstackC<2*maxStack) {
    4529                     fixThis++;
    4530 #ifdef PRINT_DEBUG
    4531                     printf("upper bound on %d decreased from %g to %g by row %d %g %g\n",kcol,colUpper[kcol],newUpper,irow,rowLower[irow],rowUpper[irow]);
    4532                     value=0.0;
    4533                     for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    4534                       int ii=column[jj];
    4535                       if (colUpper[ii]-colLower[ii]>primalTolerance_) {
    4536                         printf("(%d, %g) ",ii,rowElements[jj]);
    4537                       } else {
    4538                         value += rowElements[jj]*colLower[ii];
    4539                       }
    4540                     }
    4541                     printf(" - fixed %g\n",value);
    4542                     for (jj=rowStart[irow];jj<rowStart[irow]+rowLength[irow];jj++) {
    4543                       int ii=column[jj];
    4544                       if (colUpper[ii]-colLower[ii]<primalTolerance_) {
    4545                         printf("(%d, %g, %g) ",ii,rowElements[jj],colLower[ii]);
    4546                       }
    4547                     }
    4548                     printf("\n");
    4549 #endif
    4550                     if (!onList) {
    4551                       stackC[nstackC]=kcol;
    4552                       saveL[nstackC]=colLower[kcol];
    4553                       saveU[nstackC]=colUpper[kcol];
    4554                       assert (saveU[nstackC]>saveL[nstackC]);
    4555                       nstackC++;
    4556                       onList=true;
    4557                     }
    4558                     if (intVar[kcol])
    4559                       newUpper = CoinMin(colUpper[kcol],floor(newUpper+1.0e-4));
    4560                     colUpper[kcol]=newUpper;
    4561                     if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
    4562                       markC[kcol]=3; // say fixed
    4563                     }
    4564                     /* update immediately */
    4565                     for (jj=columnStart[kcol];jj<columnStart[kcol]+columnLength[kcol];jj++) {
    4566                       krow = row[jj];
    4567                       value = columnElements[jj];
    4568                       if (markR[krow]==-1) {
    4569                         stackR[nstackR]=krow;
    4570                         markR[krow]=nstackR;
    4571                         saveMin[nstackR]=minR[krow];
    4572                         saveMax[nstackR]=maxR[krow];
    4573                         nstackR++;
    4574                       } else if (markR[krow]==-2) {
    4575                         continue;
    4576                       }
    4577                       /* could check immediately if violation */
    4578                       /* down */
    4579                       if (value<0.0) {
    4580                         /* up does not change - down does */
    4581                         if (minR[krow]>-1.0e10)
    4582                           minR[krow] += value*moveDown;
    4583                         if (minR[krow]>rowUpper[krow]+1.0e-5) {
    4584                           colUpper[kcol]=-1.0e50; /* force infeasible */
    4585                           break;
    4586                         }
    4587                       } else {
    4588                         /* down does not change - up does */
    4589                         if (maxR[krow]<1.0e10)
    4590                           maxR[krow] += value*moveDown;
    4591                         if (maxR[krow]<rowLower[krow]-1.0e-5) {
    4592                           colUpper[kcol]=-1.0e50; /* force infeasible */
    4593                           break;
    4594                         }
    4595                       }
    4596                     }
    4597                   }
    4598                   if (colLower[kcol]>colUpper[kcol]+primalTolerance_) {
    4599                     notFeasible=1;;
    4600                     k=columnStart[jcol]+columnLength[jcol];
    4601                     istackC=nstackC+1;
    4602 #ifdef PRINT_DEBUG
    4603                     printf("** not feasible this way\n");
    4604 #endif
    4605                     break;
    4606                   }
    4607                 }
    4608               }
    4609             }
    4610           }
    4611           istackC++;
    4612         }
    4613         if (!notFeasible) {
    4614           feasible |= feas[iway];
    4615         } else if (iway==1&&feasible==0) {
    4616           /* not feasible at all */
    4617           ninfeas=1;
    4618           j=nCols-1;
    4619           iLook=nLook;
    4620           ipass=maxPass;
    4621           break;
    4622         }
    4623         if (iway==2||(iway==1&&feasible==2)) {
    4624           /* keep */
    4625           iway=3;
    4626           nfixed++;
    4627           if (mode_) {
    4628             OsiColCut cc;
    4629             int nTot=0,nFix=0,nInt=0;
    4630             bool ifCut=false;
    4631             for (istackC=0;istackC<nstackC;istackC++) {
    4632               int icol=stackC[istackC];
    4633               if (intVar[icol]) {
    4634                 if (colUpper[icol]<currentColUpper[icol]-1.0e-4) {
    4635                   element[nFix]=colUpper[icol];
    4636                   index[nFix++]=icol;
    4637                   nInt++;
    4638                   if (colsol[icol]>colUpper[icol]+primalTolerance_) {
    4639                     ifCut=true;
    4640                     anyColumnCuts=true;
    4641                   }
    4642                 }
    4643               }
    4644             }
    4645             if (nFix) {
    4646               nTot=nFix;
    4647               cc.setUbs(nFix,index,element);
    4648               nFix=0;
    4649             }
    4650             for (istackC=0;istackC<nstackC;istackC++) {
    4651               int icol=stackC[istackC];
    4652               if (intVar[icol]) {
    4653                 if (colLower[icol]>currentColLower[icol]+1.0e-4) {
    4654                   element[nFix]=colLower[icol];
    4655                   index[nFix++]=icol;
    4656                   nInt++;
    4657                   if (colsol[icol]<colLower[icol]-primalTolerance_) {
    4658                     ifCut=true;
    4659                     anyColumnCuts=true;
    4660                   }
    4661                 }
    4662               }
    4663             }
    4664             if (nFix) {
    4665               nTot+=nFix;
    4666               cc.setLbs(nFix,index,element);
    4667             }
    4668             // could tighten continuous as well
    4669             if (nInt) {
    4670               if (ifCut) {
    4671                 cc.setEffectiveness(100.0);
    4672               } else {
    4673                 cc.setEffectiveness(1.0e-5);
    4674               }
    4675 #if CGL_DEBUG > 0
    4676               CglProbingDebug::checkBounds(si,cc);
    4677 #endif
    4678               cs.insert(cc);
    4679             }
    4680           }
    4681           for (istackC=0;istackC<nstackC;istackC++) {
    4682             int icol=stackC[istackC];
    4683             if (colUpper[icol]-colLower[icol]>primalTolerance_) {
    4684               markC[icol]=0;
    4685             } else {
    4686               markC[icol]=3;
    4687             }
    4688           }
    4689           for (istackR=0;istackR<nstackR;istackR++) {
    4690             int irow=stackR[istackR];
    4691             markR[irow]=-1;
    4692           }
    4693         } else {
    4694           /* is it worth seeing if can increase coefficients
    4695              or maybe better see if it is a cut */
    4696           if (iway==0) {
    4697             nstackC0=CoinMin(nstackC,maxStack);
    4698             if (notFeasible) {
    4699               nstackC0=0;
    4700             } else {
    4701               for (istackC=0;istackC<nstackC0;istackC++) {
    4702                 int icol=stackC[istackC];
    4703                 stackC0[istackC]=icol;
    4704                 lo0[istackC]=colLower[icol];
    4705                 up0[istackC]=colUpper[icol];
    4706               }
    4707             }
    4708             /* restore all */
    4709             assert (iway==0);
    4710             for (istackC=nstackC-1;istackC>=0;istackC--) {
    4711               int icol=stackC[istackC];
    4712               double oldU=saveU[istackC];
    4713               double oldL=saveL[istackC];
    4714               colUpper[icol]=oldU;
    4715               colLower[icol]=oldL;
    4716               markC[icol]=0;
    4717             }
    4718             for (istackR=0;istackR<nstackR;istackR++) {
    4719               int irow=stackR[istackR];
    4720               // switch off strengthening if not wanted
    4721               if ((rowCuts&2)!=0) {
    4722                 bool ifCut=anyColumnCuts;
    4723                 double gap = rowUpper[irow]-maxR[irow];
    4724                 double sum=0.0;
    4725                 if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    4726                   // see if the strengthened row is a cut
    4727                   for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    4728                        kk++) {
    4729                     sum += rowElements[kk]*colsol[column[kk]];
    4730                   }
    4731                   if (sum-gap*colsol[j]>maxR[irow]+primalTolerance_||info->strengthenRow) {
    4732                     // can be a cut
    4733                     // subtract gap from upper and integer coefficient
    4734                     // saveU and saveL spare
    4735                     int * index = reinterpret_cast<int *>(saveL);
    4736                     double * element = saveU;
    4737                     int n=0;
    4738                     bool coefficientExists=false;
    4739                     for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    4740                          kk++) {
    4741                       if (column[kk]!=j) {
    4742                         index[n]=column[kk];
    4743                         element[n++]=rowElements[kk];
    4744                       } else {
    4745                         double value=rowElements[kk]-gap;
    4746                         if (fabs(value)>1.0e-12) {
    4747                           index[n]=column[kk];
    4748                           element[n++]=value;
    4749                         }
    4750                         coefficientExists=true;
    4751                       }
    4752                     }
    4753                     if (!coefficientExists) {
    4754                       index[n]=j;
    4755                       element[n++]=-gap;
    4756                     }
    4757                     OsiRowCut rc;
    4758                     rc.setLb(-DBL_MAX);
    4759                     rc.setUb(rowUpper[irow]-gap*(colLower[j]+1.0));
    4760                     // effectiveness is how far j moves
    4761                     rc.setEffectiveness((sum-gap*colsol[j]-maxR[irow])/gap);
    4762                     if (rc.effectiveness()>needEffectiveness) {
    4763                       rc.setRow(n,index,element,false);
    4764 #if CGL_DEBUG > 0
    4765                       if (debugger) assert(!debugger->invalidCut(rc));
    4766 #endif
    4767                       // If strengthenRow point to row
    4768                       //if(info->strengthenRow)
    4769                       //printf("a point to row %d\n",irow);
    4770 #ifdef STRENGTHEN_PRINT
    4771                       {
    4772                         printf("9Cut %g <= ",rc.lb());
    4773                         int k;
    4774                         for ( k=0;k<n;k++) {
    4775                           int iColumn = index[k];
    4776                           printf("%g*",element[k]);
    4777                           if (si.isInteger(iColumn))
    4778                             printf("i%d ",iColumn);
    4779                           else
    4780                             printf("x%d ",iColumn);
    4781                         }
    4782                         printf("<= %g\n",rc.ub());
    4783                         printf("Row %g <= ",rowLower[irow]);
    4784                         for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    4785                           int iColumn = column[k];
    4786                           printf("%g*",rowElements[k]);
    4787                           if (si.isInteger(iColumn))
    4788                             printf("i%d ",iColumn);
    4789                           else
    4790                             printf("x%d ",iColumn);
    4791                         }
    4792                         printf("<= %g\n",rowUpper[irow]);
    4793                       }
    4794 #endif
    4795                       rc.setWhichRow(irow) ;
    4796                       rowCut.addCutIfNotDuplicate(rc);
    4797                     }
    4798                   }
    4799                 }
    4800                 gap = minR[irow]-rowLower[irow];
    4801                 if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    4802                   // see if the strengthened row is a cut
    4803                   if (!sum) {
    4804                     for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    4805                          kk++) {
    4806                       sum += rowElements[kk]*colsol[column[kk]];
    4807                     }
    4808                   }
    4809                   if (sum+gap*colsol[j]<minR[irow]+primalTolerance_||info->strengthenRow) {
    4810                     // can be a cut
    4811                     // add gap to lower and integer coefficient
    4812                     // saveU and saveL spare
    4813                     int * index = reinterpret_cast<int *>(saveL);
    4814                     double * element = saveU;
    4815                     int n=0;
    4816                     bool coefficientExists=false;
    4817                     for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    4818                          kk++) {
    4819                       if (column[kk]!=j) {
    4820                         index[n]=column[kk];
    4821                         element[n++]=rowElements[kk];
    4822                       } else {
    4823                         double value=rowElements[kk]+gap;
    4824                         if (fabs(value)>1.0e-12) {
    4825                           index[n]=column[kk];
    4826                           element[n++]=value;
    4827                         }
    4828                         coefficientExists=true;
    4829                       }
    4830                     }
    4831                     if (!coefficientExists) {
    4832                       index[n]=j;
    4833                       element[n++]=gap;
    4834                     }
    4835                     OsiRowCut rc;
    4836                     rc.setLb(rowLower[irow]+gap*(colLower[j]+1.0));
    4837                     rc.setUb(DBL_MAX);
    4838                     // effectiveness is how far j moves
    4839                     rc.setEffectiveness((minR[irow]-sum-gap*colsol[j])/gap);
    4840                     if (rc.effectiveness()>needEffectiveness) {
    4841                       rc.setRow(n,index,element,false);
    4842 #if CGL_DEBUG > 0
    4843                       if (debugger) assert(!debugger->invalidCut(rc));
    4844 #endif
    4845                       //if(info->strengthenRow)
    4846                       //printf("b point to row %d\n",irow);
    4847 #ifdef STRENGTHEN_PRINT
    4848                       {
    4849                         printf("10Cut %g <= ",rc.lb());
    4850                         int k;
    4851                         for ( k=0;k<n;k++) {
    4852                           int iColumn = index[k];
    4853                           printf("%g*",element[k]);
    4854                           if (si.isInteger(iColumn))
    4855                             printf("i%d ",iColumn);
    4856                           else
    4857                             printf("x%d ",iColumn);
    4858                         }
    4859                         printf("<= %g\n",rc.ub());
    4860                         printf("Row %g <= ",rowLower[irow]);
    4861                         for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    4862                           int iColumn = column[k];
    4863                           printf("%g*",rowElements[k]);
    4864                           if (si.isInteger(iColumn))
    4865                             printf("i%d ",iColumn);
    4866                           else
    4867                             printf("x%d ",iColumn);
    4868                         }
    4869                         printf("<= %g\n",rowUpper[irow]);
    4870                       }
    4871 #endif
    4872                       rc.setWhichRow(irow) ;
    4873                       rowCut.addCutIfNotDuplicate(rc);
    4874                     }
    4875                   }
    4876                 }
    4877               }
    4878               minR[irow]=saveMin[istackR];
    4879               maxR[irow]=saveMax[istackR];
    4880               markR[irow]=-1;
    4881             }
    4882           } else {
    4883             if (iway==1&&feasible==3) {
    4884               iway=3;
    4885               /* point back to stack */
    4886               for (istackC=nstackC-1;istackC>=0;istackC--) {
    4887                 int icol=stackC[istackC];
    4888                 markC[icol]=istackC+1000;
    4889               }
    4890               if (mode_) {
    4891                 OsiColCut cc;
    4892                 int nTot=0,nFix=0,nInt=0;
    4893                 bool ifCut=false;
    4894                 for (istackC=0;istackC<nstackC0;istackC++) {
    4895                   int icol=stackC0[istackC];
    4896                   int istackC1=markC[icol]-1000;
    4897                   if (istackC1>=0) {
    4898                     if (CoinMin(lo0[istackC],colLower[icol])>saveL[istackC1]+1.0e-4) {
    4899                       saveL[istackC1]=CoinMin(lo0[istackC],colLower[icol]);
    4900                       if (intVar[icol]) {
    4901                         element[nFix]=saveL[istackC1];
    4902                         index[nFix++]=icol;
    4903                         nInt++;
    4904                         if (colsol[icol]<saveL[istackC1]-primalTolerance_)
    4905                           ifCut=true;
    4906                       }
    4907                       nfixed++;
    4908                     }
    4909                   }
    4910                 }
    4911                 if (nFix) {
    4912                   nTot=nFix;
    4913                   cc.setLbs(nFix,index,element);
    4914                   nFix=0;
    4915                 }
    4916                 for (istackC=0;istackC<nstackC0;istackC++) {
    4917                   int icol=stackC0[istackC];
    4918                   int istackC1=markC[icol]-1000;
    4919                   if (istackC1>=0) {
    4920                     if (CoinMax(up0[istackC],colUpper[icol])<saveU[istackC1]-1.0e-4) {
    4921                       saveU[istackC1]=CoinMax(up0[istackC],colUpper[icol]);
    4922                       if (intVar[icol]) {
    4923                         element[nFix]=saveU[istackC1];
    4924                         index[nFix++]=icol;
    4925                         nInt++;
    4926                         if (colsol[icol]>saveU[istackC1]+primalTolerance_)
    4927                           ifCut=true;
    4928                       }
    4929                       nfixed++;
    4930                     }
    4931                   }
    4932                 }
    4933                 if (nFix) {
    4934                   nTot+=nFix;
    4935                   cc.setUbs(nFix,index,element);
    4936                 }
    4937                 // could tighten continuous as well
    4938                 if (nInt) {
    4939                   if (ifCut) {
    4940                     cc.setEffectiveness(100.0);
    4941                   } else {
    4942                     cc.setEffectiveness(1.0e-5);
    4943                   }
    4944 #if CGL_DEBUG > 0
    4945                   CglProbingDebug::checkBounds(si,cc);
    4946 #endif
    4947                   cs.insert(cc);
    4948                 }
    4949               }
    4950             }
    4951             /* restore all */
    4952             for (istackC=nstackC-1;istackC>=0;istackC--) {
    4953               int icol=stackC[istackC];
    4954               double oldU=saveU[istackC];
    4955               double oldL=saveL[istackC];
    4956               colUpper[icol]=oldU;
    4957               colLower[icol]=oldL;
    4958               if (oldU>oldL+1.0e-4)
    4959                 markC[icol]=0;
    4960               else
    4961                 markC[icol]=3;
    4962             }
    4963             for (istackR=0;istackR<nstackR;istackR++) {
    4964               int irow=stackR[istackR];
    4965               // switch off strengthening if not wanted
    4966               if ((rowCuts&2)!=0) {
    4967                 bool ifCut=anyColumnCuts;
    4968                 double gap = rowUpper[irow]-maxR[irow];
    4969                 double sum=0.0;
    4970                 if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    4971                   // see if the strengthened row is a cut
    4972                   for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    4973                        kk++) {
    4974                     sum += rowElements[kk]*colsol[column[kk]];
    4975                   }
    4976                   if (sum+gap*colsol[j]>rowUpper[irow]+primalTolerance_||info->strengthenRow) {
    4977                     // can be a cut
    4978                     // add gap to integer coefficient
    4979                     // saveU and saveL spare
    4980                     int * index = reinterpret_cast<int *>(saveL);
    4981                     double * element = saveU;
    4982                     int n=0;
    4983                     bool coefficientExists=false;
    4984                     for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    4985                          kk++) {
    4986                       if (column[kk]!=j) {
    4987                         index[n]=column[kk];
    4988                         element[n++]=rowElements[kk];
    4989                         } else {
    4990                           double value=rowElements[kk]+gap;
    4991                           if (fabs(value)>1.0e-12) {
    4992                             index[n]=column[kk];
    4993                             element[n++]=value;
    4994                           }
    4995                           coefficientExists=true;
    4996                         }
    4997                       }
    4998                       if (!coefficientExists) {
    4999                         index[n]=j;
    5000                         element[n++]=gap;
    5001                       }
    5002                       OsiRowCut rc;
    5003                       rc.setLb(-DBL_MAX);
    5004                       rc.setUb(rowUpper[irow]+gap*(colUpper[j]-1.0));
    5005                       // effectiveness is how far j moves
    5006                       rc.setEffectiveness((sum+gap*colsol[j]-rowUpper[irow])/gap);
    5007                       if (rc.effectiveness()>needEffectiveness) {
    5008                         rc.setRow(n,index,element,false);
    5009 #if CGL_DEBUG > 0
    5010                         if (debugger) assert(!debugger->invalidCut(rc));
    5011 #endif
    5012                         //if(info->strengthenRow)
    5013                         //printf("c point to row %d\n",irow);
    5014 #ifdef STRENGTHEN_PRINT
    5015                       {
    5016                         printf("11Cut %g <= ",rc.lb());
    5017                         int k;
    5018                         for ( k=0;k<n;k++) {
    5019                           int iColumn = index[k];
    5020                           printf("%g*",element[k]);
    5021                           if (si.isInteger(iColumn))
    5022                             printf("i%d ",iColumn);
    5023                           else
    5024                             printf("x%d ",iColumn);
    5025                         }
    5026                         printf("<= %g\n",rc.ub());
    5027                         printf("Row %g <= ",rowLower[irow]);
    5028                         for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    5029                           int iColumn = column[k];
    5030                           printf("%g*",rowElements[k]);
    5031                           if (si.isInteger(iColumn))
    5032                             printf("i%d ",iColumn);
    5033                           else
    5034                             printf("x%d ",iColumn);
    5035                         }
    5036                         printf("<= %g\n",rowUpper[irow]);
    5037                       }
    5038 #endif
    5039                         rc.setWhichRow(irow) ;
    5040                         rowCut.addCutIfNotDuplicate(rc);
    5041                       }
    5042                     }
    5043                   }
    5044                   gap = minR[irow]-rowLower[irow];
    5045                   if (!ifCut&&(gap>primalTolerance_&&gap<1.0e8)) {
    5046                     // see if the strengthened row is a cut
    5047                     if (!sum) {
    5048                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    5049                            kk++) {
    5050                         sum += rowElements[kk]*colsol[column[kk]];
    5051                       }
    5052                     }
    5053                     if (sum-gap*colsol[j]<rowLower[irow]+primalTolerance_||info->strengthenRow) {
    5054                       // can be a cut
    5055                       // subtract gap from integer coefficient
    5056                       // saveU and saveL spare
    5057                       int * index = reinterpret_cast<int *>(saveL);
    5058                       double * element = saveU;
    5059                       int n=0;
    5060                       bool coefficientExists=false;
    5061                       for (kk=rowStart[irow];kk<rowStart[irow]+rowLength[irow];
    5062                            kk++) {
    5063                         if (column[kk]!=j) {
    5064                           index[n]=column[kk];
    5065                           element[n++]=rowElements[kk];
    5066                         } else {
    5067                           double value=rowElements[kk]-gap;
    5068                           if (fabs(value)>1.0e-12) {
    5069                             index[n]=column[kk];
    5070                             element[n++]=value;
    5071                           }
    5072                           coefficientExists=true;
    5073                         }
    5074                       }
    5075                       if (!coefficientExists) {
    5076                         index[n]=j;
    5077                         element[n++]=-gap;
    5078                       }
    5079                       OsiRowCut rc;
    5080                       rc.setLb(rowLower[irow]-gap*(colUpper[j]-1));
    5081                       rc.setUb(DBL_MAX);
    5082                       // effectiveness is how far j moves
    5083                       rc.setEffectiveness((rowLower[irow]-sum+gap*colsol[j])/gap);
    5084                       if (rc.effectiveness()>needEffectiveness) {
    5085                         rc.setRow(n,index,element,false);
    5086 #if CGL_DEBUG > 0
    5087                         if (debugger) assert(!debugger->invalidCut(rc));
    5088 #endif
    5089                         //if(info->strengthenRow)
    5090                         //printf("d point to row %d\n",irow);
    5091 #ifdef STRENGTHEN_PRINT
    5092                       {
    5093                         printf("12Cut %g <= ",rc.lb());
    5094                         int k;
    5095                         for ( k=0;k<n;k++) {
    5096                           int iColumn = index[k];
    5097                           printf("%g*",element[k]);
    5098                           if (si.isInteger(iColumn))
    5099                             printf("i%d ",iColumn);
    5100                           else
    5101                             printf("x%d ",iColumn);
    5102                         }
    5103                         printf("<= %g\n",rc.ub());
    5104                         printf("Row %g <= ",rowLower[irow]);
    5105                         for (k=rowStart[irow];k<rowStart[irow]+rowLength[irow];k++) {
    5106                           int iColumn = column[k];
    5107                           printf("%g*",rowElements[k]);
    5108                           if (si.isInteger(iColumn))
    5109                             printf("i%d ",iColumn);
    5110                           else
    5111                             printf("x%d ",iColumn);
    5112                         }
    5113                         printf("<= %g\n",rowUpper[irow]);
    5114                       }
    5115 #endif
    5116                         rc.setWhichRow(irow) ;
    5117                         rowCut.addCutIfNotDuplicate(rc);
    5118                       }
    5119                     }
    5120                   }
    5121                 }
    5122                 minR[irow]=saveMin[istackR];
    5123                 maxR[irow]=saveMax[istackR];
    5124                 markR[irow]=-1;
    5125               }
    5126             }
    5127           }
    5128         }
    5129     }
    5130   }
    5131   delete [] stackC0;
    5132   delete [] lo0;
    5133   delete [] up0;
    5134   delete [] tempL;
    5135   delete [] tempU;
    5136   delete [] markC;
    5137   delete [] stackC;
    5138   delete [] stackR;
    5139   delete [] saveL;
    5140   delete [] saveU;
    5141   delete [] saveMin;
    5142   delete [] saveMax;
    5143   delete [] index;
    5144   delete [] element;
    5145   delete [] colsol;
    5146   // Add in row cuts
    5147   if (!ninfeas) {
    5148     rowCut.addCuts(cs,info->strengthenRow,0);
    5149   }
    5150   delete [] array;
    5151   abort();
    5152   return (ninfeas);
    5153 }
    5154 
    5155 /* Set mode
     2637  Set mode
    51562638
    51572639  Take care to change only the basic mode, without affecting higher bits.
     
    57113193}
    57123194
    5713 
    5714 /* Creates cliques for use by probing.
    5715    Can also try and extend cliques as a result of probing (root node).
    5716    Returns number of cliques found.
     3195/*
     3196  Creates cliques for use by probing.
     3197  Can also try and extend cliques as a result of probing (root node).
     3198  Returns number of cliques found.
     3199
     3200  minimumSize and maximumSize limit the size of the cliques reported.
     3201
     3202  Mathematically, we look for constraints of the form
     3203      L<i> <= SUM{P}x<j> - SUM{M}x<j> <= U<i>
     3204  where all x<j> are binary, all x<j> j in P have a coefficient of +1, and
     3205  all x<j>, j in M, have a coefficient of -1. The smallest value of the lhs
     3206  is -|M|, the largest, |P|. We can make the following conclusions (numbered to
     3207  match state codes below):
     3208
     3209  3) -|M| > U<i> or |P| < L<i>  ==> infeasible
     3210
     3211  2) -|M| = U<i> ==> all x<j>, j in M, at 1, all x<j>, j in P, at 0
     3212      |P| = L<i> ==> all x<j>, j in P, at 1, all x<j>, j in M, at 0
     3213
     3214  1) U<i> = -|M|+1 ==> any x<t>, t in P, at 1 ==> all x<j>, j in M, at 1
     3215                                              ==> all x<j>, j in P\t, at 0
     3216                   ==> any x<t>, t in M, at 0 ==> all x<j>, j in P, at 0
     3217                                              ==> all x<j>, j in M\t, at 1
     3218     L<i> =  |P|-1 ==> any x<t>, j in M, at 1 ==> all x<j>, j in P, at 1
     3219                                              ==> all x<j>, j in M\t, at 0
     3220                   ==> any x<t>, j in P, at 0 ==> all x<j>, j in M, at 0
     3221                                              ==> all x<j>, j in P\t, at 1
     3222
     3223  1) is the only case where we have an interesting clique to work with.
     3224
    57173225*/
    57183226int
    5719 CglProbing::createCliques( OsiSolverInterface & si,
    5720                           int minimumSize, int maximumSize)
    5721 {
    5722   // get rid of what is there
    5723   deleteCliques();
    5724   CoinPackedMatrix matrixByRow(*si.getMatrixByRow());
    5725   int numberRows = si.getNumRows();
     3227CglProbing::createCliques (OsiSolverInterface &si,
     3228                           int minimumSize, int maximumSize)
     3229{
     3230  // Remove existing clique information
     3231  deleteCliques() ;
     3232/*
     3233  Generate cliques from the original matrix. Don't use a snapshot, even if
     3234  we have one.
     3235*/
     3236  CoinPackedMatrix matrixByRow(*si.getMatrixByRow()) ;
     3237  int numberRows = si.getNumRows() ;
    57263238  if (!rowCopy_)
    5727     numberRows_=numberRows;
    5728   numberColumns_ = si.getNumCols();
    5729 
    5730   numberCliques_=0;
    5731   int numberEntries=0;
    5732   int numberIntegers=0;
    5733   int * lookup = new int[numberColumns_];
    5734   int i;
    5735   for (i=0;i<numberColumns_;i++) {
     3239    numberRows_ = numberRows ;
     3240  numberColumns_ = si.getNumCols() ;
     3241/*
     3242  Initialise some bookkeeping. lookup will enable us to find the integer
     3243  variable's sequence number given the column number.
     3244*/
     3245  numberCliques_ = 0 ;
     3246  int numberEntries = 0 ;
     3247  int numberIntegers = 0 ;
     3248  int *lookup = new int[numberColumns_] ;
     3249  int i ;
     3250  for (i = 0 ; i < numberColumns_ ; i++) {
    57363251    if (si.isBinary(i))
    5737       lookup[i]=numberIntegers++;
     3252      lookup[i] = numberIntegers++ ;
    57383253    else
    5739       lookup[i]=-1;
    5740   }
    5741 
    5742   int * which = new int[numberColumns_];
    5743   int * whichRow = new int[numberRows];
     3254      lookup[i] = -1 ;
     3255  }
     3256  int *which = new int[numberColumns_] ;
     3257  int *whichRow = new int[numberRows] ;
     3258
    57443259  // Statistics
    5745   int totalP1=0,totalM1=0;
    5746   int numberBig=0,totalBig=0;
    5747   int numberFixed=0;
    5748 
    5749   // Row copy
    5750   const double * elementByRow = matrixByRow.getElements();
    5751   const int * column = matrixByRow.getIndices();
    5752   const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
    5753   const int * rowLength = matrixByRow.getVectorLengths();
    5754 
    5755   // Column lengths for slacks
    5756   const int * columnLength = si.getMatrixByCol()->getVectorLengths();
    5757 
    5758   const double * lower = si.getColLower();
    5759   const double * upper = si.getColUpper();
    5760   const double * rowLower = si.getRowLower();
    5761   const double * rowUpper = si.getRowUpper();
    5762   int iRow;
    5763   for (iRow=0;iRow<numberRows;iRow++) {
    5764     int numberP1=0, numberM1=0;
    5765     int j;
    5766     double upperValue=rowUpper[iRow];
    5767     double lowerValue=rowLower[iRow];
    5768     bool good=true;
    5769     int slack = -1;
    5770     for (j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    5771       int iColumn = column[j];
    5772       int iInteger=lookup[iColumn];
    5773       if (upper[iColumn]-lower[iColumn]<1.0e-8) {
    5774         // fixed
    5775         upperValue -= lower[iColumn]*elementByRow[j];
    5776         lowerValue -= lower[iColumn]*elementByRow[j];
    5777         continue;
    5778       } else if (upper[iColumn]!=1.0||lower[iColumn]!=0.0) {
    5779         good = false;
    5780         break;
    5781       } else if (iInteger<0) {
    5782         good = false;
    5783         break;
     3260  int totalP1 = 0 ;
     3261  int totalM1 = 0 ;
     3262  int numberBig = 0 ;
     3263  int totalBig = 0 ;
     3264  int numberFixed = 0 ;
     3265
     3266  // Expose the internal arrays of the matrix
     3267  const double *elementByRow = matrixByRow.getElements() ;
     3268  const int *column = matrixByRow.getIndices() ;
     3269  const CoinBigIndex *rowStart = matrixByRow.getVectorStarts() ;
     3270  const int *rowLength = matrixByRow.getVectorLengths() ;
     3271
     3272  /*
     3273    jjf: Column lengths for slacks
     3274   
     3275    See note below re. slacks. Vestigial  -- lh, 100924
     3276  */
     3277  const int *columnLength = si.getMatrixByCol()->getVectorLengths() ;
     3278
     3279  const double *lower = si.getColLower() ;
     3280  const double *upper = si.getColUpper() ;
     3281  const double *rowLower = si.getRowLower() ;
     3282  const double *rowUpper = si.getRowUpper() ;
     3283/*
     3284  Get down to the business of scanning rows.
     3285*/
     3286  int iRow ;
     3287  for (iRow = 0 ; iRow < numberRows ; iRow++) {
     3288    int numberP1 = 0 ;
     3289    int numberM1 = 0 ;
     3290    int j ;
     3291    double upperValue = rowUpper[iRow] ;
     3292    double lowerValue = rowLower[iRow] ;
     3293    bool good = true ;
     3294    int slack = -1 ;
     3295/*
     3296  Scan the columns of the current row. In order to generate cliques, the row
     3297  must contain only binary variables with coefficients of +/- 1.0. Adjust the
     3298  row bounds for fixed variables while we're here. If the row survives, the
     3299  small end of which contains the indices of the positive coefficients, the
     3300  large end the indices of the negative coefficients.
     3301
     3302  TODO: We set slack when we run across a column of length 1, but don't use it
     3303        for anything.   -- lh, 100924 --
     3304*/
     3305    for (j = rowStart[iRow] ; j < rowStart[iRow]+rowLength[iRow] ; j++) {
     3306      int iColumn = column[j] ;
     3307      int iInteger = lookup[iColumn] ;
     3308      if (upper[iColumn]-lower[iColumn] < 1.0e-8) {
     3309        upperValue -= lower[iColumn]*elementByRow[j] ;
     3310        lowerValue -= lower[iColumn]*elementByRow[j] ;
     3311        continue ;
     3312      } else if (upper[iColumn] != 1.0 || lower[iColumn] != 0.0) {
     3313        good = false ;
     3314        break ;
     3315      } else if (iInteger < 0) {
     3316        good = false ;
     3317        break ;
    57843318      } else {
    5785         if (columnLength[iColumn]==1)
    5786           slack = iInteger;
    5787       }
    5788       if (fabs(elementByRow[j])!=1.0) {
    5789         good=false;
    5790         break;
    5791       } else if (elementByRow[j]>0.0) {
    5792         which[numberP1++]=iColumn;
     3319        if (columnLength[iColumn] == 1)
     3320          slack = iInteger ;
     3321      }
     3322      if (fabs(elementByRow[j]) != 1.0) {
     3323        good = false ;
     3324        break ;
     3325      } else if (elementByRow[j] > 0.0) {
     3326        which[numberP1++] = iColumn ;
    57933327      } else {
    5794         numberM1++;
    5795         which[numberIntegers-numberM1]=iColumn;
    5796       }
    5797     }
    5798     int iUpper = static_cast<int> (floor(upperValue+1.0e-5));
    5799     int iLower = static_cast<int> (ceil(lowerValue-1.0e-5));
    5800     int state=0;
    5801     if (upperValue<1.0e6) {
    5802       if (iUpper==1-numberM1)
    5803         state=1;
    5804       else if (iUpper==-numberM1)
    5805         state=2;
    5806       else if (iUpper<-numberM1)
    5807         state=3;
    5808     }
    5809     if (!state&&lowerValue>-1.0e6) {
    5810       if (-iLower==1-numberP1)
    5811         state=-1;
    5812       else if (-iLower==-numberP1)
    5813         state=-2;
    5814       else if (-iLower<-numberP1)
    5815         state=-3;
    5816     }
    5817     if (good&&state) {
    5818       if (abs(state)==3) {
    5819         // infeasible
    5820         numberCliques_ = -99999;
    5821         break;
    5822       } else if (abs(state)==2) {
    5823         // we can fix all
    5824         numberFixed += numberP1+numberM1;
    5825         if (state>0) {
    5826           // fix all +1 at 0, -1 at 1
    5827           for (i=0;i<numberP1;i++)
    5828             si.setColUpper(which[i],0.0);
    5829           for (i=0;i<numberM1;i++)
    5830             si.setColLower(which[numberIntegers-i-1],
    5831                                  1.0);
     3328        numberM1++ ;
     3329        which[numberIntegers-numberM1] = iColumn ;
     3330      }
     3331    }
     3332/*
     3333  See if we have any strength. As explained at the head of the routine,
     3334  3: infeasible; 2: all variables fixed ; 1: clique of some sort.
     3335*/
     3336    int iUpper = static_cast<int>(floor(upperValue+1.0e-5)) ;
     3337    int iLower = static_cast<int>(ceil(lowerValue-1.0e-5)) ;
     3338    int state = 0 ;
     3339    if (upperValue < 1.0e6) {
     3340      if (iUpper == 1-numberM1)
     3341        state = 1 ;
     3342      else if (iUpper == -numberM1)
     3343        state = 2 ;
     3344      else if (iUpper < -numberM1)
     3345        state = 3 ;
     3346    }
     3347    if (!state && lowerValue > -1.0e6) {
     3348      if (-iLower == 1-numberP1)
     3349        state = -1 ;
     3350      else if (-iLower == -numberP1)
     3351        state = -2 ;
     3352      else if (-iLower < -numberP1)
     3353        state = -3 ;
     3354    }
     3355/*
     3356  +/- 3 is infeasible; we're done.
     3357*/
     3358    if (good && state) {
     3359      if (abs(state) == 3) {
     3360        numberCliques_ = -99999 ;
     3361        break ;
     3362      }
     3363/*
     3364  +2 means all x<j>, j in M, must be 1, all x<j>, j in P, must be 0.
     3365  -2 is the opposite.
     3366*/
     3367      else if (abs(state) == 2) {
     3368        numberFixed += numberP1+numberM1 ;
     3369        if (state > 0) {
     3370          for (i = 0 ; i < numberP1 ; i++)
     3371            si.setColUpper(which[i],0.0) ;
     3372          for (i = 0 ; i < numberM1 ; i++)
     3373            si.setColLower(which[numberIntegers-i-1],1.0) ;
    58323374        } else {
    5833           // fix all +1 at 1, -1 at 0
    5834           for (i=0;i<numberP1;i++)
    5835             si.setColLower(which[i],1.0);
    5836           for (i=0;i<numberM1;i++)
    5837             si.setColUpper(which[numberIntegers-i-1],
    5838                                  0.0);
     3375          for (i = 0 ; i < numberP1 ; i++)
     3376            si.setColLower(which[i],1.0) ;
     3377          for (i = 0 ; i < numberM1 ; i++)
     3378            si.setColUpper(which[numberIntegers-i-1],0.0) ;
    58393379        }
    5840       } else {
    5841         int length = numberP1+numberM1;
    5842         totalP1 += numberP1;
    5843         totalM1 += numberM1;
    5844         if (length >= minimumSize&&length<maximumSize) {
    5845           whichRow[numberCliques_++]=iRow;
    5846           numberEntries += length;
     3380      }
     3381/*
     3382  +1 means that any x<j>, j in P, at 1 implies all x<j>, j in M, at 1.
     3383  -1 means that any x<j>, j in M, at 1 implies all x<j>, j in P, at 1.
     3384*/
     3385      else {
     3386        int length = numberP1+numberM1 ;
     3387        totalP1 += numberP1 ;
     3388        totalM1 += numberM1 ;
     3389        if (length >= minimumSize && length < maximumSize) {
     3390          whichRow[numberCliques_++] = iRow ;
     3391          numberEntries += length ;
    58473392        } else if (numberP1+numberM1 >= maximumSize) {
    5848           // too big
    5849           numberBig++;
    5850           totalBig += numberP1+numberM1;
     3393          numberBig++ ;
     3394          totalBig += numberP1+numberM1 ;
    58513395        }
    58523396      }
    58533397    }
    58543398  }
    5855   if (numberCliques_<0) {
     3399/*
     3400  Done with the row scan. Print a bit of information.
     3401*/
     3402  if (numberCliques_ < 0) {
    58563403    if (logLevel_)
    5857       printf("*** Problem infeasible\n");
     3404      printf("*** Problem infeasible\n") ;
    58583405  } else {
    58593406    if (numberCliques_) {
     
    58633410               (static_cast<double>(totalP1+totalM1))/
    58643411               (static_cast<double> (numberCliques_)),
    5865                totalP1,totalM1);
     3412               totalP1,totalM1) ;
    58663413    } else {
    5867       if (logLevel_>1)
    5868         printf("No cliques found\n");
     3414      if (logLevel_ > 1)
     3415        printf("No cliques found\n") ;
    58693416    }
    58703417    if (numberBig) {
    58713418      if (logLevel_)
    58723419        printf("%d large cliques ( >= %d) found, total %d\n",
    5873              numberBig,maximumSize,totalBig);
     3420             numberBig,maximumSize,totalBig) ;
    58743421    }
    58753422    if (numberFixed) {
    58763423      if (logLevel_)
    5877         printf("%d variables fixed\n",numberFixed);
    5878     }
    5879   }
    5880   if (numberCliques_>0) {
    5881     cliqueType_ = new cliqueType [numberCliques_];
    5882     cliqueStart_ = new int [numberCliques_+1];
    5883     cliqueEntry_ = new cliqueEntry [numberEntries];
    5884     oneFixStart_ = new int [numberColumns_];
    5885     zeroFixStart_ = new int [numberColumns_];
    5886     endFixStart_ = new int [numberColumns_];
    5887     whichClique_ = new int [numberEntries];
    5888     numberEntries=0;
    5889     cliqueStart_[0]=0;
    5890     for (i=0;i<numberColumns_;i++) {
    5891       oneFixStart_[i]=-1;
    5892       zeroFixStart_[i]=-1;
    5893       endFixStart_[i]=-1;
    5894     }
    5895     int iClique;
    5896     // Possible some have been fixed
    5897     int numberCliques=numberCliques_;
    5898     numberCliques_=0;
    5899     for (iClique=0;iClique<numberCliques;iClique++) {
    5900       int iRow=whichRow[iClique];
    5901       whichRow[numberCliques_]=iRow;
    5902       int numberP1=0, numberM1=0;
    5903       int j;
    5904       double upperValue=rowUpper[iRow];
    5905       double lowerValue=rowLower[iRow];
    5906       for (j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    5907         int iColumn = column[j];
    5908         if (upper[iColumn]-lower[iColumn]<1.0e-8) {
    5909           // fixed
    5910           upperValue -= lower[iColumn]*elementByRow[j];
    5911           lowerValue -= lower[iColumn]*elementByRow[j];
    5912           continue;
     3424        printf("%d variables fixed\n",numberFixed) ;
     3425    }
     3426  }
     3427/*
     3428  Deal with any cliques we found. The only information saved from the row scan
     3429  is the row indices of the cliques. Start by setting up the bookkeeping arrays.
     3430
     3431    cliqueType_  is a bit array; if set, the clique row is an equality
     3432    cliqueStart_ is an index array; each entry points to the start of a clique
     3433                 in cliqueEntry_
     3434    cliqueEntry_ is an array with one entry per clique member; each entry
     3435                 gives the column number of the variable and a single bit
     3436                 indicating whether 1 or 0 is the strong direction.
     3437    whichClique_ crossreferences a variable back to its clique; for each
     3438                 variable there's a block of clique indices where the variable
     3439                 is strong-1, followed by a block where the variable is
     3440                 strong-0.
     3441    oneFixStart_ points to the start of the strong-1 block for the variable
     3442    zeroFixStart_ points to the start of the strong-0 block for the variable
     3443    endFixStart_ is one past the end of the entire block for the variable
     3444*/
     3445  if (numberCliques_ > 0) {
     3446    cliqueType_ = new cliqueType [numberCliques_] ;
     3447    cliqueStart_ = new int [numberCliques_+1] ;
     3448    cliqueEntry_ = new cliqueEntry [numberEntries] ;
     3449    oneFixStart_ = new int [numberColumns_] ;
     3450    zeroFixStart_ = new int [numberColumns_] ;
     3451    endFixStart_ = new int [numberColumns_] ;
     3452    whichClique_ = new int [numberEntries] ;
     3453    numberEntries = 0 ;
     3454    cliqueStart_[0] = 0 ;
     3455    for (i = 0 ; i < numberColumns_ ; i++) {
     3456      oneFixStart_[i] = -1 ;
     3457      zeroFixStart_[i] = -1 ;
     3458      endFixStart_[i] = -1 ;
     3459    }
     3460/*
     3461  Process each clique into cliqueEntry storage format. Keep in mind that some
     3462  cliques may now be vacant because their variables have been fixed, so the
     3463  first thing we need to do is evaluate each row again to see if the clique
     3464  conditions still hold, and to recover the index sets P and M. At least we
     3465  know the row has the correct form.
     3466
     3467  TODO: Is there any possibility that a type 1 row has converted to type 2? Is
     3468        it worth checking here? -- lh, 101007 --
     3469*/
     3470    int iClique ;
     3471    int numberCliques = numberCliques_ ;
     3472    numberCliques_ = 0 ;
     3473    for (iClique = 0 ; iClique < numberCliques ; iClique++) {
     3474      int iRow = whichRow[iClique] ;
     3475      whichRow[numberCliques_] = iRow ;
     3476      int numberP1 = 0 ;
     3477      int numberM1 = 0 ;
     3478      int j ;
     3479      double upperValue = rowUpper[iRow] ;
     3480      double lowerValue = rowLower[iRow] ;
     3481      for (j = rowStart[iRow] ; j < rowStart[iRow]+rowLength[iRow] ; j++) {
     3482        int iColumn = column[j] ;
     3483        if (upper[iColumn]-lower[iColumn] < 1.0e-8) {
     3484          upperValue -= lower[iColumn]*elementByRow[j] ;
     3485          lowerValue -= lower[iColumn]*elementByRow[j] ;
     3486          continue ;
    59133487        }
    5914         if (elementByRow[j]>0.0) {
    5915           which[numberP1++]=iColumn;
     3488        if (elementByRow[j] > 0.0) {
     3489          which[numberP1++] = iColumn ;
    59163490        } else {
    5917           numberM1++;
    5918           which[numberIntegers-numberM1]=iColumn;
     3491          numberM1++ ;
     3492          which[numberIntegers-numberM1] = iColumn ;
    59193493        }
    59203494      }
    5921       int iUpper = static_cast<int> (floor(upperValue+1.0e-5));
    5922       int iLower = static_cast<int> (ceil(lowerValue-1.0e-5));
    5923       int state=0;
    5924       if (upperValue<1.0e6) {
    5925         if (iUpper==1-numberM1)
    5926           state=1;
    5927       }
    5928       if (!state&&lowerValue>-1.0e6) {
    5929         state=-1;
    5930       }
    5931       if (abs(state)!=1)
    5932         continue; // must have been fixed
    5933       if (iLower==iUpper) {
    5934         cliqueType_[numberCliques_].equality=1;
     3495      int iUpper = static_cast<int>(floor(upperValue+1.0e-5)) ;
     3496      int iLower = static_cast<int>(ceil(lowerValue-1.0e-5)) ;
     3497      int state = 0 ;
     3498      if (upperValue < 1.0e6) {
     3499        if (iUpper == 1-numberM1)
     3500          state = 1 ;
     3501      }
     3502      if (!state && lowerValue > -1.0e6) {
     3503        state = -1 ;
     3504      }
     3505/*
     3506  If the clique conditions no longer hold, move on to the next candidate.
     3507*/
     3508      if (abs(state) != 1)
     3509        continue ;
     3510/*
     3511  Mark whether the clique row is an equality in cliqueType_. Record the
     3512  column index of each clique member and the strong branching direction in
     3513  cliqueEntry_.
     3514
     3515  TODO: Unless there's some indication that other information will go into
     3516        cliqueType, converting it to a simple boolean array, or even a struct
     3517        with a simple boolean array, would be better than this single-bit
     3518        setup.  -- lh, 101007 --
     3519*/
     3520      if (iLower == iUpper) {
     3521        cliqueType_[numberCliques_].equality = 1 ;
    59353522      } else {
    5936         cliqueType_[numberCliques_].equality=0;
    5937       }
    5938       if (state>0) {
    5939         for (i=0;i<numberP1;i++) {
    5940           // 1 is strong branch
    5941           int iColumn = which[i];
    5942           setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn);
    5943           setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],true);
    5944           numberEntries++;
    5945           // zero counts
    5946           oneFixStart_[iColumn]=0;
    5947           zeroFixStart_[iColumn]=0;
     3523        cliqueType_[numberCliques_].equality = 0 ;
     3524      }
     3525/*
     3526  U<i> = -|M|+1. 1 is the strong value for x<j> in P, 0 is the strong value for
     3527  x<j> in M.
     3528*/
     3529      if (state > 0) {
     3530        for (i = 0 ; i < numberP1 ; i++) {
     3531          int iColumn = which[i] ;
     3532          setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn) ;
     3533          setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],true) ;
     3534          numberEntries++ ;
     3535          oneFixStart_[iColumn] = 0 ;
     3536          zeroFixStart_[iColumn] = 0 ;
    59483537        }
    5949         for (i=0;i<numberM1;i++) {
    5950           // 0 is strong branch
    5951           int iColumn = which[numberIntegers-i-1];
    5952           setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn);
    5953           setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],false);
    5954           numberEntries++;
    5955           // zero counts
    5956           oneFixStart_[iColumn]=0;
    5957           zeroFixStart_[iColumn]=0;
     3538        for (i = 0 ; i < numberM1 ; i++) {
     3539          int iColumn = which[numberIntegers-i-1] ;
     3540          setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn) ;
     3541          setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],false) ;
     3542          numberEntries++ ;
     3543          oneFixStart_[iColumn] = 0 ;
     3544          zeroFixStart_[iColumn] = 0 ;
    59583545        }
    5959       } else {
    5960         for (i=0;i<numberP1;i++) {
    5961           // 0 is strong branch
    5962           int iColumn = which[i];
    5963           setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn);
    5964           setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],false);
    5965           numberEntries++;
    5966           // zero counts
    5967           oneFixStart_[iColumn]=0;
    5968           zeroFixStart_[iColumn]=0;
     3546      }
     3547/*
     3548  L<i> = |P|-1. 0 is the strong value for x<j> in P, 1 is the strong value for
     3549  x<j> in M.
     3550*/
     3551      else {
     3552        for (i = 0 ; i < numberP1 ; i++) {
     3553          int iColumn = which[i] ;
     3554          setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn) ;
     3555          setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],false) ;
     3556          numberEntries++ ;
     3557          oneFixStart_[iColumn] = 0 ;
     3558          zeroFixStart_[iColumn] = 0 ;
    59693559        }
    5970         for (i=0;i<numberM1;i++) {
    5971           // 1 is strong branch
    5972           int iColumn = which[numberIntegers-i-1];
    5973           setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn);
    5974           setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],true);
    5975           numberEntries++;
    5976           // zero counts
    5977           oneFixStart_[iColumn]=0;
    5978           zeroFixStart_[iColumn]=0;
     3560        for (i = 0 ; i < numberM1 ; i++) {
     3561          int iColumn = which[numberIntegers-i-1] ;
     3562          setSequenceInCliqueEntry(cliqueEntry_[numberEntries],iColumn) ;
     3563          setOneFixesInCliqueEntry(cliqueEntry_[numberEntries],true) ;
     3564          numberEntries++ ;
     3565          oneFixStart_[iColumn] = 0 ;
     3566          zeroFixStart_[iColumn] = 0 ;
    59793567        }
    59803568      }
    5981       numberCliques_++;
    5982       cliqueStart_[numberCliques_]=numberEntries;
    5983     }
    5984     // Now do column lists
    5985     // First do counts
    5986     for (iClique=0;iClique<numberCliques_;iClique++) {
    5987       for (int j=cliqueStart_[iClique];j<cliqueStart_[iClique+1];j++) {
    5988         int iColumn = sequenceInCliqueEntry(cliqueEntry_[j]);
     3569      numberCliques_++ ;
     3570      cliqueStart_[numberCliques_] = numberEntries ;
     3571    }
     3572/*
     3573  cliqueEntry now contains a block for each clique, P variables followed by
     3574  M variables. In a U<i> clique, the strong-1 entries precede the strong-0
     3575  entries. In a L<i> clique, the strong-0 entries precede the strong-1
     3576  entries. cliqueStart_ points to the start of the block for each clique.
     3577
     3578  Now count the number of times each variable occurs as strong-1 or
     3579  strong-0.  By virtue of previous initialisation, entries for variables not
     3580  in any clique will be -1.
     3581
     3582  TODO: I can't see any point in the nested loops here. Just walking
     3583        cliqueEntry_ would achieve the same result.  -- lh, 101007 --
     3584*/
     3585    for (iClique = 0 ; iClique < numberCliques_ ; iClique++) {
     3586      for (int j = cliqueStart_[iClique] ; j < cliqueStart_[iClique+1] ; j++) {
     3587        int iColumn = sequenceInCliqueEntry(cliqueEntry_[j]) ;
    59893588        if (oneFixesInCliqueEntry(cliqueEntry_[j]))
    5990           oneFixStart_[iColumn]++;
     3589          oneFixStart_[iColumn]++ ;
    59913590        else
    5992           zeroFixStart_[iColumn]++;
    5993       }
    5994     }
    5995     // now get starts and use which and end as counters
    5996     numberEntries=0;
    5997     for (int iColumn=0;iColumn<numberColumns_;iColumn++) {
    5998       if (oneFixStart_[iColumn]>=0) {
    5999         int n1=oneFixStart_[iColumn];
    6000         int n2=zeroFixStart_[iColumn];
    6001         oneFixStart_[iColumn]=numberEntries;
    6002         which[iColumn]=numberEntries;
    6003         numberEntries += n1;
    6004         zeroFixStart_[iColumn]=numberEntries;
    6005         endFixStart_[iColumn]=numberEntries;
    6006         numberEntries += n2;
    6007       }
    6008     }
    6009     // now put in
    6010     for (iClique=0;iClique<numberCliques_;iClique++) {
    6011       for (int j=cliqueStart_[iClique];j<cliqueStart_[iClique+1];j++) {
    6012         int iColumn = sequenceInCliqueEntry(cliqueEntry_[j]);
     3591          zeroFixStart_[iColumn]++ ;
     3592      }
     3593    }
     3594/*
     3595  We want to set up whichClique_ with blocks of entries for each variable,
     3596  pointing back to the cliques that use the variable.  The block of entries
     3597  for a variable will be divided into two subblocks, first the strong-1
     3598  cliques, then the strong-0 cliques.
     3599
     3600  This first loop establishes the start of each strong-1 and strong-0
     3601  subblock in oneFixStart_ and zeroFixStart_. The arrays which and
     3602  endFixStart mirror oneFixStart_ and zeroFixStart_, respectively.
     3603*/
     3604    numberEntries = 0 ;
     3605    for (int iColumn = 0 ; iColumn < numberColumns_ ; iColumn++) {
     3606      if (oneFixStart_[iColumn] >= 0) {
     3607        int n1 = oneFixStart_[iColumn] ;
     3608        int n2 = zeroFixStart_[iColumn] ;
     3609        oneFixStart_[iColumn] = numberEntries ;
     3610        which[iColumn] = numberEntries ;
     3611        numberEntries += n1 ;
     3612        zeroFixStart_[iColumn] = numberEntries ;
     3613        endFixStart_[iColumn] = numberEntries ;
     3614        numberEntries += n2 ;
     3615      }
     3616    }
     3617/*
     3618  Now drop the clique indices into whichClique. which and endFixStart advance
     3619  through each subblock as we process the cliques. When we're done, we can
     3620  discard which (it's equal to zeroFixStart). endFixStart will be one past
     3621  the end of the block for the variable. Note that endFixStart is not quite
     3622  equal to oneFixStart shifted one index, because not all variables are in use
     3623  in cliques --- all of oneFixStart, zeroFixStart, and endFixStart have entries
     3624  of -1 interspersed with meaningful entries.
     3625*/
     3626    for (iClique=0 ;iClique<numberCliques_ ;iClique++) {
     3627      for (int j=cliqueStart_[iClique] ;j<cliqueStart_[iClique+1] ;j++) {
     3628        int iColumn = sequenceInCliqueEntry(cliqueEntry_[j]) ;
    60133629        if (oneFixesInCliqueEntry(cliqueEntry_[j])) {
    6014           int put = which[iColumn];
    6015           which[iColumn]++;
    6016           whichClique_[put]=iClique;
     3630          int put = which[iColumn] ;
     3631          which[iColumn]++ ;
     3632          whichClique_[put]=iClique ;
    60173633        } else {
    6018           int put = endFixStart_[iColumn];
    6019           endFixStart_[iColumn]++;
    6020           whichClique_[put]=iClique;
     3634          int put = endFixStart_[iColumn] ;
     3635          endFixStart_[iColumn]++ ;
     3636          whichClique_[put]=iClique ;
    60213637        }
    60223638      }
    60233639    }
    60243640  }
    6025   delete [] which;
    6026   delete [] whichRow;
    6027   delete [] lookup;
    6028   return numberCliques_;
     3641  delete [] which ;
     3642  delete [] whichRow ;
     3643  delete [] lookup ;
     3644  return numberCliques_ ;
     3645}
     3646
     3647// Delete all clique information
     3648void
     3649CglProbing::deleteCliques()
     3650{
     3651  delete [] cliqueType_ ;
     3652  delete [] cliqueStart_ ;
     3653  delete [] cliqueEntry_ ;
     3654  delete [] oneFixStart_ ;
     3655  delete [] zeroFixStart_ ;
     3656  delete [] endFixStart_ ;
     3657  delete [] whichClique_ ;
     3658  delete [] cliqueRow_ ;
     3659  delete [] cliqueRowStart_ ;
     3660  cliqueType_ = NULL ;
     3661  cliqueStart_ = NULL ;
     3662  cliqueEntry_ = NULL ;
     3663  oneFixStart_ = NULL ;
     3664  zeroFixStart_ = NULL ;
     3665  endFixStart_ = NULL ;
     3666  whichClique_ = NULL ;
     3667  cliqueRow_ = NULL ;
     3668  cliqueRowStart_ = NULL ;
     3669  numberCliques_ = 0 ;
    60293670}
    60303671
     
    60733714
    60743715
    6075 // Sets up clique information for each row
     3716/*
     3717  jjf: Sets up clique information for each row
     3718
     3719  Boils down clique information into a form suitable for easy use when
     3720  calculating row bounds U<i> and L<i>. To motivate the process, consider
     3721  calculating U<i> for a row that is entangled with several cliques. If
     3722  the cliques don't intersect, we can treat each of them as a super-variable
     3723  and calculate
     3724
     3725    U<i> = (upper bound for non-clique variables) +
     3726                SUM{cliques}(upper bound clique contribution)
     3727
     3728  A nontrivial amount of thinking will convince you that all we need to know
     3729  to calculate the upper bound on the contribution of a clique is which
     3730  variables are in a given clique, and the strong direction for the variable.
     3731  See the written doc'n for details. But what if the cliques intersect? We
     3732  don't want to double count. So this method partitions the variables. The
     3733  largest clique goes in intact. Then the clique with the most unclaimed
     3734  variables goes in. And so on, until the number of unclaimed variables is one
     3735  or zero. We're not getting the full strength from cliques later in the
     3736  sequence, but that's conservative.
     3737
     3738  TODO: The current algorithm allows later cliques to poach variables that
     3739        are in the intersection with earlier cliques. It's not clear to me
     3740        that this is a good idea. It weakens the earlier (larger) cliques
     3741        with no guarantee of major benefits in terms of strengthening the
     3742        later (smaller) cliques.  -- lh, 101008 --
     3743*/
     3744
    60763745void
    6077 CglProbing::setupRowCliqueInformation(const OsiSolverInterface & si)
     3746CglProbing::setupRowCliqueInformation (const OsiSolverInterface &si)
    60783747{
    60793748  if (!numberCliques_)
    6080     return;
    6081   CoinPackedMatrix * rowCopy;
     3749    return ;
     3750/*
     3751  Use a snapshot if we have one, otherwise ask the solver for a copy of the
     3752  constraint matrix. It can happen that CglProbing managed to eliminate rows
     3753  when generating the snapshot, but the column count should be the same.
     3754*/
     3755  CoinPackedMatrix *rowCopy ;
    60823756  if (!rowCopy_) {
    6083     // create from current
    6084     numberRows_=si.getNumRows();
    6085     numberColumns_=si.getNumCols();
    6086     rowCopy = new CoinPackedMatrix(*si.getMatrixByRow());
     3757    numberRows_ = si.getNumRows() ;
     3758    numberColumns_ = si.getNumCols() ;
     3759    rowCopy = new CoinPackedMatrix(*si.getMatrixByRow()) ;
    60873760  } else {
    6088     rowCopy = rowCopy_;
    6089     assert(numberRows_<=si.getNumRows());
    6090     assert(numberColumns_==si.getNumCols());
    6091   }
    6092   assert(numberRows_&&numberColumns_);
    6093   cliqueRowStart_ = new int [numberRows_+1];
    6094   cliqueRowStart_[0]=0;
     3761    rowCopy = rowCopy_ ;
     3762    assert(numberRows_ <= si.getNumRows()) ;
     3763    assert(numberColumns_ == si.getNumCols()) ;
     3764  }
     3765  assert(numberRows_ && numberColumns_) ;
     3766  const int *column = rowCopy->getIndices() ;
     3767  const CoinBigIndex *rowStart = rowCopy->getVectorStarts() ;
     3768  const int *rowLength = rowCopy->getVectorLengths() ;
     3769
     3770  cliqueRowStart_ = new int [numberRows_+1] ;
     3771  cliqueRowStart_[0] = 0 ;
     3772
    60953773  // Temporary array while building list
    6096   cliqueEntry ** array = new cliqueEntry * [numberRows_];
    6097   // Which cliques in use
    6098   int * which = new int[numberCliques_];
    6099   int * count = new int[numberCliques_];
    6100   int * back =new int[numberColumns_];
    6101   CoinZeroN(count,numberCliques_);
    6102   CoinFillN(back,numberColumns_,-1);
    6103   const int * column = rowCopy->getIndices();
    6104   const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    6105   const int * rowLength = rowCopy->getVectorLengths();
    6106   const double * lower = si.getColLower();
    6107   const double * upper = si.getColUpper();
    6108   int iRow;
    6109   for (iRow=0;iRow<numberRows_;iRow++) {
    6110     int j;
    6111     int numberFree=0;
    6112     int numberUsed=0;
    6113     for (j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    6114       int iColumn=column[j];
    6115       if (upper[iColumn]>lower[iColumn]) {
    6116         back[iColumn]=j-rowStart[iRow];
    6117         numberFree++;
    6118         for (int k=oneFixStart_[iColumn];k<endFixStart_[iColumn];k++) {
    6119           int iClique = whichClique_[k];
     3774  cliqueEntry **array = new cliqueEntry* [numberRows_] ;
     3775
     3776  // Clique ids
     3777  int *which = new int [numberCliques_] ;
     3778  // Number of times a given clique is encountered
     3779  int *count = new int [numberCliques_] ;
     3780  // Position of column in row (unfixed variables only; -1 otherwise).
     3781  int *back = new int [numberColumns_] ;
     3782  CoinZeroN(count,numberCliques_) ;
     3783  CoinFillN(back,numberColumns_,-1) ;
     3784
     3785  const double *lower = si.getColLower() ;
     3786  const double *upper = si.getColUpper() ;
     3787
     3788/*
     3789  Open the main loop to scan rows. It won't be apparent for a while yet, but
     3790  we're only interested in rows that a) were not used to form cliques, and
     3791  b) contain variables that are members of cliques. We need to collect some
     3792  data before we can answer the question.
     3793*/
     3794  int iRow ;
     3795  for (iRow = 0 ; iRow < numberRows_ ; iRow++) {
     3796    int j ;
     3797/*
     3798  Scan this row. Count the number of free variables and record the position
     3799  of each variable in the row. Then update clique entanglement (the clique
     3800  id, and the number of variables entangled with the clique).
     3801*/
     3802    int numberFree = 0 ;
     3803    int numberUsed = 0 ;
     3804    for (j = rowStart[iRow] ; j < rowStart[iRow]+rowLength[iRow] ; j++) {
     3805      int iColumn = column[j] ;
     3806      if (upper[iColumn] > lower[iColumn]) {
     3807        back[iColumn] = j-rowStart[iRow] ;
     3808        numberFree++ ;
     3809        for (int k = oneFixStart_[iColumn] ; k < endFixStart_[iColumn] ; k++) {
     3810          int iClique = whichClique_[k] ;
    61203811          if (!count[iClique]) {
    6121             which[numberUsed++]=iClique;
     3812            which[numberUsed++] = iClique ;
    61223813          }
    6123           count[iClique]++;
     3814          count[iClique]++ ;
    61243815        }
    61253816      }
    61263817    }
    6127     // find largest cliques
    6128     bool finished=false;
    6129     int numberInThis=0;
    6130     cliqueEntry * entries = NULL;
    6131     array[iRow]=entries;
     3818/*
     3819  jjf: find largest cliques
     3820
     3821  Recast the clique entanglement data into a form that's easily useable
     3822  when calculating row bounds.
     3823*/
     3824    bool finished = false ;
     3825    int numberInThis = 0 ;
     3826    cliqueEntry *entries = NULL ;
     3827    array[iRow] = entries ;
    61323828    while (!finished) {
    6133       int largest=1;
    6134       int whichClique=-1;
    6135       for (int i=0;i<numberUsed;i++) {
    6136         int iClique = which[i];
    6137         if (count[iClique]>largest) {
    6138           largest=count[iClique];
    6139           whichClique=iClique;
     3829/*
     3830  Find the clique entangled with the largest number of variables from this
     3831  row. As the loop progresses, a more accurate description will be `Find the
     3832  clique with the largest number of entangled variables not yet claimed by
     3833  cliques already processed.'
     3834*/
     3835      int largest = 1 ;
     3836      int whichClique = -1 ;
     3837      for (int i = 0 ; i < numberUsed ; i++) {
     3838        int iClique = which[i] ;
     3839        if (count[iClique] > largest) {
     3840          largest = count[iClique] ;
     3841          whichClique = iClique ;
    61403842        }
    61413843      }
    6142       // Add in if >1 (but not if all as that means clique==row)
    6143       if (whichClique>=0&&largest<numberFree) {
     3844/*
     3845  jjf: Add in if >1 (but not if all as that means clique==row)
     3846
     3847  If we formed a clique from this row then largest == numberFree by
     3848  construction and we're done. We'll skip down to the else, set finished =
     3849  true, and move on to the next row. So the rows we're processing here are
     3850  rows that were not used to form cliques but contain 2 or more variables
     3851  that are part of the same clique.
     3852
     3853  TODO: This might be overly restrictive. You could picture a scenario where
     3854        one constraint defines a clique, and another constraint, with the same
     3855        variables but arbitrary coefficients, enforces some other restriction.
     3856        (Admittedly, no example leaps to mind.)  -- lh, 101008 --
     3857
     3858  TODO: It should be the case that if this is the row that formed the clique,
     3859        count[which[0]] == numberFree. Because we have to suck in all the
     3860        variables, and that means that the very first variable needs to be
     3861        part of a clique, and that clique id will be in which[0]. We could make
     3862        this check immediately after generating the entanglement data.
     3863        -- lh, 101008 --
     3864*/
     3865      if (whichClique >= 0 && largest < numberFree) {
    61443866        if (!numberInThis) {
    6145           int length=rowLength[iRow];
    6146           entries = new cliqueEntry [length];
    6147           array[iRow]=entries;
    6148           for (int i=0;i<length;i++) {
    6149             setOneFixesInCliqueEntry(entries[i],false);
    6150             setSequenceInCliqueEntry(entries[i],numberColumns_+1);
     3867          int length = rowLength[iRow] ;
     3868          entries = new cliqueEntry [length] ;
     3869          array[iRow] = entries ;
     3870          for (int i = 0 ; i < length ; i++) {
     3871            setOneFixesInCliqueEntry(entries[i],false) ;
     3872            setSequenceInCliqueEntry(entries[i],numberColumns_+1) ;
    61513873          }
    61523874        }
    6153         // put in (and take out all counts)
    6154         for (j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    6155           int iColumn=column[j];
    6156           if (upper[iColumn]>lower[iColumn]) {
    6157             bool found=false;
    6158             int k;
    6159             for ( k=oneFixStart_[iColumn];k<endFixStart_[iColumn];k++) {
    6160               int iClique = whichClique_[k];
    6161               if (iClique==whichClique) {
    6162                 found=true;
    6163                 break;
     3875     
     3876/*
     3877  jjf: put in (and take out all counts)
     3878
     3879  Scan the row again. For each unfixed variable, look to see if it's a member
     3880  of the clique we're processing. If so, `remove' the variable from the count
     3881  for other cliques where it's a member, and record the variable as a member
     3882  of this clique, along with its strong direction.
     3883
     3884  Presumably the rationale for scanning the block for this variable in
     3885  whichClique_, as opposed to the block for this clique in cliqueEntry_,
     3886  is that on average the number of variables in a clique is larger than the
     3887  number of cliques using a variable, so the block in whichClique_ is smaller.
     3888*/
     3889        for (j = rowStart[iRow] ; j < rowStart[iRow]+rowLength[iRow] ; j++) {
     3890          int iColumn = column[j] ;
     3891          if (upper[iColumn] > lower[iColumn]) {
     3892            bool found = false ;
     3893            int k ;
     3894            for (k = oneFixStart_[iColumn] ; k < endFixStart_[iColumn] ; k++) {
     3895              int iClique = whichClique_[k] ;
     3896              if (iClique == whichClique) {
     3897                found = true ;
     3898                break ;
    61643899              }
    61653900            }
     3901
    61663902            if (found) {
    6167               for ( k=oneFixStart_[iColumn];k<endFixStart_[iColumn];k++) {
    6168                 int iClique = whichClique_[k];
    6169