More stable ratio test (#78)

* initial work on sparse matrix representation

* store/restore functionality

* addLastRow functionality

* getRow and getColumn

* column-merging functionality

* added an interface class

* introducing also sparse vectors

* added addLastColumn functionality

* another unittest

* get sparse columns/matrices in dense form

* WIP on storing the constraint matrix inside the tableau in sparse form

* more WIP, fixed a few bugs, still have a couple of failing tests

* fixed some test issues

* initialization

* resize _a along with the rest

* sparse lu factors

* store also the transposed versions of F and V

* starting work on the sparse GE

* some work on changing the values within an existing sparse
representation, needed for sparse factorization.
still WIP

* refactoring and new functionality for CSRMatrix: any kind of
insertions and deletions

* support for empty initialization and counting elements

* sparse GE is now working. minor bug fixes elsewhere

* compute Ft and Vt as part of the G-elimination process

* tests

* basis oracles can return also sparse columns, not just dense

* sparse LU factorization class

* switch to using the sparse factorization in the engine/tableau

* bug fix in mergeColumns, and nicer printing

* bug fix

* bug fix

* bug fix: merging columns does not delete the actual column, just
leaves it empty

* configuration changes

* optimization: since the sparse columns of A are needed all the time,
just compute them once-and-for-all

* a more efficient implementation of sparse vectors

* comments and unit tests

* cleanup

* keep _A in dense column-major format, too, instead of repeatedly
invoking toDense() to gets its columns

* bad deletes

* bug fix in test

* bug fixes: relu constraint propagation, and the handling of merged
variables in the preprocessor

* new test

* compute Ft incrementally, use it whenever sparse columns are required

* did the todo

* valgrind fixes

* debugging

* un-initialized memory

* cleanup

* fixing an issue with cost function degradation

* reinstating the anti-looping trick

* debug prints

* verify invariants

* debug info for crash

* fixing a numerical stability issue
new assertions

* cleanup

* cleanup

* removing the code with the failing assertion

* weaken a couple of too-storng assertions

* bug fix and some fine-tuning of PSE

* bug fix in LU factorization
optimizations to PSE

* WIP

* cleanup

* cleanup

* fine-tuning the computation of basic costs

* more robust comparisons

src/configuration/GlobalConfiguration.cpp

@@ -17,8 +17,13 @@
 const double GlobalConfiguration::DEFAULT_EPSILON_FOR_COMPARISONS = 0.0000000001;
 const unsigned GlobalConfiguration::DEFAULT_DOUBLE_TO_STRING_PRECISION = 10;
 const unsigned GlobalConfiguration::STATISTICS_PRINTING_FREQUENCY = 100;
-const double GlobalConfiguration::BOUND_COMPARISON_TOLERANCE = 0.001;
+const double GlobalConfiguration::BOUND_COMPARISON_ADDITIVE_TOLERANCE = 0.0000001;
+const double GlobalConfiguration::BOUND_COMPARISON_MULTIPLICATIVE_TOLERANCE = 0.001 * 0.0000001;
 const double GlobalConfiguration::PIVOT_CHANGE_COLUMN_TOLERANCE = 0.000000001;
+const double GlobalConfiguration::RATIO_CONSTRAINT_ADDITIVE_TOLERANCE = 0.0000001 * 0.3;
+const double GlobalConfiguration::RATIO_CONSTRAINT_MULTIPLICATIVE_TOLERANCE = 0.001 * 0.0000001 * 0.3;
+const double GlobalConfiguration::BASIC_COSTS_ADDITIVE_TOLERANCE = 0.0000001 * 0.97;
+const double GlobalConfiguration::BASIC_COSTS_MULTIPLICATIVE_TOLERANCE = 0.001 * 0.0000001 * 0.97;
 const unsigned GlobalConfiguration::DEGRADATION_CHECKING_FREQUENCY = 100;
 const double GlobalConfiguration::DEGRADATION_THRESHOLD = 0.1;
 const double GlobalConfiguration::ACCEPTABLE_SIMPLEX_PIVOT_THRESHOLD = 0.0001;
@@ -36,6 +41,8 @@ const bool GlobalConfiguration::PREPROCESSOR_PL_CONSTRAINTS_ADD_AUX_EQUATIONS =
 const unsigned GlobalConfiguration::PSE_ITERATIONS_BEFORE_RESET = 1000;
 const double GlobalConfiguration::PSE_GAMMA_ERROR_THRESHOLD = 0.001;
 
+const double GlobalConfiguration::RELU_CONSTRAINT_COMPARISON_TOLERANCE = 0.001;
+
 const GlobalConfiguration::ExplicitBasisBoundTighteningType GlobalConfiguration::EXPLICIT_BASIS_BOUND_TIGHTENING_TYPE =
     GlobalConfiguration::COMPUTE_INVERTED_BASIS_MATRIX;
 const bool GlobalConfiguration::EXPLICIT_BOUND_TIGHTENING_UNTIL_SATURATION = false;
@@ -61,8 +68,13 @@ void GlobalConfiguration::print()
     printf( "  DEFAULT_EPSILON_FOR_COMPARISONS: %.15lf\n", DEFAULT_EPSILON_FOR_COMPARISONS );
     printf( "  DEFAULT_DOUBLE_TO_STRING_PRECISION: %u\n", DEFAULT_DOUBLE_TO_STRING_PRECISION );
     printf( "  STATISTICS_PRINTING_FREQUENCY: %u\n", STATISTICS_PRINTING_FREQUENCY );
-    printf( "  BOUND_COMPARISON_TOLERANCE: %.15lf\n", BOUND_COMPARISON_TOLERANCE );
+    printf( "  BOUND_COMPARISON_ADDITIVE_TOLERANCE: %.15lf\n", BOUND_COMPARISON_ADDITIVE_TOLERANCE );
+    printf( "  BOUND_COMPARISON_MULTIPLICATIVE_TOLERANCE: %.15lf\n", BOUND_COMPARISON_MULTIPLICATIVE_TOLERANCE );
     printf( "  PIVOT_CHANGE_COLUMN_TOLERANCE: %.15lf\n", PIVOT_CHANGE_COLUMN_TOLERANCE );
+    printf( "  RATIO_CONSTRAINT_ADDITIVE_TOLERANCE: %.15lf\n", RATIO_CONSTRAINT_ADDITIVE_TOLERANCE );
+    printf( "  RATIO_CONSTRAINT_MULTIPLICATIVE_TOLERANCE: %.15lf\n", RATIO_CONSTRAINT_MULTIPLICATIVE_TOLERANCE );
+    printf( "  BASIC_COSTS_ADDITIVE_TOLERANCE: %.15lf\n", BASIC_COSTS_ADDITIVE_TOLERANCE );
+    printf( "  BASIC_COSTS_MULTIPLICATIVE_TOLERANCE: %.15lf\n", BASIC_COSTS_MULTIPLICATIVE_TOLERANCE );
     printf( "  DEGRADATION_CHECKING_FREQUENCY: %u\n", DEGRADATION_CHECKING_FREQUENCY );
     printf( "  DEGRADATION_THRESHOLD: %.15lf\n", DEGRADATION_THRESHOLD );
     printf( "  ACCEPTABLE_SIMPLEX_PIVOT_THRESHOLD: %.15lf\n", ACCEPTABLE_SIMPLEX_PIVOT_THRESHOLD );
@@ -78,6 +90,7 @@ void GlobalConfiguration::print()
             PREPROCESSOR_PL_CONSTRAINTS_ADD_AUX_EQUATIONS ? "Yes" : "No" );
     printf( "  PSE_ITERATIONS_BEFORE_RESET: %u\n", PSE_ITERATIONS_BEFORE_RESET );
     printf( "  PSE_GAMMA_ERROR_THRESHOLD: %.15lf\n", PSE_GAMMA_ERROR_THRESHOLD );
+    printf( "  RELU_CONSTRAINT_COMPARISON_TOLERANCE: %.15lf\n", RELU_CONSTRAINT_COMPARISON_TOLERANCE );
 
     String basisBoundTighteningType;
     switch ( EXPLICIT_BASIS_BOUND_TIGHTENING_TYPE )

src/configuration/GlobalConfiguration.h

@@ -28,12 +28,21 @@ class GlobalConfiguration
     static const unsigned STATISTICS_PRINTING_FREQUENCY;
 
     // Tolerance when checking whether the value computed for a basic variable is out of bounds
-    static const double BOUND_COMPARISON_TOLERANCE;
+    static const double BOUND_COMPARISON_ADDITIVE_TOLERANCE;
+    static const double BOUND_COMPARISON_MULTIPLICATIVE_TOLERANCE;
 
     // Tolerance when checking whether a basic variable depends on a non-basic variable, by looking
     // at the change column, as part of a pivot operation.
     static const double PIVOT_CHANGE_COLUMN_TOLERANCE;
 
+    // Ratio test tolerance constants
+    static const double RATIO_CONSTRAINT_ADDITIVE_TOLERANCE;
+    static const double RATIO_CONSTRAINT_MULTIPLICATIVE_TOLERANCE;
+
+    // Cost function tolerance constants
+    static const double BASIC_COSTS_ADDITIVE_TOLERANCE;
+    static const double BASIC_COSTS_MULTIPLICATIVE_TOLERANCE;
+
     // Toggle query-preprocessing on/off.
 	static const bool PREPROCESS_INPUT_QUERY;
 
@@ -81,6 +90,8 @@ class GlobalConfiguration
     // An error threshold which, when crossed, causes projected steepest edge to reset the reference space
     static const double PSE_GAMMA_ERROR_THRESHOLD;
 
+    // The tolerance for checking whether f = Relu( b ), to determine a ReLU's statisfaction
+    static const double RELU_CONSTRAINT_COMPARISON_TOLERANCE;
     /*
       Bound tightening options
     */

src/engine/CostFunctionManager.cpp

@@ -155,14 +155,38 @@ void CostFunctionManager::computeCoreCostFunction()
 
 void CostFunctionManager::computeBasicOOBCosts()
 {
+    unsigned variable;
+    double assignment, lb, relaxedLb, ub, relaxedUb;
     for ( unsigned i = 0; i < _m; ++i )
     {
-        if ( _tableau->basicTooLow( i ) )
+        variable = _tableau->basicIndexToVariable( i );
+        assignment = _tableau->getBasicAssignment( i );
+
+        lb = _tableau->getLowerBound( variable );
+        relaxedLb =
+            lb -
+            ( GlobalConfiguration::BASIC_COSTS_ADDITIVE_TOLERANCE +
+              GlobalConfiguration::BASIC_COSTS_MULTIPLICATIVE_TOLERANCE * FloatUtils::abs( lb ) );
+
+        if ( assignment < relaxedLb )
+        {
             _basicCosts[i] = -1;
-        else if ( _tableau->basicTooHigh( i ) )
+            continue;
+        }
+
+        ub = _tableau->getUpperBound( variable );
+        relaxedUb =
+            ub +
+            ( GlobalConfiguration::BASIC_COSTS_ADDITIVE_TOLERANCE +
+              GlobalConfiguration::BASIC_COSTS_MULTIPLICATIVE_TOLERANCE * FloatUtils::abs( ub ) );
+
+        if ( assignment > relaxedUb )
+        {
             _basicCosts[i] = 1;
-        else
-            _basicCosts[i] = 0;
+            continue;
+        }
+
+        _basicCosts[i] = 0;
     }
 }
 

src/engine/ITableau.h

@@ -31,10 +31,8 @@ class ITableau
 public:
     enum BasicStatus {
         BELOW_LB = 0,
-        AT_LB = 1,
-        BETWEEN = 2,
-        AT_UB = 3,
-        ABOVE_UB = 4,
+        BETWEEN = 1,
+        ABOVE_UB = 2,
     };
 
     enum CostFunctionStatus {
@@ -162,6 +160,7 @@ class ITableau
     virtual void backwardTransformation( const double *y, double *x ) const = 0;
     virtual double getSumOfInfeasibilities() const = 0;
     virtual BasicAssignmentStatus getBasicAssignmentStatus() const = 0;
+    virtual double getBasicAssignment( unsigned basicIndex ) const = 0;
     virtual void setBasicAssignmentStatus( ITableau::BasicAssignmentStatus status ) = 0;
     virtual bool basicOutOfBounds( unsigned basic ) const = 0;
     virtual bool basicTooHigh( unsigned basic ) const = 0;

src/engine/ReluConstraint.cpp

@@ -54,7 +54,7 @@ void ReluConstraint::unregisterAsWatcher( ITableau *tableau )
 
 void ReluConstraint::notifyVariableValue( unsigned variable, double value )
 {
-    if ( FloatUtils::isZero( value, GlobalConfiguration::BOUND_COMPARISON_TOLERANCE ) )
+    if ( FloatUtils::isZero( value, GlobalConfiguration::RELU_CONSTRAINT_COMPARISON_TOLERANCE ) )
         value = 0.0;
 
     _assignment[variable] = value;
@@ -112,7 +112,7 @@ bool ReluConstraint::satisfied() const
         return false;
 
     if ( FloatUtils::isPositive( fValue ) )
-        return FloatUtils::areEqual( bValue, fValue, GlobalConfiguration::BOUND_COMPARISON_TOLERANCE );
+        return FloatUtils::areEqual( bValue, fValue, GlobalConfiguration::RELU_CONSTRAINT_COMPARISON_TOLERANCE );
     else
         return !FloatUtils::isPositive( bValue );
 }