Browse Source

worked on parametric model simplifier

main
TimQu 8 years ago
parent
commit
732bbc85d2
  1. 4
      src/storm/modelchecker/reachability/SparseDtmcEliminationModelChecker.cpp
  2. 10
      src/storm/modelchecker/region/ApproximationModel.cpp
  3. 4
      src/storm/modelchecker/region/SamplingModel.cpp
  4. 32
      src/storm/models/sparse/StandardRewardModel.cpp
  5. 9
      src/storm/models/sparse/StandardRewardModel.h
  6. 93
      src/storm/storage/FlexibleSparseMatrix.cpp
  7. 24
      src/storm/storage/FlexibleSparseMatrix.h
  8. 190
      src/storm/transformer/GoalStateMerger.cpp
  9. 161
      src/storm/transformer/GoalStateMerger.h
  10. 77
      src/storm/transformer/SparseParametricDtmcSimplifier.cpp
  11. 6
      src/storm/transformer/SparseParametricModelSimplifier.cpp

4
src/storm/modelchecker/reachability/SparseDtmcEliminationModelChecker.cpp

@ -186,9 +186,9 @@ namespace storm {
// Then, we convert the reduced matrix to a more flexible format to be able to perform state elimination more easily.
storm::storage::FlexibleSparseMatrix<ValueType> flexibleMatrix(transitionMatrix);
flexibleMatrix.createSubmatrix(maybeStates, maybeStates);
flexibleMatrix.filterEntries(maybeStates, maybeStates);
storm::storage::FlexibleSparseMatrix<ValueType> flexibleBackwardTransitions(backwardTransitions);
flexibleBackwardTransitions.createSubmatrix(maybeStates, maybeStates);
flexibleBackwardTransitions.filterEntries(maybeStates, maybeStates);
auto conversionEnd = std::chrono::high_resolution_clock::now();
std::chrono::high_resolution_clock::time_point modelCheckingStart = std::chrono::high_resolution_clock::now();

10
src/storm/modelchecker/region/ApproximationModel.cpp

@ -33,7 +33,7 @@ namespace storm {
this->computeRewards=true;
STORM_LOG_THROW(this->typeOfParametricModel==storm::models::ModelType::Dtmc, storm::exceptions::InvalidArgumentException, "Approximation with rewards is only implemented for Dtmcs");
STORM_LOG_THROW(parametricModel.hasUniqueRewardModel(), storm::exceptions::InvalidArgumentException, "The rewardmodel of the approximation model should be unique");
STORM_LOG_THROW(parametricModel.getUniqueRewardModel().hasOnlyStateRewards(), storm::exceptions::InvalidArgumentException, "The rewardmodel of the approximation model should have state rewards only");
STORM_LOG_THROW(parametricModel.getUniqueRewardModel().hasStateActionRewards() && !parametricModel.getUniqueRewardModel().hasStateRewards() && !parametricModel.getUniqueRewardModel().hasTransitionRewards(), storm::exceptions::InvalidArgumentException, "The rewardmodel of the approximation model should have state action rewards only");
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Invalid formula: " << formula << ". Approximation model only supports eventually or reachability reward formulae.");
}
@ -205,13 +205,13 @@ namespace storm {
STORM_LOG_THROW(this->vectorData.vector.size()==this->matrixData.matrix.getRowCount(), storm::exceptions::UnexpectedException, "The size of the eq-sys vector does not match to the number of rows in the eq-sys matrix");
this->vectorData.assignment.reserve(vectorData.vector.size());
// run through the state reward vector of the parametric model.
// run through the state action reward vector of the parametric model.
// Constant entries can be set directly.
// For Parametric entries we set a dummy value and insert the corresponding function and the assignment
ConstantType dummyNonZeroValue = storm::utility::one<ConstantType>();
auto vectorIt = this->vectorData.vector.begin();
for(auto oldState : this->maybeStates){
if(storm::utility::isConstant(parametricModel.getUniqueRewardModel().getStateRewardVector()[oldState])){
if(storm::utility::isConstant(parametricModel.getUniqueRewardModel().getStateActionRewardVector()[oldState])){
ConstantType reward = storm::utility::region::convertNumber<ConstantType>(storm::utility::region::getConstantPart(parametricModel.getUniqueRewardModel().getStateRewardVector()[oldState]));
//Add one of these entries for every row in the row group of oldState
for(auto matrixRow=this->matrixData.matrix.getRowGroupIndices()[oldState]; matrixRow<this->matrixData.matrix.getRowGroupIndices()[oldState+1]; ++matrixRow){
@ -220,7 +220,7 @@ namespace storm {
}
} else {
std::set<VariableType> occurringRewVariables;
storm::utility::region::gatherOccurringVariables(parametricModel.getUniqueRewardModel().getStateRewardVector()[oldState], occurringRewVariables);
storm::utility::region::gatherOccurringVariables(parametricModel.getUniqueRewardModel().getStateActionRewardVector()[oldState], occurringRewVariables);
// For each row in the row group of oldState, we get the corresponding substitution and insert the FunctionSubstitution
for(auto matrixRow=this->matrixData.matrix.getRowGroupIndices()[oldState]; matrixRow<this->matrixData.matrix.getRowGroupIndices()[oldState+1]; ++matrixRow){
//Extend the substitution for the probabilities with the reward parameters
@ -230,7 +230,7 @@ namespace storm {
substitution.insert(typename std::map<VariableType, RegionBoundary>::value_type(rewardVar, RegionBoundary::UNSPECIFIED));
}
// insert the FunctionSubstitution
auto functionsIt = this->funcSubData.functions.insert(FunctionEntry(FunctionSubstitution(parametricModel.getUniqueRewardModel().getStateRewardVector()[oldState], this->matrixData.rowSubstitutions[matrixRow]), dummyNonZeroValue)).first;
auto functionsIt = this->funcSubData.functions.insert(FunctionEntry(FunctionSubstitution(parametricModel.getUniqueRewardModel().getStateActionRewardVector()[oldState], this->matrixData.rowSubstitutions[matrixRow]), dummyNonZeroValue)).first;
//insert assignment and dummy data
this->vectorData.assignment.emplace_back(vectorIt, functionsIt->second);
*vectorIt = dummyNonZeroValue;

4
src/storm/modelchecker/region/SamplingModel.cpp

@ -43,7 +43,7 @@ namespace storm {
this->computeRewards=true;
STORM_LOG_THROW(this->typeOfParametricModel==storm::models::ModelType::Dtmc, storm::exceptions::InvalidArgumentException, "Sampling with rewards is only implemented for Dtmcs");
STORM_LOG_THROW(parametricModel.hasUniqueRewardModel(), storm::exceptions::InvalidArgumentException, "The rewardmodel of the sampling model should be unique");
STORM_LOG_THROW(parametricModel.getUniqueRewardModel().hasOnlyStateRewards(), storm::exceptions::InvalidArgumentException, "The rewardmodel of the sampling model should have state rewards only");
STORM_LOG_THROW(parametricModel.getUniqueRewardModel().hasStateActionRewards() && !parametricModel.getUniqueRewardModel().hasStateRewards() && !parametricModel.getUniqueRewardModel().hasTransitionRewards(), storm::exceptions::InvalidArgumentException, "The rewardmodel of the sempling model should have state action rewards only");
STORM_LOG_THROW(formula->getSubformula().isEventuallyFormula(), storm::exceptions::InvalidArgumentException, "The subformula should be a reachabilityreward formula");
STORM_LOG_THROW(formula->getSubformula().asEventuallyFormula().getSubformula().isInFragment(storm::logic::propositional()), storm::exceptions::InvalidArgumentException, "The subsubformula should be a propositional formula");
} else {
@ -131,7 +131,7 @@ namespace storm {
std::vector<ConstantType> b;
if(this->computeRewards){
b.resize(submatrix.getRowCount());
storm::utility::vector::selectVectorValues(b, this->maybeStates, instantiatedModel.getUniqueRewardModel().getStateRewardVector());
storm::utility::vector::selectVectorValues(b, this->maybeStates, instantiatedModel.getUniqueRewardModel().getStateActionRewardVector());
} else {
b = instantiatedModel.getTransitionMatrix().getConstrainedRowSumVector(this->maybeStates, this->targetStates);
}

32
src/storm/models/sparse/StandardRewardModel.cpp

@ -232,6 +232,34 @@ namespace storm {
return result;
}
template<typename ValueType>
template<typename MatrixValueType>
storm::storage::BitVector StandardRewardModel<ValueType>::getStatesWithZeroReward(storm::storage::SparseMatrix<MatrixValueType> const& transitionMatrix) const {
storm::storage::BitVector result = this->hasStateRewards() ? storm::utility::vector::filterZero(this->getStateRewardVector()) : storm::storage::BitVector(transitionMatrix.getRowGroupCount(), true);
if (this->hasStateActionRewards()) {
for (uint_fast64_t state = 0; state < transitionMatrix.getRowGroupCount(); ++state) {
for (uint_fast64_t row = transitionMatrix.getRowGroupIndices()[state]; row < transitionMatrix.getRowGroupIndices()[state+1]; ++row) {
if(!storm::utility::isZero(this->getStateActionRewardVector()[row])) {
result.set(state, false);
break;
}
}
}
}
if (this->hasTransitionRewards()) {
for (uint_fast64_t state = 0; state < transitionMatrix.getRowGroupCount(); ++state) {
for (auto const& rewardMatrixEntry : this->getTransitionRewardMatrix().getRowGroup(state)) {
if(!storm::utility::isZero(rewardMatrixEntry.getValue())) {
result.set(state, false);
break;
}
}
}
}
return result;
}
template<typename ValueType>
void StandardRewardModel<ValueType>::setStateActionRewardValue(uint_fast64_t row, ValueType const& value) {
this->optionalStateActionRewardVector.get()[row] = value;
@ -304,6 +332,7 @@ namespace storm {
template std::vector<double> StandardRewardModel<double>::getTotalRewardVector(uint_fast64_t numberOfRows, storm::storage::SparseMatrix<double> const& transitionMatrix, storm::storage::BitVector const& filter) const;
template std::vector<double> StandardRewardModel<double>::getTotalRewardVector(storm::storage::SparseMatrix<double> const& transitionMatrix, std::vector<double> const& weights, bool scaleTransAndActions) const;
template std::vector<double> StandardRewardModel<double>::getTotalActionRewardVector(storm::storage::SparseMatrix<double> const& transitionMatrix, std::vector<double> const& stateRewardWeights) const;
template storm::storage::BitVector StandardRewardModel<double>::getStatesWithZeroReward(storm::storage::SparseMatrix<double> const& transitionMatrix) const;
template void StandardRewardModel<double>::reduceToStateBasedRewards(storm::storage::SparseMatrix<double> const& transitionMatrix, bool reduceToStateRewards);
template void StandardRewardModel<double>::setStateActionReward(uint_fast64_t choiceIndex, double const & newValue);
template void StandardRewardModel<double>::setStateReward(uint_fast64_t state, double const & newValue);
@ -325,6 +354,7 @@ namespace storm {
template std::vector<storm::RationalNumber> StandardRewardModel<storm::RationalNumber>::getTotalRewardVector(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix) const;
template std::vector<storm::RationalNumber> StandardRewardModel<storm::RationalNumber>::getTotalRewardVector(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix, std::vector<storm::RationalNumber> const& weights, bool scaleTransAndActions) const;
template std::vector<storm::RationalNumber> StandardRewardModel<storm::RationalNumber>::getTotalActionRewardVector(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix, std::vector<storm::RationalNumber> const& stateRewardWeights) const;
template storm::storage::BitVector StandardRewardModel<storm::RationalNumber>::getStatesWithZeroReward(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix) const;
template void StandardRewardModel<storm::RationalNumber>::reduceToStateBasedRewards(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix, bool reduceToStateRewards);
template void StandardRewardModel<storm::RationalNumber>::setStateActionReward(uint_fast64_t choiceIndex, storm::RationalNumber const & newValue);
template void StandardRewardModel<storm::RationalNumber>::setStateReward(uint_fast64_t state, storm::RationalNumber const & newValue);
@ -334,6 +364,8 @@ namespace storm {
template std::vector<storm::RationalFunction> StandardRewardModel<storm::RationalFunction>::getTotalRewardVector(uint_fast64_t numberOfRows, storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix, storm::storage::BitVector const& filter) const;
template std::vector<storm::RationalFunction> StandardRewardModel<storm::RationalFunction>::getTotalRewardVector(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix) const;
template std::vector<storm::RationalFunction> StandardRewardModel<storm::RationalFunction>::getTotalRewardVector(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix, std::vector<storm::RationalFunction> const& weights, bool scaleTransAndActions) const;
template storm::storage::BitVector StandardRewardModel<storm::RationalFunction>::getStatesWithZeroReward(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix) const;
template std::vector<storm::RationalFunction> StandardRewardModel<storm::RationalFunction>::getTotalActionRewardVector(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix, std::vector<storm::RationalFunction> const& stateRewardWeights) const;
template void StandardRewardModel<storm::RationalFunction>::reduceToStateBasedRewards(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix, bool reduceToStateRewards);
template void StandardRewardModel<storm::RationalFunction>::setStateActionReward(uint_fast64_t choiceIndex, storm::RationalFunction const & newValue);

9
src/storm/models/sparse/StandardRewardModel.h

@ -230,6 +230,15 @@ namespace storm {
template<typename MatrixValueType>
std::vector<ValueType> getTotalActionRewardVector(storm::storage::SparseMatrix<MatrixValueType> const& transitionMatrix, std::vector<MatrixValueType> const& stateRewardWeights) const;
/*!
* Returns the set of states at which a all rewards (state-, action- and transition-rewards) are zero.
*
* @param transitionMatrix the transition matrix of the model (used to determine the actions and transitions that belong to a state)
* @ return a vector representing all states at which the reward is zero
*/
template<typename MatrixValueType>
storm::storage::BitVector getStatesWithZeroReward(storm::storage::SparseMatrix<MatrixValueType> const& transitionMatrix) const;
/*!
* Sets the given value in the state-action reward vector at the given row. This assumes that the reward
* model has state-action rewards.

93
src/storm/storage/FlexibleSparseMatrix.cpp

@ -109,16 +109,12 @@ namespace storm {
template<typename ValueType>
typename FlexibleSparseMatrix<ValueType>::index_type FlexibleSparseMatrix<ValueType>::getRowGroupCount() const {
return rowGroupIndices.size();
return rowGroupIndices.size() - 1;
}
template<typename ValueType>
typename FlexibleSparseMatrix<ValueType>::index_type FlexibleSparseMatrix<ValueType>::getRowGroupSize(index_type group) const {
if (group == getRowGroupCount() - 1) {
return getRowCount() - rowGroupIndices[group];
} else {
return rowGroupIndices[group + 1] - rowGroupIndices[group];
}
return rowGroupIndices[group + 1] - rowGroupIndices[group];
}
template<typename ValueType>
@ -157,9 +153,9 @@ namespace storm {
bool FlexibleSparseMatrix<ValueType>::hasTrivialRowGrouping() const {
return trivialRowGrouping;
}
template<typename ValueType>
void FlexibleSparseMatrix<ValueType>::createSubmatrix(storm::storage::BitVector const& rowConstraint, storm::storage::BitVector const& columnConstraint) {
void FlexibleSparseMatrix<ValueType>::filterEntries(storm::storage::BitVector const& rowConstraint, storm::storage::BitVector const& columnConstraint) {
for (uint_fast64_t rowIndex = 0; rowIndex < this->data.size(); ++rowIndex) {
auto& row = this->data[rowIndex];
if (!rowConstraint.get(rowIndex)) {
@ -179,12 +175,83 @@ namespace storm {
template<typename ValueType>
storm::storage::SparseMatrix<ValueType> FlexibleSparseMatrix<ValueType>::createSparseMatrix() {
storm::storage::SparseMatrixBuilder<ValueType> matrixBuilder(getRowCount(), getColumnCount());
for (uint_fast64_t rowIndex = 0; rowIndex < getRowCount(); ++rowIndex) {
auto& row = this->data[rowIndex];
uint_fast64_t numEntries = 0;
for (auto const& row : this->data) {
numEntries += row.size();
}
storm::storage::SparseMatrixBuilder<ValueType> matrixBuilder(getRowCount(), getColumnCount(), numEntries, hasTrivialRowGrouping(), hasTrivialRowGrouping() ? 0 : getRowGroupCount());
uint_fast64_t currRowIndex = 0;
auto rowGroupIndexIt = getRowGroupIndices().begin();
for (auto const& row : this->data) {
if(!hasTrivialRowGrouping()) {
while (currRowIndex >= *rowGroupIndexIt) {
matrixBuilder.newRowGroup(currRowIndex);
++rowGroupIndexIt;
}
}
for (auto const& entry : row) {
matrixBuilder.addNextValue(rowIndex, entry.getColumn(), entry.getValue());
matrixBuilder.addNextValue(currRowIndex, entry.getColumn(), entry.getValue());
}
++currRowIndex;
}
// The matrix might end with one or more empty row groups
if(!hasTrivialRowGrouping()) {
while (currRowIndex >= *rowGroupIndexIt) {
matrixBuilder.newRowGroup(currRowIndex);
++rowGroupIndexIt;
}
}
return matrixBuilder.build();
}
template<typename ValueType>
storm::storage::SparseMatrix<ValueType> FlexibleSparseMatrix<ValueType>::createSparseMatrix(storm::storage::BitVector const& rowConstraint, storm::storage::BitVector const& columnConstraint) {
uint_fast64_t numEntries = 0;
for (auto const& rowIndex : rowConstraint) {
auto const& row = data[rowIndex];
for(auto const& entry : row) {
if (columnConstraint.get(entry.getColumn())) {
++numEntries;
}
}
}
uint_fast64_t numRowGroups = 0;
if (!hasTrivialRowGrouping()) {
auto lastRowGroupIndexIt = getRowGroupIndices().end() - 1;
auto rowGroupIndexIt = getRowGroupIndices().begin();
while (rowGroupIndexIt != lastRowGroupIndexIt) {
// Check whether the rowGroup will be nonempty
if(rowConstraint.getNextSetIndex(*rowGroupIndexIt) < *(++rowGroupIndexIt)) {
++numRowGroups;
}
}
}
std::vector<uint_fast64_t> oldToNewColumnIndexMapping(getColumnCount(), getColumnCount());
uint_fast64_t newColumnIndex = 0;
for (auto const& oldColumnIndex : columnConstraint) {
oldToNewColumnIndexMapping[oldColumnIndex] = newColumnIndex++;
}
storm::storage::SparseMatrixBuilder<ValueType> matrixBuilder(rowConstraint.getNumberOfSetBits(), newColumnIndex, numEntries, hasTrivialRowGrouping(), numRowGroups);
uint_fast64_t currRowIndex = 0;
auto rowGroupIndexIt = getRowGroupIndices().begin();
for (auto const& oldRowIndex : rowConstraint) {
if(!hasTrivialRowGrouping() && oldRowIndex >= *rowGroupIndexIt) {
matrixBuilder.newRowGroup(currRowIndex);
// Skip empty row groups
do {
++rowGroupIndexIt;
} while (currRowIndex >= *rowGroupIndexIt);
}
auto const& row = data[oldRowIndex];
for (auto const& entry : row) {
if(columnConstraint.get(entry.getColumn())) {
matrixBuilder.addNextValue(currRowIndex, oldToNewColumnIndexMapping[entry.getColumn()], entry.getValue());
}
}
++currRowIndex;
}
return matrixBuilder.build();
}
@ -237,7 +304,7 @@ namespace storm {
FlexibleIndex rowGroupCount = matrix.getRowGroupCount();
for (FlexibleIndex rowGroup = 0; rowGroup < rowGroupCount; ++rowGroup) {
out << "\t---- group " << rowGroup << "/" << (rowGroupCount - 1) << " ---- " << std::endl;
FlexibleIndex endRow = rowGroup+1 < rowGroupCount ? matrix.rowGroupIndices[rowGroup + 1] : matrix.getRowCount();
FlexibleIndex endRow = matrix.rowGroupIndices[rowGroup + 1];
// Iterate over all rows.
for (FlexibleIndex row = matrix.rowGroupIndices[rowGroup]; row < endRow; ++row) {
// Print the actual row.

24
src/storm/storage/FlexibleSparseMatrix.h

@ -158,19 +158,31 @@ namespace storm {
bool hasTrivialRowGrouping() const;
/*!
* Creates a submatrix of the current matrix in place by dropping all rows and columns whose bits are not
* set to one in the given bit vector.
* Erases all entries whose row and column does not satisfy the given rowConstraint and the given columnConstraint
*
* @param rowConstraint A bit vector indicating which rows to keep.
* @param columnConstraint A bit vector indicating which columns to keep.
* @param rowConstraint A bit vector indicating which row entries to keep.
* @param columnConstraint A bit vector indicating which column entries to keep.
*/
void createSubmatrix(storm::storage::BitVector const& rowConstraint, storm::storage::BitVector const& columnConstraint);
void filterEntries(storm::storage::BitVector const& rowConstraint, storm::storage::BitVector const& columnConstraint);
/*!
* Creates a sparse matrix from the flexible sparse matrix.
* @return The sparse matrix.
*/
storm::storage::SparseMatrix<ValueType> createSparseMatrix();
/*!
* Creates a sparse matrix from the flexible sparse matrix.
* Only the selected rows and columns will be considered.
* Empty rowGroups will be ignored
*
* @param rowConstraint A bit vector indicating which rows to keep.
* @param columnConstraint A bit vector indicating which columns to keep.
*
* @return The sparse matrix.
*/
storm::storage::SparseMatrix<ValueType> createSparseMatrix(storm::storage::BitVector const& rowConstraint, storm::storage::BitVector const& columnConstraint);
/*!
* Checks whether the given state has a self-loop with an arbitrary probability in the probability matrix.

190
src/storm/transformer/GoalStateMerger.cpp

@ -0,0 +1,190 @@
#include "storm/transformer/GoalStateMerger.h"
#include <limits>
#include <memory>
#include "storm/utility/constants.h"
#include "storm/utility/macros.h"
#include "storm/utility/vector.h"
#include "storm/models/sparse/Dtmc.h"
#include "storm/models/sparse/Mdp.h"
#include "storm/models/sparse/StandardRewardModel.h"
#include "storm/exceptions/InvalidArgumentException.h"
namespace storm {
namespace transformer {
template <typename SparseModelType>
GoalStateMerger<SparseModelType>::GoalStateMerger(SparseModelType const& model) : originalModel(model) {
// Intentionally left empty
}
template <typename SparseModelType>
std::shared_ptr<SparseModelType> GoalStateMerger<SparseModelType>::mergeTargetAndSinkStates(storm::storage::BitVector const& maybeStates, storm::storage::BitVector& targetStates, storm::storage::BitVector& sinkStates, std::vector<std::string> const& selectedRewardModels) {
STORM_LOG_THROW(maybeStates.isDisjointFrom(targetStates) && targetStates.isDisjointFrom(sinkStates) && sinkStates.isDisjointFrom(maybeStates), storm::exceptions::InvalidArgumentException, "maybestates, targetstates, and sinkstates are assumed to be disjoint when creating the submodel. However, this is not the case.");
boost::optional<uint_fast64_t> targetState, sinkState;
auto builder = initializeTransitionMatrixBuilder(maybeStates, targetStates, sinkStates, targetState, sinkState);
auto transitionMatrix = buildTransitionMatrix(maybeStates, targetStates, sinkStates, targetState, sinkState, builder);
uint_fast64_t resNumStates = transitionMatrix.getRowGroupCount();
// Get the labeling for the initial states
storm::storage::BitVector initialStates = originalModel.getInitialStates() % maybeStates;
initialStates.resize(resNumStates, false);
if(!originalModel.getInitialStates().isDisjointFrom(targetStates)) {
initialStates.set(*targetState, true);
}
if(!originalModel.getInitialStates().isDisjointFrom(sinkStates)) {
initialStates.set(*sinkState, true);
}
storm::models::sparse::StateLabeling labeling(resNumStates);
labeling.addLabel("init", std::move(initialStates));
// Get the reward models
std::unordered_map<std::string, typename SparseModelType::RewardModelType> rewardModels;
for (auto rewardModelName : selectedRewardModels) {
auto origTotalRewards = originalModel.getRewardModel(rewardModelName).getTotalRewardVector(originalModel.getTransitionMatrix());
auto resTotalRewards = storm::utility::vector::filterVector(origTotalRewards, maybeStates);
resTotalRewards.resize(resNumStates, storm::utility::zero<typename SparseModelType::RewardModelType::ValueType>());
rewardModels.insert(std::make_pair(rewardModelName, typename SparseModelType::RewardModelType(boost::none, resTotalRewards)));
}
// modify the given target and sink states
targetStates = storm::storage::BitVector(resNumStates, false);
if(targetState) {
targetStates.set(*targetState, true);
}
sinkStates = storm::storage::BitVector(resNumStates, false);
if(sinkState) {
sinkStates.set(*sinkState, true);
}
// Return the result
return std::make_shared<SparseModelType>(std::move(transitionMatrix), std::move(labeling), std::move(rewardModels));
}
template <typename SparseModelType>
storm::storage::SparseMatrixBuilder<typename SparseModelType::ValueType> GoalStateMerger<SparseModelType>::initializeTransitionMatrixBuilder(storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& targetStates, storm::storage::BitVector const& sinkStates, boost::optional<uint_fast64_t>& newTargetState, boost::optional<uint_fast64_t>& newSinkState) {
storm::storage::SparseMatrix<typename SparseModelType::ValueType> const& origMatrix = originalModel.getTransitionMatrix();
// Get the number of rows, cols and entries that the resulting transition matrix will have.
uint_fast64_t resNumStates(maybeStates.getNumberOfSetBits()), resNumActions(0), resNumTransitions(0);
bool targetStateRequired = !originalModel.getInitialStates().isDisjointFrom(targetStates);
bool sinkStateRequired = !originalModel.getInitialStates().isDisjointFrom(targetStates);
for( auto state : maybeStates) {
resNumActions += origMatrix.getRowGroupSize(state);
bool hasTransitionToTarget(false), hasTransitionToSink(false);
auto const& endOfRowGroup = origMatrix.getRowGroupIndices()[state+1];
for (uint_fast64_t row = origMatrix.getRowGroupIndices()[state]; row < endOfRowGroup; ++row) {
for (auto const& entry : origMatrix.getRow(row)) {
if(maybeStates.get(entry.getColumn())) {
++resNumTransitions;
} else if (targetStates.get(entry.getColumn())) {
hasTransitionToTarget = true;
} else if (sinkStates.get(entry.getColumn())) {
hasTransitionToSink = true;
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "There is a transition originating from a maybestate that does not lead to a maybe-, target-, or sinkstate.");
}
}
if(hasTransitionToTarget) {
++resNumTransitions;
targetStateRequired = true;
}
if(hasTransitionToSink) {
++resNumTransitions;
sinkStateRequired = true;
}
}
}
// Get the index of the target/ sink state in the resulting model (if these states will exist)
if(targetStateRequired) {
newTargetState = resNumStates;
++resNumStates;
++resNumActions;
++resNumTransitions;
}
if(sinkStateRequired) {
newSinkState = resNumStates;
++resNumStates;
++resNumActions;
++resNumTransitions;
}
return storm::storage::SparseMatrixBuilder<typename SparseModelType::ValueType>(resNumActions, resNumStates, resNumTransitions, true, true, resNumStates);
}
template <typename SparseModelType>
storm::storage::SparseMatrix<typename SparseModelType::ValueType> GoalStateMerger<SparseModelType>::buildTransitionMatrix(storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& targetStates, storm::storage::BitVector const& sinkStates, boost::optional<uint_fast64_t> const& newTargetState, boost::optional<uint_fast64_t>const& newSinkState, storm::storage::SparseMatrixBuilder<typename SparseModelType::ValueType>& builder) {
// Get a Mapping that yields for each column in the old matrix the corresponding column in the new matrix
std::vector<uint_fast64_t> oldToNewIndexMap(maybeStates.size(), std::numeric_limits<uint_fast64_t>::max()); // init with some invalid state
uint_fast64_t newStateIndex = 0;
for (auto maybeState : maybeStates) {
oldToNewIndexMap[maybeState] = newStateIndex;
++newStateIndex;
}
// Build the transition matrix
storm::storage::SparseMatrix<typename SparseModelType::ValueType> const& origMatrix = originalModel.getTransitionMatrix();
uint_fast64_t currRow = 0;
for (auto state : maybeStates) {
builder.newRowGroup(currRow);
boost::optional<typename SparseModelType::ValueType> targetProbability, sinkProbability;
auto const& endOfRowGroup = origMatrix.getRowGroupIndices()[state+1];
for (uint_fast64_t row = origMatrix.getRowGroupIndices()[state]; row < endOfRowGroup; ++row) {
for (auto const& entry : origMatrix.getRow(row)) {
if(maybeStates.get(entry.getColumn())) {
builder.addNextValue(currRow, oldToNewIndexMap[entry.getColumn()], entry.getValue());
} else if (targetStates.get(entry.getColumn())) {
targetProbability = targetProbability.is_initialized() ? *targetProbability + entry.getValue() : entry.getValue();
} else if (sinkStates.get(entry.getColumn())) {
sinkProbability = sinkProbability.is_initialized() ? *sinkProbability + entry.getValue() : entry.getValue();
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "There is a transition originating from a maybestate that does not lead to a maybe-, target-, or sinkstate.");
}
}
if(targetProbability) {
assert(newTargetState);
builder.addNextValue(currRow, *newTargetState, storm::utility::simplify(*targetProbability));
}
if(sinkProbability) {
assert(newSinkState);
builder.addNextValue(currRow, *newSinkState, storm::utility::simplify(*sinkProbability));
}
++currRow;
}
}
// Add the selfloops at target and sink
if(newTargetState) {
builder.newRowGroup(currRow);
builder.addNextValue(currRow, *newTargetState, storm::utility::one<typename SparseModelType::ValueType>());
++currRow;
}
if(newSinkState) {
builder.newRowGroup(currRow);
builder.addNextValue(currRow, *newSinkState, storm::utility::one<typename SparseModelType::ValueType>());
++currRow;
}
return builder.build();
}
template class GoalStateMerger<storm::models::sparse::Dtmc<double>>;
template class GoalStateMerger<storm::models::sparse::Mdp<double>>;
template class GoalStateMerger<storm::models::sparse::Dtmc<storm::RationalNumber>>;
template class GoalStateMerger<storm::models::sparse::Mdp<storm::RationalNumber>>;
template class GoalStateMerger<storm::models::sparse::Dtmc<storm::RationalFunction>>;
template class GoalStateMerger<storm::models::sparse::Mdp<storm::RationalFunction>>;
}
}

161
src/storm/transformer/GoalStateMerger.h

@ -1,15 +1,12 @@
#pragma once
#include <limits>
#include <memory>
#include <string>
#include <vector>
#include <boost/optional.hpp>
#include "storm/utility/constants.h"
#include "storm/utility/macros.h"
#include "storm/models/sparse/Dtmc.h"
#include "storm/models/sparse/Mdp.h"
#include "storm/exceptions/InvalidArgumentException.h"
#include "storm/storage/BitVector.h"
#include "storm/storage/SparseMatrix.h"
namespace storm {
@ -21,7 +18,9 @@ namespace storm {
template <typename SparseModelType>
class GoalStateMerger {
public:
GoalStateMerger(SparseModelType const& model);
/* Computes a submodel of the specified model that only considers the states given by maybeStates as well as
* * one target state to which all transitions to a state selected by targetStates are redirected and
* * one sink state to which all transitions to a state selected by sinkStates are redirected.
@ -30,133 +29,31 @@ namespace storm {
* The two given bitvectors targetStates and sinkStates are modified such that they represent the corresponding state in the obtained submodel.
*
* Notes:
* * The resulting model will not have any rewardmodels, labels (other then "init") etc.
* * The resulting model will not have any labels (other then "init") and only the selected reward models.
* * Rewards are reduced to stateActionRewards.
* * The target and sink states will not get any reward
* * It is assumed that the given set of maybeStates can only be left via either a target or a sink state. Otherwise an exception is thrown.
* * It is assumed that maybeStates, targetStates, and sinkStates are all disjoint. Otherwise an exception is thrown.
*/
static std::shared_ptr<SparseModelType> mergeTargetAndSinkStates(SparseModelType const& model, storm::storage::BitVector const& maybeStates, storm::storage::BitVector& targetStates, storm::storage::BitVector& sinkStates) {
STORM_LOG_THROW(maybeStates.isDisjointFrom(targetStates) && targetStates.isDisjointFrom(sinkStates) && sinkStates.isDisjointFrom(maybeStates), storm::exceptions::InvalidArgumentException, "maybestates, targetstates, and sinkstates are assumed to be disjoint when creating the submodel. However, this is not the case.");
storm::storage::SparseMatrix<typename SparseModelType::ValueType> const& origMatrix = model.getTransitionMatrix();
// Get the number of rows, cols and entries that the resulting transition matrix will have.
uint_fast64_t resNumStates(maybeStates.getNumberOfSetBits()), resNumActions(0), resNumTransitions(0);
bool targetStateRequired = !model.getInitialStates().isDisjointFrom(targetStates);
bool sinkStateRequired = !model.getInitialStates().isDisjointFrom(targetStates);
for( auto state : maybeStates) {
resNumActions += origMatrix.getRowGroupSize(state);
bool hasTransitionToTarget(false), hasTransitionToSink(false);
auto const& endOfRowGroup = origMatrix.getRowGroupIndices()[state+1];
for (uint_fast64_t row = origMatrix.getRowGroupIndices()[state]; row < endOfRowGroup; ++row) {
for (auto const& entry : origMatrix.getRow(row)) {
if(maybeStates.get(entry.getColumn())) {
++resNumTransitions;
} else if (targetStates.get(entry.getColumn())) {
hasTransitionToTarget = true;
} else if (sinkStates.get(entry.getColumn())) {
hasTransitionToSink = true;
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "There is a transition originating from a maybestate that does not lead to a maybe-, target-, or sinkstate.");
}
}
if(hasTransitionToTarget) {
++resNumTransitions;
targetStateRequired = true;
}
if(hasTransitionToSink) {
++resNumTransitions;
sinkStateRequired = true;
}
}
}
// Get the index of the target/ sink state in the resulting model (if these states will exist)
uint_fast64_t targetState(std::numeric_limits<uint_fast64_t>::max()), sinkState(std::numeric_limits<uint_fast64_t>::max()); // init with some invalid state
if(targetStateRequired) {
targetState = resNumStates;
++resNumStates;
++resNumActions;
++resNumTransitions;
}
if(sinkStateRequired) {
sinkState = resNumStates;
++resNumStates;
++resNumActions;
++resNumTransitions;
}
// Get a Mapping that yields for each column in the old matrix the corresponding column in the new matrix
std::vector<uint_fast64_t> oldToNewIndexMap(maybeStates.size(), std::numeric_limits<uint_fast64_t>::max()); // init with some invalid state
uint_fast64_t newStateIndex = 0;
for (auto maybeState : maybeStates) {
oldToNewIndexMap[maybeState] = newStateIndex;
++newStateIndex;
}
// Build the transition matrix
storm::storage::SparseMatrixBuilder<typename SparseModelType::ValueType> builder(resNumActions, resNumStates, resNumTransitions, true, true, resNumStates);
uint_fast64_t currRow = 0;
for (auto state : maybeStates) {
builder.newRowGroup(currRow);
boost::optional<typename SparseModelType::ValueType> targetProbability, sinkProbability;
auto const& endOfRowGroup = origMatrix.getRowGroupIndices()[state+1];
for (uint_fast64_t row = origMatrix.getRowGroupIndices()[state]; row < endOfRowGroup; ++row) {
for (auto const& entry : origMatrix.getRow(row)) {
if(maybeStates.get(entry.getColumn())) {
builder.addNextValue(currRow, oldToNewIndexMap[entry.getColumn()], entry.getValue());
} else if (targetStates.get(entry.getColumn())) {
targetProbability = targetProbability.is_initialized() ? *targetProbability + entry.getValue() : entry.getValue();
} else if (sinkStates.get(entry.getColumn())) {
sinkProbability = sinkProbability.is_initialized() ? *sinkProbability + entry.getValue() : entry.getValue();
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "There is a transition originating from a maybestate that does not lead to a maybe-, target-, or sinkstate.");
}
}
if(targetProbability) {
builder.addNextValue(currRow, targetState, storm::utility::simplify(*targetProbability));
}
if(sinkProbability) {
builder.addNextValue(currRow, sinkState, storm::utility::simplify(*sinkProbability));
}
++currRow;
}
}
// Add the selfloops at target and sink
if(targetStateRequired) {
builder.newRowGroup(currRow);
builder.addNextValue(currRow, targetState, storm::utility::one<typename SparseModelType::ValueType>());
++currRow;
}
if(sinkStateRequired) {
builder.newRowGroup(currRow);
builder.addNextValue(currRow, sinkState, storm::utility::one<typename SparseModelType::ValueType>());
++currRow;
}
// Get the labeling for the initial states
storm::storage::BitVector initialStates = model.getInitialStates() % maybeStates;
initialStates.resize(resNumStates, false);
if(!model.getInitialStates().isDisjointFrom(targetStates)) {
initialStates.set(targetState, true);
}
if(!model.getInitialStates().isDisjointFrom(sinkStates)) {
initialStates.set(sinkState, true);
}
storm::models::sparse::StateLabeling labeling(resNumStates);
labeling.addLabel("init", std::move(initialStates));
// modify the given target and sink states
targetStates = storm::storage::BitVector(resNumStates, false);
if(targetStateRequired) {
targetStates.set(targetState, true);
}
sinkStates = storm::storage::BitVector(resNumStates, false);
if(sinkStateRequired) {
sinkStates.set(sinkState, true);
}
// Return the result
return std::make_shared<SparseModelType>(builder.build(), std::move(labeling));
}
std::shared_ptr<SparseModelType> mergeTargetAndSinkStates(storm::storage::BitVector const& maybeStates, storm::storage::BitVector& targetStates, storm::storage::BitVector& sinkStates, std::vector<std::string> const& selectedRewardModels = std::vector<std::string>());
private:
SparseModelType const& originalModel;
/*!
* Initializes the matrix builder for the transition matrix of the resulting model
*
* @param newTargetState Will be set to the index of the target state of the resulting model (Only if it has a target state)
* @param newSinkState Will be set to the index of the sink state of the resulting model (Only if it has a sink state)
*/
storm::storage::SparseMatrixBuilder<typename SparseModelType::ValueType> initializeTransitionMatrixBuilder(storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& targetStates, storm::storage::BitVector const& sinkStates, boost::optional<uint_fast64_t>& newTargetState, boost::optional<uint_fast64_t>& newSinkState);
/*!
* Builds the transition matrix of the resulting model
*/
storm::storage::SparseMatrix<typename SparseModelType::ValueType> buildTransitionMatrix(storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& targetStates, storm::storage::BitVector const& sinkStates, boost::optional<uint_fast64_t> const& newTargetState, boost::optional<uint_fast64_t> const& newSinkState, storm::storage::SparseMatrixBuilder<typename SparseModelType::ValueType>& builder);
};
}
}

77
src/storm/transformer/SparseParametricDtmcSimplifier.cpp

@ -9,6 +9,8 @@
#include "storm/transformer/GoalStateMerger.h"
#include "storm/utility/graph.h"
#include <storm/exceptions/NotSupportedException.h>
#include <storm/exceptions/UnexpectedException.h>
namespace storm {
namespace transformer {
@ -30,11 +32,12 @@ namespace storm {
storm::storage::BitVector psiStates = std::move(propositionalChecker.check(formula.getSubformula().asUntilFormula().getRightSubformula())->asExplicitQualitativeCheckResult().getTruthValuesVector());
std::pair<storm::storage::BitVector, storm::storage::BitVector> statesWithProbability01 = storm::utility::graph::performProb01(this->originalModel.getBackwardTransitions(), phiStates, psiStates);
// Only consider the maybestates that are reachable from one initial state without hopping over a target (i.e., prob1) state
storm::storage::BitVector reachableGreater0States = storm::utility::graph::getReachableStates(this->originalModel.getTransitionMatrix(), this->originalModel.getInitialStates(), ~statesWithProbability01.first, statesWithProbability01.second);
storm::storage::BitVector reachableGreater0States = storm::utility::graph::getReachableStates(this->originalModel.getTransitionMatrix(), this->originalModel.getInitialStates() & ~statesWithProbability01.first, ~statesWithProbability01.first, statesWithProbability01.second);
storm::storage::BitVector maybeStates = reachableGreater0States & ~statesWithProbability01.second;
// obtain the resulting subsystem
this->simplifiedModel = storm::transformer::GoalStateMerger<SparseModelType>::mergeTargetAndSinkStates(this->originalModel, maybeStates, statesWithProbability01.second, statesWithProbability01.first);
storm::transformer::GoalStateMerger<SparseModelType> goalStateMerger(this->originalModel);
this->simplifiedModel = goalStateMerger.mergeTargetAndSinkStates(maybeStates, statesWithProbability01.second, statesWithProbability01.first);
this->simplifiedModel->getStateLabeling().addLabel("target", statesWithProbability01.second);
this->simplifiedModel->getStateLabeling().addLabel("sink", statesWithProbability01.first);
@ -51,14 +54,76 @@ namespace storm {
template<typename SparseModelType>
bool SparseParametricDtmcSimplifier<SparseModelType>::simplifyForBoundedUntilProbabilities(storm::logic::ProbabilityOperatorFormula const& formula) {
// If this method was not overriden by any subclass, simplification is not possible
return false;
STORM_LOG_THROW(!formula.getSubformula().asBoundedUntilFormula().hasLowerBound(), storm::exceptions::NotSupportedException, "Lower step bounds are not supported.");
STORM_LOG_THROW(formula.getSubformula().asBoundedUntilFormula().hasUpperBound(), storm::exceptions::UnexpectedException, "Expected a bounded until formula with an upper bound.");
STORM_LOG_THROW(formula.getSubformula().asBoundedUntilFormula().hasIntegerUpperBound(), storm::exceptions::UnexpectedException, "Expected a bounded until formula with an integral upper bound.");
// Get the prob0, target, and the maybeStates
storm::modelchecker::SparsePropositionalModelChecker<SparseModelType> propositionalChecker(this->originalModel);
if(!propositionalChecker.canHandle(formula.getSubformula().asBoundedUntilFormula().getLeftSubformula()) || !propositionalChecker.canHandle(formula.getSubformula().asBoundedUntilFormula().getRightSubformula())) {
STORM_LOG_DEBUG("Can not simplify when Until-formula has non-propositional subformula(s). Formula: " << formula);
return false;
}
storm::storage::BitVector phiStates = std::move(propositionalChecker.check(formula.getSubformula().asBoundedUntilFormula().getLeftSubformula())->asExplicitQualitativeCheckResult().getTruthValuesVector());
storm::storage::BitVector psiStates = std::move(propositionalChecker.check(formula.getSubformula().asBoundedUntilFormula().getRightSubformula())->asExplicitQualitativeCheckResult().getTruthValuesVector());
storm::storage::BitVector probGreater0States = storm::utility::graph::performProbGreater0(this->originalModel.getBackwardTransitions(), phiStates, psiStates, true, formula.getSubformula().asBoundedUntilFormula().getUpperBound().evaluateAsInt());
storm::storage::BitVector prob0States = ~probGreater0States;
// Only consider the maybestates that are reachable from one initial probGreater0 state within the given amount of steps and without hopping over a target state
storm::storage::BitVector reachableGreater0States = storm::utility::graph::getReachableStates(this->originalModel.getTransitionMatrix(), this->originalModel.getInitialStates() & probGreater0States, probGreater0States, psiStates, true, formula.getSubformula().asBoundedUntilFormula().getUpperBound().evaluateAsInt());
storm::storage::BitVector maybeStates = reachableGreater0States & ~psiStates;
// obtain the resulting subsystem
storm::transformer::GoalStateMerger<SparseModelType> goalStateMerger(this->originalModel);
this->simplifiedModel = goalStateMerger.mergeTargetAndSinkStates(maybeStates, prob0States, psiStates);
this->simplifiedModel->getStateLabeling().addLabel("target", psiStates);
this->simplifiedModel->getStateLabeling().addLabel("sink", prob0States);
// obtain the simplified formula for the simplified model
auto labelFormula = std::make_shared<storm::logic::AtomicLabelFormula const> ("target");
auto boundedUntilFormula = std::make_shared<storm::logic::BoundedUntilFormula const>(storm::logic::Formula::getTrueFormula(), labelFormula, boost::none, storm::logic::);
this->simplifiedFormula = std::make_shared<storm::logic::ProbabilityOperatorFormula const>(eventuallyFormula, formula.getOperatorInformation());
return true;
}
template<typename SparseModelType>
bool SparseParametricDtmcSimplifier<SparseModelType>::simplifyForReachabilityRewards(storm::logic::RewardOperatorFormula const& formula) {
// If this method was not overriden by any subclass, simplification is not possible
return false;
typename SparseModelType::RewardModelType const& originalRewardModel = formula.hasRewardModelName() ? this->originalModel.getRewardModel(formula.getRewardModelName()) : this->originalModel.getUniqueRewardModel();
// Get the prob1 and the maybeStates
storm::modelchecker::SparsePropositionalModelChecker<SparseModelType> propositionalChecker(this->originalModel);
if(!propositionalChecker.canHandle(formula.getSubformula().asEventuallyFormula().getSubformula())) {
STORM_LOG_DEBUG("Can not simplify when reachability reward formula has non-propositional subformula(s). Formula: " << formula);
return false;
}
storm::storage::BitVector targetStates = std::move(propositionalChecker.check(formula.getSubformula().asEventuallyFormula().getSubformula())->asExplicitQualitativeCheckResult().getTruthValuesVector());
// The set of target states can be extended by the states that reach target with probability 1 without collecting any reward
targetStates = storm::utility::graph::performProb1(this->originalModel.getBackwardTransitions(), originalRewardModel.getStatesWithZeroReward(this->originalModel.getTransitionMatrix()), targetStates);
storm::storage::BitVector statesWithProb1 = storm::utility::graph::performProb1(this->originalModel.getBackwardTransitions(), storm::storage::BitVector(this->originalModel.getNumberOfStates(), true), targetStates);
storm::storage::BitVector infinityStates = ~statesWithProb1;
// Only consider the states that are reachable from an initial state without hopping over a target state
storm::storage::BitVector reachableStates = storm::utility::graph::getReachableStates(this->originalModel.getTransitionMatrix(), this->originalModel.getInitialStates() & statesWithProb1, statesWithProb1, targetStates);
storm::storage::BitVector maybeStates = reachableStates & ~targetStates;
targetStates = targetStates & reachableStates;
// obtain the resulting subsystem
std::vector<std::string> rewardModelNameAsVector(1, formula.hasRewardModelName() ? formula.getRewardModelName() : this->originalModel.getRewardModels().begin()->first);
storm::transformer::GoalStateMerger<SparseModelType> goalStateMerger(this->originalModel);
this->simplifiedModel = goalStateMerger.mergeTargetAndSinkStates(maybeStates, targetStates, infinityStates, rewardModelNameAsVector);
this->simplifiedModel->getStateLabeling().addLabel("target", targetStates);
this->simplifiedModel->getStateLabeling().addLabel("sink", infinityStates);
// obtain the simplified formula for the simplified model
auto labelFormula = std::make_shared<storm::logic::AtomicLabelFormula const> ("target");
auto eventuallyFormula = std::make_shared<storm::logic::EventuallyFormula const>(labelFormula, storm::logic::FormulaContext::Reward);
this->simplifiedFormula = std::make_shared<storm::logic::RewardOperatorFormula const>(eventuallyFormula, rewardModelNameAsVector.front(), formula.getOperatorInformation(), storm::logic::RewardMeasureType::Expectation);
// Eliminate all states for which all outgoing transitions are constant
this->simplifiedModel = this->eliminateConstantDeterministicStates(*this->simplifiedModel);
return true;
}
template<typename SparseModelType>

6
src/storm/transformer/SparseParametricModelSimplifier.cpp

@ -10,8 +10,10 @@
#include "storm/solver/stateelimination/PrioritizedStateEliminator.h"
#include "storm/storage/FlexibleSparseMatrix.h"
#include "storm/utility/vector.h"
#include "storm/exceptions/InvalidStateException.h"
#include "storm/exceptions/NotImplementedException.h"
#include "storm/exceptions/InvalidPropertyException.h"
namespace storm {
namespace transformer {
@ -37,6 +39,7 @@ namespace storm {
}
} else if (formula.isRewardOperatorFormula()) {
storm::logic::RewardOperatorFormula rewOpForm = formula.asRewardOperatorFormula();
STORM_LOG_THROW((rewOpForm.hasRewardModelName() && originalModel.hasRewardModel(rewOpForm.getRewardModelName())) || (!rewOpForm.hasRewardModelName() && originalModel.hasUniqueRewardModel()), storm::exceptions::InvalidPropertyException, "The reward model specified by formula " << formula << " is not available in the given model.");
if (rewOpForm.getSubformula().isReachabilityRewardFormula()) {
return simplifyForReachabilityRewards(rewOpForm);
} else if (rewOpForm.getSubformula().isCumulativeRewardFormula()) {
@ -141,7 +144,6 @@ namespace storm {
storm::solver::stateelimination::PrioritizedStateEliminator<typename SparseModelType::ValueType> stateEliminator(flexibleMatrix, flexibleBackwardTransitions, statesToEliminate, stateValues);
stateEliminator.eliminateAll();
flexibleMatrix.createSubmatrix(keptRows, keptStates);
stateValues = storm::utility::vector::filterVector(stateValues, keptRows);
std::unordered_map<std::string, typename SparseModelType::RewardModelType> rewardModels;
@ -149,7 +151,7 @@ namespace storm {
rewardModels.insert(std::make_pair(model.getRewardModels().begin()->first, typename SparseModelType::RewardModelType(boost::none, stateValues)));
}
return std::make_shared<SparseModelType>(flexibleMatrix.createSparseMatrix(), model.getStateLabeling().getSubLabeling(keptStates), std::move(rewardModels));
return std::make_shared<SparseModelType>(flexibleMatrix.createSparseMatrix(keptRows, keptStates), model.getStateLabeling().getSubLabeling(keptStates), std::move(rewardModels));
}

Loading…
Cancel
Save