Merge branch 'master' into dft

resources/3rdparty/gmm-5.2/include/gmm/gmm_blas.h

@@ -44,7 +44,7 @@
 
 #ifdef STORM_HAVE_INTELTBB
 #include <new> // This fixes a potential dependency ordering problem between GMM and TBB
-#include "tbb/tbb.h"
+#include "storm/adapters/IntelTbbAdapter.h"
 #include <iterator>
 #endif
 

src/storm-cli-utilities/cli.cpp

@@ -17,9 +17,8 @@
 
 
 // Includes for the linked libraries and versions header.
-#ifdef STORM_HAVE_INTELTBB
-#	include "tbb/tbb_stddef.h"
-#endif
+#include "storm/adapters/IntelTbbAdapter.h"
+
 #ifdef STORM_HAVE_GLPK
 #	include "glpk.h"
 #endif

src/storm/adapters/IntelTbbAdapter.h

@@ -0,0 +1,32 @@
+#pragma once
+
+// To detect whether the usage of TBB is possible, this include is neccessary
+#include "storm-config.h"
+
+
+/* On macOS, TBB includes the header
+ * /Library/Developer/CommandLineTools/SDKs/MacOSX10.14.sdk/usr/include/mach/boolean.h
+ * which #defines macros TRUE and FALSE. Since these also occur as identifiers, we #undef them after including TBB.
+ * We still issue a warning in case these macros have been defined before.
+*/
+
+#ifdef TRUE
+#warning "Undefining previously defined macro 'TRUE'."
+#endif
+
+#ifdef FALSE
+#warning "Undefining previously defined macro 'FALSE'."
+#endif
+
+#ifdef STORM_HAVE_INTELTBB
+#include "tbb/tbb.h"
+#include "tbb/tbb_stddef.h"
+#endif
+
+#ifdef TRUE
+#undef TRUE
+#endif
+
+#ifdef FALSE
+#undef FALSE
+#endif

src/storm/analysis/GraphConditions.h

@@ -9,6 +9,7 @@
 namespace storm {
     namespace analysis {
 
+        
         template <typename ValueType, typename Enable=void>
         struct ConstraintType {
             typedef void* val;
@@ -16,7 +17,7 @@ namespace storm {
 
         template<typename ValueType>
         struct ConstraintType<ValueType, typename std::enable_if<std::is_same<storm::RationalFunction, ValueType>::value>::type> {
-            typedef carl::Formula<typename ValueType::PolyType::PolyType> val;
+            typedef carl::Formula<typename storm::RationalFunction::PolyType::PolyType> val;
         };
 
         /**

src/storm/environment/solver/MinMaxSolverEnvironment.h

@@ -25,8 +25,6 @@ namespace storm {
         void setRelativeTerminationCriterion(bool value);
         storm::solver::MultiplicationStyle const& getMultiplicationStyle() const;
         void setMultiplicationStyle(storm::solver::MultiplicationStyle value);
-        bool isForceBoundsSet() const;
-        void setForceBounds(bool value);
         bool isSymmetricUpdatesSet() const;
         void setSymmetricUpdates(bool value);
         
@@ -37,7 +35,6 @@ namespace storm {
         storm::RationalNumber precision;
         bool considerRelativeTerminationCriterion;
         storm::solver::MultiplicationStyle multiplicationStyle;
-        bool forceBounds;
         bool symmetricUpdates;
     };
 }

src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp

@@ -783,7 +783,7 @@ namespace storm {
                             newRew0AStates.set(ecElimResult.oldToNewStateMapping[oldRew0AState]);
                         }
                         
-                        return computeReachabilityRewardsHelper(env, std::move(goal), ecElimResult.matrix, ecElimResult.matrix.transpose(true),
+                        MDPSparseModelCheckingHelperReturnType<ValueType> result = computeReachabilityRewardsHelper(env, std::move(goal), ecElimResult.matrix, ecElimResult.matrix.transpose(true),
                                                                 [&] (uint_fast64_t rowCount, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates) {
                                                                     std::vector<ValueType> result;
                                                                     std::vector<ValueType> oldChoiceRewards = rewardModel.getTotalRewardVector(transitionMatrix);
@@ -813,6 +813,11 @@ namespace storm {
                                                                     }
                                                                     return newChoicesWithoutReward;
                                                                 });
+                        
+                        std::vector<ValueType> resultInEcQuotient = std::move(result.values);
+                        result.values.resize(ecElimResult.oldToNewStateMapping.size());
+                        storm::utility::vector::selectVectorValues(result.values, ecElimResult.oldToNewStateMapping, resultInEcQuotient);
+                        return result;
                     }
                 }
             }

src/storm/models/sparse/MarkovAutomaton.cpp

@@ -5,6 +5,7 @@
 #include "storm/solver/stateelimination/StateEliminator.h"
 #include "storm/storage/FlexibleSparseMatrix.h"
 #include "storm/utility/constants.h"
+#include "storm/utility/ConstantsComparator.h"
 #include "storm/utility/vector.h"
 #include "storm/utility/macros.h"
 #include "storm/utility/graph.h"
@@ -137,6 +138,7 @@ namespace storm {
                     this->exitRates.reserve(this->getNumberOfStates());
                 }
                 
+                storm::utility::ConstantsComparator<ValueType> comparator;
                 for (uint_fast64_t state = 0; state< this->getNumberOfStates(); ++state) {
                     uint_fast64_t row = this->getTransitionMatrix().getRowGroupIndices()[state];
                     if (this->markovianStates.get(state)) {
@@ -151,13 +153,13 @@ namespace storm {
                         ++row;
                     } else {
                         if (assertRates) {
-                            STORM_LOG_THROW(storm::utility::isZero<ValueType>(this->exitRates[state]), storm::exceptions::InvalidArgumentException, "The specified exit rate for (non-Markovian) choice should be 0.");
+                            STORM_LOG_THROW(comparator.isZero(this->exitRates[state]), storm::exceptions::InvalidArgumentException, "The specified exit rate for (non-Markovian) choice should be 0.");
                         } else {
                             this->exitRates.push_back(storm::utility::zero<ValueType>());
                         }
                     }
                     for (; row < this->getTransitionMatrix().getRowGroupIndices()[state+1]; ++row) {
-                        STORM_LOG_THROW(storm::utility::isOne(this->getTransitionMatrix().getRowSum(row)), storm::exceptions::InvalidArgumentException, "Entries of transition matrix do not sum up to one for (non-Markovian) choice " << row << " of state " << state << " (sum is " << this->getTransitionMatrix().getRowSum(row) << ").");
+                        STORM_LOG_THROW(comparator.isOne(this->getTransitionMatrix().getRowSum(row)), storm::exceptions::InvalidArgumentException, "Entries of transition matrix do not sum up to one for (non-Markovian) choice " << row << " of state " << state << " (sum is " << this->getTransitionMatrix().getRowSum(row) << ").");
                     }
                 }
             }

src/storm/settings/ArgumentValidators.cpp

@@ -129,6 +129,10 @@ namespace storm {
             return createRangeValidatorExcluding<uint64_t>(lowerBound, upperBound);
         }
         
+        std::shared_ptr<ArgumentValidator<uint64_t>> ArgumentValidatorFactory::createUnsignedRangeValidatorIncluding(uint64_t lowerBound, uint64_t upperBound) {
+            return createRangeValidatorIncluding<uint64_t>(lowerBound, upperBound);
+        }
+        
         std::shared_ptr<ArgumentValidator<double>> ArgumentValidatorFactory::createDoubleRangeValidatorExcluding(double lowerBound, double upperBound) {
             return createRangeValidatorExcluding<double>(lowerBound, upperBound);
         }

src/storm/settings/ArgumentValidators.h

@@ -70,6 +70,7 @@ namespace storm {
         public:
             static std::shared_ptr<ArgumentValidator<int64_t>> createIntegerRangeValidatorExcluding(int_fast64_t lowerBound, int_fast64_t upperBound);
             static std::shared_ptr<ArgumentValidator<uint64_t>> createUnsignedRangeValidatorExcluding(uint64_t lowerBound, uint64_t upperBound);
+            static std::shared_ptr<ArgumentValidator<uint64_t>> createUnsignedRangeValidatorIncluding(uint64_t lowerBound, uint64_t upperBound);
             static std::shared_ptr<ArgumentValidator<double>> createDoubleRangeValidatorExcluding(double lowerBound, double upperBound);
             
             static std::shared_ptr<ArgumentValidator<double>> createDoubleRangeValidatorIncluding(double lowerBound, double upperBound);

src/storm/settings/modules/CoreSettings.cpp

@@ -9,6 +9,7 @@
 #include "storm/solver/SolverSelectionOptions.h"
 
 #include "storm/storage/dd/DdType.h"
+#include "storm/adapters/IntelTbbAdapter.h"
 
 #include "storm/utility/macros.h"
 #include "storm/exceptions/IllegalArgumentValueException.h"

src/storm/settings/modules/GurobiSettings.cpp

@@ -15,6 +15,8 @@ namespace storm {
             const std::string GurobiSettings::integerToleranceOption = "inttol";
             const std::string GurobiSettings::threadsOption = "threads";
             const std::string GurobiSettings::outputOption = "output";
+            const std::string GurobiSettings::mipFocusOption = "mipfocus";
+            const std::string GurobiSettings::concurrentMipThreadsOption = "concurrentmip";
             
             GurobiSettings::GurobiSettings() : ModuleSettings(moduleName) {
                 this->addOption(storm::settings::OptionBuilder(moduleName, threadsOption, true, "The number of threads that may be used by Gurobi.").addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("count", "The number of threads.").setDefaultValueUnsignedInteger(1).build()).build());
@@ -22,6 +24,10 @@ namespace storm {
                 this->addOption(storm::settings::OptionBuilder(moduleName, outputOption, true, "If set, the Gurobi output will be printed to the command line.").build());
                 
                 this->addOption(storm::settings::OptionBuilder(moduleName, integerToleranceOption, true, "Sets Gurobi's precision for integer variables.").addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("value", "The precision to achieve.").setDefaultValueDouble(1e-06).addValidatorDouble(ArgumentValidatorFactory::createDoubleRangeValidatorExcluding(0.0, 1.0)).build()).build());
+                
+                this->addOption(storm::settings::OptionBuilder(moduleName, mipFocusOption, true, "The high level solution strategy used to solve MILPs.").addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("value", "The number of the strategy.").setDefaultValueUnsignedInteger(0).addValidatorUnsignedInteger(ArgumentValidatorFactory::createUnsignedRangeValidatorIncluding(0, 3)).build()).build());
+                
+                this->addOption(storm::settings::OptionBuilder(moduleName, concurrentMipThreadsOption, true, "The number of MIP solvers Gurobi spawns in parallel. Shall not be larger then the number of threads.").addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("value", "The number of parallel solvers.").setDefaultValueUnsignedInteger(1).addValidatorUnsignedInteger(ArgumentValidatorFactory::createUnsignedGreaterEqualValidator(1)).build()).build());
             }
             
             bool GurobiSettings::isIntegerToleranceSet() const {
@@ -44,11 +50,19 @@ namespace storm {
                 return this->getOption(outputOption).getHasOptionBeenSet();
             }
             
+            uint_fast64_t GurobiSettings::getMIPFocus() const {
+                return this->getOption(mipFocusOption).getArgumentByName("value").getValueAsUnsignedInteger();
+            }
+
+            uint_fast64_t GurobiSettings::getNumberOfConcurrentMipThreads() const {
+                return this->getOption(concurrentMipThreadsOption).getArgumentByName("value").getValueAsUnsignedInteger();
+            }
+
             bool GurobiSettings::check() const {
                 if (isOutputSet() || isIntegerToleranceSet() || isNumberOfThreadsSet()) {
                     STORM_LOG_WARN_COND(storm::settings::getModule<storm::settings::modules::CoreSettings>().getLpSolver() == storm::solver::LpSolverType::Gurobi, "Gurobi is not selected as the preferred LP solver, so setting options for Gurobi might have no effect.");
                 }
-                
+                STORM_LOG_WARN_COND(getNumberOfConcurrentMipThreads() <= getNumberOfThreads(), "Requested more concurrent MILP solvers then the available threads for Gurobi." );
                 return true;
             }
             

src/storm/settings/modules/GurobiSettings.h

@@ -45,6 +45,21 @@ namespace storm {
                  */
                 uint_fast64_t getNumberOfThreads() const;
                 
+                /*!
+                 * Retrieves the selected high-level solution strategy for MILPs.
+                 *
+                 * @return The high-level solution strategy.
+                 */
+                uint_fast64_t getMIPFocus() const;
+                
+                
+                /*!
+                 * Retrieves the number of MIP solvers, Gurobi spawns in parallel.
+                 *
+                 * @return The number of MIP solvers Gurobi spawns in parallel..
+                 */
+                uint_fast64_t getNumberOfConcurrentMipThreads() const;
+                
                 /*!
                  * Retrieves whether the output option was set.
                  *
@@ -62,6 +77,8 @@ namespace storm {
                 static const std::string integerToleranceOption;
                 static const std::string threadsOption;
                 static const std::string outputOption;
+                static const std::string mipFocusOption;
+                static const std::string concurrentMipThreadsOption;
             };
             
         } // namespace modules

src/storm/solver/GmmxxMultiplier.cpp

@@ -2,6 +2,7 @@
 
 #include "storm/adapters/RationalNumberAdapter.h"
 #include "storm/adapters/RationalFunctionAdapter.h"
+#include "storm/adapters/IntelTbbAdapter.h"
 #include "storm/storage/SparseMatrix.h"
 #include "storm/settings/SettingsManager.h"
 #include "storm/settings/modules/CoreSettings.h"

src/storm/solver/GurobiLpSolver.cpp

@@ -77,6 +77,12 @@ namespace storm {
             error = GRBsetintparam(env, "Threads", storm::settings::getModule<storm::settings::modules::GurobiSettings>().getNumberOfThreads());
             STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi Parameter Threads (" << GRBgeterrormsg(env) << ", error code " << error << ").");
             
+            error = GRBsetintparam(env, "MIPFocus", storm::settings::getModule<storm::settings::modules::GurobiSettings>().getMIPFocus());
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi Parameter MIPFocus (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            
+            error = GRBsetintparam(env, "ConcurrentMIP", storm::settings::getModule<storm::settings::modules::GurobiSettings>().getNumberOfConcurrentMipThreads());
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi Parameter ConcurrentMIP (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            
             // Enable the following line to force Gurobi to be as precise about the binary variables as required by the given precision option.
             error = GRBsetdblparam(env, "IntFeasTol", storm::settings::getModule<storm::settings::modules::GurobiSettings>().getIntegerTolerance());
             STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi Parameter IntFeasTol (" << GRBgeterrormsg(env) << ", error code " << error << ").");
@@ -175,6 +181,7 @@ namespace storm {
         
         template<typename ValueType>
         void GurobiLpSolver<ValueType>::addConstraint(std::string const& name, storm::expressions::Expression const& constraint) {
+            STORM_LOG_TRACE("Adding constraint " << name << " to GurobiLpSolver:" << std::endl << "\t" << constraint);
             STORM_LOG_THROW(constraint.isRelationalExpression(), storm::exceptions::InvalidArgumentException, "Illegal constraint is not a relational expression.");
             STORM_LOG_THROW(constraint.getOperator() != storm::expressions::OperatorType::NotEqual, storm::exceptions::InvalidArgumentException, "Illegal constraint uses inequality operator.");
             
@@ -436,6 +443,132 @@ namespace storm {
             }
         }
         
+        template<typename ValueType>
+        void GurobiLpSolver<ValueType>::setMaximalSolutionCount(uint64_t value) {
+            int error = GRBsetintparam(GRBgetenv(model), "PoolSolutions", value);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi Parameter PoolSolutions (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+        }
+        
+        template<typename ValueType>
+        uint64_t GurobiLpSolver<ValueType>::getSolutionCount() const {
+            STORM_LOG_THROW(this->currentModelHasBeenOptimized, storm::exceptions::InvalidStateException, "Illegal call to GurobiLpSolver<ValueType>::getSolutionCount: model has not been optimized.");
+            int result = -1;
+            int error = GRBgetintattr(model, "SolCount", &result);
+            STORM_LOG_THROW(error == 0 && result >= 0, storm::exceptions::InvalidStateException, "Unable to get solution count or invalid number of solutions.");
+            return result;
+        }
+        
+        template<typename ValueType>
+        ValueType GurobiLpSolver<ValueType>::getContinuousValue(storm::expressions::Variable const& variable, uint64_t const& solutionIndex) const {
+            if (!this->isOptimal()) {
+                STORM_LOG_THROW(!this->isInfeasible(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(!this->isUnbounded(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(false, storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+            }
+            STORM_LOG_ASSERT(solutionIndex < getSolutionCount(), "Invalid solution index.");
+            
+            auto variableIndexPair = this->variableToIndexMap.find(variable);
+            STORM_LOG_THROW(variableIndexPair != this->variableToIndexMap.end(), storm::exceptions::InvalidAccessException, "Accessing value of unknown variable '" << variable.getName() << "'.");
+            
+            double value = 0;
+            int error = GRBsetintparam(GRBgetenv(model), GRB_INT_PAR_SOLUTIONNUMBER, solutionIndex);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi solution index (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            error = GRBgetdblattrelement(model, GRB_DBL_ATTR_Xn, variableIndexPair->second, &value);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to get Gurobi solution (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            
+            return storm::utility::convertNumber<ValueType>(value);
+        }
+        
+        template<typename ValueType>
+        int_fast64_t GurobiLpSolver<ValueType>::getIntegerValue(storm::expressions::Variable const& variable, uint64_t const& solutionIndex) const {
+            if (!this->isOptimal()) {
+                STORM_LOG_THROW(!this->isInfeasible(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(!this->isUnbounded(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(false, storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+            }
+            STORM_LOG_ASSERT(solutionIndex < getSolutionCount(), "Invalid solution index.");
+            
+            auto variableIndexPair = this->variableToIndexMap.find(variable);
+            STORM_LOG_THROW(variableIndexPair != this->variableToIndexMap.end(), storm::exceptions::InvalidAccessException, "Accessing value of unknown variable '" << variable.getName() << "'.");
+            
+            double value = 0;
+            int error = GRBsetintparam(GRBgetenv(model), GRB_INT_PAR_SOLUTIONNUMBER, solutionIndex);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi solution index (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            error = GRBgetdblattrelement(model, GRB_DBL_ATTR_Xn, variableIndexPair->second, &value);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to get Gurobi solution (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            STORM_LOG_THROW(std::abs(static_cast<int>(value) - value) <= storm::settings::getModule<storm::settings::modules::GurobiSettings>().getIntegerTolerance(), storm::exceptions::InvalidStateException, "Illegal value for integer variable in Gurobi solution (" << value << ").");
+            
+            return static_cast<int_fast64_t>(value);
+        }
+        
+        template<typename ValueType>
+        bool GurobiLpSolver<ValueType>::getBinaryValue(storm::expressions::Variable const& variable, uint64_t const& solutionIndex) const {
+            if (!this->isOptimal()) {
+                STORM_LOG_THROW(!this->isInfeasible(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(!this->isUnbounded(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(false, storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+            }
+            STORM_LOG_ASSERT(solutionIndex < getSolutionCount(), "Invalid solution index.");
+
+            auto variableIndexPair = this->variableToIndexMap.find(variable);
+            STORM_LOG_THROW(variableIndexPair != this->variableToIndexMap.end(), storm::exceptions::InvalidAccessException, "Accessing value of unknown variable '" << variable.getName() << "'.");
+            
+            double value = 0;
+            int error = GRBsetintparam(GRBgetenv(model), GRB_INT_PAR_SOLUTIONNUMBER, solutionIndex);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi solution index (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            error = GRBgetdblattrelement(model, GRB_DBL_ATTR_Xn, variableIndexPair->second, &value);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to get Gurobi solution (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+
+            if (value > 0.5) {
+                STORM_LOG_THROW(std::abs(static_cast<int>(value) - 1) <= storm::settings::getModule<storm::settings::modules::GurobiSettings>().getIntegerTolerance(), storm::exceptions::InvalidStateException, "Illegal value for integer variable in Gurobi solution (" << value << ").");
+            } else {
+                STORM_LOG_THROW(value <= storm::settings::getModule<storm::settings::modules::GurobiSettings>().getIntegerTolerance(), storm::exceptions::InvalidStateException, "Illegal value for integer variable in Gurobi solution (" << value << ").");
+            }
+            
+            return static_cast<bool>(value);
+        }
+        
+        template<typename ValueType>
+        ValueType GurobiLpSolver<ValueType>::getObjectiveValue(uint64_t const& solutionIndex) const {
+            if (!this->isOptimal()) {
+                STORM_LOG_THROW(!this->isInfeasible(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(!this->isUnbounded(), storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                STORM_LOG_THROW(false, storm::exceptions::InvalidAccessException, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+            }
+            STORM_LOG_ASSERT(solutionIndex < getSolutionCount(), "Invalid solution index.");
+
+            double value = 0;
+            int error = GRBsetintparam(GRBgetenv(model), GRB_INT_PAR_SOLUTIONNUMBER, solutionIndex);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi solution index (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            error = GRBgetdblattr(model, GRB_DBL_ATTR_POOLOBJVAL, &value);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to get Gurobi solution (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            
+            return storm::utility::convertNumber<ValueType>(value);
+        }
+        
+        template<typename ValueType>
+        void GurobiLpSolver<ValueType>::setMaximalMILPGap(ValueType const& gap, bool relative) {
+            int error = -1;
+            if (relative) {
+                error = GRBsetdblparam(GRBgetenv(model), GRB_DBL_PAR_MIPGAP, storm::utility::convertNumber<double>(gap));
+            } else {
+                error = GRBsetdblparam(GRBgetenv(model), GRB_DBL_PAR_MIPGAPABS, storm::utility::convertNumber<double>(gap));
+            }
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to set Gurobi MILP GAP (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+        }
+        
+        template<typename ValueType>
+        ValueType GurobiLpSolver<ValueType>::getMILPGap(bool relative) const {
+            double relativeGap;
+            int error = GRBgetdblattr(model, GRB_DBL_ATTR_MIPGAP, &relativeGap);
+            STORM_LOG_THROW(error == 0, storm::exceptions::InvalidStateException, "Unable to get Gurobi MILP GAP (" << GRBgeterrormsg(env) << ", error code " << error << ").");
+            auto result = storm::utility::convertNumber<ValueType>(relativeGap);
+            if (relative) {
+                return result;
+            } else {
+                return result * getObjectiveValue();
+            }
+        }
         
 #else
         template<typename ValueType>
@@ -577,6 +710,49 @@ namespace storm {
             throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
         }
 
+        template<typename ValueType>
+        void GurobiLpSolver<ValueType>::setMaximalSolutionCount(uint64_t) {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        template<typename ValueType>
+        uint64_t GurobiLpSolver<ValueType>::getSolutionCount() const {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        template<typename ValueType>
+        ValueType GurobiLpSolver<ValueType>::getContinuousValue(storm::expressions::Variable const&, uint64_t const&) const {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        template<typename ValueType>
+        int_fast64_t GurobiLpSolver<ValueType>::getIntegerValue(storm::expressions::Variable const&, uint64_t const&) const {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        template<typename ValueType>
+        bool GurobiLpSolver<ValueType>::getBinaryValue(storm::expressions::Variable const&, uint64_t const&) const {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        template<typename ValueType>
+        ValueType GurobiLpSolver<ValueType>::getObjectiveValue(uint64_t const&) const {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        
+        template<typename ValueType>
+        void GurobiLpSolver<ValueType>::setMaximalMILPGap(ValueType const&, bool) {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+        template<typename ValueType>
+        ValueType GurobiLpSolver<ValueType>::getMILPGap(bool) const {
+            throw storm::exceptions::NotImplementedException() << "This version of storm was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of storm with Gurobi support.";
+        }
+        
+
+        
 #endif
         template class GurobiLpSolver<double>;
         template class GurobiLpSolver<storm::RationalNumber>;

src/storm/solver/GurobiLpSolver.h

@@ -99,7 +99,6 @@ namespace storm {
             virtual int_fast64_t getIntegerValue(storm::expressions::Variable const& name) const override;
             virtual bool getBinaryValue(storm::expressions::Variable const& name) const override;
             virtual ValueType getObjectiveValue() const override;
-            
             // Methods to print the LP problem to a file.
             virtual void writeModelToFile(std::string const& filename) const override;
 
@@ -107,6 +106,28 @@ namespace storm {
             
             virtual void push() override;
             virtual void pop() override;
+            
+            // Methods to retrieve values of sub-optimal solutions found along the way.
+            void setMaximalSolutionCount(uint64_t value); // How many solutions will be stored (at max)
+            uint64_t getSolutionCount() const; // How many solutions have been found
+            ValueType getContinuousValue(storm::expressions::Variable const& name, uint64_t const& solutionIndex) const;
+            int_fast64_t getIntegerValue(storm::expressions::Variable const& name, uint64_t const& solutionIndex) const;
+            bool getBinaryValue(storm::expressions::Variable const& name, uint64_t const& solutionIndex) const;
+            ValueType getObjectiveValue(uint64_t const& solutionIndex) const;
+            
+            /*!
+             * Specifies the maximum difference between lower- and upper objective bounds that triggers termination.
+             * That means a solution is considered optimal if
+             * upperBound - lowerBound < gap * (relative ? |upperBound| : 1).
+             * Only relevant for programs with integer/boolean variables.
+             */
+            void setMaximalMILPGap(ValueType const& gap, bool relative);
+            
+            /*!
+             * Returns the obtained gap after a call to optimize()
+             */
+            ValueType getMILPGap(bool relative) const;
+            
         private:
             /*!
              * Sets some properties of the Gurobi environment according to parameters given by the options.

src/storm/solver/NativeMultiplier.cpp

@@ -10,6 +10,7 @@
 
 #include "storm/adapters/RationalNumberAdapter.h"
 #include "storm/adapters/RationalFunctionAdapter.h"
+#include "storm/adapters/IntelTbbAdapter.h"
 
 #include "storm/utility/macros.h"
 

src/storm/storage/SparseMatrix.cpp

@@ -1,12 +1,5 @@
 #include <boost/functional/hash.hpp>
 
-// To detect whether the usage of TBB is possible, this include is neccessary
-#include "storm-config.h"
-
-#ifdef STORM_HAVE_INTELTBB
-#include "tbb/tbb.h"
-#endif
-
 #include "storm/storage/sparse/StateType.h"
 #include "storm/storage/SparseMatrix.h"
 #include "storm/adapters/RationalFunctionAdapter.h"

src/storm/storage/SparseMatrix.h

@@ -16,6 +16,7 @@
 #include "storm/utility/macros.h"
 #include "storm/utility/constants.h"
 #include "storm/adapters/RationalFunctionAdapter.h"
+#include "storm/adapters/IntelTbbAdapter.h"
 
 // Forward declaration for adapter classes.
 namespace storm {

src/storm/storage/geometry/NativePolytope.cpp

@@ -285,7 +285,6 @@ namespace storm {
             template <typename ValueType>
             std::shared_ptr<Polytope<ValueType>> NativePolytope<ValueType>::affineTransformation(std::vector<Point> const& matrix, Point const& vector) const {
                 STORM_LOG_THROW(!matrix.empty(), storm::exceptions::InvalidArgumentException, "Invoked affine transformation with a matrix without rows.");
-                
                 StormEigen::Index rows = matrix.size();
                 StormEigen::Index columns = matrix.front().size();
                 EigenMatrix eigenMatrix(rows, columns);
@@ -296,6 +295,9 @@ namespace storm {
 
                 StormEigen::FullPivLU<EigenMatrix> luMatrix( eigenMatrix );
                 STORM_LOG_THROW(luMatrix.isInvertible(), storm::exceptions::NotImplementedException, "Affine Transformation of native polytope only implemented if the transformation matrix is invertable");
+                if (isUniversal()) {
+                    return std::make_shared<NativePolytope<ValueType>>(std::vector<Halfspace<ValueType>>());
+                }
                 EigenMatrix newA = A * luMatrix.inverse();
                 EigenVector newb = b + (newA * eigenVector);
                 return std::make_shared<NativePolytope<ValueType>>(emptyStatus, std::move(newA), std::move(newb));
@@ -398,6 +400,9 @@ namespace storm {
             
             template <typename ValueType>
             std::shared_ptr<Polytope<ValueType>> NativePolytope<ValueType>::clean() {
+                if (isEmpty()) {
+                    return create(boost::none, {});
+                }
                 return create(boost::none, getVertices());
             }
 

src/storm/utility/VectorHelper.cpp

@@ -5,8 +5,7 @@
 
 #include "storm/adapters/RationalNumberAdapter.h"
 #include "storm/adapters/RationalFunctionAdapter.h"
-
-#include "storm-config.h"
+#include "storm/adapters/IntelTbbAdapter.h"
 
 #include "storm/utility/vector.h"
 

src/storm/utility/constants.cpp

@@ -52,7 +52,7 @@ namespace storm {
         
         template<>
         bool isNan(double const& value) {
-            return isnan(value);
+            return std::isnan(value);
         }
         
         bool isAlmostZero(double const& a) {

src/storm/utility/vector.h

@@ -1,16 +1,12 @@
 #ifndef STORM_UTILITY_VECTOR_H_
 #define STORM_UTILITY_VECTOR_H_
 
-#include "storm-config.h"
-#ifdef STORM_HAVE_INTELTBB
-#include "tbb/tbb.h"
-#endif
-
 #include <iostream>
 #include <algorithm>
 #include <functional>
 #include <numeric>
 #include <storm/adapters/RationalFunctionAdapter.h>
+#include <storm/adapters/IntelTbbAdapter.h>
 
 #include <boost/optional.hpp>
 
@@ -1156,6 +1152,16 @@ namespace storm {
                 return std::any_of(v.begin(), v.end(), [](T value){return !storm::utility::isZero(value);});
             }
             
+            template<typename T>
+            bool hasZeroEntry(std::vector<T> const& v){
+                return std::any_of(v.begin(), v.end(), [](T value){return storm::utility::isZero(value);});
+            }
+            
+            template<typename T>
+            bool hasInfinityEntry(std::vector<T> const& v){
+                return std::any_of(v.begin(), v.end(), [](T value){return storm::utility::isInfinity(value);});
+            }
+
             inline std::set<storm::RationalFunctionVariable> getVariables(std::vector<storm::RationalFunction> const& vector) {
                 std::set<storm::RationalFunctionVariable> result;
                 for(auto const& entry : vector) {