Further steps towards CTMC model checking.

Former-commit-id: f057eeb17e

src/modelchecker/AbstractModelChecker.cpp

@@ -32,6 +32,8 @@ namespace storm {
                 return this->computeGloballyProbabilities(pathFormula.asGloballyFormula(), qualitative, optimalityType);
             } else if (pathFormula.isUntilFormula()) {
                 return this->computeUntilProbabilities(pathFormula.asUntilFormula(), qualitative, optimalityType);
+            } else if (pathFormula.isNextFormula()) {
+                return this->computeNextProbabilities(pathFormula.asNextFormula(), qualitative, optimalityType);
             }
             STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "The given formula '" << pathFormula << "' is invalid.");
         }

src/modelchecker/csl/SparseCtmcCslModelChecker.cpp

@@ -6,6 +6,7 @@
 #include "src/utility/vector.h"
 #include "src/utility/graph.h"
 #include "src/utility/solver.h"
+#include "src/utility/numerical.h"
 
 #include "src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h"
 #include "src/modelchecker/results/ExplicitQualitativeCheckResult.h"
@@ -42,7 +43,7 @@ namespace storm {
             ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();;
             ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();
             std::pair<double, double> const& intervalBounds =  pathFormula.getIntervalBounds();
-            std::unique_ptr<CheckResult> result = std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(this->computeBoundedUntilProbabilitiesHelper(leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, intervalBounds.first, intervalBounds.second)));
+            std::unique_ptr<CheckResult> result = std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(this->computeBoundedUntilProbabilitiesHelper(leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), this->getModel().getExitRateVector(), qualitative, intervalBounds.first, intervalBounds.second)));
             
             return result;
         }
@@ -51,7 +52,7 @@ namespace storm {
         std::unique_ptr<CheckResult> SparseCtmcCslModelChecker<ValueType>::computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative, boost::optional<storm::logic::OptimalityType> const& optimalityType) {
             std::unique_ptr<CheckResult> subResultPointer = this->check(pathFormula.getSubformula());
             ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult();
-            std::vector<ValueType> result = SparseDtmcPrctlModelChecker<ValueType>::computeNextProbabilitiesHelper(this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector(), *this->linearEquationSolver);
+            std::vector<ValueType> result = SparseDtmcPrctlModelChecker<ValueType>::computeNextProbabilitiesHelper(this->computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()), subResult.getTruthValuesVector(), *this->linearEquationSolver);
             return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(result)));
         }
         
@@ -61,7 +62,7 @@ namespace storm {
             std::unique_ptr<CheckResult> rightResultPointer = this->check(pathFormula.getRightSubformula());
             ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();
             ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();
-            return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(this->computeUntilProbabilitiesHelper(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *this->linearEquationSolver)));
+            return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(this->computeUntilProbabilitiesHelper(this->computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *this->linearEquationSolver)));
         }
         
         template<class ValueType>
@@ -74,7 +75,7 @@ namespace storm {
             // If the time bounds are [0, inf], we rather call untimed reachability.
             storm::utility::ConstantsComparator<ValueType> comparator;
             if (comparator.isZero(lowerBound) && comparator.isInfinity(upperBound)) {
-                return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(this->computeUntilProbabilitiesHelper(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), phiStates, psiStates, qualitative, *this->linearEquationSolver)));
+                return this->computeUntilProbabilitiesHelper(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), phiStates, psiStates, qualitative, *this->linearEquationSolver);
             }
 
             // From this point on, we know that we have to solve a more complicated problem [t, t'] with either t != 0
@@ -85,7 +86,8 @@ namespace storm {
             
             // If we identify the states that have probability 0 of reaching the target states, we can exclude them from the
             // further computations.
-            storm::storage::BitVector statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0(this->getModel(), this->getModel().getBackwardTransitions(), phiStates, psiStates);
+            storm::storage::SparseMatrix<ValueType> backwardTransitions = this->getModel().getBackwardTransitions();
+            storm::storage::BitVector statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0(backwardTransitions, phiStates, psiStates);
             STORM_LOG_INFO("Found " << statesWithProbabilityGreater0.getNumberOfSetBits() << " states with probability greater 0.");
             storm::storage::BitVector statesWithProbabilityGreater0NonPsi = statesWithProbabilityGreater0 & ~psiStates;
             STORM_LOG_INFO("Found " << statesWithProbabilityGreater0NonPsi.getNumberOfSetBits() << " 'maybe' states.");
@@ -95,26 +97,72 @@ namespace storm {
                     // In this case, the interval is of the form [0, 0].
                     storm::utility::vector::setVectorValues<ValueType>(result, psiStates, storm::utility::one<ValueType>());
                 } else {
-                    // Find the maximal rate of all 'maybe' states to take it as the uniformization rate.
-                    ValueType uniformizationRate = 0;
-                    for (auto const& state : statesWithProbabilityGreater0NonPsi) {
-                        uniformizationRate = std::max(uniformizationRate, exitRates[state]);
-                    }
-                    STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive.");
-                    
                     if (comparator.isZero(lowerBound)) {
                         // In this case, the interval is of the form [0, t].
                         // Note that this excludes [0, inf] since this is untimed reachability and we considered this case earlier.
 
-                        // Compute the uniformized matrix.
-                        storm::storage::SparseMatrix<ValueType> uniformizedMatrix = this->computeUniformizedMatrix(this->getModel().getTransitionMatrix(), statesWithProbabilityGreater0, psiStates, uniformizationRate, exitRates);
+                        // Find the maximal rate of all 'maybe' states to take it as the uniformization rate.
+                        ValueType uniformizationRate = 0;
+                        for (auto const& state : statesWithProbabilityGreater0NonPsi) {
+                            uniformizationRate = std::max(uniformizationRate, exitRates[state]);
+                        }
+                        STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive.");
                         
+                        // Compute the uniformized matrix.
+                        storm::storage::SparseMatrix<ValueType> uniformizedMatrix = this->computeUniformizedMatrix(this->getModel().getTransitionMatrix(), storm::storage::BitVector(this->getModel().getNumberOfStates(), true), psiStates, uniformizationRate, exitRates);
+
+                        // Finally compute the transient probabilities.
+                        std::vector<ValueType> psiStateValues(statesWithProbabilityGreater0.getNumberOfSetBits(), storm::utility::zero<ValueType>());
+                        storm::utility::vector::setVectorValues(psiStateValues, psiStates % statesWithProbabilityGreater0, storm::utility::one<ValueType>());
+                        std::vector<ValueType> subresult = this->computeTransientProbabilities(uniformizedMatrix, uniformizationRate * upperBound, psiStateValues, *this->linearEquationSolver);
+                        storm::utility::vector::setVectorValues(result, statesWithProbabilityGreater0, subresult);
                     } else if (comparator.isInfinity(upperBound)) {
                         // In this case, the interval is of the form [t, inf] with t != 0.
                         
+                        // Start by computing the (unbounded) reachability probabilities of reaching psi states while
+                        // staying in phi states.
+                        result = this->computeUntilProbabilitiesHelper(this->getModel().getTransitionMatrix(), backwardTransitions, phiStates, psiStates, qualitative, *this->linearEquationSolver);
+                        
+                        ValueType uniformizationRate = 0;
+                        for (auto const& state : statesWithProbabilityGreater0) {
+                            uniformizationRate = std::max(uniformizationRate, exitRates[state]);
+                        }
+                        
+                        // Compute the uniformized matrix.
+                        storm::storage::SparseMatrix<ValueType> uniformizedMatrix = this->computeUniformizedMatrix(this->getModel().getTransitionMatrix(), statesWithProbabilityGreater0, storm::storage::BitVector(statesWithProbabilityGreater0.getNumberOfSetBits()), uniformizationRate, exitRates);
+                        
+                        // Finally compute the transient probabilities.
+                        result = this->computeTransientProbabilities(uniformizedMatrix, uniformizationRate * lowerBound, result, *this->linearEquationSolver);
                     } else {
                         // In this case, the interval is of the form [t, t'] with t != 0 and t' != inf.
                         
+                        // Find the maximal rate of all 'maybe' states to take it as the uniformization rate.
+                        ValueType uniformizationRate = 0;
+                        for (auto const& state : statesWithProbabilityGreater0NonPsi) {
+                            uniformizationRate = std::max(uniformizationRate, exitRates[state]);
+                        }
+                        STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive.");
+                        
+                        // Compute the (first) uniformized matrix.
+                        storm::storage::SparseMatrix<ValueType> uniformizedMatrix = this->computeUniformizedMatrix(this->getModel().getTransitionMatrix(), statesWithProbabilityGreater0, psiStates, uniformizationRate, exitRates);
+                        
+                        // Start by computing the transient probabilities of reaching a psi state in time t' - t.
+                        std::vector<ValueType> psiStateValues(statesWithProbabilityGreater0.getNumberOfSetBits(), storm::utility::zero<ValueType>());
+                        storm::utility::vector::setVectorValues(psiStateValues, psiStates % statesWithProbabilityGreater0, storm::utility::one<ValueType>());
+                        std::vector<ValueType> subresult = this->computeTransientProbabilities(uniformizedMatrix, uniformizationRate * (upperBound - lowerBound), psiStateValues, *this->linearEquationSolver);
+                        storm::utility::vector::setVectorValues(result, statesWithProbabilityGreater0, subresult);
+
+                        // Then compute the transient probabilities of being in such a state after t time units. For this,
+                        // we must re-uniformize the CTMC, so we need to compute the second uniformized matrix.
+                        storm::storage::BitVector absorbingStates = ~phiStates;
+                        uniformizationRate = 0;
+                        for (auto const& state : ~absorbingStates) {
+                            uniformizationRate = std::max(uniformizationRate, exitRates[state]);
+                        }
+                        
+                        // Finally, we compute the second set of transient probabilities.
+                        uniformizedMatrix = this->computeUniformizedMatrix(this->getModel().getTransitionMatrix(), storm::storage::BitVector(this->getModel().getNumberOfStates(), true), absorbingStates, uniformizationRate, exitRates);
+                        result = this->computeTransientProbabilities(uniformizedMatrix, uniformizationRate * lowerBound, result, *this->linearEquationSolver);
                     }
                 }
             }
@@ -151,10 +199,64 @@ namespace storm {
             }
         }
         
+        template<class ValueType>
+        std::vector<ValueType> SparseCtmcCslModelChecker<ValueType>::computeTransientProbabilities(storm::storage::SparseMatrix<ValueType> const& uniformizedMatrix, ValueType const& lambda, std::vector<ValueType> values, storm::solver::LinearEquationSolver<ValueType> const& linearEquationSolver) const {
+            // Use Fox-Glynn to get the truncation points and the weights.
+            std::tuple<uint_fast64_t, uint_fast64_t, ValueType, std::vector<ValueType>> foxGlynnResult = storm::utility::numerical::getFoxGlynnCutoff(lambda, 1e-300, 1e+300, storm::settings::generalSettings().getPrecision() / 8.0);
+            
+            // Scale the weights so they add up to one.
+            for (auto& element : std::get<3>(foxGlynnResult)) {
+                element /= std::get<2>(foxGlynnResult);
+            }
+            
+            // Initialize result.
+            std::vector<ValueType> result;
+            if (std::get<0>(foxGlynnResult) == 0) {
+                result = values;
+                storm::utility::vector::scaleVectorInPlace(result, std::get<3>(foxGlynnResult)[0]);
+            } else {
+                result = std::vector<ValueType>(values.size());
+            }
+            std::vector<ValueType> multiplicationResult(result.size());
+            
+            // Perform the matrix-vector multiplications (without adding)
+            if (std::get<0>(foxGlynnResult) > 1) {
+                linearEquationSolver.performMatrixVectorMultiplication(uniformizedMatrix, values, &result, std::get<0>(foxGlynnResult) - 1, &multiplicationResult);
+            }
+            
+            // For the indices that fall in between the truncation points, we need to perform the matrix-vector
+            // multiplication, scale and add the result.
+            ValueType weight = 0;
+            std::function<ValueType(ValueType const&, ValueType const&)> addAndScale = [&weight] (ValueType const& a, ValueType const& b) { return a + weight * b; };
+            for (uint_fast64_t index = std::get<0>(foxGlynnResult); index <= std::get<1>(foxGlynnResult); ++index) {
+                linearEquationSolver.performMatrixVectorMultiplication(uniformizedMatrix, values, nullptr, 1, &multiplicationResult);
+                
+                ValueType weight = std::get<3>(foxGlynnResult)[index - std::get<0>(foxGlynnResult)];
+                storm::utility::vector::applyPointwiseInPlace(result, values, addAndScale);
+            }
+            
+            return result;
+        }
+        
         template<class ValueType>
         std::vector<ValueType> SparseCtmcCslModelChecker<ValueType>::computeUntilProbabilitiesHelper(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::solver::LinearEquationSolver<ValueType> const& linearEquationSolver) {
             return SparseDtmcPrctlModelChecker<ValueType>::computeUntilProbabilitiesHelper(transitionMatrix, backwardTransitions, phiStates, psiStates, qualitative, linearEquationSolver);
         }
         
+        template<class ValueType>
+        storm::storage::SparseMatrix<ValueType> SparseCtmcCslModelChecker<ValueType>::computeProbabilityMatrix(storm::storage::SparseMatrix<ValueType> const& rateMatrix, std::vector<ValueType> const& exitRates) {
+            // Turn the rates into probabilities by scaling each row with the exit rate of the state.
+            storm::storage::SparseMatrix<ValueType> result(rateMatrix);
+            for (uint_fast64_t row = 0; row < result.getRowCount(); ++row) {
+                for (auto& entry : result.getRow(row)) {
+                    entry.setValue(entry.getValue() / exitRates[row]);
+                }
+            }
+            return result;
+        }
+        
+        // Explicitly instantiate the model checker.
+        template class SparseCtmcCslModelChecker<double>;
+        
     } // namespace modelchecker
 } // namespace storm

src/modelchecker/csl/SparseCtmcCslModelChecker.h

@@ -36,8 +36,28 @@ namespace storm {
              * @param absorbingStates The states that need to be made absorbing.
              * @param uniformizationRate The rate to be used for uniformization.
              * @param exitRates The exit rates of all states.
+             * @return The uniformized matrix.
              */
-            storm::storage::SparseMatrix<ValueType> computeUniformizedMatrix(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& absorbingStates, ValueType uniformizationRate, std::vector<ValueType> const& exitRates);
+            static storm::storage::SparseMatrix<ValueType> computeUniformizedMatrix(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& absorbingStates, ValueType uniformizationRate, std::vector<ValueType> const& exitRates);
+            
+            /*!
+             * Computes the transient probabilities for lambda time steps.
+             *
+             * @param uniformizedMatrix The uniformized transition matrix.
+             * @param lambda The number of time steps.
+             * @param values A vector mapping each state to an initial probability.
+             * @param linearEquationSolver The linear equation solver to use.
+             * @return The vector of transient probabilities.
+             */
+            std::vector<ValueType> computeTransientProbabilities(storm::storage::SparseMatrix<ValueType> const& uniformizedMatrix, ValueType const& lambda, std::vector<ValueType> values, storm::solver::LinearEquationSolver<ValueType> const& linearEquationSolver) const;
+            
+            /*!
+             * Converts the given rate-matrix into a time-abstract probability matrix.
+             *
+             * @param rateMatrix The rate matrix.
+             * @param exitRates The exit rate vector.
+             */
+            static storm::storage::SparseMatrix<ValueType> computeProbabilityMatrix(storm::storage::SparseMatrix<ValueType> const& rateMatrix, std::vector<ValueType> const& exitRates);
             
             // An object that is used for solving linear equations and performing matrix-vector multiplication.
             std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> linearEquationSolver;

src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h

@@ -14,6 +14,8 @@ namespace storm {
         template<class ValueType>
         class SparseDtmcPrctlModelChecker : public SparsePropositionalModelChecker<ValueType> {
         public:
+            friend class SparseCtmcCslModelChecker<ValueType>;
+            
             explicit SparseDtmcPrctlModelChecker(storm::models::sparse::Dtmc<ValueType> const& model);
             explicit SparseDtmcPrctlModelChecker(storm::models::sparse::Dtmc<ValueType> const& model, std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>>&& linearEquationSolver);
             

src/modelchecker/propositional/SparsePropositionalModelChecker.cpp

@@ -3,6 +3,7 @@
 #include "src/adapters/CarlAdapter.h"
 
 #include "src/models/sparse/Dtmc.h"
+#include "src/models/sparse/Ctmc.h"
 #include "src/models/sparse/Mdp.h"
 
 #include "src/modelchecker/results/ExplicitQualitativeCheckResult.h"
@@ -50,6 +51,7 @@ namespace storm {
         
         // Explicitly instantiate the template class.
         template storm::models::sparse::Dtmc<double> const& SparsePropositionalModelChecker<double>::getModelAs() const;
+        template storm::models::sparse::Ctmc<double> const& SparsePropositionalModelChecker<double>::getModelAs() const;
         template storm::models::sparse::Mdp<double> const& SparsePropositionalModelChecker<double>::getModelAs() const;
         template class SparsePropositionalModelChecker<double>;
         

src/models/sparse/Ctmc.cpp

@@ -12,7 +12,7 @@ namespace storm {
                                   boost::optional<storm::storage::SparseMatrix<ValueType>> const& optionalTransitionRewardMatrix,
                                   boost::optional<std::vector<boost::container::flat_set<uint_fast64_t>>> const& optionalChoiceLabeling)
             : DeterministicModel<ValueType>(storm::models::ModelType::Ctmc, rateMatrix, stateLabeling, optionalStateRewardVector, optionalTransitionRewardMatrix, optionalChoiceLabeling) {
-                // Intentionally left empty.
+                exitRates = createExitRateVector(this->getTransitionMatrix());
             }
             
             template<typename ValueType>
@@ -21,7 +21,22 @@ namespace storm {
                                   boost::optional<storm::storage::SparseMatrix<ValueType>>&& optionalTransitionRewardMatrix,
                                   boost::optional<std::vector<boost::container::flat_set<uint_fast64_t>>>&& optionalChoiceLabeling)
             : DeterministicModel<ValueType>(storm::models::ModelType::Ctmc, std::move(rateMatrix), std::move(stateLabeling), std::move(optionalStateRewardVector), std::move(optionalTransitionRewardMatrix), std::move(optionalChoiceLabeling)) {
-                // Intentionally left empty.
+                // It is important to refer to the transition matrix here, because the given rate matrix has been move elsewhere.
+                exitRates = createExitRateVector(this->getTransitionMatrix());
+            }
+            
+            template<typename ValueType>
+            std::vector<ValueType> const& Ctmc<ValueType>::getExitRateVector() const {
+                return exitRates;
+            }
+            
+            template<typename ValueType>
+            std::vector<ValueType> Ctmc<ValueType>::createExitRateVector(storm::storage::SparseMatrix<ValueType> const& rateMatrix) {
+                std::vector<ValueType> exitRates(rateMatrix.getRowCount());
+                for (uint_fast64_t row = 0; row < rateMatrix.getRowCount(); ++row) {
+                    exitRates[row] = rateMatrix.getRowSum(row);
+                }
+                return exitRates;
             }
             
             template class Ctmc<double>;

src/models/sparse/Ctmc.h

@@ -52,7 +52,24 @@ namespace storm {
                 Ctmc(Ctmc<ValueType>&& ctmc) = default;
                 Ctmc& operator=(Ctmc<ValueType>&& ctmc) = default;
 #endif
-
+                /*!
+                 * Retrieves the vector of exit rates of the model.
+                 *
+                 * @return The exit rate vector.
+                 */
+                std::vector<ValueType> const& getExitRateVector() const;
+                
+            private:
+                /*!
+                 * Computes the exit rate vector based on the given rate matrix.
+                 *
+                 * @param rateMatrix The rate matrix.
+                 * @return The exit rate vector.
+                 */
+                static std::vector<ValueType> createExitRateVector(storm::storage::SparseMatrix<ValueType> const& rateMatrix);
+                
+                // A vector containing the exit rates of all states.
+                std::vector<ValueType> exitRates;
             };
             
         } // namespace sparse

src/utility/cli.h

@@ -67,6 +67,7 @@ log4cplus::Logger printer;
 #include "src/modelchecker/reachability/SparseDtmcEliminationModelChecker.h"
 #include "src/modelchecker/prctl/SparseMdpPrctlModelChecker.h"
 #include "src/modelchecker/prctl/TopologicalValueIterationMdpPrctlModelChecker.h"
+#include "src/modelchecker/csl/SparseCtmcCslModelChecker.h"
 
 // Headers for counterexample generation.
 #include "src/counterexamples/MILPMinimalLabelSetGenerator.h"
@@ -441,6 +442,11 @@ namespace storm {
                         storm::modelchecker::SparseMdpPrctlModelChecker<ValueType> modelchecker(*mdp);
                         result = modelchecker.check(*formula.get());
 #endif
+                    } else if (model->getType() == storm::models::ModelType::Ctmc) {
+                        std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> ctmc = sparseModel->template as<storm::models::sparse::Ctmc<ValueType>>();
+
+                        storm::modelchecker::SparseCtmcCslModelChecker<ValueType> modelchecker(*ctmc);
+                        result = modelchecker.check(*formula.get());
                     }
                     
                     if (result) {

src/utility/constants.cpp

@@ -110,6 +110,10 @@ namespace storm {
             return std::abs(value) <= precision;
         }
         
+        bool ConstantsComparator<double>::isInfinity(double const& value) const {
+            return value == infinity<double>();
+        }
+        
         bool ConstantsComparator<double>::isEqual(double const& value1, double const& value2) const {
             return std::abs(value1 - value2) <= precision;
         }

src/utility/constants.h

@@ -93,6 +93,8 @@ namespace storm {
             
             bool isZero(double const& value) const;
             
+            bool isInfinity(double const& value) const;
+            
             bool isEqual(double const& value1, double const& value2) const;
             
             bool isConstant(double const& value) const;