diff --git a/src/modelchecker/AbstractModelChecker.cpp b/src/modelchecker/AbstractModelChecker.cpp index 7ba1de099..2b6d2ef3a 100644 --- a/src/modelchecker/AbstractModelChecker.cpp +++ b/src/modelchecker/AbstractModelChecker.cpp @@ -216,8 +216,8 @@ namespace storm { } } - std::unique_ptr AbstractModelChecker::checkExpectedTimeOperatorFormula(storm::logic::ExpectedTimeOperatorFormula const& stateFormula) { - STORM_LOG_THROW(stateFormula.getSubformula().isStateFormula(), storm::exceptions::InvalidArgumentException, "The given formula is invalid."); + std::unique_ptr AbstractModelChecker::checkExpectedTimeOperatorFormula(storm::logic::ExpectedTimeOperatorFormula const& stateFormula) { + STORM_LOG_THROW(stateFormula.getSubformula().isEventuallyFormula(), storm::exceptions::InvalidArgumentException, "The given formula is invalid."); // If the reward bound is 0, is suffices to do qualitative model checking. bool qualitative = false; @@ -248,8 +248,8 @@ namespace storm { } } - std::unique_ptr AbstractModelChecker::checkLongRunAverageOperatorFormula(storm::logic::LongRunAverageOperatorFormula const& stateFormula) { - STORM_LOG_THROW(stateFormula.getSubformula().isEventuallyFormula(), storm::exceptions::InvalidArgumentException, "The given formula is invalid."); + std::unique_ptr AbstractModelChecker::checkLongRunAverageOperatorFormula(storm::logic::LongRunAverageOperatorFormula const& stateFormula) { + STORM_LOG_THROW(stateFormula.getSubformula().isStateFormula(), storm::exceptions::InvalidArgumentException, "The given formula is invalid."); std::unique_ptr result; if (stateFormula.hasOptimalityType()) { diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp index 02a2b5758..9c70023a7 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp @@ -340,7 +340,6 @@ namespace storm { storm::storage::BitVector statesNotInBsccs = ~statesInBsccs; // calculate steady state distribution for all BSCCs by calculating an eigenvector for the eigenvalue 1 of the transposed transition matrix for the bsccs - std::unique_ptr> solver = storm::utility::solver::getLinearEquationSolver(); storm::storage::SparseMatrix bsccEquationSystem = transitionMatrix.getSubmatrix(false, statesInBsccs, statesInBsccs, true).transpose(); ValueType one = storm::utility::one(); @@ -356,7 +355,10 @@ namespace storm { std::vector bsccEquationSystemRightSide(bsccEquationSystem.getColumnCount(), zero); std::vector bsccEquationSystemSolution(bsccEquationSystem.getColumnCount(), one); - solver->solveEquationSystem(bsccEquationSystem, bsccEquationSystemSolution, bsccEquationSystemRightSide); + { + auto solver = this->linearEquationSolverFactory->create(bsccEquationSystem); + solver->solveEquationSystem(bsccEquationSystemSolution, bsccEquationSystemRightSide); + } //calculate LRA Value for each BSCC from steady state distribution in BSCCs // we have to do some scaling, as the probabilities for each BSCC have to sum up to one, which they don't necessarily do in the solution of the equation system @@ -395,7 +397,10 @@ namespace storm { std::vector rewardSolution(rewardEquationSystemMatrix.getColumnCount(), one); - solver->solveEquationSystem(rewardEquationSystemMatrix, rewardSolution, rewardRightSide); + { + auto solver = this->linearEquationSolverFactory->create(rewardEquationSystemMatrix); + solver->solveEquationSystem(rewardSolution, rewardRightSide); + } // now fill the result vector std::vector result(numOfStates); @@ -453,7 +458,7 @@ namespace storm { template - ValueType SparseDtmcPrctlModelChecker::computeLraForBSCC(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, storm::storage::StronglyConnectedComponent const& bscc) { + ValueType SparseDtmcPrctlModelChecker::computeLraForBSCC(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, storm::storage::StronglyConnectedComponent const& bscc, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { //if size is one we can immediately derive the result if (bscc.size() == 1){ if (psiStates.get(*(bscc.begin()))) { @@ -462,7 +467,6 @@ namespace storm { return storm::utility::zero(); } } - std::unique_ptr> solver = storm::utility::solver::getLinearEquationSolver(); storm::storage::BitVector subsystem = storm::storage::BitVector(transitionMatrix.getRowCount()); subsystem.set(bscc.begin(), bscc.end()); @@ -502,7 +506,11 @@ namespace storm { std::vector b(subsystemMatrix.getRowCount(), zero); b[subsystemMatrix.getRowCount() - 1] = one; - solver->solveEquationSystem(subsystemMatrix, steadyStateDist, b); + + { + auto solver = linearEquationSolverFactory.create(subsystemMatrix); + solver->solveEquationSystem(steadyStateDist, b); + } //remove the last entry of the vector as it was just there to enforce the unique solution steadyStateDist.pop_back(); diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h index a84413495..b55f00faa 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h @@ -43,7 +43,7 @@ namespace storm { static std::vector computeReachabilityRewardsHelper(storm::storage::SparseMatrix const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative); std::vector computeLongRunAverageHelper(storm::storage::BitVector const& psiStates, bool qualitative) const; - static ValueType computeLraForBSCC(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& goalStates, storm::storage::StronglyConnectedComponent const& bscc); + static ValueType computeLraForBSCC(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& goalStates, storm::storage::StronglyConnectedComponent const& bscc, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); // An object that is used for retrieving linear equation solvers. std::unique_ptr> linearEquationSolverFactory; diff --git a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp index ba88c4ea7..9fa78f041 100644 --- a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp @@ -366,14 +366,12 @@ namespace storm { //This transitions have the LRA of the MEC as reward. //The expected reward corresponds to sum of LRAs in MEC weighted by the reachability probability of the MEC - std::unique_ptr> solver = storm::utility::solver::getLinearEquationSolver(); - //we now build the submatrix of the transition matrix of the system with the auxiliary state, that only contains the states from //the original state, i.e. all "maybe-states" storm::storage::SparseMatrixBuilder rewardEquationSystemMatrixBuilder(transitionMatrix.getRowCount() + statesInMecs.getNumberOfSetBits(), transitionMatrix.getColumnCount(), transitionMatrix.getEntryCount(), - true, + false, true, transitionMatrix.getRowGroupCount()); @@ -384,19 +382,32 @@ namespace storm { for (uint_fast64_t rowGroupIndex = 0; rowGroupIndex < transitionMatrix.getRowGroupCount(); ++rowGroupIndex) { rewardEquationSystemMatrixBuilder.newRowGroup(rowIndex); for (uint_fast64_t i = 0; i < transitionMatrix.getRowGroupSize(rowGroupIndex); ++i) { + //we have to make sure that an entry exists for all diagonal elements, even if it is zero. Other wise the call to convertToEquationSystem will produce wrong results or fail. + bool foundDiagonal = false; for (auto entry : transitionMatrix.getRow(oldRowIndex)) { + if (!foundDiagonal) { + if (entry.getColumn() > rowGroupIndex) { + foundDiagonal = true; + rewardEquationSystemMatrixBuilder.addNextValue(rowIndex, rowGroupIndex, zero); + } else if (entry.getColumn() == rowGroupIndex) { + foundDiagonal = true; + } + } //copy over values from transition matrix of the actual system rewardEquationSystemMatrixBuilder.addNextValue(rowIndex, entry.getColumn(), entry.getValue()); } + if (!foundDiagonal) { + rewardEquationSystemMatrixBuilder.addNextValue(rowIndex, rowGroupIndex, zero); + } ++oldRowIndex; ++rowIndex; } - ++rowGroupIndex; if (statesInMecs.get(rowGroupIndex)) { + //put the transition-reward on the right side + rewardRightSide[rowIndex] = mecLra[stateToMecIndexMap[rowGroupIndex]]; //add the choice where we go to the auxiliary state, which is a row with all zeros in the submatrix we build + rewardEquationSystemMatrixBuilder.addNextValue(rowIndex, rowGroupIndex, zero); ++rowIndex; - //put the transition-reward on the right side - rewardRightSide[rowIndex] = mecLra[rowGroupIndex]; } } @@ -405,7 +416,11 @@ namespace storm { std::vector result(rewardEquationSystemMatrix.getColumnCount(), one); - solver->solveEquationSystem(rewardEquationSystemMatrix, result, rewardRightSide); + { + auto solver = this->MinMaxLinearEquationSolverFactory->create(rewardEquationSystemMatrix); + solver->solveEquationSystem(minimize, result, rewardRightSide); + } + return result; } diff --git a/test/functional/modelchecker/SparseMdpPrctlModelCheckerTest.cpp b/test/functional/modelchecker/SparseMdpPrctlModelCheckerTest.cpp index 9a1ab6b1d..3ea26261f 100644 --- a/test/functional/modelchecker/SparseMdpPrctlModelCheckerTest.cpp +++ b/test/functional/modelchecker/SparseMdpPrctlModelCheckerTest.cpp @@ -197,7 +197,7 @@ TEST(SparseMdpPrctlModelCheckerTest, AsynchronousLeader) { TEST(SparseMdpPrctlModelCheckerTest, LRA) { storm::storage::SparseMatrixBuilder matrixBuilder; - std::shared_ptr> mdp; + std::shared_ptr> mdp; { matrixBuilder = storm::storage::SparseMatrixBuilder(2, 2, 2); @@ -205,16 +205,16 @@ TEST(SparseMdpPrctlModelCheckerTest, LRA) { matrixBuilder.addNextValue(1, 0, 1.); storm::storage::SparseMatrix transitionMatrix = matrixBuilder.build(); - storm::models::AtomicPropositionsLabeling ap(2, 1); - ap.addAtomicProposition("a"); - ap.addAtomicPropositionToState("a", 1); + storm::models::sparse::StateLabeling ap(2); + ap.addLabel("a"); + ap.addLabelToState("a", 1); - mdp.reset(new storm::models::Mdp(transitionMatrix, ap, boost::none, boost::none, boost::none)); + mdp.reset(new storm::models::sparse::Mdp(transitionMatrix, ap, boost::none, boost::none, boost::none)); - storm::modelchecker::SparseMdpPrctlModelChecker checker(*mdp, std::shared_ptr>(new storm::solver::NativeNondeterministicLinearEquationSolver())); + storm::modelchecker::SparseMdpPrctlModelChecker checker(*mdp, std::unique_ptr>(new storm::utility::solver::NativeMinMaxLinearEquationSolverFactory())); auto labelFormula = std::make_shared("a"); - auto lraFormula = std::make_shared(storm::logic::OptimalityType::Minimize, labelFormula); + auto lraFormula = std::make_shared(storm::logic::OptimalityType::Maximize, labelFormula); std::unique_ptr result = checker.check(*lraFormula); storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult = result->asExplicitQuantitativeCheckResult();