Handle MalformedBasis exceptions when adding a new equation (#50)

* a more efficient computation of the explicit basis:
1. Compute from left to write instead of from right to left. This
makes multipliciation by eta matrices easier
2. Compute invLP once-and-for-all, when the basis is factorized

* pivoting stats

* bug fix in store/restore

* supprt for a "basis column oracle", which allows the basis
factorization classes to query the tableau for constraint matrix columns

* cleanup in FT factorization

* grab a fresh basis from the oracle instead of condensing th etas

* bug fix in test

* started separating out the gaussian elimination functionality from the
basis factorization classes

* error when adding an equation

* handle the initial assignment of basic variables as part of Tableau::initializeTableau(0

* make RowBoundTightener remember the tableau, instead of passing the
tableau again and again every time

* Tableau's addRow() informs objects registered as "resizeWatchers" that
the tableau dimensions have changed.
Tableau handles the case where the basis factorization fails after a
row is added, by choosing a new set of basic variables using gaussian
elimination of the constraint matrix

* oops

* oops

src/engine/AutoRowBoundTightener.h

@@ -19,9 +19,9 @@
 class AutoRowBoundTightener
 {
 public:
-	AutoRowBoundTightener()
+	AutoRowBoundTightener( const ITableau &tableau )
 	{
-		_rowBoundTightener = T::createRowBoundTightener();
+		_rowBoundTightener = T::createRowBoundTightener( tableau );
 	}
 
 	~AutoRowBoundTightener()

src/engine/Engine.cpp

@@ -25,7 +25,8 @@
 #include "TimeUtils.h"
 
 Engine::Engine()
-    : _smtCore( this )
+    : _rowBoundTightener( *_tableau )
+    , _smtCore( this )
     , _numPlConstraintsDisabledByValidSplits( 0 )
     , _preprocessingEnabled( false )
     , _work( NULL )
@@ -407,7 +408,7 @@ void Engine::performSimplexStep()
     if ( !fakePivot )
     {
         _tableau->computePivotRow();
-        _rowBoundTightener->examinePivotRow( _tableau );
+        _rowBoundTightener->examinePivotRow();
     }
 
     // Perform the actual pivot
@@ -562,37 +563,16 @@ bool Engine::processInputQuery( InputQuery &inputQuery, bool preprocess )
             ++equationIndex;
         }
 
-        // Placeholder: better constraint matrix analysis as part
-        // of the preprocessing phase.
-
-        AutoConstraintMatrixAnalyzer analyzer;
-        analyzer->analyze( _tableau->getA(), _tableau->getM(), _tableau->getN() );
-
-        if ( analyzer->getRank() != _tableau->getM() )
-        {
-            printf( "Warning!! Contraint matrix rank is %u (out of %u)\n",
-                    analyzer->getRank(), _tableau->getM() );
-        }
-
-        List<unsigned> independentColumns = analyzer->getIndependentColumns();
-
-        unsigned assigned = 0;
-        for( unsigned basicVar : independentColumns )
-        {
-            _tableau->markAsBasic( basicVar );
-            _tableau->assignIndexToBasicVariable( basicVar, assigned );
-            assigned++;
-        }
-
         for ( unsigned i = 0; i < n; ++i )
         {
             _tableau->setLowerBound( i, _preprocessedQuery.getLowerBound( i ) );
             _tableau->setUpperBound( i, _preprocessedQuery.getUpperBound( i ) );
         }
 
         _tableau->registerToWatchAllVariables( _rowBoundTightener );
+        _tableau->registerResizeWatcher( _rowBoundTightener );
 
-        _rowBoundTightener->initialize( _tableau );
+        _rowBoundTightener->setDimensions();
 
         _plConstraints = _preprocessedQuery.getPiecewiseLinearConstraints();
         for ( const auto &constraint : _plConstraints )
@@ -601,7 +581,20 @@ bool Engine::processInputQuery( InputQuery &inputQuery, bool preprocess )
             constraint->setStatistics( &_statistics );
         }
 
-        _tableau->initializeTableau();
+        // Placeholder: better constraint matrix analysis as part
+        // of the preprocessing phase.
+
+        AutoConstraintMatrixAnalyzer analyzer;
+        analyzer->analyze( _tableau->getA(), _tableau->getM(), _tableau->getN() );
+
+        if ( analyzer->getRank() != _tableau->getM() )
+        {
+            printf( "Warning!! Contraint matrix rank is %u (out of %u)\n",
+                    analyzer->getRank(), _tableau->getM() );
+        }
+        List<unsigned> independentColumns = analyzer->getIndependentColumns();
+        _tableau->initializeTableau( independentColumns );
+
         _costFunctionManager->initialize();
         _tableau->registerCostFunctionManager( _costFunctionManager );
         _activeEntryStrategy->initialize( _tableau );
@@ -711,7 +704,7 @@ void Engine::restoreState( const EngineState &state )
     _numPlConstraintsDisabledByValidSplits = state._numPlConstraintsDisabledByValidSplits;
 
     // Make sure the data structures are initialized to the correct size
-    _rowBoundTightener->initialize( _tableau );
+    _rowBoundTightener->setDimensions();
     adjustWorkMemorySize();
     _activeEntryStrategy->resizeHook( _tableau );
     _costFunctionManager->initialize();
@@ -755,7 +748,7 @@ void Engine::applySplit( const PiecewiseLinearCaseSplit &split )
 
     adjustWorkMemorySize();
 
-    _rowBoundTightener->initialize( _tableau );
+    _rowBoundTightener->resetBounds();
 
     for ( auto &bound : bounds )
     {
@@ -878,7 +871,7 @@ void Engine::tightenBoundsOnConstraintMatrix()
     if ( _statistics.getNumMainLoopIterations() %
          GlobalConfiguration::BOUND_TIGHTING_ON_CONSTRAINT_MATRIX_FREQUENCY == 0 )
     {
-        _rowBoundTightener->examineConstraintMatrix( _tableau, true );
+        _rowBoundTightener->examineConstraintMatrix( true );
         _statistics.incNumBoundTighteningOnConstraintMatrix();
     }
 
@@ -897,15 +890,15 @@ void Engine::explicitBasisBoundTightening()
     switch ( GlobalConfiguration::EXPLICIT_BASIS_BOUND_TIGHTENING_TYPE )
     {
     case GlobalConfiguration::USE_BASIS_MATRIX:
-        _rowBoundTightener->examineBasisMatrix( _tableau, saturation );
+        _rowBoundTightener->examineBasisMatrix( saturation );
         break;
 
     case GlobalConfiguration::COMPUTE_INVERTED_BASIS_MATRIX:
-        _rowBoundTightener->examineInvertedBasisMatrix( _tableau, saturation );
+        _rowBoundTightener->examineInvertedBasisMatrix( saturation );
         break;
 
     case GlobalConfiguration::USE_IMPLICIT_INVERTED_BASIS_MATRIX:
-        _rowBoundTightener->examineImplicitInvertedBasisMatrix( _tableau, saturation );
+        _rowBoundTightener->examineImplicitInvertedBasisMatrix( saturation );
         break;
     }
 
@@ -926,7 +919,7 @@ void Engine::performPrecisionRestoration( PrecisionRestorer::RestoreBasics resto
     _statistics.addTimeForPrecisionRestoration( TimeUtils::timePassed( start, end ) );
 
     _statistics.incNumPrecisionRestorations();
-    _rowBoundTightener->clear( _tableau );
+    _rowBoundTightener->clear();
 
     // debug
     double after = _degradationChecker.computeDegradation( *_tableau );
@@ -946,7 +939,7 @@ void Engine::performPrecisionRestoration( PrecisionRestorer::RestoreBasics resto
         _statistics.addTimeForPrecisionRestoration( TimeUtils::timePassed( start, end ) );
         _statistics.incNumPrecisionRestorations();
 
-        _rowBoundTightener->clear( _tableau );
+        _rowBoundTightener->clear();
 
         // debug
         double afterSecond = _degradationChecker.computeDegradation( *_tableau );

src/engine/IRowBoundTightener.h

@@ -16,37 +16,41 @@
 #include "ITableau.h"
 #include "Tightening.h"
 
-class IRowBoundTightener : public ITableau::VariableWatcher
+class IRowBoundTightener : public ITableau::VariableWatcher, public ITableau::ResizeWatcher
 {
 public:
     virtual ~IRowBoundTightener() {};
 
     /*
-      Allocate internal work memory according to the tableau size and
-      initialize tightest lower/upper bounds using the talbeau.
+      Allocate internal work memory according to the tableau size.
     */
-    virtual void initialize( const ITableau &tableau ) = 0;
+    virtual void setDimensions() = 0;
+
+    /*
+      Initialize tightest lower/upper bounds using the talbeau.
+    */
+    virtual void resetBounds() = 0;
 
     /*
       Clear all learned bounds, without reallocating memory.
     */
-    virtual void clear( const ITableau &tableau ) = 0;
+    virtual void clear() = 0;
 
     /*
       Derive and enqueue new bounds for all varaibles, using the
       explicit basis matrix B0 that should be available through the
       tableau. Can also do this until saturation, meaning that we
       continue until no new bounds are learned.
      */
-    virtual void examineBasisMatrix( const ITableau &tableau, bool untilSaturation ) = 0;
+    virtual void examineBasisMatrix( bool untilSaturation ) = 0;
 
     /*
       Derive and enqueue new bounds for all varaibles, using the
       inverse of the explicit basis matrix, inv(B0), which should be available
       through the tableau. Can also do this until saturation, meaning that we
       continue until no new bounds are learned.
      */
-    virtual void examineInvertedBasisMatrix( const ITableau &tableau, bool untilSaturation ) = 0;
+    virtual void examineInvertedBasisMatrix( bool untilSaturation ) = 0;
 
     /*
       Derive and enqueue new bounds for all varaibles, implicitly using the
@@ -55,22 +59,22 @@ class IRowBoundTightener : public ITableau::VariableWatcher
       is performed via FTRANs. Can also do this until saturation, meaning that we
       continue until no new bounds are learned.
     */
-    virtual void examineImplicitInvertedBasisMatrix( const ITableau &tableau, bool untilSaturation ) = 0;
+    virtual void examineImplicitInvertedBasisMatrix( bool untilSaturation ) = 0;
 
     /*
       Derive and enqueue new bounds for all varaibles, using the
       original constraint matrix A and right hands side vector b. Can
       also do this until saturation, meaning that we continue until no
       new bounds are learned.
     */
-    virtual void examineConstraintMatrix( const ITableau &tableau, bool untilSaturation ) = 0;
+    virtual void examineConstraintMatrix( bool untilSaturation ) = 0;
 
     /*
       Derive and enqueue new bounds immedaitely following a pivot
       operation in the given tableau. The tightening is performed for
       the entering variable (which is now basic).
     */
-    virtual void examinePivotRow( ITableau &tableau ) = 0;
+    virtual void examinePivotRow() = 0;
 
     /*
       Get the tightenings entailed by the constraint.

src/engine/ITableau.h

@@ -69,10 +69,22 @@ class ITableau
         virtual void notifyUpperBound( unsigned /* variable */, double /* bound */ ) {}
     };
 
+    class ResizeWatcher
+    {
+    public:
+        /*
+          This callback will be invoked when the tableau size changes,
+          typically when new variables are added.
+        */
+        virtual void notifyDimensionChange( unsigned /* m */, unsigned /* n */ ) {}
+    };
+
     virtual void registerToWatchAllVariables( VariableWatcher *watcher ) = 0;
     virtual void registerToWatchVariable( VariableWatcher *watcher, unsigned variable ) = 0;
     virtual void unregisterToWatchVariable( VariableWatcher *watcher, unsigned variable ) = 0;
 
+    virtual void registerResizeWatcher( ResizeWatcher *watcher ) = 0;
+
     virtual void registerCostFunctionManager( ICostFunctionManager *costFunctionManager ) = 0;
 
     virtual ~ITableau() {};
@@ -82,7 +94,7 @@ class ITableau
     virtual void setRightHandSide( const double *b ) = 0;
     virtual void setRightHandSide( unsigned index, double value ) = 0;
     virtual void markAsBasic( unsigned variable ) = 0;
-    virtual void initializeTableau() = 0;
+    virtual void initializeTableau( const List<unsigned> &initialBasicVariables ) = 0;
     virtual double getValue( unsigned variable ) = 0;
     virtual bool allBoundsValid() const = 0;
     virtual double getLowerBound( unsigned variable ) const = 0;