diff --git a/src/builder/DftExplorationHeuristic.cpp b/src/builder/DftExplorationHeuristic.cpp
index 75556a06f..db26fc60d 100644
--- a/src/builder/DftExplorationHeuristic.cpp
+++ b/src/builder/DftExplorationHeuristic.cpp
@@ -8,16 +8,22 @@ namespace storm {
     namespace builder {
 
         template<typename ValueType>
-        DFTExplorationHeuristic<ValueType>::DFTExplorationHeuristic(size_t id) : id(id), expand(false), depth(0), probability(storm::utility::zero<ValueType>()) {
+        DFTExplorationHeuristic<ValueType>::DFTExplorationHeuristic(size_t id) : id(id), expand(false), lowerBound(storm::utility::zero<ValueType>()), upperBound(storm::utility::infinity<ValueType>()), depth(0), probability(storm::utility::one<ValueType>()) {
             // Intentionally left empty
         }
 
         template<typename ValueType>
-        DFTExplorationHeuristic<ValueType>::DFTExplorationHeuristic(size_t id, DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate) : id(id), expand(false), depth(predecessor.depth + 1) {
+        DFTExplorationHeuristic<ValueType>::DFTExplorationHeuristic(size_t id, DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate, ValueType lowerBound, ValueType upperBound) : id(id), expand(false), lowerBound(lowerBound), upperBound(upperBound), depth(predecessor.depth + 1) {
             STORM_LOG_ASSERT(storm::utility::zero<ValueType>() < exitRate, "Exit rate is 0");
             probability = predecessor.probability * rate/exitRate;
         }
 
+        template<typename ValueType>
+        void DFTExplorationHeuristic<ValueType>::setBounds(ValueType lowerBound, ValueType upperBound) {
+            this->lowerBound = lowerBound;
+            this->upperBound = upperBound;
+        }
+
         template<>
         bool DFTExplorationHeuristicProbability<double>::updateHeuristicValues(DFTExplorationHeuristic<double> const& predecessor, double rate, double exitRate) {
             STORM_LOG_ASSERT(exitRate > 0, "Exit rate is 0");
@@ -41,18 +47,29 @@ namespace storm {
             STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Heuristic rate ration does not work for rational functions.");
         }
 
+        template<>
+        double DFTExplorationHeuristicBoundDifference<double>::getPriority() const {
+            return upperBound / lowerBound;
+        }
+
+        template<>
+        double DFTExplorationHeuristicBoundDifference<storm::RationalFunction>::getPriority() const {
+            STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Heuristic bound difference does not work for rational functions.");
+        }
 
         // Instantiate templates.
         template class DFTExplorationHeuristic<double>;
         template class DFTExplorationHeuristicNone<double>;
         template class DFTExplorationHeuristicDepth<double>;
         template class DFTExplorationHeuristicProbability<double>;
+        template class DFTExplorationHeuristicBoundDifference<double>;
 
 #ifdef STORM_HAVE_CARL
         template class DFTExplorationHeuristic<storm::RationalFunction>;
         template class DFTExplorationHeuristicNone<storm::RationalFunction>;
         template class DFTExplorationHeuristicDepth<storm::RationalFunction>;
         template class DFTExplorationHeuristicProbability<storm::RationalFunction>;
+        template class DFTExplorationHeuristicBoundDifference<storm::RationalFunction>;
 #endif
     }
 }
diff --git a/src/builder/DftExplorationHeuristic.h b/src/builder/DftExplorationHeuristic.h
index f24d7defb..12f178f25 100644
--- a/src/builder/DftExplorationHeuristic.h
+++ b/src/builder/DftExplorationHeuristic.h
@@ -22,7 +22,9 @@ namespace storm {
         public:
             DFTExplorationHeuristic(size_t id);
 
-            DFTExplorationHeuristic(size_t id, DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate);
+            DFTExplorationHeuristic(size_t id, DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate, ValueType lowerBound, ValueType upperBound);
+
+            void setBounds(ValueType lowerBound, ValueType upperBound);
 
             virtual bool updateHeuristicValues(DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate) = 0;
 
@@ -38,6 +40,10 @@ namespace storm {
                 return id;
             }
 
+            bool isExpand() {
+                return expand;
+            }
+
             size_t getDepth() const {
                 return depth;
             }
@@ -46,9 +52,19 @@ namespace storm {
                 return probability;
             }
 
+            ValueType getLowerBound() const {
+                return lowerBound;
+            }
+
+            ValueType getUpperBound() const {
+                return upperBound;
+            }
+
         protected:
             size_t id;
             bool expand;
+            ValueType lowerBound;
+            ValueType upperBound;
             size_t depth;
             ValueType probability;
         };
@@ -60,7 +76,7 @@ namespace storm {
                 // Intentionally left empty
             }
 
-            DFTExplorationHeuristicNone(size_t id, DFTExplorationHeuristicNone<ValueType> const& predecessor, ValueType rate, ValueType exitRate) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) {
+            DFTExplorationHeuristicNone(size_t id, DFTExplorationHeuristicNone<ValueType> const& predecessor, ValueType rate, ValueType exitRate, ValueType lowerBound, ValueType upperBound) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate, lowerBound, upperBound) {
                 // Intentionally left empty
             }
 
@@ -88,7 +104,7 @@ namespace storm {
                 // Intentionally left empty
             }
 
-            DFTExplorationHeuristicDepth(size_t id, DFTExplorationHeuristicDepth<ValueType> const& predecessor, ValueType rate, ValueType exitRate) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) {
+            DFTExplorationHeuristicDepth(size_t id, DFTExplorationHeuristicDepth<ValueType> const& predecessor, ValueType rate, ValueType exitRate, ValueType lowerBound, ValueType upperBound) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate, lowerBound, upperBound) {
                 // Intentionally left empty
             }
 
@@ -121,7 +137,7 @@ namespace storm {
                 // Intentionally left empty
             }
 
-            DFTExplorationHeuristicProbability(size_t id, DFTExplorationHeuristicProbability<ValueType> const& predecessor, ValueType rate, ValueType exitRate) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) {
+            DFTExplorationHeuristicProbability(size_t id, DFTExplorationHeuristicProbability<ValueType> const& predecessor, ValueType rate, ValueType exitRate, ValueType lowerBound, ValueType upperBound) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate, lowerBound, upperBound) {
                 // Intentionally left empty
             }
 
@@ -138,6 +154,37 @@ namespace storm {
             }
         };
 
+        template<typename ValueType>
+        class DFTExplorationHeuristicBoundDifference : public DFTExplorationHeuristic<ValueType> {
+        public:
+            DFTExplorationHeuristicBoundDifference(size_t id) : DFTExplorationHeuristic<ValueType>(id) {
+                // Intentionally left empty
+            }
+
+            DFTExplorationHeuristicBoundDifference(size_t id, DFTExplorationHeuristicBoundDifference<ValueType> const& predecessor, ValueType rate, ValueType exitRate, ValueType lowerBound, ValueType upperBound) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate, lowerBound, upperBound) {
+                // Intentionally left empty
+            }
+
+            bool updateHeuristicValues(DFTExplorationHeuristic<ValueType> const& predecessor, ValueType rate, ValueType exitRate) override {
+                return false;
+            }
+
+            double getPriority() const override;
+
+            bool isSkip(double approximationThreshold) const override {
+                return !this->expand && this->getPriority() < approximationThreshold;
+            }
+
+            bool operator<(DFTExplorationHeuristicBoundDifference<ValueType> const& other) const {
+                return this->getPriority() < other.getPriority();
+            }
+
+        private:
+            ValueType lowerBound;
+            ValueType upperBound;
+        };
+
+
     }
 }
 
diff --git a/src/builder/ExplicitDFTModelBuilderApprox.cpp b/src/builder/ExplicitDFTModelBuilderApprox.cpp
index 7de1e0bb4..0328cf055 100644
--- a/src/builder/ExplicitDFTModelBuilderApprox.cpp
+++ b/src/builder/ExplicitDFTModelBuilderApprox.cpp
@@ -3,6 +3,7 @@
 #include <src/models/sparse/Ctmc.h>
 #include <src/utility/constants.h>
 #include <src/utility/vector.h>
+#include "src/utility/bitoperations.h"
 #include <src/exceptions/UnexpectedException.h>
 #include "src/settings/modules/DFTSettings.h"
 #include "src/settings/SettingsManager.h"
@@ -32,12 +33,12 @@ namespace storm {
                 matrixBuilder(!generator.isDeterministicModel()),
                 stateStorage(((dft.stateVectorSize() / 64) + 1) * 64),
                 // TODO Matthias: make choosable
-                explorationQueue(dft.nrElements()+1, 0, 1)
-                //explorationQueue(1001, 0, 0.001)
+                //explorationQueue(dft.nrElements()+1, 0, 1)
+                explorationQueue(1001, 0, 0.001)
         {
             // Intentionally left empty.
             // TODO Matthias: remove again
-            heuristic = storm::builder::ApproximationHeuristic::NONE;
+            heuristic = storm::builder::ApproximationHeuristic::PROBABILITY;
         }
 
         template<typename ValueType, typename StateType>
@@ -256,6 +257,9 @@ namespace storm {
                 if (currentState->isPseudoState()) {
                     // Create concrete state from pseudo state
                     currentState->construct();
+                    ValueType lowerBound = getLowerBound(currentState);
+                    ValueType upperBound = getUpperBound(currentState);
+                    currentExplorationHeuristic->setBounds(lowerBound, upperBound);
                 }
                 STORM_LOG_ASSERT(!currentState->isPseudoState(), "State is pseudo state.");
 
@@ -267,8 +271,8 @@ namespace storm {
                 // Try to explore the next state
                 generator.load(currentState);
 
-                //if (nrExpandedStates > approximationThreshold) {
-                if (currentExplorationHeuristic->isSkip(approximationThreshold)) {
+                if (nrExpandedStates > approximationThreshold && !currentExplorationHeuristic->isExpand()) {
+                //if (currentExplorationHeuristic->isSkip(approximationThreshold)) {
                     // Skip the current state
                     ++nrSkippedStates;
                     STORM_LOG_TRACE("Skip expansion of state: " << dft.getStateString(currentState));
@@ -299,11 +303,24 @@ namespace storm {
                                 // Update heuristic values
                                 DFTStatePointer state = iter->second.first;
                                 if (!iter->second.second) {
+                                    ValueType lowerBound;
+                                    ValueType upperBound;
+                                    if (state->isPseudoState()) {
+                                        lowerBound = storm::utility::infinity<ValueType>();
+                                        upperBound = storm::utility::infinity<ValueType>();
+                                    } else {
+                                        lowerBound = getLowerBound(state);
+                                        upperBound = getUpperBound(state);
+                                    }
                                     // Initialize heuristic values
-                                    ExplorationHeuristicPointer heuristic = std::make_shared<ExplorationHeuristic>(stateProbabilityPair.first, *currentExplorationHeuristic, stateProbabilityPair.second, choice.getTotalMass());
+                                    ExplorationHeuristicPointer heuristic = std::make_shared<ExplorationHeuristic>(stateProbabilityPair.first, *currentExplorationHeuristic, stateProbabilityPair.second, choice.getTotalMass(), lowerBound, upperBound);
                                     iter->second.second = heuristic;
+                                    if (state->hasFailed(dft.getTopLevelIndex()) || state->isFailsafe(dft.getTopLevelIndex()) || state->nrFailableDependencies() > 0 || (state->nrFailableDependencies() == 0 && state->nrFailableBEs() == 0)) {
+                                        // Do not skip absorbing state or if reached by dependencies
+                                        iter->second.second->markExpand();
+                                    }
                                     explorationQueue.push(heuristic);
-                                } else {
+                                } else if (!iter->second.second->isExpand()) {
                                     double oldPriority = iter->second.second->getPriority();
                                     if (iter->second.second->updateHeuristicValues(*currentExplorationHeuristic, stateProbabilityPair.second, choice.getTotalMass())) {
                                         // Update priority queue
@@ -311,11 +328,6 @@ namespace storm {
                                         explorationQueue.update(iter->second.second, oldPriority);
                                     }
                                 }
-                                if (state->hasFailed(dft.getTopLevelIndex()) || state->isFailsafe(dft.getTopLevelIndex()) || state->nrFailableDependencies() > 0 || (state->nrFailableDependencies() == 0 && state->nrFailableBEs() == 0)) {
-                                    // Do not skip absorbing state or if reached by dependencies
-                                    iter->second.second->markExpand();
-                                    // TODO Matthias give highest priority to ensure expanding before all skipped
-                                }
                             }
                         }
                         matrixBuilder.finishRow();
@@ -383,11 +395,7 @@ namespace storm {
         template<typename ValueType, typename StateType>
         std::shared_ptr<storm::models::sparse::Model<ValueType>> ExplicitDFTModelBuilderApprox<ValueType, StateType>::getModelApproximation(bool lowerBound) {
             // TODO Matthias: handle case with no skipped states
-            if (lowerBound) {
-                changeMatrixLowerBound(modelComponents.transitionMatrix);
-            } else {
-                changeMatrixUpperBound(modelComponents.transitionMatrix);
-            }
+            changeMatrixBound(modelComponents.transitionMatrix, lowerBound);
             return createModel(true);
         }
 
@@ -443,44 +451,72 @@ namespace storm {
         }
 
         template<typename ValueType, typename StateType>
-        void ExplicitDFTModelBuilderApprox<ValueType, StateType>::changeMatrixLowerBound(storm::storage::SparseMatrix<ValueType> & matrix) const {
+        void ExplicitDFTModelBuilderApprox<ValueType, StateType>::changeMatrixBound(storm::storage::SparseMatrix<ValueType> & matrix, bool lowerBound) const {
             // Set lower bound for skipped states
             for (auto it = skippedStates.begin(); it != skippedStates.end(); ++it) {
                 auto matrixEntry = matrix.getRow(it->first, 0).begin();
                 STORM_LOG_ASSERT(matrixEntry->getColumn() == failedStateId, "Transition has wrong target state.");
-                // Get the lower bound by considering the failure of all possible BEs
-                DFTStatePointer state = it->second.first;
-                ValueType newRate = storm::utility::zero<ValueType>();
-                for (size_t index = 0; index < state->nrFailableBEs(); ++index) {
-                    newRate += state->getFailableBERate(index);
-                }
-                for (size_t index = 0; index < state->nrNotFailableBEs(); ++index) {
-                    newRate += state->getNotFailableBERate(index);
+                STORM_LOG_ASSERT(!it->second.first->isPseudoState(), "State is still pseudo state.");
+                if (lowerBound) {
+                    matrixEntry->setValue(it->second.second->getLowerBound());
+                } else {
+                    matrixEntry->setValue(it->second.second->getUpperBound());
                 }
-                matrixEntry->setValue(newRate);
             }
         }
 
         template<typename ValueType, typename StateType>
-        void ExplicitDFTModelBuilderApprox<ValueType, StateType>::changeMatrixUpperBound(storm::storage::SparseMatrix<ValueType> & matrix) const {
-            // Set uppper bound for skipped states
-            for (auto it = skippedStates.begin(); it != skippedStates.end(); ++it) {
-                auto matrixEntry = matrix.getRow(it->first, 0).begin();
-                STORM_LOG_ASSERT(matrixEntry->getColumn() == failedStateId, "Transition has wrong target state.");
-                // Get the upper bound by considering the failure of all BEs
-                // The used formula for the rate is 1/( 1/a + 1/b + ...)
-                // TODO Matthias: improve by using closed formula for AND of all BEs
-                DFTStatePointer state = it->second.first;
-                ValueType newRate = storm::utility::zero<ValueType>();
-                for (size_t index = 0; index < state->nrFailableBEs(); ++index) {
-                    newRate += storm::utility::one<ValueType>() / state->getFailableBERate(index);
-                }
-                for (size_t index = 0; index < state->nrNotFailableBEs(); ++index) {
-                    newRate += storm::utility::one<ValueType>() / state->getNotFailableBERate(index);
+        ValueType ExplicitDFTModelBuilderApprox<ValueType, StateType>::getLowerBound(DFTStatePointer const& state) const {
+            // Get the lower bound by considering the failure of all possible BEs
+            ValueType lowerBound = storm::utility::zero<ValueType>();
+            for (size_t index = 0; index < state->nrFailableBEs(); ++index) {
+                lowerBound += state->getFailableBERate(index);
+            }
+            for (size_t index = 0; index < state->nrNotFailableBEs(); ++index) {
+                lowerBound += state->getNotFailableBERate(index);
+            }
+            return lowerBound;
+        }
+
+        template<typename ValueType, typename StateType>
+        ValueType ExplicitDFTModelBuilderApprox<ValueType, StateType>::getUpperBound(DFTStatePointer const& state) const {
+            // Get the upper bound by considering the failure of all BEs
+            // The used formula for the rate is 1/( 1/a + 1/b + ...)
+            // TODO Matthias: improve by using closed formula for AND of all BEs
+            ValueType upperBound = storm::utility::zero<ValueType>();
+
+            // Get all possible rates
+            std::vector<ValueType> rates(state->nrFailableBEs() + state->nrNotFailableBEs());
+            for (size_t index = 0; index < state->nrFailableBEs(); ++index) {
+                rates[index] = state->getFailableBERate(index);
+            }
+            for (size_t index = 0; index < state->nrNotFailableBEs(); ++index) {
+                rates[index + state->nrFailableBEs()] = state->getNotFailableBERate(index);
+            }
+
+            // TODO Matthias: works only for <64 BEs!
+            for (size_t i = 1; i < 4 && i <= rates.size(); ++i) {
+                size_t permutation = smallestIntWithNBitsSet(static_cast<size_t>(i));
+                ValueType sum = storm::utility::zero<ValueType>();
+                do {
+                    ValueType permResult = storm::utility::zero<ValueType>();
+                    for(size_t j = 0; j < rates.size(); ++j) {
+                        if(permutation & (1 << j)) {
+                            permResult += rates[j];
+                        }
+                    }
+                    permutation = nextBitPermutation(permutation);
+                    STORM_LOG_ASSERT(!storm::utility::isZero(permResult), "PermResult is 0");
+                    sum += storm::utility::one<ValueType>() / permResult;
+                } while(permutation < (1 << rates.size()) && permutation != 0);
+                if (i % 2 == 0) {
+                    upperBound -= sum;
+                } else {
+                    upperBound += sum;
                 }
-                newRate = storm::utility::one<ValueType>() / newRate;
-                matrixEntry->setValue(newRate);
             }
+            STORM_LOG_ASSERT(!storm::utility::isZero(upperBound), "UpperBound is 0");
+            return storm::utility::one<ValueType>() / upperBound;
         }
 
         template<typename ValueType, typename StateType>
diff --git a/src/builder/ExplicitDFTModelBuilderApprox.h b/src/builder/ExplicitDFTModelBuilderApprox.h
index 1f779e44f..446bed6ad 100644
--- a/src/builder/ExplicitDFTModelBuilderApprox.h
+++ b/src/builder/ExplicitDFTModelBuilderApprox.h
@@ -28,7 +28,7 @@ namespace storm {
 
             using DFTStatePointer = std::shared_ptr<storm::storage::DFTState<ValueType>>;
             // TODO Matthias: make choosable
-            using ExplorationHeuristic = DFTExplorationHeuristicDepth<ValueType>;
+            using ExplorationHeuristic = DFTExplorationHeuristicProbability<ValueType>;
             using ExplorationHeuristicPointer = std::shared_ptr<ExplorationHeuristic>;
 
 
@@ -217,18 +217,26 @@ namespace storm {
             void setMarkovian(bool markovian);
 
             /**
-             * Change matrix to reflect the lower approximation bound.
+             * Change matrix to reflect the lower or upper approximation bound.
              *
-             * @param matrix Matrix to change. The change are reflected here.
+             * @param matrix     Matrix to change. The change are reflected here.
+             * @param lowerBound Flag indicating if the lower bound should be used. Otherwise the upper bound is used.
              */
-            void changeMatrixLowerBound(storm::storage::SparseMatrix<ValueType> & matrix) const;
+            void changeMatrixBound(storm::storage::SparseMatrix<ValueType> & matrix, bool lowerBound) const;
 
             /*!
-             * Change matrix to reflect the upper approximation bound.
+             * Get lower bound approximation for state.
              *
-             * @param matrix Matrix to change. The change are reflected here.
+             * @return Lower bound approximation.
              */
-            void changeMatrixUpperBound(storm::storage::SparseMatrix<ValueType> & matrix) const;
+            ValueType getLowerBound(DFTStatePointer const& state) const;
+
+            /*!
+             * Get upper bound approximation for state.
+             *
+             * @return Upper bound approximation.
+             */
+            ValueType getUpperBound(DFTStatePointer const& state) const;
 
             /*!
              * Compares the priority of two states.
diff --git a/src/modelchecker/dft/DFTModelChecker.cpp b/src/modelchecker/dft/DFTModelChecker.cpp
index f19943f2f..59d85f5c2 100644
--- a/src/modelchecker/dft/DFTModelChecker.cpp
+++ b/src/modelchecker/dft/DFTModelChecker.cpp
@@ -255,6 +255,7 @@ namespace storm {
 
         template<>
         bool DFTModelChecker<double>::isApproximationSufficient(double lowerBound, double upperBound, double approximationError) {
+            STORM_LOG_THROW(!std::isnan(lowerBound) && !std::isnan(upperBound), storm::exceptions::NotSupportedException, "Approximation does not work if result is NaN.");
             return upperBound - lowerBound <= approximationError * (lowerBound + upperBound) / 2;
         }
 
diff --git a/src/storage/BucketPriorityQueue.cpp b/src/storage/BucketPriorityQueue.cpp
index fe1ca95b5..eb174cce0 100644
--- a/src/storage/BucketPriorityQueue.cpp
+++ b/src/storage/BucketPriorityQueue.cpp
@@ -16,19 +16,19 @@ namespace storm {
         void BucketPriorityQueue<ValueType>::fix() {
             if (currentBucket < buckets.size() && nrUnsortedItems > 0) {
                 // Fix current bucket
-                //std::make_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare);
+                std::make_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare);
                 nrUnsortedItems = 0;
             }
         }
 
         template<typename ValueType>
         bool BucketPriorityQueue<ValueType>::empty() const {
-            return currentBucket == buckets.size();
+            return currentBucket == buckets.size() && immediateBucket.empty();
         }
 
         template<typename ValueType>
         std::size_t BucketPriorityQueue<ValueType>::size() const {
-            size_t size = 0;
+            size_t size = immediateBucket.size();
             for (size_t i = currentBucket; currentBucket < buckets.size(); ++i) {
                 size += buckets[i].size();
             }
@@ -37,6 +37,9 @@ namespace storm {
 
         template<typename ValueType>
         typename BucketPriorityQueue<ValueType>::HeuristicPointer const& BucketPriorityQueue<ValueType>::top() const {
+            if (!immediateBucket.empty()) {
+                return immediateBucket.back();
+            }
             STORM_LOG_ASSERT(!empty(), "BucketPriorityQueue is empty");
             STORM_LOG_ASSERT(nrUnsortedItems == 0, "First bucket is not sorted");
             return buckets[currentBucket].front();
@@ -44,6 +47,10 @@ namespace storm {
 
         template<typename ValueType>
         void BucketPriorityQueue<ValueType>::push(HeuristicPointer const& item) {
+            if (item->isExpand()) {
+                immediateBucket.push_back(item);
+                return;
+            }
             size_t bucket = getBucket(item->getPriority());
             if (bucket < currentBucket) {
                 currentBucket = bucket;
@@ -62,6 +69,7 @@ namespace storm {
 
         template<typename ValueType>
         void BucketPriorityQueue<ValueType>::update(HeuristicPointer const& item, double oldPriority) {
+            STORM_LOG_ASSERT(!item->isExpand(), "Item is marked for expansion");
             size_t newBucket = getBucket(item->getPriority());
             size_t oldBucket = getBucket(oldPriority);
 
@@ -113,6 +121,10 @@ namespace storm {
 
         template<typename ValueType>
         void BucketPriorityQueue<ValueType>::pop() {
+            if (!immediateBucket.empty()) {
+                immediateBucket.pop_back();
+                return;
+            }
             STORM_LOG_ASSERT(!empty(), "BucketPriorityQueue is empty");
             STORM_LOG_ASSERT(nrUnsortedItems == 0, "First bucket is not sorted");
             std::pop_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare);
@@ -153,6 +165,11 @@ namespace storm {
         template<typename ValueType>
         void BucketPriorityQueue<ValueType>::print(std::ostream& out) const {
             out << "Bucket priority queue with size " << buckets.size() << ", lower value: " << lowerValue << " and step per bucket: " << stepPerBucket << std::endl;
+            out << "Immediate bucket: ";
+            for (HeuristicPointer heuristic : immediateBucket) {
+                out << heuristic->getId() << ", ";
+            }
+            out << std::endl;
             out << "Current bucket (" << currentBucket << ") has " << nrUnsortedItems  << " unsorted items" << std::endl;
             for (size_t bucket = 0; bucket < buckets.size(); ++bucket) {
                 if (!buckets[bucket].empty()) {
diff --git a/src/storage/BucketPriorityQueue.h b/src/storage/BucketPriorityQueue.h
index 3d3c2d35d..253b34e19 100644
--- a/src/storage/BucketPriorityQueue.h
+++ b/src/storage/BucketPriorityQueue.h
@@ -13,7 +13,7 @@ namespace storm {
         template<typename ValueType>
         class BucketPriorityQueue {
 
-            using HeuristicPointer = std::shared_ptr<storm::builder::DFTExplorationHeuristicDepth<ValueType>>;
+            using HeuristicPointer = std::shared_ptr<storm::builder::DFTExplorationHeuristicProbability<ValueType>>;
 
         public:
             explicit BucketPriorityQueue(size_t nrBuckets, double lowerValue, double stepPerBucket);
@@ -43,6 +43,9 @@ namespace storm {
             // List of buckets
             std::vector<std::vector<HeuristicPointer>> buckets;
 
+            // Bucket containing all items which should be considered immediately
+            std::vector<HeuristicPointer> immediateBucket;
+
             // Index of first bucket which contains items
             size_t currentBucket;
 
diff --git a/src/storage/dft/DFTState.cpp b/src/storage/dft/DFTState.cpp
index 369bed3f7..b9381598f 100644
--- a/src/storage/dft/DFTState.cpp
+++ b/src/storage/dft/DFTState.cpp
@@ -18,7 +18,9 @@ namespace storm {
             }
 
             for (auto elem : mDft.getBasicElements()) {
-                mCurrentlyNotFailableBE.push_back(elem->id());
+                if (!storm::utility::isZero(elem->activeFailureRate())) {
+                    mCurrentlyNotFailableBE.push_back(elem->id());
+                }
             }
 
             // Initialize activation
@@ -242,8 +244,10 @@ namespace storm {
         ValueType DFTState<ValueType>::getBERate(size_t id, bool considerPassive) const {
             STORM_LOG_ASSERT(mDft.isBasicElement(id), "Element is no BE.");
             if (considerPassive && mDft.hasRepresentant(id) && !isActive(mDft.getRepresentant(id))) {
+                STORM_LOG_ASSERT(!storm::utility::isZero(mDft.getBasicElement(id)->passiveFailureRate()), "Failure rate of BE " << mDft.getBasicElement(id)->id() << " is 0 in state " << mDft.getStateString(mStatus, mStateGenerationInfo, mId));
                 return mDft.getBasicElement(id)->passiveFailureRate();
             } else {
+                STORM_LOG_ASSERT(!storm::utility::isZero(mDft.getBasicElement(id)->activeFailureRate()), "Failure rate of BE " << mDft.getBasicElement(id)->id() << " is 0 in state " << mDft.getStateString(mStatus, mStateGenerationInfo, mId));
                 return mDft.getBasicElement(id)->activeFailureRate();
             }
         }
@@ -257,7 +261,7 @@ namespace storm {
         template<typename ValueType>
         ValueType DFTState<ValueType>::getNotFailableBERate(size_t index) const {
             STORM_LOG_ASSERT(index < nrNotFailableBEs(), "Index invalid.");
-            STORM_LOG_ASSERT(storm::utility::isZero<ValueType>(mDft.getBasicElement(mCurrentlyNotFailableBE[index])->activeFailureRate()) ||
+            STORM_LOG_ASSERT(isPseudoState() || storm::utility::isZero<ValueType>(mDft.getBasicElement(mCurrentlyNotFailableBE[index])->activeFailureRate()) ||
                              (mDft.hasRepresentant(mCurrentlyNotFailableBE[index]) && !isActive(mDft.getRepresentant(mCurrentlyNotFailableBE[index]))), "BE " << mCurrentlyNotFailableBE[index] << " can fail in state: " << mDft.getStateString(mStatus, mStateGenerationInfo, mId));
             // Use active failure rate as passive failure rate is 0.
             return getBERate(mCurrentlyNotFailableBE[index], false);