From f5b1259f3c72ef90c1dd300cd78c39103e404231 Mon Sep 17 00:00:00 2001 From: dehnert Date: Tue, 28 Nov 2017 11:06:05 +0100 Subject: [PATCH] fixed issue related to Markov automata without proababilistic states --- .../helper/SparseMarkovAutomatonCslHelper.cpp | 67 ++++++++++++------- 1 file changed, 44 insertions(+), 23 deletions(-) diff --git a/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp index bdf030f31..ad3db92f0 100644 --- a/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp @@ -44,9 +44,16 @@ namespace storm { // * a matrix aProbabilistic with all (non-discretized) transitions from probabilistic non-goal states to other probabilistic non-goal states. // * a matrix aProbabilisticToMarkovian with all (non-discretized) transitions from probabilistic non-goal states to all Markovian non-goal states. typename storm::storage::SparseMatrix aMarkovian = transitionMatrix.getSubmatrix(true, markovianNonGoalStates, markovianNonGoalStates, true); - typename storm::storage::SparseMatrix aMarkovianToProbabilistic = transitionMatrix.getSubmatrix(true, markovianNonGoalStates, probabilisticNonGoalStates); - typename storm::storage::SparseMatrix aProbabilistic = transitionMatrix.getSubmatrix(true, probabilisticNonGoalStates, probabilisticNonGoalStates); - typename storm::storage::SparseMatrix aProbabilisticToMarkovian = transitionMatrix.getSubmatrix(true, probabilisticNonGoalStates, markovianNonGoalStates); + + bool existProbabilisticStates = !probabilisticNonGoalStates.empty(); + typename storm::storage::SparseMatrix aMarkovianToProbabilistic; + typename storm::storage::SparseMatrix aProbabilistic; + typename storm::storage::SparseMatrix aProbabilisticToMarkovian; + if (existProbabilisticStates) { + aMarkovianToProbabilistic = transitionMatrix.getSubmatrix(true, markovianNonGoalStates, probabilisticNonGoalStates); + aProbabilistic = transitionMatrix.getSubmatrix(true, probabilisticNonGoalStates, probabilisticNonGoalStates); + aProbabilisticToMarkovian = transitionMatrix.getSubmatrix(true, probabilisticNonGoalStates, markovianNonGoalStates); + } // The matrices with transitions from Markovian states need to be digitized. // Digitize aMarkovian. Based on whether the transition is a self-loop or not, we apply the two digitization rules. @@ -64,20 +71,25 @@ namespace storm { } // Digitize aMarkovianToProbabilistic. As there are no self-loops in this case, we only need to apply the digitization formula for regular successors. - rowIndex = 0; - for (auto state : markovianNonGoalStates) { - for (auto& element : aMarkovianToProbabilistic.getRow(rowIndex)) { - element.setValue((1 - std::exp(-exitRates[state] * delta)) * element.getValue()); + if (existProbabilisticStates) { + rowIndex = 0; + for (auto state : markovianNonGoalStates) { + for (auto& element : aMarkovianToProbabilistic.getRow(rowIndex)) { + element.setValue((1 - std::exp(-exitRates[state] * delta)) * element.getValue()); + } + ++rowIndex; } - ++rowIndex; } // Initialize the two vectors that hold the variable one-step probabilities to all target states for probabilistic and Markovian (non-goal) states. - std::vector bProbabilistic(aProbabilistic.getRowCount()); + std::vector bProbabilistic(existProbabilisticStates ? aProbabilistic.getRowCount() : 0); std::vector bMarkovian(markovianNonGoalStates.getNumberOfSetBits()); // Compute the two fixed right-hand side vectors, one for Markovian states and one for the probabilistic ones. - std::vector bProbabilisticFixed = transitionMatrix.getConstrainedRowGroupSumVector(probabilisticNonGoalStates, goalStates); + std::vector bProbabilisticFixed; + if (existProbabilisticStates) { + bProbabilisticFixed = transitionMatrix.getConstrainedRowGroupSumVector(probabilisticNonGoalStates, goalStates); + } std::vector bMarkovianFixed; bMarkovianFixed.reserve(markovianNonGoalStates.getNumberOfSetBits()); for (auto state : markovianNonGoalStates) { @@ -110,25 +122,34 @@ namespace storm { // * perform one timed-step using v_MS := A_MSwG * v_MS + A_MStoPS * v_PS + (A * 1_G)|MS std::vector markovianNonGoalValuesSwap(markovianNonGoalValues); for (uint64_t currentStep = 0; currentStep < numberOfSteps; ++currentStep) { - // Start by (re-)computing bProbabilistic = bProbabilisticFixed + aProbabilisticToMarkovian * vMarkovian. - aProbabilisticToMarkovian.multiplyWithVector(markovianNonGoalValues, bProbabilistic); - storm::utility::vector::addVectors(bProbabilistic, bProbabilisticFixed, bProbabilistic); + if (existProbabilisticStates) { + // Start by (re-)computing bProbabilistic = bProbabilisticFixed + aProbabilisticToMarkovian * vMarkovian. + aProbabilisticToMarkovian.multiplyWithVector(markovianNonGoalValues, bProbabilistic); + storm::utility::vector::addVectors(bProbabilistic, bProbabilisticFixed, bProbabilistic); - // Now perform the inner value iteration for probabilistic states. - solver->solveEquations(env, dir, probabilisticNonGoalValues, bProbabilistic); + // Now perform the inner value iteration for probabilistic states. + solver->solveEquations(env, dir, probabilisticNonGoalValues, bProbabilistic); + + // (Re-)compute bMarkovian = bMarkovianFixed + aMarkovianToProbabilistic * vProbabilistic. + aMarkovianToProbabilistic.multiplyWithVector(probabilisticNonGoalValues, bMarkovian); + storm::utility::vector::addVectors(bMarkovian, bMarkovianFixed, bMarkovian); + } - // (Re-)compute bMarkovian = bMarkovianFixed + aMarkovianToProbabilistic * vProbabilistic. - aMarkovianToProbabilistic.multiplyWithVector(probabilisticNonGoalValues, bMarkovian); - storm::utility::vector::addVectors(bMarkovian, bMarkovianFixed, bMarkovian); aMarkovian.multiplyWithVector(markovianNonGoalValues, markovianNonGoalValuesSwap); std::swap(markovianNonGoalValues, markovianNonGoalValuesSwap); - storm::utility::vector::addVectors(markovianNonGoalValues, bMarkovian, markovianNonGoalValues); + if (existProbabilisticStates) { + storm::utility::vector::addVectors(markovianNonGoalValues, bMarkovian, markovianNonGoalValues); + } else { + storm::utility::vector::addVectors(markovianNonGoalValues, bMarkovianFixed, markovianNonGoalValues); + } } - // After the loop, perform one more step of the value iteration for PS states. - aProbabilisticToMarkovian.multiplyWithVector(markovianNonGoalValues, bProbabilistic); - storm::utility::vector::addVectors(bProbabilistic, bProbabilisticFixed, bProbabilistic); - solver->solveEquations(env, dir, probabilisticNonGoalValues, bProbabilistic); + if (existProbabilisticStates) { + // After the loop, perform one more step of the value iteration for PS states. + aProbabilisticToMarkovian.multiplyWithVector(markovianNonGoalValues, bProbabilistic); + storm::utility::vector::addVectors(bProbabilistic, bProbabilisticFixed, bProbabilistic); + solver->solveEquations(env, dir, probabilisticNonGoalValues, bProbabilistic); + } } template ::SupportsExponential, int>::type>