diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp index 634eaf806..2f247b8bd 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp @@ -52,7 +52,7 @@ namespace storm { storm::storage::BitVector statesWithProbabilityGreater0NonPsi = statesWithProbabilityGreater0 & ~psiStates; STORM_LOG_INFO("Found " << statesWithProbabilityGreater0NonPsi.getNumberOfSetBits() << " 'maybe' states."); - if (!statesWithProbabilityGreater0NonPsi.empty()) { + if (!statesWithProbabilityGreater0.empty()) { if (storm::utility::isZero(upperBound)) { // In this case, the interval is of the form [0, 0]. result = std::vector(numberOfStates, storm::utility::zero()); @@ -62,30 +62,31 @@ namespace storm { // 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. - // 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]); - } - uniformizationRate *= 1.02; - STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive."); - - // Compute the uniformized matrix. - storm::storage::SparseMatrix uniformizedMatrix = computeUniformizedMatrix(rateMatrix, statesWithProbabilityGreater0NonPsi, uniformizationRate, exitRates); - - // Compute the vector that is to be added as a compensation for removing the absorbing states. - std::vector b = rateMatrix.getConstrainedRowSumVector(statesWithProbabilityGreater0NonPsi, psiStates); - for (auto& element : b) { - element /= uniformizationRate; - } - - // Finally compute the transient probabilities. - std::vector values(statesWithProbabilityGreater0NonPsi.getNumberOfSetBits(), storm::utility::zero()); - std::vector subresult = computeTransientProbabilities(env, uniformizedMatrix, &b, upperBound, uniformizationRate, values); result = std::vector(numberOfStates, storm::utility::zero()); - - storm::utility::vector::setVectorValues(result, statesWithProbabilityGreater0NonPsi, subresult); - storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); + storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); + if (!statesWithProbabilityGreater0NonPsi.empty()) { + // 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]); + } + uniformizationRate *= 1.02; + STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive."); + + // Compute the uniformized matrix. + storm::storage::SparseMatrix uniformizedMatrix = computeUniformizedMatrix(rateMatrix, statesWithProbabilityGreater0NonPsi, uniformizationRate, exitRates); + + // Compute the vector that is to be added as a compensation for removing the absorbing states. + std::vector b = rateMatrix.getConstrainedRowSumVector(statesWithProbabilityGreater0NonPsi, psiStates); + for (auto& element : b) { + element /= uniformizationRate; + } + + // Finally compute the transient probabilities. + std::vector values(statesWithProbabilityGreater0NonPsi.getNumberOfSetBits(), storm::utility::zero()); + std::vector subresult = computeTransientProbabilities(env, uniformizedMatrix, &b, upperBound, uniformizationRate, values); + storm::utility::vector::setVectorValues(result, statesWithProbabilityGreater0NonPsi, subresult); + } } else if (upperBound == storm::utility::infinity()) { // In this case, the interval is of the form [t, inf] with t != 0. @@ -120,35 +121,38 @@ namespace storm { if (lowerBound != upperBound) { // In this case, the interval is of the form [t, t'] with t != 0, t' != inf and t != t'. - // Find the maximal rate of all 'maybe' states to take it as the uniformization rate. - ValueType uniformizationRate = storm::utility::zero(); - for (auto const& state : statesWithProbabilityGreater0NonPsi) { - uniformizationRate = std::max(uniformizationRate, exitRates[state]); - } - uniformizationRate *= 1.02; - STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive."); - - // Compute the (first) uniformized matrix. - storm::storage::SparseMatrix uniformizedMatrix = computeUniformizedMatrix(rateMatrix, statesWithProbabilityGreater0NonPsi, uniformizationRate, exitRates); - - // Compute the vector that is to be added as a compensation for removing the absorbing states. - std::vector b = rateMatrix.getConstrainedRowSumVector(statesWithProbabilityGreater0NonPsi, psiStates); - for (auto& element : b) { - element /= uniformizationRate; - } - - // Start by computing the transient probabilities of reaching a psi state in time t' - t. - std::vector values(statesWithProbabilityGreater0NonPsi.getNumberOfSetBits(), storm::utility::zero()); - std::vector subresult = computeTransientProbabilities(env, uniformizedMatrix, &b, upperBound - lowerBound, uniformizationRate, values); storm::storage::BitVector relevantStates = statesWithProbabilityGreater0 & phiStates; - std::vector newSubresult = std::vector(relevantStates.getNumberOfSetBits()); - storm::utility::vector::setVectorValues(newSubresult, statesWithProbabilityGreater0NonPsi % relevantStates, subresult); + std::vector newSubresult(relevantStates.getNumberOfSetBits(), storm::utility::zero()); storm::utility::vector::setVectorValues(newSubresult, psiStates % relevantStates, storm::utility::one()); + if (!statesWithProbabilityGreater0NonPsi.empty()) { + // Find the maximal rate of all 'maybe' states to take it as the uniformization rate. + ValueType uniformizationRate = storm::utility::zero(); + for (auto const& state : statesWithProbabilityGreater0NonPsi) { + uniformizationRate = std::max(uniformizationRate, exitRates[state]); + } + uniformizationRate *= 1.02; + STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive."); + + // Compute the (first) uniformized matrix. + storm::storage::SparseMatrix uniformizedMatrix = computeUniformizedMatrix(rateMatrix, statesWithProbabilityGreater0NonPsi, uniformizationRate, exitRates); + + // Compute the vector that is to be added as a compensation for removing the absorbing states. + std::vector b = rateMatrix.getConstrainedRowSumVector(statesWithProbabilityGreater0NonPsi, psiStates); + for (auto& element : b) { + element /= uniformizationRate; + } + + // Start by computing the transient probabilities of reaching a psi state in time t' - t. + std::vector values(statesWithProbabilityGreater0NonPsi.getNumberOfSetBits(), storm::utility::zero()); + std::vector subresult = computeTransientProbabilities(env, uniformizedMatrix, &b, upperBound - lowerBound, uniformizationRate, values); + newSubresult = std::vector(relevantStates.getNumberOfSetBits()); + storm::utility::vector::setVectorValues(newSubresult, statesWithProbabilityGreater0NonPsi % relevantStates, 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. - uniformizationRate = storm::utility::zero(); + ValueType uniformizationRate = storm::utility::zero(); for (auto const& state : relevantStates) { uniformizationRate = std::max(uniformizationRate, exitRates[state]); } @@ -156,7 +160,7 @@ namespace storm { STORM_LOG_THROW(uniformizationRate > 0, storm::exceptions::InvalidStateException, "The uniformization rate must be positive."); // Finally, we compute the second set of transient probabilities. - uniformizedMatrix = computeUniformizedMatrix(rateMatrix, relevantStates, uniformizationRate, exitRates); + storm::storage::SparseMatrix uniformizedMatrix = computeUniformizedMatrix(rateMatrix, relevantStates, uniformizationRate, exitRates); newSubresult = computeTransientProbabilities(env, uniformizedMatrix, nullptr, lowerBound, uniformizationRate, newSubresult); // Fill in the correct values.