Bound tightenings using INVPROP and Partial Multi-Neuron Relaxation algorithms.

CHANGELOG.md

@@ -5,9 +5,12 @@
   - Dropped support for parsing Tensorflow network format. Newest Marabou version that supports Tensorflow is at commit 190555573e4702.
   - Fixed bug in the parsing of `transpose` nodes in command line C++ parser.
   - Implemented forward-backward abstract interpretation, symbolic bound tightening, interval arithmetic and simulations for all activation functions.
+  - Implemented backward analysis using preimage-approximation algorithm for `Relu`, `LeakyRelu`, `Sign` and `Bilinear` Layers.
   - Added the BaBSR heuristic as a new branching strategy for ReLU Splitting
   - Support Sub of two variables, "Mul" of two constants, Slice, and ConstantOfShape in the python onnx parser
   - Renamed SmtCore module to SearchTreeHandler
+  - Implemented backward analysis using INVPROP algorithm with added support for all activation functions.
+  - Implemented backward analysis using partial multi-neuron relaxation with BBPS-based heuristic for neuron selection.
 
 ## Version 2.0.0
 

src/configuration/GlobalConfiguration.cpp

@@ -68,6 +68,24 @@ const double GlobalConfiguration::COST_FUNCTION_ERROR_THRESHOLD = 0.0000000001;
 
 const unsigned GlobalConfiguration::SIMULATION_RANDOM_SEED = 1;
 
+const unsigned GlobalConfiguration::VOLUME_ESTIMATION_RANDOM_SEED = 1;
+const unsigned GlobalConfiguration::VOLUME_ESTIMATION_ITERATIONS = 25000;
+const unsigned GlobalConfiguration::PREIMAGE_APPROXIMATION_OPTIMIZATION_RANDOM_SEED = 1;
+const unsigned GlobalConfiguration::PREIMAGE_APPROXIMATION_OPTIMIZATION_MAX_ITERATIONS = 25;
+const double GlobalConfiguration::PREIMAGE_APPROXIMATION_OPTIMIZATION_STEP_SIZE = 0.025;
+const double GlobalConfiguration::PREIMAGE_APPROXIMATION_OPTIMIZATION_LEARNING_RATE = 0.25;
+const double GlobalConfiguration::PREIMAGE_APPROXIMATION_OPTIMIZATION_WEIGHT_DECAY = 0;
+
+const unsigned GlobalConfiguration::INVPROP_MAX_ITERATIONS = 10;
+const double GlobalConfiguration::INVPROP_STEP_SIZE = 0.025;
+const double GlobalConfiguration::INVPROP_LEARNING_RATE = 0.005;
+const double GlobalConfiguration::INVPROP_WEIGHT_DECAY = 0.5;
+const double GlobalConfiguration::INVPROP_INITIAL_GAMMA = 0.025;
+
+const unsigned GlobalConfiguration::PMNR_RANDOM_SEED = 1;
+const unsigned GlobalConfiguration::PMNR_SELECTED_NEURONS = 2;
+const unsigned GlobalConfiguration::PMNR_BBPS_BRANCHING_CANDIDATES = 100;
+
 const bool GlobalConfiguration::USE_HARRIS_RATIO_TEST = true;
 
 const double GlobalConfiguration::SYMBOLIC_TIGHTENING_ROUNDING_CONSTANT = 0.00000000001;
@@ -118,6 +136,7 @@ const bool GlobalConfiguration::WRITE_JSON_PROOF = false;
 
 const unsigned GlobalConfiguration::BACKWARD_BOUND_PROPAGATION_DEPTH = 3;
 const unsigned GlobalConfiguration::MAX_ROUNDS_OF_BACKWARD_ANALYSIS = 10;
+const unsigned GlobalConfiguration::MAX_ROUNDS_OF_PMNR_BACKWARD_ANALYSIS = 10;
 
 const bool GlobalConfiguration::ANALYZE_PROOF_DEPENDENCIES = true;
 const bool GlobalConfiguration::MINIMIZE_PROOF_DEPENDENCIES = true;

src/configuration/GlobalConfiguration.h

@@ -151,6 +151,51 @@ class GlobalConfiguration
     // Random seed for generating simulation values.
     static const unsigned SIMULATION_RANDOM_SEED;
 
+    // Random seed for EstimateVolume procedure (PreimageApproximation).
+    static const unsigned VOLUME_ESTIMATION_RANDOM_SEED;
+
+    // Number of iterations for EstimateVolume procedure (PreimageApproximation).
+    static const unsigned VOLUME_ESTIMATION_ITERATIONS;
+
+    // Random seed for PreimageApproximation optimization.
+    static const unsigned PREIMAGE_APPROXIMATION_OPTIMIZATION_RANDOM_SEED;
+
+    // Maximum iterations for PreimageApproximation optimization.
+    static const unsigned PREIMAGE_APPROXIMATION_OPTIMIZATION_MAX_ITERATIONS;
+
+    // Step size for PreimageApproximation optimization.
+    static const double PREIMAGE_APPROXIMATION_OPTIMIZATION_STEP_SIZE;
+
+    // Learning rate for PreimageApproximation optimization.
+    static const double PREIMAGE_APPROXIMATION_OPTIMIZATION_LEARNING_RATE;
+
+    // Weight decay for PreimageApproximation optimization.
+    static const double PREIMAGE_APPROXIMATION_OPTIMIZATION_WEIGHT_DECAY;
+
+    // Maximum iterations for INVPROP optimization.
+    static const unsigned INVPROP_MAX_ITERATIONS;
+
+    // Step size for INVPROP optimization.
+    static const double INVPROP_STEP_SIZE;
+
+    // Learning rate for INVPROP optimization.
+    static const double INVPROP_LEARNING_RATE;
+
+    // Weight decay for INVPROP optimization.
+    static const double INVPROP_WEIGHT_DECAY;
+
+    // Initial gamma values for INVPROP optimization.
+    static const double INVPROP_INITIAL_GAMMA;
+
+    // Random seed for PMNR (with randomized hyperplanes).
+    static const unsigned PMNR_RANDOM_SEED;
+
+    // Number of selected neurons for PMNR (with heuristically selected hyperplanes).
+    static const unsigned PMNR_SELECTED_NEURONS;
+
+    // Number of candidates for PMNR-BBPS branching points.
+    static const unsigned PMNR_BBPS_BRANCHING_CANDIDATES;
+
     // How often should projected steepest edge reset the reference space?
     static const unsigned PSE_ITERATIONS_BEFORE_RESET;
 
@@ -263,6 +308,10 @@ class GlobalConfiguration
      */
     static const unsigned MAX_ROUNDS_OF_BACKWARD_ANALYSIS;
 
+    /* How many rounds of backward analysis to perform for PMNR algorithm?
+     */
+    static const unsigned MAX_ROUNDS_OF_PMNR_BACKWARD_ANALYSIS;
+
     /* Analyze lemma dependencies when producing proofs
      */
     static const bool ANALYZE_PROOF_DEPENDENCIES;
@@ -271,7 +320,6 @@ class GlobalConfiguration
      */
     static const bool MINIMIZE_PROOF_DEPENDENCIES;
 
-
 #ifdef ENABLE_GUROBI
     /*
       The number of threads Gurobi spawns

src/configuration/OptionParser.cpp

@@ -267,7 +267,8 @@ void OptionParser::initialize()
                 &( ( *_stringOptions )[Options::MILP_SOLVER_BOUND_TIGHTENING_TYPE] ) )
                 ->default_value( ( *_stringOptions )[Options::MILP_SOLVER_BOUND_TIGHTENING_TYPE] ),
             "The MILP solver bound tightening type: "
-            "lp/backward-once/backward-converge/lp-inc/milp/milp-inc/iter-prop/none." )
+            "lp/backward-once/backward-converge/backward-preimage-approx/backward-pmnr/lp-inc/milp/"
+            "milp-inc/iter-prop/none." )
 #endif
         ;
 

src/configuration/Options.cpp

@@ -209,6 +209,10 @@ MILPSolverBoundTighteningType Options::getMILPSolverBoundTighteningType() const
             return MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_ONCE;
         if ( strategyString == "backward-converge" )
             return MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_CONVERGE;
+        if ( strategyString == "backward-preimage-approx" )
+            return MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PREIMAGE_APPROX;
+        if ( strategyString == "backward-pmnr" )
+            return MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PMNR;
         else if ( strategyString == "milp" )
             return MILPSolverBoundTighteningType::MILP_ENCODING;
         else if ( strategyString == "milp-inc" )

src/engine/Engine.cpp

@@ -1598,6 +1598,8 @@ void Engine::performMILPSolverBoundedTightening( Query *inputQuery )
         case MILPSolverBoundTighteningType::LP_RELAXATION_INCREMENTAL:
         case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_ONCE:
         case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_CONVERGE:
+        case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PREIMAGE_APPROX:
+        case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PMNR:
             _networkLevelReasoner->lpRelaxationPropagation();
             break;
         case MILPSolverBoundTighteningType::MILP_ENCODING:
@@ -1662,6 +1664,31 @@ void Engine::performAdditionalBackwardAnalysisIfNeeded()
                 printf( "Backward analysis tightened %u bounds\n", tightened );
         }
     }
+
+    if ( _milpSolverBoundTighteningType ==
+         MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PREIMAGE_APPROX )
+    {
+        performMILPSolverBoundedTightening( &( *_preprocessedQuery ) );
+        unsigned tightened = performSymbolicBoundTightening( &( *_preprocessedQuery ) );
+        if ( _verbosity > 0 )
+            printf( "Backward analysis tightened %u bounds\n", tightened );
+    }
+
+    if ( _milpSolverBoundTighteningType == MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PMNR )
+    {
+        unsigned iter = 1;
+        unsigned tightened = performSymbolicBoundTightening( &( *_preprocessedQuery ) );
+        if ( _verbosity > 0 )
+            printf( "Backward analysis tightened %u bounds\n", tightened );
+        while ( tightened && iter < GlobalConfiguration::MAX_ROUNDS_OF_PMNR_BACKWARD_ANALYSIS )
+        {
+            performMILPSolverBoundedTightening( &( *_preprocessedQuery ) );
+            tightened = performSymbolicBoundTightening( &( *_preprocessedQuery ) );
+            if ( _verbosity > 0 )
+                printf( "Backward analysis tightened %u bounds\n", tightened );
+            ++iter;
+        }
+    }
 }
 
 void Engine::performMILPSolverBoundedTighteningForSingleLayer( unsigned targetIndex )
@@ -1689,6 +1716,8 @@ void Engine::performMILPSolverBoundedTighteningForSingleLayer( unsigned targetIn
             return;
         case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_ONCE:
         case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_CONVERGE:
+        case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PREIMAGE_APPROX:
+        case MILPSolverBoundTighteningType::BACKWARD_ANALYSIS_PMNR:
         case MILPSolverBoundTighteningType::ITERATIVE_PROPAGATION:
         case MILPSolverBoundTighteningType::NONE:
             return;

src/engine/MILPEncoder.cpp

@@ -16,9 +16,9 @@
 
 #include "MILPEncoder.h"
 
-#include "DeepPolySoftmaxElement.h"
 #include "FloatUtils.h"
 #include "GurobiWrapper.h"
+#include "Layer.h"
 #include "TimeUtils.h"
 
 MILPEncoder::MILPEncoder( const ITableau &tableau )
@@ -629,14 +629,14 @@ void MILPEncoder::encodeSoftmaxConstraint( GurobiWrapper &gurobi, SoftmaxConstra
             }
             if ( !useLSE2 )
             {
-                symbolicLowerBias = NLR::DeepPolySoftmaxElement::LSELowerBound(
-                    sourceMids, sourceLbs, sourceUbs, i );
+                symbolicLowerBias =
+                    NLR::Layer::LSELowerBound( sourceMids, sourceLbs, sourceUbs, i );
                 if ( !FloatUtils::wellFormed( symbolicLowerBias ) )
                     wellFormed = false;
                 for ( unsigned j = 0; j < size; ++j )
                 {
-                    double dldj = NLR::DeepPolySoftmaxElement::dLSELowerBound(
-                        sourceMids, sourceLbs, sourceUbs, i, j );
+                    double dldj =
+                        NLR::Layer::dLSELowerBound( sourceMids, sourceLbs, sourceUbs, i, j );
                     if ( !FloatUtils::wellFormed( dldj ) )
                         wellFormed = false;
                     terms.append(
@@ -646,14 +646,14 @@ void MILPEncoder::encodeSoftmaxConstraint( GurobiWrapper &gurobi, SoftmaxConstra
             }
             else
             {
-                symbolicLowerBias = NLR::DeepPolySoftmaxElement::LSELowerBound2(
-                    sourceMids, sourceLbs, sourceUbs, i );
+                symbolicLowerBias =
+                    NLR::Layer::LSELowerBound2( sourceMids, sourceLbs, sourceUbs, i );
                 if ( !FloatUtils::wellFormed( symbolicLowerBias ) )
                     wellFormed = false;
                 for ( unsigned j = 0; j < size; ++j )
                 {
-                    double dldj = NLR::DeepPolySoftmaxElement::dLSELowerBound2(
-                        sourceMids, sourceLbs, sourceUbs, i, j );
+                    double dldj =
+                        NLR::Layer::dLSELowerBound2( sourceMids, sourceLbs, sourceUbs, i, j );
                     if ( !FloatUtils::wellFormed( dldj ) )
                         wellFormed = false;
                     terms.append(
@@ -667,15 +667,14 @@ void MILPEncoder::encodeSoftmaxConstraint( GurobiWrapper &gurobi, SoftmaxConstra
             // Upper-bound
             wellFormed = true;
             double symbolicUpperBias =
-                NLR::DeepPolySoftmaxElement::LSEUpperBound( sourceMids, targetLbs, targetUbs, i );
+                NLR::Layer::LSEUpperBound( sourceMids, targetLbs, targetUbs, i );
             if ( !FloatUtils::wellFormed( symbolicUpperBias ) )
                 wellFormed = false;
             terms.clear();
             terms.append( GurobiWrapper::Term( 1, Stringf( "x%u", targetVariables[i] ) ) );
             for ( unsigned j = 0; j < size; ++j )
             {
-                double dudj = NLR::DeepPolySoftmaxElement::dLSEUpperbound(
-                    sourceMids, targetLbs, targetUbs, i, j );
+                double dudj = NLR::Layer::dLSEUpperbound( sourceMids, targetLbs, targetUbs, i, j );
                 if ( !FloatUtils::wellFormed( dudj ) )
                     wellFormed = false;
                 terms.append( GurobiWrapper::Term( -dudj, Stringf( "x%u", sourceVariables[j] ) ) );

src/engine/MILPSolverBoundTighteningType.h

@@ -34,8 +34,14 @@ enum class MILPSolverBoundTighteningType {
     // Perform backward analysis
     BACKWARD_ANALYSIS_ONCE = 5,
     BACKWARD_ANALYSIS_CONVERGE = 6,
+    // Perform backward analysis using the PreimageApproximation Algorithm (arXiv:2305.03686v4
+    // [cs.SE])
+    BACKWARD_ANALYSIS_PREIMAGE_APPROX = 7,
+    // Perform backward analysis using PMNR with INVPROP and BBPS-based neuron selection
+    // (arXiv:2302.01404v4 [cs.LG], arXiv:2405.21063v3 [cs.LG]).
+    BACKWARD_ANALYSIS_PMNR = 8,
     // Option to have no MILP bound tightening performed
-    NONE = 10,
+    NONE = 9,
 };
 
 #endif // __MILPSolverBoundTighteningType_h__

src/engine/PolygonalTightening.h

@@ -0,0 +1,124 @@
+/*********************                                                        */
+/*! \file PolygonalTightening.h
+ ** \verbatim
+ ** Top contributors (to current version):
+ **   Duligur Ibeling, Guy Katz, Ido Shmuel
+ ** This file is part of the Marabou project.
+ ** Copyright (c) 2017-2024 by the authors listed in the file AUTHORS
+ ** in the top-level source directory) and their institutional affiliations.
+ ** All rights reserved. See the file COPYING in the top-level source
+ ** directory for licensing information.\endverbatim
+ **
+ ** [[ Add lengthier description here ]]
+ **/
+
+#ifndef __PolygonalTightening_h__
+#define __PolygonalTightening_h__
+
+#include "FloatUtils.h"
+#include "MStringf.h"
+#include "Map.h"
+#include "NeuronIndex.h"
+
+#include <cstdio>
+
+class PolygonalTightening
+{
+public:
+    enum PolygonalBoundType {
+        LB = 0,
+        UB = 1,
+    };
+
+    PolygonalTightening( Map<NLR::NeuronIndex, double> neuronToCoefficient,
+                         double value,
+                         PolygonalBoundType type )
+        : _neuronToCoefficient( neuronToCoefficient )
+        , _value( value )
+        , _type( type )
+    {
+    }
+
+    /*
+      The coefficient of each neuron.
+    */
+    Map<NLR::NeuronIndex, double> _neuronToCoefficient;
+
+    /*
+      Its new value.
+    */
+    double _value;
+
+    /*
+      Whether the tightening tightens the
+      lower bound or the upper bound.
+    */
+    PolygonalBoundType _type;
+
+    /*
+      Equality operator.
+    */
+    bool operator==( const PolygonalTightening &other ) const
+    {
+        bool allFound = true;
+        for ( const auto &pair : _neuronToCoefficient )
+        {
+            bool currentFound = false;
+            for ( const auto &otherPair : other._neuronToCoefficient )
+            {
+                currentFound |= ( pair.first._layer == otherPair.first._layer &&
+                                  pair.first._neuron == otherPair.first._neuron &&
+                                  pair.second == otherPair.second );
+            }
+            allFound &= currentFound;
+        }
+        bool result = allFound && _value == other._value && _type == other._type &&
+                      _neuronToCoefficient.size() == other._neuronToCoefficient.size();
+        return result;
+    }
+
+    /*
+       Get matching coefficent for a NeuronIndex,
+       return 0 if the NeuronIndex is not present.
+    */
+    double getCoeff( NLR::NeuronIndex index ) const
+    {
+        if ( _neuronToCoefficient.exists( index ) )
+            return _neuronToCoefficient[index];
+        return 0;
+    }
+
+    void dump() const
+    {
+        String output = "PolygonalTightening: ";
+        unsigned count = 0;
+        for ( const auto &pair : _neuronToCoefficient )
+        {
+            double coeff = pair.second;
+            if ( FloatUtils::isZero( coeff ) )
+                continue;
+
+            if ( count )
+            {
+                output += Stringf( "%s %.2lf neuron%u_%u ",
+                                   FloatUtils::isPositive( coeff ) ? "+" : "-",
+                                   FloatUtils::abs( coeff ),
+                                   pair.first._layer,
+                                   pair.first._neuron );
+            }
+            else
+            {
+                output +=
+                    Stringf( "%.2lf neuron%u_%u ", coeff, pair.first._layer, pair.first._neuron );
+            }
+            ++count;
+        }
+        if ( count == 0 )
+        {
+            output += Stringf( "%.2lf ", 0 );
+        }
+        output += Stringf( "%s %.2lf", _type == LB ? ">=" : "<=", _value );
+        printf( "%s\n", output.ascii() );
+    }
+};
+#endif // __PolygonalTightening_h