diff --git a/examples/dft/and_approx.dft b/examples/dft/and_approx.dft new file mode 100644 index 000000000..ff6393dd5 --- /dev/null +++ b/examples/dft/and_approx.dft @@ -0,0 +1,5 @@ +toplevel "A"; +"A" and "B" "C" "D"; +"B" lambda=1 dorm=0; +"C" lambda=100 dorm=0; +"D" lambda=50 dorm=0; diff --git a/examples/dft/and_approx_param.dft b/examples/dft/and_approx_param.dft new file mode 100644 index 000000000..bbd00bd80 --- /dev/null +++ b/examples/dft/and_approx_param.dft @@ -0,0 +1,6 @@ +param x; +toplevel "A"; +"A" and "B" "C" "D"; +"B" lambda=1 dorm=0; +"C" lambda=100 dorm=0; +"D" lambda=100*x dorm=0; diff --git a/examples/dft/approx.dft b/examples/dft/approx.dft new file mode 100644 index 000000000..38c7faba0 --- /dev/null +++ b/examples/dft/approx.dft @@ -0,0 +1,5 @@ +toplevel "A"; +"A" and "B" "C" "D"; +"B" lambda=0.1 dorm=0; +"C" lambda=10 dorm=0; +"D" lambda=0.1 dorm=0; diff --git a/examples/dft/fdep4.dft b/examples/dft/fdep4.dft new file mode 100644 index 000000000..7e3c5642a --- /dev/null +++ b/examples/dft/fdep4.dft @@ -0,0 +1,7 @@ +toplevel "A"; +"A" or "F" "B"; +"F" fdep "E" "C" "D"; +"B" wsp "C" "D"; +"C" lambda=1 dorm=0; +"D" lambda=1 dorm=0.5; +"E" lambda=0.5 dorm=0; diff --git a/examples/dft/mas.dft b/examples/dft/mas.dft new file mode 100644 index 000000000..bed57d59c --- /dev/null +++ b/examples/dft/mas.dft @@ -0,0 +1,58 @@ +toplevel "MAS"; +"MAS" or "CPU" "DB" "MB" "VMB" "MEM" "VMS"; +"CPU" or "CW" "SO1" "SO2" "PG" "SM"; +"CW" and "CWA" "CWB"; +"SO1" and "SO1A" "SO1B"; +"SO2" and "SO2A" "SO2B"; +"PG" and "PGA" "PGB"; +"SM" and "SMA" "SMB"; +"CWA" csp "CWAev" "S1" "S2"; +"CWB" csp "CWBev" "S1" "S2"; +"SO1A" csp "SO1Aev" "S1" "S2"; +"SO1B" csp "SO1Bev" "S1" "S2"; +"SO2A" csp "SO2Aev" "S1" "S2"; +"SO2B" csp "SO2Bev" "S1" "S2"; +"PGA" csp "PGAev" "S1" "S2"; +"PGB" csp "PGBev" "S1" "S2"; +"SMA" csp "SMAev" "S1" "S2"; +"SMB" csp "SMBev" "S1" "S2"; +"CWAev" lambda=1.0e-6 dorm=0; +"CWBev" lambda=1.0e-6 dorm=0; +"SO1Aev" lambda=1.0e-6 dorm=0; +"SO1Bev" lambda=1.0e-6 dorm=0; +"SO2Aev" lambda=1.0e-6 dorm=0; +"SO2Bev" lambda=1.0e-6 dorm=0; +"PGAev" lambda=1.0e-6 dorm=0; +"PGBev" lambda=1.0e-6 dorm=0; +"SMAev" lambda=1.0e-6 dorm=0; +"SMBev" lambda=1.0e-6 dorm=0; +"S1" lambda=1.0e-6 dorm=0; +"S2" lambda=1.0e-6 dorm=0; +"DB" and "DB1" "DB2" "DB3"; +"DB1" lambda=5.0e-6 dorm=0; +"DB2" lambda=5.0e-6 dorm=0; +"DB3" lambda=5.0e-6 dorm=0; +"MB" and "MB1" "MB2" "MB3"; +"MB1" lambda=5.0e-6 dorm=0; +"MB2" lambda=5.0e-6 dorm=0; +"MB3" lambda=5.0e-6 dorm=0; +"VMB" and "VMB1" "VMB2"; +"VMB1" lambda=5.0e-6 dorm=0; +"VMB2" lambda=5.0e-6 dorm=0; +"MEM" and "MEM1" "MEM2"; +"MEM1" lambda=1.0e-5 dorm=0; +"MEM2" lambda=1.0e-5 dorm=0; +"VMS" or "VM1" "VM2"; +"VM1" and "VM1A" "VM1B"; +"VM2" and "VM2A" "VM2B"; +"VM1A" csp "VM1Aev" "VMS1" "VMS2"; +"VM1B" csp "VM1Bev" "VMS1" "VMS2"; +"VM2A" csp "VM2Aev" "VMS1" "VMS2"; +"VM2B" csp "VM2Bev" "VMS1" "VMS2"; +"VM1Aev" lambda=1.0e-6 dorm=0; +"VM1Bev" lambda=1.0e-6 dorm=0; +"VM2Aev" lambda=1.0e-6 dorm=0; +"VM2Bev" lambda=1.0e-6 dorm=0; +"VMS1" lambda=1.0e-6 dorm=0; +"VMS2" lambda=1.0e-6 dorm=0; + diff --git a/examples/dft/nonmonoton.dft b/examples/dft/nonmonoton.dft new file mode 100644 index 000000000..0575028e8 --- /dev/null +++ b/examples/dft/nonmonoton.dft @@ -0,0 +1,6 @@ +toplevel "A"; +"A" or "B" "Z"; +"B" pand "D" "S"; +"Z" lambda=1 dorm=0; +"D" lambda=100 dorm=0; +"S" lambda=50 dorm=0; diff --git a/examples/dft/nonmonoton_param.dft b/examples/dft/nonmonoton_param.dft index e83f3b012..fab011e64 100644 --- a/examples/dft/nonmonoton_param.dft +++ b/examples/dft/nonmonoton_param.dft @@ -2,7 +2,7 @@ param x; param y; toplevel "A"; "A" or "B" "Z"; -"B" pand "D" "S"; -"Z" lambda=y dorm=0; -"D" lambda=100 dorm=0; -"S" lambda=100*x dorm=0; +"Z" pand "C" "D"; +"B" lambda=y dorm=0; +"C" lambda=100 dorm=0; +"D" lambda=100*x dorm=0; diff --git a/resources/3rdparty/CMakeLists.txt b/resources/3rdparty/CMakeLists.txt index ffd5fc54c..57d36770f 100644 --- a/resources/3rdparty/CMakeLists.txt +++ b/resources/3rdparty/CMakeLists.txt @@ -15,7 +15,6 @@ mark_as_advanced(AUTORECONF) ## ############################################################# - # Do not take a branch, needs internet connection. ExternalProject_Add( l3pp diff --git a/src/builder/DftExplorationHeuristic.cpp b/src/builder/DftExplorationHeuristic.cpp new file mode 100644 index 000000000..527513134 --- /dev/null +++ b/src/builder/DftExplorationHeuristic.cpp @@ -0,0 +1,95 @@ +#include "src/builder/DftExplorationHeuristic.h" +#include "src/adapters/CarlAdapter.h" +#include "src/utility/macros.h" +#include "src/utility/constants.h" +#include "src/exceptions/NotImplementedException.h" + +namespace storm { + namespace builder { + + template<typename ValueType> + DFTExplorationHeuristic<ValueType>::DFTExplorationHeuristic(size_t id) : id(id), expand(true), 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), lowerBound(storm::utility::zero<ValueType>()), upperBound(storm::utility::zero<ValueType>()), 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"); + probability += predecessor.getProbability() * rate/exitRate; + return true; + } + + template<> + bool DFTExplorationHeuristicProbability<storm::RationalFunction>::updateHeuristicValues(DFTExplorationHeuristic<storm::RationalFunction> const& predecessor, storm::RationalFunction rate, storm::RationalFunction exitRate) { + STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Heuristic rate ration does not work for rational functions."); + return false; + } + + template<> + double DFTExplorationHeuristicProbability<double>::getPriority() const { + return probability; + } + + template<> + double DFTExplorationHeuristicProbability<storm::RationalFunction>::getPriority() const { + STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Heuristic rate ration does not work for rational functions."); + } + + template<> + bool DFTExplorationHeuristicBoundDifference<double>::updateHeuristicValues(DFTExplorationHeuristic<double> const& predecessor, double rate, double exitRate) { + STORM_LOG_ASSERT(exitRate > 0, "Exit rate is 0"); + probability += predecessor.getProbability() * rate/exitRate; + return true; + } + + template<> + bool DFTExplorationHeuristicBoundDifference<storm::RationalFunction>::updateHeuristicValues(DFTExplorationHeuristic<storm::RationalFunction> const& predecessor, storm::RationalFunction rate, storm::RationalFunction exitRate) { + STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Heuristic rate ration does not work for rational functions."); + return false; + } + + template<typename ValueType> + void DFTExplorationHeuristicBoundDifference<ValueType>::setBounds(ValueType lowerBound, ValueType upperBound) { + this->lowerBound = lowerBound; + this->upperBound = upperBound; + difference = (storm::utility::one<ValueType>() / upperBound) - (storm::utility::one<ValueType>() / lowerBound); + } + + template<> + double DFTExplorationHeuristicBoundDifference<double>::getPriority() const { + return probability * difference; + } + + 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 new file mode 100644 index 000000000..44931dfbf --- /dev/null +++ b/src/builder/DftExplorationHeuristic.h @@ -0,0 +1,190 @@ +#ifndef STORM_BUILDER_DFTEXPLORATIONHEURISTIC_H_ +#define STORM_BUILDER_DFTEXPLORATIONHEURISTIC_H_ + +#include <memory> + +namespace storm { + + namespace builder { + + /*! + * Enum representing the heuristic used for deciding which states to expand. + */ + enum class ApproximationHeuristic { NONE, DEPTH, PROBABILITY, BOUNDDIFFERENCE }; + + + /*! + * General super class for appoximation heuristics. + */ + template<typename ValueType> + class DFTExplorationHeuristic { + + public: + DFTExplorationHeuristic(size_t id); + + DFTExplorationHeuristic(size_t id, DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate); + + void setBounds(ValueType lowerBound, ValueType upperBound); + + virtual bool updateHeuristicValues(DFTExplorationHeuristic const& predecessor, ValueType rate, ValueType exitRate) = 0; + + virtual double getPriority() const = 0; + + virtual bool isSkip(double approximationThreshold) const = 0; + + void markExpand() { + expand = true; + } + + size_t getId() const { + return id; + } + + bool isExpand() { + return expand; + } + + size_t getDepth() const { + return depth; + } + + ValueType getProbability() const { + 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; + }; + + template<typename ValueType> + class DFTExplorationHeuristicNone : public DFTExplorationHeuristic<ValueType> { + public: + DFTExplorationHeuristicNone(size_t id) : DFTExplorationHeuristic<ValueType>(id) { + // Intentionally left empty + } + + DFTExplorationHeuristicNone(size_t id, DFTExplorationHeuristicNone<ValueType> const& predecessor, ValueType rate, ValueType exitRate) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) { + // Intentionally left empty + } + + bool updateHeuristicValues(DFTExplorationHeuristic<ValueType> const& predecessor, ValueType rate, ValueType exitRate) override { + return false; + } + + double getPriority() const override { + return 0; + } + + bool isSkip(double approximationThreshold) const override { + return false; + } + + bool operator<(DFTExplorationHeuristicNone<ValueType> const& other) const { + return this->id > other.id; + } + }; + + template<typename ValueType> + class DFTExplorationHeuristicDepth : public DFTExplorationHeuristic<ValueType> { + public: + DFTExplorationHeuristicDepth(size_t id) : DFTExplorationHeuristic<ValueType>(id) { + // Intentionally left empty + } + + DFTExplorationHeuristicDepth(size_t id, DFTExplorationHeuristicDepth<ValueType> const& predecessor, ValueType rate, ValueType exitRate) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) { + // Intentionally left empty + } + + bool updateHeuristicValues(DFTExplorationHeuristic<ValueType> const& predecessor, ValueType rate, ValueType exitRate) override { + if (predecessor.getDepth() < this->depth) { + this->depth = predecessor.getDepth(); + return true; + } + return false; + } + + + double getPriority() const override { + return this->depth; + } + + bool isSkip(double approximationThreshold) const override { + return !this->expand && this->depth > approximationThreshold; + } + + bool operator<(DFTExplorationHeuristicDepth<ValueType> const& other) const { + return this->depth > other.depth; + } + }; + + template<typename ValueType> + class DFTExplorationHeuristicProbability : public DFTExplorationHeuristic<ValueType> { + public: + DFTExplorationHeuristicProbability(size_t id) : DFTExplorationHeuristic<ValueType>(id) { + // Intentionally left empty + } + + DFTExplorationHeuristicProbability(size_t id, DFTExplorationHeuristicProbability<ValueType> const& predecessor, ValueType rate, ValueType exitRate) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) { + // Intentionally left empty + } + + bool updateHeuristicValues(DFTExplorationHeuristic<ValueType> const& predecessor, ValueType rate, ValueType exitRate) override; + + double getPriority() const override; + + bool isSkip(double approximationThreshold) const override { + return !this->expand && this->getPriority() < approximationThreshold; + } + + bool operator<(DFTExplorationHeuristicProbability<ValueType> const& other) const { + return this->getPriority() < other.getPriority(); + } + }; + + 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) : DFTExplorationHeuristic<ValueType>(id, predecessor, rate, exitRate) { + // Intentionally left empty + } + + void setBounds(ValueType lowerBound, ValueType upperBound); + + bool updateHeuristicValues(DFTExplorationHeuristic<ValueType> const& predecessor, ValueType rate, ValueType exitRate) override; + + 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 difference; + }; + + + } +} + +#endif /* STORM_BUILDER_DFTEXPLORATIONHEURISTIC_H_ */ diff --git a/src/builder/DftSmtBuilder.cpp b/src/builder/DftSmtBuilder.cpp new file mode 100644 index 000000000..c1226531c --- /dev/null +++ b/src/builder/DftSmtBuilder.cpp @@ -0,0 +1,120 @@ +#include "src/builder/DftSmtBuilder.h" +#include "src/exceptions/NotImplementedException.h" + +namespace storm { + namespace builder { + + template <typename ValueType> + DFTSMTBuilder<ValueType>::DFTSMTBuilder() : manager(std::make_shared<storm::expressions::ExpressionManager>()) { + solver = storm::utility::solver::SmtSolverFactory().create(*manager); + } + + template <typename ValueType> + void DFTSMTBuilder<ValueType>::convertToSMT(storm::storage::DFT<ValueType> const& dft) { + std::cout << "Convert DFT to SMT" << std::endl; + timeMax = manager->integer(dft.nrBasicElements()); + timeFailSafe = manager->integer(dft.nrBasicElements() + 1); + timeZero = manager->integer(0); + + // Convert all elements + for (size_t i = 0; i < dft.nrElements(); ++i) { + std::shared_ptr<storm::storage::DFTElement<ValueType> const> element = dft.getElement(i); + std::cout << "Consider " << element->toString() << std::endl; + if (element->isBasicElement()) { + storm::expressions::Variable varBE = convert(std::static_pointer_cast<storm::storage::DFTBE<ValueType> const>(element)); + varsBE.push_back(varBE); + } else if (element->isGate()) { + storm::expressions::Variable varGate = convert(std::static_pointer_cast<storm::storage::DFTGate<ValueType> const>(element)); + if (dft.getTopLevelIndex() == i) { + topLevel = varGate; + } + } else if (element->isDependency()) { + convert(std::static_pointer_cast<storm::storage::DFTDependency<ValueType> const>(element)); + } else if (element->isRestriction()) { + convert(std::static_pointer_cast<storm::storage::DFTRestriction<ValueType> const>(element)); + } + } + + // No simultaneous fails can occur + for (size_t i = 0; i < varsBE.size() - 1; ++i) { + storm::expressions::Expression be = varsBE[i]; + for (size_t j = i + 1; j < varsBE.size(); ++j) { + storm::expressions::Expression assertion = be != varsBE[j]; + solver->add(assertion); + } + } + + // For every time-point one BE must fail + for (size_t time = 1; time <= dft.nrBasicElements(); ++time) { + storm::expressions::Expression exprTime = manager->integer(time); + storm::expressions::Expression assertion = varsBE[0] == exprTime; + for (size_t i = 1; i < varsBE.size(); ++i) { + assertion = assertion || varsBE[i] == exprTime; + } + assertion = assertion || topLevel <= exprTime; + solver->add(assertion); + } + + } + + template <typename ValueType> + storm::expressions::Variable DFTSMTBuilder<ValueType>::convert(std::shared_ptr<storm::storage::DFTBE<ValueType> const> const& be) { + storm::expressions::Variable var = manager->declareIntegerVariable(be->name()); + storm::expressions::Expression assertion = timeZero < var && var < timeFailSafe; + solver->add(assertion); + return var; + } + + template <typename ValueType> + storm::expressions::Variable DFTSMTBuilder<ValueType>::convert(std::shared_ptr<storm::storage::DFTGate<ValueType> const> const& gate) { + storm::expressions::Variable var = manager->declareIntegerVariable(gate->name()); + storm::expressions::Expression assertion = timeZero < var && var <= timeFailSafe; + solver->add(assertion); + return var; + } + + template <typename ValueType> + storm::expressions::Variable DFTSMTBuilder<ValueType>::convert(std::shared_ptr<storm::storage::DFTDependency<ValueType> const> const& dependency) { + STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "The dependency cannot be converted into SMT."); + } + + template <typename ValueType> + storm::expressions::Variable DFTSMTBuilder<ValueType>::convert(std::shared_ptr<storm::storage::DFTRestriction<ValueType> const> const& restriction) { + STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "The restriction cannot be converted into SMT."); + } + + template <typename ValueType> + bool DFTSMTBuilder<ValueType>::check() const { + std::cout << "Check" << std::endl; + storm::solver::SmtSolver::CheckResult result = solver->check(); + switch (result) { + case solver::SmtSolver::CheckResult::Sat: + { + std::cout << "SAT with model:" << std::endl; + std::shared_ptr<storm::solver::SmtSolver::ModelReference> model = solver->getModel(); + for (auto const& pair : *manager) { + std::cout << pair.first.getName() << "->" << model->getIntegerValue(pair.first) << std::endl; + } + return true; + } + case solver::SmtSolver::CheckResult::Unsat: + return false; + case solver::SmtSolver::CheckResult::Unknown: + default: + STORM_LOG_ASSERT(false, "Result is unknown."); + return false; + } + } + + + // Explicitly instantiate the class. + template class DFTSMTBuilder<double>; + +#ifdef STORM_HAVE_CARL + template class DFTSMTBuilder<storm::RationalFunction>; +#endif + + } // namespace builder +} // namespace storm + + diff --git a/src/builder/DftSmtBuilder.h b/src/builder/DftSmtBuilder.h new file mode 100644 index 000000000..596e8106f --- /dev/null +++ b/src/builder/DftSmtBuilder.h @@ -0,0 +1,49 @@ +#ifndef DFTSMTBUILDER_H +#define DFTSMTBUILDER_H + +#include "src/solver/SmtSolver.h" +#include "src/utility/solver.h" +#include "src/storage/dft/DFT.h" + +namespace storm { + namespace builder { + + template<typename ValueType> + class DFTSMTBuilder { + + public: + DFTSMTBuilder(); + + void convertToSMT(storm::storage::DFT<ValueType> const& dft); + + bool check() const; + + private: + + std::shared_ptr<storm::solver::SmtSolver> solver; + + std::shared_ptr<storm::expressions::ExpressionManager> manager; + + storm::expressions::Expression timeMax; + + storm::expressions::Expression timeFailSafe; + + storm::expressions::Expression timeZero; + + storm::expressions::Expression topLevel; + + std::vector<storm::expressions::Variable> varsBE; + + storm::expressions::Variable convert(std::shared_ptr<storm::storage::DFTBE<ValueType> const> const& be); + + storm::expressions::Variable convert(std::shared_ptr<storm::storage::DFTGate<ValueType> const> const& gate); + + storm::expressions::Variable convert(std::shared_ptr<storm::storage::DFTDependency<ValueType> const> const& dependency); + + storm::expressions::Variable convert(std::shared_ptr<storm::storage::DFTRestriction<ValueType> const> const& restriction); + + }; + } +} + +#endif /* DFTSMTBUILDER_H */ diff --git a/src/builder/ExplicitDFTModelBuilder.cpp b/src/builder/ExplicitDFTModelBuilder.cpp index adc44b6ec..4090e0e12 100644 --- a/src/builder/ExplicitDFTModelBuilder.cpp +++ b/src/builder/ExplicitDFTModelBuilder.cpp @@ -4,6 +4,8 @@ #include <src/utility/constants.h> #include "src/utility/vector.h" #include <src/exceptions/UnexpectedException.h> +#include "src/settings/modules/DFTSettings.h" +#include "src/settings/SettingsManager.h" #include <map> namespace storm { @@ -56,6 +58,7 @@ namespace storm { STORM_LOG_ASSERT(mStates.getValue(pseudoStatePair.second) >= OFFSET_PSEUDO_STATE, "State is no pseudo state."); STORM_LOG_TRACE("Create pseudo state from bit vector " << pseudoStatePair.second); DFTStatePointer pseudoState = std::make_shared<storm::storage::DFTState<ValueType>>(pseudoStatePair.second, mDft, *mStateGenerationInfo, newIndex); + pseudoState->construct(); STORM_LOG_ASSERT(pseudoStatePair.second == pseudoState->status(), "Pseudo states do not coincide."); STORM_LOG_TRACE("Explore pseudo state " << mDft.getStateString(pseudoState) << " with id " << pseudoState->getId()); auto exploreResult = exploreStates(pseudoState, rowOffset, transitionMatrixBuilder, tmpMarkovianStates, modelComponents.exitRates); @@ -157,7 +160,7 @@ namespace storm { return model; } - + template <typename ValueType> std::pair<uint_fast64_t, bool> ExplicitDFTModelBuilder<ValueType>::exploreStates(DFTStatePointer const& state, size_t& rowOffset, storm::storage::SparseMatrixBuilder<ValueType>& transitionMatrixBuilder, std::vector<uint_fast64_t>& markovianStates, std::vector<ValueType>& exitRates) { STORM_LOG_TRACE("Explore state: " << mDft.getStateString(state)); @@ -288,7 +291,7 @@ namespace storm { if (mDft.hasRepresentant(nextBE->id())) { // Active must be checked for the state we are coming from as this state is responsible for the // rate and not the new state we are going to - isActive = state->isActive(mDft.getRepresentant(nextBE->id())->id()); + isActive = state->isActive(mDft.getRepresentant(nextBE->id())); } ValueType rate = isActive ? nextBE->activeFailureRate() : nextBE->passiveFailureRate(); STORM_LOG_ASSERT(!storm::utility::isZero(rate), "Rate is 0."); diff --git a/src/builder/ExplicitDFTModelBuilder.h b/src/builder/ExplicitDFTModelBuilder.h index 3e78180d1..b9c54d2b6 100644 --- a/src/builder/ExplicitDFTModelBuilder.h +++ b/src/builder/ExplicitDFTModelBuilder.h @@ -60,7 +60,7 @@ namespace storm { size_t initialStateIndex = 0; public: - struct LabelOptions { + struct LabelOptions { bool buildFailLabel = true; bool buildFailSafeLabel = false; std::set<std::string> beLabels = {}; diff --git a/src/builder/ExplicitDFTModelBuilderApprox.cpp b/src/builder/ExplicitDFTModelBuilderApprox.cpp new file mode 100644 index 000000000..2f302b0db --- /dev/null +++ b/src/builder/ExplicitDFTModelBuilderApprox.cpp @@ -0,0 +1,731 @@ +#include "src/builder/ExplicitDFTModelBuilderApprox.h" +#include <src/models/sparse/MarkovAutomaton.h> +#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" +#include <map> + +namespace storm { + namespace builder { + + template<typename ValueType, typename StateType> + ExplicitDFTModelBuilderApprox<ValueType, StateType>::ModelComponents::ModelComponents() : transitionMatrix(), stateLabeling(), markovianStates(), exitRates(), choiceLabeling() { + // Intentionally left empty. + } + + template<typename ValueType, typename StateType> + ExplicitDFTModelBuilderApprox<ValueType, StateType>::MatrixBuilder::MatrixBuilder(bool canHaveNondeterminism) : mappingOffset(0), stateRemapping(), currentRowGroup(0), currentRow(0), canHaveNondeterminism((canHaveNondeterminism)) { + // Create matrix builder + builder = storm::storage::SparseMatrixBuilder<ValueType>(0, 0, 0, false, canHaveNondeterminism, 0); + } + + template<typename ValueType, typename StateType> + ExplicitDFTModelBuilderApprox<ValueType, StateType>::ExplicitDFTModelBuilderApprox(storm::storage::DFT<ValueType> const& dft, storm::storage::DFTIndependentSymmetries const& symmetries, bool enableDC) : + dft(dft), + stateGenerationInfo(std::make_shared<storm::storage::DFTStateGenerationInfo>(dft.buildStateGenerationInfo(symmetries))), + enableDC(enableDC), + usedHeuristic(storm::settings::getModule<storm::settings::modules::DFTSettings>().getApproximationHeuristic()), + generator(dft, *stateGenerationInfo, enableDC, mergeFailedStates), + matrixBuilder(!generator.isDeterministicModel()), + stateStorage(((dft.stateVectorSize() / 64) + 1) * 64), + // TODO Matthias: make choosable + //explorationQueue(dft.nrElements()+1, 0, 1) + explorationQueue(200, 0, 0.9) + { + // Intentionally left empty. + // TODO Matthias: remove again + usedHeuristic = storm::builder::ApproximationHeuristic::BOUNDDIFFERENCE; + + // Compute independent subtrees + if (dft.topLevelType() == storm::storage::DFTElementType::OR) { + // We only support this for approximation with top level element OR + for (auto const& child : dft.getGate(dft.getTopLevelIndex())->children()) { + // Consider all children of the top level gate + std::vector<size_t> isubdft; + if (child->nrParents() > 1 || child->hasOutgoingDependencies()) { + STORM_LOG_TRACE("child " << child->name() << "does not allow modularisation."); + isubdft.clear(); + } else if (dft.isGate(child->id())) { + isubdft = dft.getGate(child->id())->independentSubDft(false); + } else { + STORM_LOG_ASSERT(dft.isBasicElement(child->id()), "Child is no BE."); + if(dft.getBasicElement(child->id())->hasIngoingDependencies()) { + STORM_LOG_TRACE("child " << child->name() << "does not allow modularisation."); + isubdft.clear(); + } else { + isubdft = {child->id()}; + } + } + if(isubdft.empty()) { + subtreeBEs.clear(); + break; + } else { + // Add new independent subtree + std::vector<size_t> subtree; + for (size_t id : isubdft) { + if (dft.isBasicElement(id)) { + subtree.push_back(id); + } + } + subtreeBEs.push_back(subtree); + } + } + } + if (subtreeBEs.empty()) { + // Consider complete tree + std::vector<size_t> subtree; + for (size_t id = 0; id < dft.nrElements(); ++id) { + if (dft.isBasicElement(id)) { + subtree.push_back(id); + } + } + subtreeBEs.push_back(subtree); + } + + for (auto tree : subtreeBEs) { + std::stringstream stream; + stream << "Subtree: "; + for (size_t id : tree) { + stream << id << ", "; + } + STORM_LOG_TRACE(stream.str()); + } + } + + template<typename ValueType, typename StateType> + void ExplicitDFTModelBuilderApprox<ValueType, StateType>::buildModel(LabelOptions const& labelOpts, size_t iteration, double approximationThreshold) { + STORM_LOG_TRACE("Generating DFT state space"); + + if (iteration < 1) { + // Initialize + modelComponents.markovianStates = storm::storage::BitVector(INITIAL_BITVECTOR_SIZE); + + if(mergeFailedStates) { + // Introduce explicit fail state + storm::generator::StateBehavior<ValueType, StateType> behavior = generator.createMergeFailedState([this] (DFTStatePointer const& state) { + this->failedStateId = newIndex++; + matrixBuilder.stateRemapping.push_back(0); + return this->failedStateId; + } ); + + matrixBuilder.setRemapping(failedStateId); + STORM_LOG_ASSERT(!behavior.empty(), "Behavior is empty."); + matrixBuilder.newRowGroup(); + setMarkovian(behavior.begin()->isMarkovian()); + + // Now add self loop. + // TODO Matthias: maybe use general method. + STORM_LOG_ASSERT(behavior.getNumberOfChoices() == 1, "Wrong number of choices for failed state."); + STORM_LOG_ASSERT(behavior.begin()->size() == 1, "Wrong number of transitions for failed state."); + std::pair<StateType, ValueType> stateProbabilityPair = *(behavior.begin()->begin()); + STORM_LOG_ASSERT(stateProbabilityPair.first == failedStateId, "No self loop for failed state."); + STORM_LOG_ASSERT(storm::utility::isOne<ValueType>(stateProbabilityPair.second), "Probability for failed state != 1."); + matrixBuilder.addTransition(stateProbabilityPair.first, stateProbabilityPair.second); + matrixBuilder.finishRow(); + } + + // Build initial state + this->stateStorage.initialStateIndices = generator.getInitialStates(std::bind(&ExplicitDFTModelBuilderApprox::getOrAddStateIndex, this, std::placeholders::_1)); + STORM_LOG_ASSERT(stateStorage.initialStateIndices.size() == 1, "Only one initial state assumed."); + initialStateIndex = stateStorage.initialStateIndices[0]; + STORM_LOG_TRACE("Initial state: " << initialStateIndex); + // Initialize heuristic values for inital state + STORM_LOG_ASSERT(!statesNotExplored.at(initialStateIndex).second, "Heuristic for initial state is already initialized"); + ExplorationHeuristicPointer heuristic = std::make_shared<ExplorationHeuristic>(initialStateIndex); + heuristic->markExpand(); + statesNotExplored[initialStateIndex].second = heuristic; + explorationQueue.push(heuristic); + } else { + initializeNextIteration(); + } + + if (approximationThreshold > 0) { + switch (usedHeuristic) { + case storm::builder::ApproximationHeuristic::NONE: + // Do not change anything + approximationThreshold = dft.nrElements()+10; + break; + case storm::builder::ApproximationHeuristic::DEPTH: + approximationThreshold = iteration; + break; + case storm::builder::ApproximationHeuristic::PROBABILITY: + case storm::builder::ApproximationHeuristic::BOUNDDIFFERENCE: + approximationThreshold = 10 * std::pow(2, iteration); + break; + } + } + exploreStateSpace(approximationThreshold); + + size_t stateSize = stateStorage.getNumberOfStates() + (mergeFailedStates ? 1 : 0); + modelComponents.markovianStates.resize(stateSize); + modelComponents.deterministicModel = generator.isDeterministicModel(); + + // Fix the entries in the transition matrix according to the mapping of ids to row group indices + STORM_LOG_ASSERT(matrixBuilder.stateRemapping[initialStateIndex] == initialStateIndex, "Initial state should not be remapped."); + // TODO Matthias: do not consider all rows? + STORM_LOG_TRACE("Remap matrix: " << matrixBuilder.stateRemapping << ", offset: " << matrixBuilder.mappingOffset); + matrixBuilder.remap(); + + STORM_LOG_TRACE("State remapping: " << matrixBuilder.stateRemapping); + STORM_LOG_TRACE("Markovian states: " << modelComponents.markovianStates); + STORM_LOG_DEBUG("Model has " << stateSize << " states"); + STORM_LOG_DEBUG("Model is " << (generator.isDeterministicModel() ? "deterministic" : "non-deterministic")); + + // Build transition matrix + modelComponents.transitionMatrix = matrixBuilder.builder.build(stateSize, stateSize); + if (stateSize <= 15) { + STORM_LOG_TRACE("Transition matrix: " << std::endl << modelComponents.transitionMatrix); + } else { + STORM_LOG_TRACE("Transition matrix: too big to print"); + } + + buildLabeling(labelOpts); + } + + template<typename ValueType, typename StateType> + void ExplicitDFTModelBuilderApprox<ValueType, StateType>::initializeNextIteration() { + STORM_LOG_TRACE("Refining DFT state space"); + + // TODO Matthias: should be easier now as skipped states all are at the end of the matrix + // Push skipped states to explore queue + // TODO Matthias: remove + for (auto const& skippedState : skippedStates) { + statesNotExplored[skippedState.second.first->getId()] = skippedState.second; + explorationQueue.push(skippedState.second.second); + } + + // Initialize matrix builder again + // TODO Matthias: avoid copy + std::vector<uint_fast64_t> copyRemapping = matrixBuilder.stateRemapping; + matrixBuilder = MatrixBuilder(!generator.isDeterministicModel()); + matrixBuilder.stateRemapping = copyRemapping; + StateType nrStates = modelComponents.transitionMatrix.getRowGroupCount(); + STORM_LOG_ASSERT(nrStates == matrixBuilder.stateRemapping.size(), "No. of states does not coincide with mapping size."); + + // Start by creating a remapping from the old indices to the new indices + std::vector<StateType> indexRemapping = std::vector<StateType>(nrStates, 0); + auto iterSkipped = skippedStates.begin(); + size_t skippedBefore = 0; + for (size_t i = 0; i < indexRemapping.size(); ++i) { + while (iterSkipped != skippedStates.end() && iterSkipped->first <= i) { + STORM_LOG_ASSERT(iterSkipped->first < indexRemapping.size(), "Skipped is too high."); + ++skippedBefore; + ++iterSkipped; + } + indexRemapping[i] = i - skippedBefore; + } + + // Set remapping + size_t nrExpandedStates = nrStates - skippedBefore; + matrixBuilder.mappingOffset = nrStates; + STORM_LOG_TRACE("# expanded states: " << nrExpandedStates); + StateType skippedIndex = nrExpandedStates; + std::map<StateType, std::pair<DFTStatePointer, ExplorationHeuristicPointer>> skippedStatesNew; + for (size_t id = 0; id < matrixBuilder.stateRemapping.size(); ++id) { + StateType index = matrixBuilder.stateRemapping[id]; + auto itFind = skippedStates.find(index); + if (itFind != skippedStates.end()) { + // Set new mapping for skipped state + matrixBuilder.stateRemapping[id] = skippedIndex; + skippedStatesNew[skippedIndex] = itFind->second; + indexRemapping[index] = skippedIndex; + ++skippedIndex; + } else { + // Set new mapping for expanded state + matrixBuilder.stateRemapping[id] = indexRemapping[index]; + } + } + STORM_LOG_TRACE("New state remapping: " << matrixBuilder.stateRemapping); + std::stringstream ss; + ss << "Index remapping:" << std::endl; + for (auto tmp : indexRemapping) { + ss << tmp << " "; + } + STORM_LOG_TRACE(ss.str()); + + // Remap markovian states + storm::storage::BitVector markovianStatesNew = storm::storage::BitVector(modelComponents.markovianStates.size(), true); + // Iterate over all not set bits + modelComponents.markovianStates.complement(); + size_t index = modelComponents.markovianStates.getNextSetIndex(0); + while (index < modelComponents.markovianStates.size()) { + markovianStatesNew.set(indexRemapping[index], false); + index = modelComponents.markovianStates.getNextSetIndex(index+1); + } + STORM_LOG_ASSERT(modelComponents.markovianStates.size() - modelComponents.markovianStates.getNumberOfSetBits() == markovianStatesNew.getNumberOfSetBits(), "Remapping of markovian states is wrong."); + STORM_LOG_ASSERT(markovianStatesNew.size() == nrStates, "No. of states does not coincide with markovian size."); + modelComponents.markovianStates = markovianStatesNew; + + // Build submatrix for expanded states + // TODO Matthias: only use row groups when necessary + for (StateType oldRowGroup = 0; oldRowGroup < modelComponents.transitionMatrix.getRowGroupCount(); ++oldRowGroup) { + if (indexRemapping[oldRowGroup] < nrExpandedStates) { + // State is expanded -> copy to new matrix + matrixBuilder.newRowGroup(); + for (StateType oldRow = modelComponents.transitionMatrix.getRowGroupIndices()[oldRowGroup]; oldRow < modelComponents.transitionMatrix.getRowGroupIndices()[oldRowGroup+1]; ++oldRow) { + for (typename storm::storage::SparseMatrix<ValueType>::const_iterator itEntry = modelComponents.transitionMatrix.begin(oldRow); itEntry != modelComponents.transitionMatrix.end(oldRow); ++itEntry) { + auto itFind = skippedStates.find(itEntry->getColumn()); + if (itFind != skippedStates.end()) { + // Set id for skipped states as we remap it later + matrixBuilder.addTransition(matrixBuilder.mappingOffset + itFind->second.first->getId(), itEntry->getValue()); + } else { + // Set newly remapped index for expanded states + matrixBuilder.addTransition(indexRemapping[itEntry->getColumn()], itEntry->getValue()); + } + } + matrixBuilder.finishRow(); + } + } + } + + skippedStates = skippedStatesNew; + + STORM_LOG_ASSERT(matrixBuilder.getCurrentRowGroup() == nrExpandedStates, "Row group size does not match."); + skippedStates.clear(); + } + + template<typename ValueType, typename StateType> + void ExplicitDFTModelBuilderApprox<ValueType, StateType>::exploreStateSpace(double approximationThreshold) { + size_t nrExpandedStates = 0; + size_t nrSkippedStates = 0; + // TODO Matthias: do not empty queue every time but break before + while (!explorationQueue.empty()) { + // Get the first state in the queue + ExplorationHeuristicPointer currentExplorationHeuristic = explorationQueue.popTop(); + StateType currentId = currentExplorationHeuristic->getId(); + auto itFind = statesNotExplored.find(currentId); + STORM_LOG_ASSERT(itFind != statesNotExplored.end(), "Id " << currentId << " not found"); + DFTStatePointer currentState = itFind->second.first; + STORM_LOG_ASSERT(currentExplorationHeuristic == itFind->second.second, "Exploration heuristics do not match"); + STORM_LOG_ASSERT(currentState->getId() == currentId, "Ids do not match"); + // Remove it from the list of not explored states + statesNotExplored.erase(itFind); + STORM_LOG_ASSERT(stateStorage.stateToId.contains(currentState->status()), "State is not contained in state storage."); + STORM_LOG_ASSERT(stateStorage.stateToId.getValue(currentState->status()) == currentId, "Ids of states do not coincide."); + + // Get concrete state if necessary + if (currentState->isPseudoState()) { + // Create concrete state from pseudo state + currentState->construct(); + } + STORM_LOG_ASSERT(!currentState->isPseudoState(), "State is pseudo state."); + + // Remember that the current row group was actually filled with the transitions of a different state + matrixBuilder.setRemapping(currentId); + + matrixBuilder.newRowGroup(); + + // Try to explore the next state + generator.load(currentState); + + if (approximationThreshold > 0.0 && nrExpandedStates > approximationThreshold && !currentExplorationHeuristic->isExpand()) { + //if (currentExplorationHeuristic->isSkip(approximationThreshold)) { + // Skip the current state + ++nrSkippedStates; + STORM_LOG_TRACE("Skip expansion of state: " << dft.getStateString(currentState)); + setMarkovian(true); + // Add transition to target state with temporary value 0 + // TODO Matthias: what to do when there is no unique target state? + matrixBuilder.addTransition(failedStateId, storm::utility::zero<ValueType>()); + // Remember skipped state + skippedStates[matrixBuilder.getCurrentRowGroup() - 1] = std::make_pair(currentState, currentExplorationHeuristic); + matrixBuilder.finishRow(); + } else { + // Explore the current state + ++nrExpandedStates; + storm::generator::StateBehavior<ValueType, StateType> behavior = generator.expand(std::bind(&ExplicitDFTModelBuilderApprox::getOrAddStateIndex, this, std::placeholders::_1)); + STORM_LOG_ASSERT(!behavior.empty(), "Behavior is empty."); + setMarkovian(behavior.begin()->isMarkovian()); + + // Now add all choices. + for (auto const& choice : behavior) { + // Add the probabilistic behavior to the matrix. + for (auto const& stateProbabilityPair : choice) { + STORM_LOG_ASSERT(!storm::utility::isZero(stateProbabilityPair.second), "Probability zero."); + // Set transition to state id + offset. This helps in only remapping all previously skipped states. + matrixBuilder.addTransition(matrixBuilder.mappingOffset + stateProbabilityPair.first, stateProbabilityPair.second); + // Set heuristic values for reached states + auto iter = statesNotExplored.find(stateProbabilityPair.first); + if (iter != statesNotExplored.end()) { + // Update heuristic values + DFTStatePointer state = iter->second.first; + if (!iter->second.second) { + // Initialize heuristic values + ExplorationHeuristicPointer heuristic = std::make_shared<ExplorationHeuristic>(stateProbabilityPair.first, *currentExplorationHeuristic, stateProbabilityPair.second, choice.getTotalMass()); + 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(); + } + if (usedHeuristic == storm::builder::ApproximationHeuristic::BOUNDDIFFERENCE) { + // Compute bounds for heuristic now + if (state->isPseudoState()) { + // Create concrete state from pseudo state + state->construct(); + } + STORM_LOG_ASSERT(!currentState->isPseudoState(), "State is pseudo state."); + + // Initialize bounds + ValueType lowerBound = getLowerBound(state); + ValueType upperBound = getUpperBound(state); + heuristic->setBounds(lowerBound, upperBound); + } + + explorationQueue.push(heuristic); + } 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 + explorationQueue.update(iter->second.second, oldPriority); + } + } + } + } + matrixBuilder.finishRow(); + } + } + } // end exploration + + STORM_LOG_INFO("Expanded " << nrExpandedStates << " states"); + STORM_LOG_INFO("Skipped " << nrSkippedStates << " states"); + STORM_LOG_ASSERT(nrSkippedStates == skippedStates.size(), "Nr skipped states is wrong"); + } + + template<typename ValueType, typename StateType> + void ExplicitDFTModelBuilderApprox<ValueType, StateType>::buildLabeling(LabelOptions const& labelOpts) { + // Build state labeling + modelComponents.stateLabeling = storm::models::sparse::StateLabeling(modelComponents.transitionMatrix.getRowGroupCount()); + // Initial state + modelComponents.stateLabeling.addLabel("init"); + modelComponents.stateLabeling.addLabelToState("init", initialStateIndex); + // Label all states corresponding to their status (failed, failsafe, failed BE) + if(labelOpts.buildFailLabel) { + modelComponents.stateLabeling.addLabel("failed"); + } + if(labelOpts.buildFailSafeLabel) { + modelComponents.stateLabeling.addLabel("failsafe"); + } + + // Collect labels for all BE + std::vector<std::shared_ptr<storage::DFTBE<ValueType>>> basicElements = dft.getBasicElements(); + for (std::shared_ptr<storage::DFTBE<ValueType>> elem : basicElements) { + if(labelOpts.beLabels.count(elem->name()) > 0) { + modelComponents.stateLabeling.addLabel(elem->name() + "_fail"); + } + } + + // Set labels to states + if(mergeFailedStates) { + modelComponents.stateLabeling.addLabelToState("failed", failedStateId); + } + for (auto const& stateIdPair : stateStorage.stateToId) { + storm::storage::BitVector state = stateIdPair.first; + size_t stateId = stateIdPair.second; + if (!mergeFailedStates && labelOpts.buildFailLabel && dft.hasFailed(state, *stateGenerationInfo)) { + modelComponents.stateLabeling.addLabelToState("failed", stateId); + } + if (labelOpts.buildFailSafeLabel && dft.isFailsafe(state, *stateGenerationInfo)) { + modelComponents.stateLabeling.addLabelToState("failsafe", stateId); + }; + // Set fail status for each BE + for (std::shared_ptr<storage::DFTBE<ValueType>> elem : basicElements) { + if (labelOpts.beLabels.count(elem->name()) > 0 && storm::storage::DFTState<ValueType>::hasFailed(state, stateGenerationInfo->getStateIndex(elem->id())) ) { + modelComponents.stateLabeling.addLabelToState(elem->name() + "_fail", stateId); + } + } + } + } + + template<typename ValueType, typename StateType> + std::shared_ptr<storm::models::sparse::Model<ValueType>> ExplicitDFTModelBuilderApprox<ValueType, StateType>::getModel() { + STORM_LOG_ASSERT(skippedStates.size() == 0, "Concrete model has skipped states"); + return createModel(false); + } + + 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 + changeMatrixBound(modelComponents.transitionMatrix, lowerBound); + return createModel(true); + } + + template<typename ValueType, typename StateType> + std::shared_ptr<storm::models::sparse::Model<ValueType>> ExplicitDFTModelBuilderApprox<ValueType, StateType>::createModel(bool copy) { + std::shared_ptr<storm::models::sparse::Model<ValueType>> model; + + if (modelComponents.deterministicModel) { + // Build CTMC + if (copy) { + model = std::make_shared<storm::models::sparse::Ctmc<ValueType>>(modelComponents.transitionMatrix, modelComponents.stateLabeling); + } else { + model = std::make_shared<storm::models::sparse::Ctmc<ValueType>>(std::move(modelComponents.transitionMatrix), std::move(modelComponents.stateLabeling)); + } + } else { + // Build MA + // Compute exit rates + // TODO Matthias: avoid computing multiple times + modelComponents.exitRates = std::vector<ValueType>(modelComponents.markovianStates.size()); + std::vector<typename storm::storage::SparseMatrix<ValueType>::index_type> indices = modelComponents.transitionMatrix.getRowGroupIndices(); + for (StateType stateIndex = 0; stateIndex < modelComponents.markovianStates.size(); ++stateIndex) { + if (modelComponents.markovianStates[stateIndex]) { + modelComponents.exitRates[stateIndex] = modelComponents.transitionMatrix.getRowSum(indices[stateIndex]); + } else { + modelComponents.exitRates[stateIndex] = storm::utility::zero<ValueType>(); + } + } + STORM_LOG_TRACE("Exit rates: " << modelComponents.exitRates); + + std::shared_ptr<storm::models::sparse::MarkovAutomaton<ValueType>> ma; + if (copy) { + ma = std::make_shared<storm::models::sparse::MarkovAutomaton<ValueType>>(modelComponents.transitionMatrix, modelComponents.stateLabeling, modelComponents.markovianStates, modelComponents.exitRates); + } else { + ma = std::make_shared<storm::models::sparse::MarkovAutomaton<ValueType>>(std::move(modelComponents.transitionMatrix), std::move(modelComponents.stateLabeling), std::move(modelComponents.markovianStates), std::move(modelComponents.exitRates)); + } + if (ma->hasOnlyTrivialNondeterminism()) { + // Markov automaton can be converted into CTMC + // TODO Matthias: change components which were not moved accordingly + model = ma->convertToCTMC(); + } else { + model = ma; + } + } + + STORM_LOG_DEBUG("No. states: " << model->getNumberOfStates()); + STORM_LOG_DEBUG("No. transitions: " << model->getNumberOfTransitions()); + if (model->getNumberOfStates() <= 15) { + STORM_LOG_TRACE("Transition matrix: " << std::endl << model->getTransitionMatrix()); + } else { + STORM_LOG_TRACE("Transition matrix: too big to print"); + } + return model; + } + + template<typename ValueType, typename StateType> + 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."); + STORM_LOG_ASSERT(!it->second.first->isPseudoState(), "State is still pseudo state."); + + ExplorationHeuristicPointer heuristic = it->second.second; + if (storm::utility::isZero(heuristic->getUpperBound())) { + // Initialize bounds + ValueType lowerBound = getLowerBound(it->second.first); + ValueType upperBound = getUpperBound(it->second.first); + heuristic->setBounds(lowerBound, upperBound); + } + + // Change bound + if (lowerBound) { + matrixEntry->setValue(it->second.second->getLowerBound()); + } else { + matrixEntry->setValue(it->second.second->getUpperBound()); + } + } + } + + template<typename ValueType, typename StateType> + 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); + } + 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 + ValueType upperBound = storm::utility::one<ValueType>(); + ValueType rateSum = storm::utility::zero<ValueType>(); + + // Compute for each independent subtree + for (std::vector<size_t> const& subtree : subtreeBEs) { + // Get all possible rates + std::vector<ValueType> rates; + storm::storage::BitVector coldBEs(subtree.size(), false); + for (size_t i = 0; i < subtree.size(); ++i) { + size_t id = subtree[i]; + if (state->isOperational(id)) { + // Get BE rate + ValueType rate = state->getBERate(id); + if (storm::utility::isZero<ValueType>(rate)) { + // Get active failure rate for cold BE + rate = dft.getBasicElement(id)->activeFailureRate(); + // Mark BE as cold + coldBEs.set(i, true); + } + rates.push_back(rate); + rateSum += rate; + } + } + + STORM_LOG_ASSERT(rates.size() > 0, "No rates failable"); + + // Sort rates + std::sort(rates.begin(), rates.end()); + std::vector<size_t> rateCount(rates.size(), 0); + size_t lastIndex = 0; + for (size_t i = 0; i < rates.size(); ++i) { + if (rates[lastIndex] != rates[i]) { + lastIndex = i; + } + ++rateCount[lastIndex]; + } + + for (size_t i = 0; i < rates.size(); ++i) { + // Cold BEs cannot fail in the first step + if (!coldBEs.get(i) && rateCount[i] > 0) { + std::iter_swap(rates.begin() + i, rates.end() - 1); + // Compute AND MTTF of subtree without current rate and scale with current rate + upperBound += rates.back() * rateCount[i] * computeMTTFAnd(rates, rates.size() - 1); + // Swap back + // TODO try to avoid swapping back + std::iter_swap(rates.begin() + i, rates.end() - 1); + } + } + } + + STORM_LOG_TRACE("Upper bound is " << (rateSum / upperBound) << " for state " << state->getId()); + STORM_LOG_ASSERT(!storm::utility::isZero(upperBound), "Upper bound is 0"); + STORM_LOG_ASSERT(!storm::utility::isZero(rateSum), "State is absorbing"); + return rateSum / upperBound; + } + + template<typename ValueType, typename StateType> + ValueType ExplicitDFTModelBuilderApprox<ValueType, StateType>::computeMTTFAnd(std::vector<ValueType> const& rates, size_t size) const { + ValueType result = storm::utility::zero<ValueType>(); + if (size == 0) { + return result; + } + + // TODO Matthias: works only for <64 BEs! + // WARNING: this code produces wrong results for more than 32 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 & static_cast<size_t>(1 << static_cast<size_t>(j))) { + // WARNING: if the first bit is set, it also recognizes the 32nd bit as set + // TODO: fix + permResult += rates[j]; + } + } + permutation = nextBitPermutation(permutation); + STORM_LOG_ASSERT(!storm::utility::isZero(permResult), "PermResult is 0"); + sum += storm::utility::one<ValueType>() / permResult; + } while(permutation < (static_cast<size_t>(1) << rates.size()) && permutation != 0); + if (i % 2 == 0) { + result -= sum; + } else { + result += sum; + } + }*/ + + // Compute result with permutations of size <= 3 + ValueType one = storm::utility::one<ValueType>(); + for (size_t i1 = 0; i1 < size; ++i1) { + // + 1/a + ValueType sum = rates[i1]; + result += one / sum; + for (size_t i2 = 0; i2 < i1; ++i2) { + // - 1/(a+b) + ValueType sum2 = sum + rates[i2]; + result -= one / sum2; + for (size_t i3 = 0; i3 < i2; ++i3) { + // + 1/(a+b+c) + result += one / (sum2 + rates[i3]); + } + } + } + + STORM_LOG_ASSERT(!storm::utility::isZero(result), "UpperBound is 0"); + return result; + } + + template<typename ValueType, typename StateType> + StateType ExplicitDFTModelBuilderApprox<ValueType, StateType>::getOrAddStateIndex(DFTStatePointer const& state) { + StateType stateId; + bool changed = false; + + if (stateGenerationInfo->hasSymmetries()) { + // Order state by symmetry + STORM_LOG_TRACE("Check for symmetry: " << dft.getStateString(state)); + changed = state->orderBySymmetry(); + STORM_LOG_TRACE("State " << (changed ? "changed to " : "did not change") << (changed ? dft.getStateString(state) : "")); + } + + if (stateStorage.stateToId.contains(state->status())) { + // State already exists + stateId = stateStorage.stateToId.getValue(state->status()); + STORM_LOG_TRACE("State " << dft.getStateString(state) << " with id " << stateId << " already exists"); + if (!changed) { + // Check if state is pseudo state + // If state is explored already the possible pseudo state was already constructed + auto iter = statesNotExplored.find(stateId); + if (iter != statesNotExplored.end() && iter->second.first->isPseudoState()) { + // Create pseudo state now + STORM_LOG_ASSERT(iter->second.first->getId() == stateId, "Ids do not match."); + STORM_LOG_ASSERT(iter->second.first->status() == state->status(), "Pseudo states do not coincide."); + state->setId(stateId); + // Update mapping to map to concrete state now + // TODO Matthias: just change pointer? + statesNotExplored[stateId] = std::make_pair(state, iter->second.second); + // We do not push the new state on the exploration queue as the pseudo state was already pushed + STORM_LOG_TRACE("Created pseudo state " << dft.getStateString(state)); + } + } + } else { + // State does not exist yet + STORM_LOG_ASSERT(state->isPseudoState() == changed, "State type (pseudo/concrete) wrong."); + // Create new state + state->setId(newIndex++); + stateId = stateStorage.stateToId.findOrAdd(state->status(), state->getId()); + STORM_LOG_ASSERT(stateId == state->getId(), "Ids do not match."); + // Insert state as not yet explored + ExplorationHeuristicPointer nullHeuristic; + statesNotExplored[stateId] = std::make_pair(state, nullHeuristic); + // Reserve one slot for the new state in the remapping + matrixBuilder.stateRemapping.push_back(0); + STORM_LOG_TRACE("New " << (state->isPseudoState() ? "pseudo" : "concrete") << " state: " << dft.getStateString(state)); + } + return stateId; + } + + template<typename ValueType, typename StateType> + void ExplicitDFTModelBuilderApprox<ValueType, StateType>::setMarkovian(bool markovian) { + if (matrixBuilder.getCurrentRowGroup() > modelComponents.markovianStates.size()) { + // Resize BitVector + modelComponents.markovianStates.resize(modelComponents.markovianStates.size() + INITIAL_BITVECTOR_SIZE); + } + modelComponents.markovianStates.set(matrixBuilder.getCurrentRowGroup() - 1, markovian); + } + + template<typename ValueType, typename StateType> + void ExplicitDFTModelBuilderApprox<ValueType, StateType>::printNotExplored() const { + std::cout << "states not explored:" << std::endl; + for (auto it : statesNotExplored) { + std::cout << it.first << " -> " << dft.getStateString(it.second.first) << std::endl; + } + } + + + // Explicitly instantiate the class. + template class ExplicitDFTModelBuilderApprox<double>; + +#ifdef STORM_HAVE_CARL + template class ExplicitDFTModelBuilderApprox<storm::RationalFunction>; +#endif + + } // namespace builder +} // namespace storm + + diff --git a/src/builder/ExplicitDFTModelBuilderApprox.h b/src/builder/ExplicitDFTModelBuilderApprox.h new file mode 100644 index 000000000..0e9b9dc57 --- /dev/null +++ b/src/builder/ExplicitDFTModelBuilderApprox.h @@ -0,0 +1,335 @@ +#ifndef EXPLICITDFTMODELBUILDERAPPROX_H +#define EXPLICITDFTMODELBUILDERAPPROX_H + +#include "src/builder/DftExplorationHeuristic.h" +#include "src/models/sparse/StateLabeling.h" +#include "src/models/sparse/StandardRewardModel.h" +#include "src/models/sparse/Model.h" +#include "src/generator/DftNextStateGenerator.h" +#include "src/storage/SparseMatrix.h" +#include "src/storage/sparse/StateStorage.h" +#include "src/storage/dft/DFT.h" +#include "src/storage/dft/SymmetricUnits.h" +#include "src/storage/BucketPriorityQueue.h" +#include <boost/container/flat_set.hpp> +#include <boost/optional/optional.hpp> +#include <stack> +#include <unordered_set> +#include <limits> + +namespace storm { + namespace builder { + + /*! + * Build a Markov chain from DFT. + */ + template<typename ValueType, typename StateType = uint32_t> + class ExplicitDFTModelBuilderApprox { + + using DFTStatePointer = std::shared_ptr<storm::storage::DFTState<ValueType>>; + // TODO Matthias: make choosable + using ExplorationHeuristic = DFTExplorationHeuristicBoundDifference<ValueType>; + using ExplorationHeuristicPointer = std::shared_ptr<ExplorationHeuristic>; + + + // A structure holding the individual components of a model. + struct ModelComponents { + // Constructor + ModelComponents(); + + // The transition matrix. + storm::storage::SparseMatrix<ValueType> transitionMatrix; + + // The state labeling. + storm::models::sparse::StateLabeling stateLabeling; + + // The Markovian states. + storm::storage::BitVector markovianStates; + + // The exit rates. + std::vector<ValueType> exitRates; + + // A vector that stores a labeling for each choice. + boost::optional<std::vector<boost::container::flat_set<StateType>>> choiceLabeling; + + // A flag indicating if the model is deterministic. + bool deterministicModel; + }; + + // A class holding the information for building the transition matrix. + class MatrixBuilder { + public: + // Constructor + MatrixBuilder(bool canHaveNondeterminism); + + /*! + * Set a mapping from a state id to the index in the matrix. + * + * @param id Id of the state. + */ + void setRemapping(StateType id) { + STORM_LOG_ASSERT(id < stateRemapping.size(), "Invalid index for remapping."); + stateRemapping[id] = currentRowGroup; + } + + /*! + * Create a new row group if the model is nondeterministic. + */ + void newRowGroup() { + if (canHaveNondeterminism) { + builder.newRowGroup(currentRow); + } + ++currentRowGroup; + } + + /*! + * Add a transition from the current row. + * + * @param index Target index + * @param value Value of transition + */ + void addTransition(StateType index, ValueType value) { + builder.addNextValue(currentRow, index, value); + } + + /*! + * Finish the current row. + */ + void finishRow() { + ++currentRow; + } + + /*! + * Remap the columns in the matrix. + */ + void remap() { + builder.replaceColumns(stateRemapping, mappingOffset); + } + + /*! + * Get the current row group. + * + * @return The current row group. + */ + StateType getCurrentRowGroup() { + return currentRowGroup; + } + + // Matrix builder. + storm::storage::SparseMatrixBuilder<ValueType> builder; + + // Offset to distinguish states which will not be remapped anymore and those which will. + size_t mappingOffset; + + // A mapping from state ids to the row group indices in which they actually reside. + // TODO Matthias: avoid hack with fixed int type + std::vector<uint_fast64_t> stateRemapping; + + private: + + // Index of the current row group. + StateType currentRowGroup; + + // Index of the current row. + StateType currentRow; + + // Flag indicating if row groups are needed. + bool canHaveNondeterminism; + }; + + public: + // A structure holding the labeling options. + struct LabelOptions { + bool buildFailLabel = true; + bool buildFailSafeLabel = false; + std::set<std::string> beLabels = {}; + }; + + /*! + * Constructor. + * + * @param dft DFT. + * @param symmetries Symmetries in the dft. + * @param enableDC Flag indicating if dont care propagation should be used. + */ + ExplicitDFTModelBuilderApprox(storm::storage::DFT<ValueType> const& dft, storm::storage::DFTIndependentSymmetries const& symmetries, bool enableDC); + + /*! + * Build model from DFT. + * + * @param labelOpts Options for labeling. + * @param iteration Current number of iteration. + * @param approximationThreshold Threshold determining when to skip exploring states. + */ + void buildModel(LabelOptions const& labelOpts, size_t iteration, double approximationThreshold = 0.0); + + /*! + * Get the built model. + * + * @return The model built from the DFT. + */ + std::shared_ptr<storm::models::sparse::Model<ValueType>> getModel(); + + /*! + * Get the built approximation model for either the lower or upper bound. + * + * @param lowerBound If true, the lower bound model is returned, else the upper bound model + * + * @return The model built from the DFT. + */ + std::shared_ptr<storm::models::sparse::Model<ValueType>> getModelApproximation(bool lowerBound = true); + + private: + + /*! + * Explore state space of DFT. + * + * @param approximationThreshold Threshold to determine when to skip states. + */ + void exploreStateSpace(double approximationThreshold); + + /*! + * Initialize the matrix for a refinement iteration. + */ + void initializeNextIteration(); + + /*! + * Build the labeling. + * + * @param labelOpts Options for labeling. + */ + void buildLabeling(LabelOptions const& labelOpts); + + /*! + * Add a state to the explored states (if not already there). It also handles pseudo states. + * + * @param state The state to add. + * + * @return Id of state. + */ + StateType getOrAddStateIndex(DFTStatePointer const& state); + + /*! + * Set markovian flag for the current state. + * + * @param markovian Flag indicating if the state is markovian. + */ + void setMarkovian(bool markovian); + + /** + * Change matrix to reflect the lower or upper approximation bound. + * + * @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 changeMatrixBound(storm::storage::SparseMatrix<ValueType> & matrix, bool lowerBound) const; + + /*! + * Get lower bound approximation for state. + * + * @return Lower bound approximation. + */ + ValueType getLowerBound(DFTStatePointer const& state) const; + + /*! + * Get upper bound approximation for state. + * + * @return Upper bound approximation. + */ + ValueType getUpperBound(DFTStatePointer const& state) const; + + /*! + * Compute the MTTF of an AND gate via a closed formula. + * The used formula is 1/( 1/a + 1/b + 1/c + ... - 1/(a+b) - 1/(a+c) - ... + 1/(a+b+c) + ... - ...) + * + * @param rates List of rates of children of AND. + * @param size Only indices < size are considered in the vector. + * @return MTTF. + */ + ValueType computeMTTFAnd(std::vector<ValueType> const& rates, size_t size) const; + + /*! + * Compares the priority of two states. + * + * @param idA Id of first state + * @param idB Id of second state + * + * @return True if the priority of the first state is greater then the priority of the second one. + */ + bool isPriorityGreater(StateType idA, StateType idB) const; + + void printNotExplored() const; + + /*! + * Create the model model from the model components. + * + * @param copy If true, all elements of the model component are copied (used for approximation). If false + * they are moved to reduce the memory overhead. + * + * @return The model built from the model components. + */ + std::shared_ptr<storm::models::sparse::Model<ValueType>> createModel(bool copy); + + // Initial size of the bitvector. + const size_t INITIAL_BITVECTOR_SIZE = 20000; + // Offset used for pseudo states. + const StateType OFFSET_PSEUDO_STATE = std::numeric_limits<StateType>::max() / 2; + + // Dft + storm::storage::DFT<ValueType> const& dft; + + // General information for state generation + // TODO Matthias: use const reference + std::shared_ptr<storm::storage::DFTStateGenerationInfo> stateGenerationInfo; + + // Flag indication if dont care propagation should be used. + bool enableDC = true; + + //TODO Matthias: make changeable + const bool mergeFailedStates = true; + + // Heuristic used for approximation + storm::builder::ApproximationHeuristic usedHeuristic; + + // Current id for new state + size_t newIndex = 0; + + // Id of failed state + size_t failedStateId = 0; + + // Id of initial state + size_t initialStateIndex = 0; + + // Next state generator for exploring the state space + storm::generator::DftNextStateGenerator<ValueType, StateType> generator; + + // Structure for the components of the model. + ModelComponents modelComponents; + + // Structure for the transition matrix builder. + MatrixBuilder matrixBuilder; + + // Internal information about the states that were explored. + storm::storage::sparse::StateStorage<StateType> stateStorage; + + // A priority queue of states that still need to be explored. + storm::storage::BucketPriorityQueue<ValueType> explorationQueue; + //storm::storage::DynamicPriorityQueue<ExplorationHeuristicPointer, std::vector<ExplorationHeuristicPointer>, std::function<bool(ExplorationHeuristicPointer, ExplorationHeuristicPointer)>> explorationQueue; + + // A mapping of not yet explored states from the id to the tuple (state object, heuristic values). + std::map<StateType, std::pair<DFTStatePointer, ExplorationHeuristicPointer>> statesNotExplored; + + // Holds all skipped states which were not yet expanded. More concretely it is a mapping from matrix indices + // to the corresponding skipped states. + // Notice that we need an ordered map here to easily iterate in increasing order over state ids. + // TODO remove again + std::map<StateType, std::pair<DFTStatePointer, ExplorationHeuristicPointer>> skippedStates; + + // List of independent subtrees and the BEs contained in them. + std::vector<std::vector<size_t>> subtreeBEs; + }; + + } +} + +#endif /* EXPLICITDFTMODELBUILDERAPPROX_H */ diff --git a/src/builder/ParallelCompositionBuilder.cpp b/src/builder/ParallelCompositionBuilder.cpp new file mode 100644 index 000000000..5d78d74fa --- /dev/null +++ b/src/builder/ParallelCompositionBuilder.cpp @@ -0,0 +1,172 @@ +#include "src/builder/ParallelCompositionBuilder.h" +#include "src/models/sparse/StandardRewardModel.h" +#include <src/utility/constants.h> + +namespace storm { + namespace builder { + + template<typename ValueType> + std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> ParallelCompositionBuilder<ValueType>::compose(std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> const& ctmcA, std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> const& ctmcB, bool labelAnd) { + STORM_LOG_TRACE("Parallel composition"); + + storm::storage::SparseMatrix<ValueType> matrixA = ctmcA->getTransitionMatrix(); + storm::storage::SparseMatrix<ValueType> matrixB = ctmcB->getTransitionMatrix(); + storm::models::sparse::StateLabeling labelingA = ctmcA->getStateLabeling(); + storm::models::sparse::StateLabeling labelingB = ctmcB->getStateLabeling(); + size_t sizeA = ctmcA->getNumberOfStates(); + size_t sizeB = ctmcB->getNumberOfStates(); + size_t size = sizeA * sizeB; + size_t rowIndex = 0; + + // Build matrix + storm::storage::SparseMatrixBuilder<ValueType> builder(size, size, 0, true, false, 0); + + for (size_t stateA = 0; stateA < sizeA; ++stateA) { + for (size_t stateB = 0; stateB < sizeB; ++stateB) { + STORM_LOG_ASSERT(rowIndex == stateA * sizeB + stateB, "Row " << rowIndex << " is not correct"); + + auto rowA = matrixA.getRow(stateA); + auto itA = rowA.begin(); + auto rowB = matrixB.getRow(stateB); + auto itB = rowB.begin(); + + // First consider all target states of A < the current stateA + while (itA != rowA.end() && itA->getColumn() < stateA) { + builder.addNextValue(rowIndex, itA->getColumn() * sizeB + stateB, itA->getValue()); + ++itA; + } + + // Then consider all target states of B + while (itB != rowB.end()) { + builder.addNextValue(rowIndex, stateA * sizeB + itB->getColumn(), itB->getValue()); + ++itB; + } + + // Last consider all remaining target states of A > the current stateA + while (itA != rowA.end()) { + builder.addNextValue(rowIndex, itA->getColumn() * sizeB + stateB, itA->getValue()); + ++itA; + } + + ++rowIndex; + } + } + + storm::storage::SparseMatrix<ValueType> matrixComposed = builder.build(); + STORM_LOG_ASSERT(matrixComposed.getRowCount() == size, "Row count is not correct"); + STORM_LOG_ASSERT(matrixComposed.getColumnCount() == size, "Column count is not correct"); + + // Build labeling + storm::models::sparse::StateLabeling labeling = storm::models::sparse::StateLabeling(sizeA * sizeB); + + if (labelAnd) { + for (std::string const& label : labelingA.getLabels()) { + if (labelingB.containsLabel(label)) { + // Only consider labels contained in both CTMCs + storm::storage::BitVector labelStates(size, false); + for (auto entryA : labelingA.getStates(label)) { + for (auto entryB : labelingB.getStates(label)) { + labelStates.set(entryA * sizeB + entryB); + } + } + labeling.addLabel(label, labelStates); + } + } + } else { + // Set labels from A + for (std::string const& label : labelingA.getLabels()) { + if (label == "init") { + // Initial states must be initial in both CTMCs + STORM_LOG_ASSERT(labelingB.containsLabel(label), "B does not have init."); + storm::storage::BitVector labelStates(size, false); + for (auto entryA : labelingA.getStates(label)) { + for (auto entryB : labelingB.getStates(label)) { + labelStates.set(entryA * sizeB + entryB); + } + } + labeling.addLabel(label, labelStates); + } else { + storm::storage::BitVector labelStates(size, false); + for (auto entry : labelingA.getStates(label)) { + for (size_t index = entry * sizeB; index < entry * sizeB + sizeB; ++index) { + labelStates.set(index, true); + } + } + labeling.addLabel(label, labelStates); + } + } + // Set labels from B + for (std::string const& label : labelingB.getLabels()) { + if (label == "init") { + continue; + } + if (labeling.containsLabel(label)) { + // Label is already there from A + for (auto entry : labelingB.getStates(label)) { + for (size_t index = 0; index < sizeA; ++index) { + labeling.addLabelToState(label, index * sizeB + entry); + } + } + } else { + storm::storage::BitVector labelStates(size, false); + for (auto entry : labelingB.getStates(label)) { + for (size_t index = 0; index < sizeA; ++index) { + labelStates.set(index * sizeB + entry, true); + } + } + labeling.addLabel(label, labelStates); + } + } + } + + + // Build CTMC + std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> composedCtmc = std::make_shared<storm::models::sparse::Ctmc<ValueType>>(matrixComposed, labeling); + + // Print for debugging + /*std::cout << "Matrix A:" << std::endl; + std::cout << matrixA << std::endl; + std::cout << "Matrix B:" << std::endl; + std::cout << matrixB << std::endl; + std::cout << "Composed matrix: " << std::endl << matrixComposed; + std::cout << "Labeling A" << std::endl; + labelingA.printLabelingInformationToStream(std::cout); + for (size_t stateA = 0; stateA < sizeA; ++stateA) { + std::cout << "State " << stateA << ": "; + for (auto label : labelingA.getLabelsOfState(stateA)) { + std::cout << label << ", "; + } + std::cout << std::endl; + } + std::cout << "Labeling B" << std::endl; + labelingB.printLabelingInformationToStream(std::cout); + for (size_t stateB = 0; stateB < sizeB; ++stateB) { + std::cout << "State " << stateB << ": "; + for (auto label : labelingB.getLabelsOfState(stateB)) { + std::cout << label << ", "; + } + std::cout << std::endl; + } + std::cout << "Labeling Composed" << std::endl; + labeling.printLabelingInformationToStream(std::cout); + for (size_t state = 0; state < size; ++state) { + std::cout << "State " << state << ": "; + for (auto label : labeling.getLabelsOfState(state)) { + std::cout << label << ", "; + } + std::cout << std::endl; + }*/ + + return composedCtmc; + } + + + // Explicitly instantiate the class. + template class ParallelCompositionBuilder<double>; + +#ifdef STORM_HAVE_CARL + template class ParallelCompositionBuilder<storm::RationalFunction>; +#endif + + } // namespace builder +} // namespace storm diff --git a/src/builder/ParallelCompositionBuilder.h b/src/builder/ParallelCompositionBuilder.h new file mode 100644 index 000000000..3b0fab790 --- /dev/null +++ b/src/builder/ParallelCompositionBuilder.h @@ -0,0 +1,24 @@ +#ifndef PARALLELCOMPOSITIONBUILDER_H +#define PARALLELCOMPOSITIONBUILDER_H + +#include <src/models/sparse/Ctmc.h> + +namespace storm { + namespace builder { + + /*! + * Build a parallel composition of Markov chains. + */ + template<typename ValueType> + class ParallelCompositionBuilder { + + public: + + static std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> compose(std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> const& ctmcA, std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> const& ctmcB, bool labelAnd); + + }; + + } +} + +#endif /* EXPLICITDFTMODELBUPARALLELCOMPOSITIONBUILDER_HILDERAPPROX_H */ diff --git a/src/generator/DftNextStateGenerator.cpp b/src/generator/DftNextStateGenerator.cpp new file mode 100644 index 000000000..8fac28f59 --- /dev/null +++ b/src/generator/DftNextStateGenerator.cpp @@ -0,0 +1,209 @@ +#include "src/generator/DftNextStateGenerator.h" + +#include "src/utility/constants.h" +#include "src/utility/macros.h" +#include "src/exceptions/NotImplementedException.h" + +namespace storm { + namespace generator { + + template<typename ValueType, typename StateType> + DftNextStateGenerator<ValueType, StateType>::DftNextStateGenerator(storm::storage::DFT<ValueType> const& dft, storm::storage::DFTStateGenerationInfo const& stateGenerationInfo, bool enableDC, bool mergeFailedStates) : mDft(dft), mStateGenerationInfo(stateGenerationInfo), state(nullptr), enableDC(enableDC), mergeFailedStates(mergeFailedStates) { + deterministicModel = !mDft.canHaveNondeterminism(); + } + + template<typename ValueType, typename StateType> + bool DftNextStateGenerator<ValueType, StateType>::isDeterministicModel() const { + return deterministicModel; + } + + template<typename ValueType, typename StateType> + std::vector<StateType> DftNextStateGenerator<ValueType, StateType>::getInitialStates(StateToIdCallback const& stateToIdCallback) { + DFTStatePointer initialState = std::make_shared<storm::storage::DFTState<ValueType>>(mDft, mStateGenerationInfo, 0); + + // Register initial state + StateType id = stateToIdCallback(initialState); + + initialState->setId(id); + return {id}; + } + + template<typename ValueType, typename StateType> + void DftNextStateGenerator<ValueType, StateType>::load(storm::storage::BitVector const& state) { + // Load the state from bitvector + size_t id = 0; //TODO Matthias: set correct id + this->state = std::make_shared<storm::storage::DFTState<ValueType>>(state, mDft, mStateGenerationInfo, id); + } + + template<typename ValueType, typename StateType> + void DftNextStateGenerator<ValueType, StateType>::load(DFTStatePointer const& state) { + // Store a pointer to the state itself, because we need to be able to access it when expanding it. + this->state = state; + } + + template<typename ValueType, typename StateType> + StateBehavior<ValueType, StateType> DftNextStateGenerator<ValueType, StateType>::expand(StateToIdCallback const& stateToIdCallback) { + STORM_LOG_TRACE("Explore state: " << mDft.getStateString(state)); + + // Prepare the result, in case we return early. + StateBehavior<ValueType, StateType> result; + + // Initialization + bool hasDependencies = state->nrFailableDependencies() > 0; + size_t failableCount = hasDependencies ? state->nrFailableDependencies() : state->nrFailableBEs(); + size_t currentFailable = 0; + Choice<ValueType, StateType> choice(0, !hasDependencies); + + // Check for absorbing state + if (mDft.hasFailed(state) || mDft.isFailsafe(state) || failableCount == 0) { + // Add self loop + choice.addProbability(state->getId(), storm::utility::one<ValueType>()); + STORM_LOG_TRACE("Added self loop for " << state->getId()); + // No further exploration required + result.addChoice(std::move(choice)); + result.setExpanded(); + return result; + } + + // Let BE fail + while (currentFailable < failableCount) { + STORM_LOG_ASSERT(!mDft.hasFailed(state), "Dft has failed."); + + // Construct new state as copy from original one + DFTStatePointer newState = state->copy(); + std::pair<std::shared_ptr<storm::storage::DFTBE<ValueType> const>, bool> nextBEPair = newState->letNextBEFail(currentFailable); + std::shared_ptr<storm::storage::DFTBE<ValueType> const>& nextBE = nextBEPair.first; + STORM_LOG_ASSERT(nextBE, "NextBE is null."); + STORM_LOG_ASSERT(nextBEPair.second == hasDependencies, "Failure due to dependencies does not match."); + STORM_LOG_TRACE("With the failure of: " << nextBE->name() << " [" << nextBE->id() << "] in " << mDft.getStateString(state)); + + // Propagate + storm::storage::DFTStateSpaceGenerationQueues<ValueType> queues; + + // Propagate failure + for (DFTGatePointer parent : nextBE->parents()) { + if (newState->isOperational(parent->id())) { + queues.propagateFailure(parent); + } + } + // Propagate failures + while (!queues.failurePropagationDone()) { + DFTGatePointer next = queues.nextFailurePropagation(); + next->checkFails(*newState, queues); + newState->updateFailableDependencies(next->id()); + } + + // Check restrictions + for (DFTRestrictionPointer restr : nextBE->restrictions()) { + queues.checkRestrictionLater(restr); + } + // Check restrictions + while(!queues.restrictionChecksDone()) { + DFTRestrictionPointer next = queues.nextRestrictionCheck(); + next->checkFails(*newState, queues); + newState->updateFailableDependencies(next->id()); + } + + if(newState->isInvalid()) { + // Continue with next possible state + ++currentFailable; + continue; + } + + // Get the id of the successor state + StateType newStateId; + if (newState->hasFailed(mDft.getTopLevelIndex()) && mergeFailedStates) { + // Use unique failed state + newStateId = mergeFailedStateId; + } else { + // Propagate failsafe + while (!queues.failsafePropagationDone()) { + DFTGatePointer next = queues.nextFailsafePropagation(); + next->checkFailsafe(*newState, queues); + } + + // Propagate dont cares + while (enableDC && !queues.dontCarePropagationDone()) { + DFTElementPointer next = queues.nextDontCarePropagation(); + next->checkDontCareAnymore(*newState, queues); + } + + // Update failable dependencies + newState->updateFailableDependencies(nextBE->id()); + newState->updateDontCareDependencies(nextBE->id()); + + // Add new state + newStateId = stateToIdCallback(newState); + } + + // Set transitions + if (hasDependencies) { + // Failure is due to dependency -> add non-deterministic choice + ValueType probability = mDft.getDependency(state->getDependencyId(currentFailable))->probability(); + choice.addProbability(newStateId, probability); + STORM_LOG_TRACE("Added transition to " << newStateId << " with probability " << probability); + + if (!storm::utility::isOne(probability)) { + // Add transition to state where dependency was unsuccessful + DFTStatePointer unsuccessfulState = state->copy(); + unsuccessfulState->letDependencyBeUnsuccessful(currentFailable); + // Add state + StateType unsuccessfulStateId = stateToIdCallback(unsuccessfulState); + ValueType remainingProbability = storm::utility::one<ValueType>() - probability; + choice.addProbability(unsuccessfulStateId, remainingProbability); + STORM_LOG_TRACE("Added transition to " << unsuccessfulStateId << " with remaining probability " << remainingProbability); + } + result.addChoice(std::move(choice)); + } else { + // Failure is due to "normal" BE failure + // Set failure rate according to activation + bool isActive = true; + if (mDft.hasRepresentant(nextBE->id())) { + // Active must be checked for the state we are coming from as this state is responsible for the rate + isActive = state->isActive(mDft.getRepresentant(nextBE->id())); + } + ValueType rate = isActive ? nextBE->activeFailureRate() : nextBE->passiveFailureRate(); + STORM_LOG_ASSERT(!storm::utility::isZero(rate), "Rate is 0."); + choice.addProbability(newStateId, rate); + STORM_LOG_TRACE("Added transition to " << newStateId << " with " << (isActive ? "active" : "passive") << " failure rate " << rate); + } + + ++currentFailable; + } // end while failing BE + + if (!hasDependencies) { + // Add all rates as one choice + result.addChoice(std::move(choice)); + } + + STORM_LOG_TRACE("Finished exploring state: " << mDft.getStateString(state)); + result.setExpanded(); + return result; + } + + template<typename ValueType, typename StateType> + StateBehavior<ValueType, StateType> DftNextStateGenerator<ValueType, StateType>::createMergeFailedState(StateToIdCallback const& stateToIdCallback) { + STORM_LOG_ASSERT(mergeFailedStates, "No unique failed state used."); + // Introduce explicit fail state + DFTStatePointer failedState = std::make_shared<storm::storage::DFTState<ValueType>>(mDft, mStateGenerationInfo, 0); + mergeFailedStateId = stateToIdCallback(failedState); + STORM_LOG_TRACE("Introduce fail state with id: " << mergeFailedStateId); + + // Add self loop + Choice<ValueType, StateType> choice(0, true); + choice.addProbability(mergeFailedStateId, storm::utility::one<ValueType>()); + + // No further exploration required + StateBehavior<ValueType, StateType> result; + result.addChoice(std::move(choice)); + result.setExpanded(); + return result; + } + + template class DftNextStateGenerator<double>; + +#ifdef STORM_HAVE_CARL + template class DftNextStateGenerator<storm::RationalFunction>; +#endif + } +} diff --git a/src/generator/DftNextStateGenerator.h b/src/generator/DftNextStateGenerator.h new file mode 100644 index 000000000..899bcb2c6 --- /dev/null +++ b/src/generator/DftNextStateGenerator.h @@ -0,0 +1,73 @@ +#ifndef STORM_GENERATOR_DFTNEXTSTATEGENERATOR_H_ +#define STORM_GENERATOR_DFTNEXTSTATEGENERATOR_H_ + +#include "src/generator/NextStateGenerator.h" +#include "src/storage/dft/DFT.h" + +#include "src/utility/ConstantsComparator.h" + +namespace storm { + namespace generator { + + /*! + * Next state generator for DFTs. + */ + template<typename ValueType, typename StateType = uint32_t> + class DftNextStateGenerator { + // TODO: inherit from NextStateGenerator + + using DFTStatePointer = std::shared_ptr<storm::storage::DFTState<ValueType>>; + using DFTElementPointer = std::shared_ptr<storm::storage::DFTElement<ValueType>>; + using DFTGatePointer = std::shared_ptr<storm::storage::DFTGate<ValueType>>; + using DFTRestrictionPointer = std::shared_ptr<storm::storage::DFTRestriction<ValueType>>; + + public: + typedef std::function<StateType (DFTStatePointer const&)> StateToIdCallback; + + DftNextStateGenerator(storm::storage::DFT<ValueType> const& dft, storm::storage::DFTStateGenerationInfo const& stateGenerationInfo, bool enableDC, bool mergeFailedStates); + + bool isDeterministicModel() const; + std::vector<StateType> getInitialStates(StateToIdCallback const& stateToIdCallback); + + void load(storm::storage::BitVector const& state); + void load(DFTStatePointer const& state); + StateBehavior<ValueType, StateType> expand(StateToIdCallback const& stateToIdCallback); + + /*! + * Create unique failed state. + * + * @param stateToIdCallback Callback for state. The callback should just return the id and not use the state. + * + * @return Behavior of state. + */ + StateBehavior<ValueType, StateType> createMergeFailedState(StateToIdCallback const& stateToIdCallback); + + private: + + // The dft used for the generation of next states. + storm::storage::DFT<ValueType> const& mDft; + + // General information for the state generation. + storm::storage::DFTStateGenerationInfo const& mStateGenerationInfo; + + // Current state + DFTStatePointer state; + + // Flag indicating if dont care propagation is enabled. + bool enableDC; + + // Flag indication if all failed states should be merged into one. + bool mergeFailedStates = true; + + // Id of the merged failed state + StateType mergeFailedStateId = 0; + + // Flag indicating if the model is deterministic. + bool deterministicModel = false; + + }; + + } +} + +#endif /* STORM_GENERATOR_DFTNEXTSTATEGENERATOR_H_ */ diff --git a/src/generator/StateBehavior.cpp b/src/generator/StateBehavior.cpp index ecd770941..c8ba1f10d 100644 --- a/src/generator/StateBehavior.cpp +++ b/src/generator/StateBehavior.cpp @@ -69,7 +69,8 @@ namespace storm { std::size_t StateBehavior<ValueType, StateType>::getNumberOfChoices() const { return choices.size(); } - + + template class StateBehavior<double>; #ifdef STORM_HAVE_CARL diff --git a/src/modelchecker/DFTAnalyser.h b/src/modelchecker/DFTAnalyser.h deleted file mode 100644 index 24c1068ba..000000000 --- a/src/modelchecker/DFTAnalyser.h +++ /dev/null @@ -1,191 +0,0 @@ -#pragma once - -#include "logic/Formula.h" -#include "parser/DFTGalileoParser.h" -#include "builder/ExplicitDFTModelBuilder.h" -#include "modelchecker/results/CheckResult.h" -#include "utility/storm.h" -#include "storage/dft/DFTIsomorphism.h" -#include "utility/bitoperations.h" - - -#include <chrono> - -template<typename ValueType> -class DFTAnalyser { - - std::chrono::duration<double> buildingTime = std::chrono::duration<double>::zero(); - std::chrono::duration<double> explorationTime = std::chrono::duration<double>::zero(); - std::chrono::duration<double> bisimulationTime = std::chrono::duration<double>::zero(); - std::chrono::duration<double> modelCheckingTime = std::chrono::duration<double>::zero(); - std::chrono::duration<double> totalTime = std::chrono::duration<double>::zero(); - ValueType checkResult = storm::utility::zero<ValueType>(); -public: - void check(storm::storage::DFT<ValueType> const& origDft , std::shared_ptr<const storm::logic::Formula> const& formula, bool symred = true, bool allowModularisation = true, bool enableDC = true) { - - // Building DFT - std::chrono::high_resolution_clock::time_point totalStart = std::chrono::high_resolution_clock::now(); - - - // Optimizing DFT - storm::storage::DFT<ValueType> dft = origDft.optimize(); - checkResult = checkHelper(dft, formula, symred, allowModularisation, enableDC); - totalTime = std::chrono::high_resolution_clock::now() - totalStart; - - - - } - -private: - ValueType checkHelper(storm::storage::DFT<ValueType> const& dft , std::shared_ptr<const storm::logic::Formula> const& formula, bool symred = true, bool allowModularisation = true, bool enableDC = true) { - STORM_LOG_TRACE("Check helper called"); - bool invResults = false; - std::vector<storm::storage::DFT<ValueType>> dfts = {dft}; - std::vector<ValueType> res; - size_t nrK = 0; // K out of M - size_t nrM = 0; // K out of M - - if(allowModularisation) { - if(dft.topLevelType() == storm::storage::DFTElementType::AND) { - STORM_LOG_TRACE("top modularisation called AND"); - dfts = dft.topModularisation(); - STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); - nrK = dfts.size(); - nrM = dfts.size(); - } - if(dft.topLevelType() == storm::storage::DFTElementType::OR) { - STORM_LOG_TRACE("top modularisation called OR"); - dfts = dft.topModularisation(); - STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); - nrK = 0; - nrM = dfts.size(); - invResults = true; - } - if(dft.topLevelType() == storm::storage::DFTElementType::VOT) { - STORM_LOG_TRACE("top modularisation called VOT"); - dfts = dft.topModularisation(); - STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); - nrK = std::static_pointer_cast<storm::storage::DFTVot<ValueType> const>(dft.getTopLevelGate())->threshold(); - nrM = dfts.size(); - if(nrK <= nrM/2) { - nrK -= 1; - invResults = true; - } - } - if(dfts.size() > 1) { - STORM_LOG_TRACE("Recursive CHECK Call"); - for(auto const ft : dfts) { - res.push_back(checkHelper(ft, formula, symred, allowModularisation)); - } - } - } - if(res.empty()) { - // Model based modularisation. - auto const& models = buildMarkovModels(dfts, formula, symred, enableDC); - - - for (auto const& model : models) { - // Model checking - STORM_LOG_INFO("Model checking..."); - std::chrono::high_resolution_clock::time_point modelCheckingStart = std::chrono::high_resolution_clock::now(); - std::unique_ptr<storm::modelchecker::CheckResult> result(storm::verifySparseModel(model, formula)); - STORM_LOG_INFO("Model checking done."); - STORM_LOG_ASSERT(result, "Result does not exist."); - result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); - modelCheckingTime += std::chrono::high_resolution_clock::now() - modelCheckingStart; - res.push_back(result->asExplicitQuantitativeCheckResult<ValueType>().getValueMap().begin()->second); - } - } - if(nrM <= 1) { - // No modularisation done. - STORM_LOG_ASSERT(res.size()==1, "Result size not 1."); - return res[0]; - } - - STORM_LOG_TRACE("Combining all results... K=" << nrK << "; M=" << nrM << "; invResults=" << (invResults?"On":"Off")); - ValueType result = storm::utility::zero<ValueType>(); - int limK = invResults ? -1 : nrM+1; - int chK = invResults ? -1 : 1; - for(int cK = nrK; cK != limK; cK += chK ) { - STORM_LOG_ASSERT(cK >= 0, "ck negative."); - size_t permutation = smallestIntWithNBitsSet(static_cast<size_t>(cK)); - do { - STORM_LOG_TRACE("Permutation="<<permutation); - ValueType permResult = storm::utility::one<ValueType>(); - for(size_t i = 0; i < res.size(); ++i) { - if(permutation & (1 << i)) { - permResult *= res[i]; - } else { - permResult *= storm::utility::one<ValueType>() - res[i]; - } - } - STORM_LOG_TRACE("Result for permutation:"<<permResult); - permutation = nextBitPermutation(permutation); - result += permResult; - } while(permutation < (1 << nrM) && permutation != 0); - } - if(invResults) { - return storm::utility::one<ValueType>() - result; - } - return result; - } - - - std::vector<std::shared_ptr<storm::models::sparse::Model<ValueType>>> buildMarkovModels(std::vector<storm::storage::DFT<ValueType>> const& dfts, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool enableDC) { - std::vector<std::shared_ptr<storm::models::sparse::Model<ValueType>>> models; - for(auto& dft : dfts) { - std::chrono::high_resolution_clock::time_point buildingStart = std::chrono::high_resolution_clock::now(); - - std::map<size_t, std::vector<std::vector<size_t>>> emptySymmetry; - storm::storage::DFTIndependentSymmetries symmetries(emptySymmetry); - if(symred) { - auto colouring = dft.colourDFT(); - symmetries = dft.findSymmetries(colouring); - STORM_LOG_INFO("Found " << symmetries.groups.size() << " symmetries."); - STORM_LOG_TRACE("Symmetries: " << std::endl << symmetries); - } - std::chrono::high_resolution_clock::time_point buildingEnd = std::chrono::high_resolution_clock::now(); - buildingTime += buildingEnd - buildingStart; - - // Building Markov Automaton - STORM_LOG_INFO("Building Model..."); - storm::builder::ExplicitDFTModelBuilder<ValueType> builder(dft, symmetries, enableDC); - typename storm::builder::ExplicitDFTModelBuilder<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula - std::shared_ptr<storm::models::sparse::Model<ValueType>> model = builder.buildModel(labeloptions); - //model->printModelInformationToStream(std::cout); - STORM_LOG_INFO("No. states (Explored): " << model->getNumberOfStates()); - STORM_LOG_INFO("No. transitions (Explored): " << model->getNumberOfTransitions()); - std::chrono::high_resolution_clock::time_point explorationEnd = std::chrono::high_resolution_clock::now(); - explorationTime += explorationEnd -buildingEnd; - - // Bisimulation - if (model->isOfType(storm::models::ModelType::Ctmc) && storm::settings::getModule<storm::settings::modules::GeneralSettings>().isBisimulationSet()) { - STORM_LOG_INFO("Bisimulation..."); - model = storm::performDeterministicSparseBisimulationMinimization<storm::models::sparse::Ctmc<ValueType>>(model->template as<storm::models::sparse::Ctmc<ValueType>>(), {formula}, storm::storage::BisimulationType::Weak)->template as<storm::models::sparse::Ctmc<ValueType>>(); - //model->printModelInformationToStream(std::cout); - } - STORM_LOG_INFO("No. states (Bisimulation): " << model->getNumberOfStates()); - STORM_LOG_INFO("No. transitions (Bisimulation): " << model->getNumberOfTransitions()); - bisimulationTime += std::chrono::high_resolution_clock::now() - explorationEnd; - models.push_back(model); - } - return models; - } - -public: - void printTimings(std::ostream& os = std::cout) { - os << "Times:" << std::endl; - os << "Building:\t" << buildingTime.count() << std::endl; - os << "Exploration:\t" << explorationTime.count() << std::endl; - os << "Bisimulation:\t" << bisimulationTime.count() << std::endl; - os << "Modelchecking:\t" << modelCheckingTime.count() << std::endl; - os << "Total:\t\t" << totalTime.count() << std::endl; - } - - void printResult(std::ostream& os = std::cout) { - - os << "Result: ["; - os << checkResult << "]" << std::endl; - } - -}; diff --git a/src/modelchecker/dft/DFTASFChecker.cpp b/src/modelchecker/dft/DFTASFChecker.cpp new file mode 100644 index 000000000..ac28ba5cf --- /dev/null +++ b/src/modelchecker/dft/DFTASFChecker.cpp @@ -0,0 +1,413 @@ +#include "DFTASFChecker.h" +#include <string> + +namespace storm { + + namespace modelchecker { + + /* + * Variable[VarIndex] is the maximum of the others + */ + class IsMaximum : public DFTConstraint { + public: + IsMaximum(uint64_t varIndex, std::vector<uint64_t> const& varIndices) : varIndex(varIndex), varIndices(varIndices) { + + } + + virtual ~IsMaximum() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(and "; + // assert it is largereq than all values. + for (auto const& ovi : varIndices) { + sstr << "(>= " << varNames.at(varIndex) << " " << varNames.at(ovi) << ") "; + } + // assert it is one of the values. + sstr << "(or "; + for (auto const& ovi : varIndices) { + sstr << "(= " << varNames.at(varIndex) << " " << varNames.at(ovi) << ") "; + } + sstr << ")"; // end of the or + sstr << ")"; // end outer and. + return sstr.str(); + + } + + + private: + uint64_t varIndex; + std::vector<uint64_t> varIndices; + }; + /* + * First is the minimum of the others + */ + class IsMinimum : public DFTConstraint { + public: + IsMinimum(uint64_t varIndex, std::vector<uint64_t> const& varIndices) : varIndex(varIndex), varIndices(varIndices) { + + } + + + virtual ~IsMinimum() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(and "; + // assert it is smallereq than all values. + for (auto const& ovi : varIndices) { + sstr << "(<= " << varNames.at(varIndex) << " " << varNames.at(ovi) << ") "; + } + // assert it is one of the values. + sstr << "(or "; + for (auto const& ovi : varIndices) { + sstr << "(= " << varNames.at(varIndex) << " " << varNames.at(ovi) << ") "; + } + sstr << ")"; // end of the or + sstr << ")"; // end outer and. + return sstr.str(); + + } + + + private: + uint64_t varIndex; + std::vector<uint64_t> varIndices; + }; + + class BetweenValues : public DFTConstraint { + public: + BetweenValues(uint64_t varIndex, uint64_t lower, uint64_t upper) : varIndex(varIndex), upperBound(upper) , lowerBound(lower) { + + } + virtual ~BetweenValues() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(and "; + sstr << "(>= " << varNames.at(varIndex) << " " << lowerBound << ")"; + sstr << "(<= " << varNames.at(varIndex) << " " << upperBound << ")"; + sstr << ")"; + return sstr.str(); + } + + private: + uint64_t varIndex; + uint64_t upperBound; + uint64_t lowerBound; + }; + + class And : public DFTConstraint { + public: + And(std::vector<std::shared_ptr<DFTConstraint>> const& constraints) : constraints(constraints) {} + + virtual ~And() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(and"; + for(auto const& c : constraints) { + sstr << " " << c->toSmtlib2(varNames); + } + sstr << ")"; + return sstr.str(); + } + private: + std::vector<std::shared_ptr<DFTConstraint>> constraints; + + }; + + class Iff : public DFTConstraint { + public: + Iff(std::shared_ptr<DFTConstraint> l, std::shared_ptr<DFTConstraint> r) : lhs(l), rhs(r) { + + } + + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(= " << lhs->toSmtlib2(varNames) << " " << rhs->toSmtlib2(varNames) << ")"; + return sstr.str(); + } + + private: + std::shared_ptr<DFTConstraint> lhs; + std::shared_ptr<DFTConstraint> rhs; + }; + + class IsConstantValue : public DFTConstraint { + public: + IsConstantValue(uint64_t varIndex, uint64_t val) : varIndex(varIndex), value(val) { + + } + + virtual ~IsConstantValue() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + assert(varIndex < varNames.size()); + sstr << "(= " << varNames.at(varIndex) << " " << value << ")"; + return sstr.str(); + } + + private: + uint64_t varIndex; + uint64_t value; + }; + + class IsEqual : public DFTConstraint { + public: + IsEqual(uint64_t varIndex1, uint64_t varIndex2) :var1Index(varIndex1), var2Index(varIndex2) { + + } + + virtual ~IsEqual() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + return "(= " + varNames.at(var1Index) + " " + varNames.at(var2Index) + " )"; + } + + private: + uint64_t var1Index; + uint64_t var2Index; + }; + + class IsLEqual : public DFTConstraint { + public: + IsLEqual(uint64_t varIndex1, uint64_t varIndex2) :var1Index(varIndex1), var2Index(varIndex2) { + + } + + virtual ~IsLEqual() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + return "(<= " + varNames.at(var1Index) + " " + varNames.at(var2Index) + " )"; + } + + private: + uint64_t var1Index; + uint64_t var2Index; + }; + + class PairwiseDifferent : public DFTConstraint { + public: + PairwiseDifferent(std::vector<uint64_t> const& indices) : varIndices(indices) { + + } + virtual ~PairwiseDifferent() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(distinct"; + // for(uint64_t i = 0; i < varIndices.size(); ++i) { + // for(uint64_t j = i + 1; j < varIndices.size(); ++j) { + // sstr << "()"; + // } + // } + for (auto const& varIndex : varIndices) { + sstr << " " << varNames.at(varIndex); + } + sstr << ")"; + return sstr.str(); + } + + private: + std::vector<uint64_t> varIndices; + }; + + class Sorted : public DFTConstraint { + public: + Sorted(std::vector<uint64_t> varIndices) : varIndices(varIndices) { + + } + + virtual ~Sorted() { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(and "; + for(uint64_t i = 1; i < varIndices.size(); ++i) { + sstr << "(<= " << varNames.at(varIndices.at(i-1)) << " " << varNames.at(varIndices.at(i)) << ")"; + } + sstr << ") "; + return sstr.str(); + } + + + private: + std::vector<uint64_t> varIndices; + }; + + class IfThenElse : public DFTConstraint { + public: + IfThenElse(std::shared_ptr<DFTConstraint> ifC, std::shared_ptr<DFTConstraint> thenC, std::shared_ptr<DFTConstraint> elseC) : ifConstraint(ifC), thenConstraint(thenC), elseConstraint(elseC) { + + } + + std::string toSmtlib2(std::vector<std::string> const& varNames) const override { + std::stringstream sstr; + sstr << "(ite " << ifConstraint->toSmtlib2(varNames) << " " << thenConstraint->toSmtlib2(varNames) << " " << elseConstraint->toSmtlib2(varNames) << ")"; + return sstr.str(); + } + + private: + std::shared_ptr<DFTConstraint> ifConstraint; + std::shared_ptr<DFTConstraint> thenConstraint; + std::shared_ptr<DFTConstraint> elseConstraint; + }; + + + DFTASFChecker::DFTASFChecker(storm::storage::DFT<double> const& dft) : dft(dft) { + // Intentionally left empty. + } + + uint64_t DFTASFChecker::getClaimVariableIndex(uint64_t spare, uint64_t child) const { + return claimVariables.at(SpareAndChildPair(spare, child)); + } + + void DFTASFChecker::convert() { + + std::vector<uint64_t> beVariables; + // Convert all elements + for (size_t i = 0; i < dft.nrElements(); ++i) { + std::shared_ptr<storm::storage::DFTElement<ValueType> const> element = dft.getElement(i); + varNames.push_back("t_" + element->name()); + timePointVariables.emplace(i, varNames.size() - 1); + switch (element->type()) { + case storm::storage::DFTElementType::BE: + beVariables.push_back(varNames.size() - 1); + break; + case storm::storage::DFTElementType::SPARE: + { + auto spare = std::static_pointer_cast<storm::storage::DFTSpare<double> const>(element); + for( auto const& spareChild : spare->children()) { + varNames.push_back("c_" + element->name() + "_" + spareChild->name()); + claimVariables.emplace(SpareAndChildPair(element->id(), spareChild->id()), varNames.size() - 1); + } + break; + } + default: + break; + } + } + + // BE + constraints.push_back(std::make_shared<PairwiseDifferent>(beVariables)); + constraints.back()->setDescription("No two BEs fail at the same time"); + for( auto const& beV : beVariables) { + constraints.push_back(std::make_shared<BetweenValues>(beV, 1, dft.nrBasicElements())); + } + + // Claim variables + for (auto const& csvV : claimVariables) { + constraints.push_back(std::make_shared<BetweenValues>(csvV.second, 0, dft.nrBasicElements() + 1)); + } + + for (size_t i = 0; i < dft.nrElements(); ++i) { + std::shared_ptr<storm::storage::DFTElement<ValueType> const> element = dft.getElement(i); + if(element->isSpareGate()) { + auto spare = std::static_pointer_cast<storm::storage::DFTSpare<double> const>(element); + auto const& spareChild = spare->children().front(); + constraints.push_back(std::make_shared<IsConstantValue>(getClaimVariableIndex(spare->id(), spareChild->id()), 0)); + constraints.back()->setDescription("Spare " + spare->name() + " claims first child"); + + } + } + + + for (size_t i = 0; i < dft.nrElements(); ++i) { + std::vector<uint64_t> childVarIndices; + if (dft.isGate(i)) { + std::shared_ptr<storm::storage::DFTGate<ValueType> const> gate = dft.getGate(i); + for (auto const& child : gate->children()) { + childVarIndices.push_back(timePointVariables.at(child->id())); + } + } + + std::shared_ptr<storm::storage::DFTElement<ValueType> const> element = dft.getElement(i); + switch (element->type()) { + case storm::storage::DFTElementType::AND: + constraints.push_back(std::make_shared<IsMaximum>(timePointVariables.at(i), childVarIndices)); + constraints.back()->setDescription("And gate"); + break; + case storm::storage::DFTElementType::OR: + constraints.push_back(std::make_shared<IsMinimum>(timePointVariables.at(i), childVarIndices)); + constraints.back()->setDescription("Or gate"); + break; + case storm::storage::DFTElementType::PAND: + constraints.push_back(std::make_shared<IfThenElse>(std::make_shared<Sorted>(childVarIndices), std::make_shared<IsEqual>(timePointVariables.at(i), timePointVariables.at(childVarIndices.back())), std::make_shared<IsConstantValue>(timePointVariables.at(i), dft.nrBasicElements()+1))); + constraints.back()->setDescription("Pand gate"); + case storm::storage::DFTElementType::SPARE: + { + auto spare = std::static_pointer_cast<storm::storage::DFTSpare<double> const>(element); + auto const& children = spare->children(); + + constraints.push_back(std::make_shared<Iff>(std::make_shared<IsLEqual>(getClaimVariableIndex(spare->id(), children.back()->id()), childVarIndices.back()), std::make_shared<IsEqual>(timePointVariables.at(i), childVarIndices.back()))); + constraints.back()->setDescription("Last child & claimed -> spare fails"); + std::vector<std::shared_ptr<DFTConstraint>> ifcs; + uint64_t xv = childVarIndices.at(childVarIndices.size()-2); // if v is the child, xv is the time x fails. + uint64_t csn = getClaimVariableIndex(spare->id(), children.back()->id()); // csn is the moment the spare claims the next child + uint64_t xn = childVarIndices.back(); // xn is the moment the next child fails + ifcs.push_back(std::make_shared<IsLEqual>(xv, xn)); + for(auto const& otherSpare : children.back()->parents()) { + if(otherSpare->id() == i) { continue; }// not a OTHER spare. + ifcs.push_back(std::make_shared<IsConstantValue>(getClaimVariableIndex(otherSpare->id(), children.back()->id()), dft.nrBasicElements()+1)); + } + std::shared_ptr<DFTConstraint> ite = std::make_shared<IfThenElse>(std::make_shared<And>(ifcs), std::make_shared<IsEqual>(csn, xv), std::make_shared<IsEqual>(timePointVariables.at(i), xv)); + constraints.push_back(std::make_shared<Iff>(std::make_shared<IsLEqual>(getClaimVariableIndex(spare->id(), children.at(children.size()-2)->id()), xv), ite)); + constraints.back()->setDescription("1 but last child & claimed"); + + for( uint64_t j = 0; j < children.size() - 2; ++j) { + uint64_t currChild = children.size() - 1 - j; + + } + break; + } + default: + break; + } + } + + constraints.push_back(std::make_shared<IsConstantValue>(dft.getTopLevelIndex(), dft.nrBasicElements()+1)); + } + + void DFTASFChecker::toFile(std::string const& filename) { + std::ofstream ofs; + std::cout << "Writing to " << filename << std::endl; + ofs.open(filename); + ofs << "; time point variables" << std::endl; + for (auto const& timeVarEntry : timePointVariables) { + ofs << "(declare-fun " << varNames[timeVarEntry.second] << "() Int)" << std::endl; + } + ofs << "; claim variables" << std::endl; + for (auto const& claimVarEntry : claimVariables) { + ofs << "(declare-fun " << varNames[claimVarEntry.second] << "() Int)" << std::endl; + } + for (auto const& constraint : constraints) { + ofs << "; " << constraint->description() << std::endl; + ofs << "(assert " << constraint->toSmtlib2(varNames) << ")" << std::endl; + } + ofs << "(check-sat)" << std::endl; + ofs.close(); + + } + } +} \ No newline at end of file diff --git a/src/modelchecker/dft/DFTASFChecker.h b/src/modelchecker/dft/DFTASFChecker.h new file mode 100644 index 000000000..73a323b03 --- /dev/null +++ b/src/modelchecker/dft/DFTASFChecker.h @@ -0,0 +1,61 @@ +#pragma once + +#include <string> +#include <vector> +#include <unordered_map> +#include "src/storage/dft/DFT.h" + + +namespace storm { + namespace modelchecker { + class DFTConstraint { + public: + virtual ~DFTConstraint() { + + } + + virtual std::string toSmtlib2(std::vector<std::string> const& varNames) const = 0; + virtual std::string description() const { return descript; } + + void setDescription(std::string const& descr) { + descript = descr; + } + + private: + std::string descript; + }; + + class SpareAndChildPair { + + public: + SpareAndChildPair(uint64_t spareIndex, uint64_t childIndex) : spareIndex(spareIndex), childIndex(childIndex) { + + } + + uint64_t spareIndex; + uint64_t childIndex; + }; + + bool operator<(SpareAndChildPair const& p1, SpareAndChildPair const& p2) { + return p1.spareIndex < p2.spareIndex || (p1.spareIndex == p2.spareIndex && p1.childIndex < p2.childIndex); + } + + class DFTASFChecker { + using ValueType = double; + + public: + DFTASFChecker(storm::storage::DFT<ValueType> const&); + void convert(); + void toFile(std::string const&); + + private: + uint64_t getClaimVariableIndex(uint64_t spareIndex, uint64_t childIndex) const; + + storm::storage::DFT<ValueType> const& dft; + std::vector<std::string> varNames; + std::unordered_map<uint64_t, uint64_t> timePointVariables; + std::vector<std::shared_ptr<DFTConstraint>> constraints; + std::map<SpareAndChildPair, uint64_t> claimVariables; + }; + } +} \ No newline at end of file diff --git a/src/modelchecker/dft/DFTModelChecker.cpp b/src/modelchecker/dft/DFTModelChecker.cpp new file mode 100644 index 000000000..748abe188 --- /dev/null +++ b/src/modelchecker/dft/DFTModelChecker.cpp @@ -0,0 +1,441 @@ +#include "src/modelchecker/dft/DFTModelChecker.h" + +#include "src/builder/ExplicitDFTModelBuilder.h" +#include "src/builder/ExplicitDFTModelBuilderApprox.h" +#include "src/builder/ParallelCompositionBuilder.h" +#include "src/storage/dft/DFTIsomorphism.h" +#include "src/settings/modules/DFTSettings.h" +#include "src/utility/bitoperations.h" + +namespace storm { + namespace modelchecker { + + template<typename ValueType> + DFTModelChecker<ValueType>::DFTModelChecker() { + checkResult = storm::utility::zero<ValueType>(); + } + + template<typename ValueType> + void DFTModelChecker<ValueType>::check(storm::storage::DFT<ValueType> const& origDft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool allowModularisation, bool enableDC, double approximationError) { + // Initialize + this->explorationTime = std::chrono::duration<double>::zero(); + this->buildingTime = std::chrono::duration<double>::zero(); + this->bisimulationTime = std::chrono::duration<double>::zero(); + this->modelCheckingTime = std::chrono::duration<double>::zero(); + this->totalTime = std::chrono::duration<double>::zero(); + this->approximationError = approximationError; + totalStart = std::chrono::high_resolution_clock::now(); + + // Optimizing DFT + storm::storage::DFT<ValueType> dft = origDft.optimize(); + + // TODO Matthias: check that all paths reach the target state for approximation + + // Checking DFT + if (formula->isProbabilityOperatorFormula() || !allowModularisation) { + checkResult = checkHelper(dft, formula, symred, allowModularisation, enableDC, approximationError); + } else { + std::shared_ptr<storm::models::sparse::Model<ValueType>> model = buildModelComposition(dft, formula, symred, allowModularisation, enableDC); + // Model checking + std::unique_ptr<storm::modelchecker::CheckResult> result = checkModel(model, formula); + result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); + checkResult = result->asExplicitQuantitativeCheckResult<ValueType>().getValueMap().begin()->second; + + } + this->totalTime = std::chrono::high_resolution_clock::now() - totalStart; + } + + template<typename ValueType> + typename DFTModelChecker<ValueType>::dft_result DFTModelChecker<ValueType>::checkHelper(storm::storage::DFT<ValueType> const& dft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool allowModularisation, bool enableDC, double approximationError) { + STORM_LOG_TRACE("Check helper called"); + bool modularisationPossible = allowModularisation; + + // Try modularisation + if(modularisationPossible) { + bool invResults = false; + std::vector<storm::storage::DFT<ValueType>> dfts; + size_t nrK = 0; // K out of M + size_t nrM = 0; // K out of M + + switch (dft.topLevelType()) { + case storm::storage::DFTElementType::AND: + STORM_LOG_TRACE("top modularisation called AND"); + dfts = dft.topModularisation(); + STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); + nrK = dfts.size(); + nrM = dfts.size(); + modularisationPossible = dfts.size() > 1; + break; + case storm::storage::DFTElementType::OR: + STORM_LOG_TRACE("top modularisation called OR"); + dfts = dft.topModularisation(); + STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); + nrK = 0; + nrM = dfts.size(); + invResults = true; + modularisationPossible = dfts.size() > 1; + break; + case storm::storage::DFTElementType::VOT: + STORM_LOG_TRACE("top modularisation called VOT"); + dfts = dft.topModularisation(); + STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); + nrK = std::static_pointer_cast<storm::storage::DFTVot<ValueType> const>(dft.getTopLevelGate())->threshold(); + nrM = dfts.size(); + if(nrK <= nrM/2) { + nrK -= 1; + invResults = true; + } + modularisationPossible = dfts.size() > 1; + break; + default: + // No static gate -> no modularisation applicable + modularisationPossible = false; + break; + } + + if(modularisationPossible) { + STORM_LOG_TRACE("Recursive CHECK Call"); + if (formula->isProbabilityOperatorFormula()) { + // Recursively call model checking + std::vector<ValueType> res; + for(auto const ft : dfts) { + dft_result ftResult = checkHelper(ft, formula, symred, true, enableDC, 0.0); + res.push_back(boost::get<ValueType>(ftResult)); + } + + // Combine modularisation results + STORM_LOG_TRACE("Combining all results... K=" << nrK << "; M=" << nrM << "; invResults=" << (invResults?"On":"Off")); + ValueType result = storm::utility::zero<ValueType>(); + int limK = invResults ? -1 : nrM+1; + int chK = invResults ? -1 : 1; + // WARNING: there is a bug for computing permutations with more than 32 elements + STORM_LOG_ASSERT(res.size() < 32, "Permutations work only for < 32 elements"); + for(int cK = nrK; cK != limK; cK += chK ) { + STORM_LOG_ASSERT(cK >= 0, "ck negative."); + size_t permutation = smallestIntWithNBitsSet(static_cast<size_t>(cK)); + do { + STORM_LOG_TRACE("Permutation="<<permutation); + ValueType permResult = storm::utility::one<ValueType>(); + for(size_t i = 0; i < res.size(); ++i) { + if(permutation & (1 << i)) { + permResult *= res[i]; + } else { + permResult *= storm::utility::one<ValueType>() - res[i]; + } + } + STORM_LOG_TRACE("Result for permutation:"<<permResult); + permutation = nextBitPermutation(permutation); + result += permResult; + } while(permutation < (1 << nrM) && permutation != 0); + } + if(invResults) { + result = storm::utility::one<ValueType>() - result; + } + return result; + } + } + } + + // If we are here, no modularisation was possible + STORM_LOG_ASSERT(!modularisationPossible, "Modularisation should not be possible."); + return checkDFT(dft, formula, symred, enableDC, approximationError); + } + + template<typename ValueType> + std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> DFTModelChecker<ValueType>::buildModelComposition(storm::storage::DFT<ValueType> const& dft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool allowModularisation, bool enableDC) { + STORM_LOG_TRACE("Build model via composition"); + // Use parallel composition for CTMCs for expected time + STORM_LOG_ASSERT(formula->isTimeOperatorFormula(), "Formula is not a time operator formula"); + bool modularisationPossible = allowModularisation; + + // Try modularisation + if(modularisationPossible) { + std::vector<storm::storage::DFT<ValueType>> dfts; + bool isAnd = true; + + switch (dft.topLevelType()) { + case storm::storage::DFTElementType::AND: + STORM_LOG_TRACE("top modularisation called AND"); + dfts = dft.topModularisation(); + STORM_LOG_TRACE("Modularisation into " << dfts.size() << " submodules."); + modularisationPossible = dfts.size() > 1; + isAnd = true; + break; + case storm::storage::DFTElementType::OR: + STORM_LOG_TRACE("top modularisation called OR"); + dfts = dft.topModularisation(); + STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); + modularisationPossible = dfts.size() > 1; + isAnd = false; + break; + case storm::storage::DFTElementType::VOT: + /*STORM_LOG_TRACE("top modularisation called VOT"); + dfts = dft.topModularisation(); + STORM_LOG_TRACE("Modularsation into " << dfts.size() << " submodules."); + std::static_pointer_cast<storm::storage::DFTVot<ValueType> const>(dft.getTopLevelGate())->threshold(); + */ + // TODO enable modularisation for voting gate + modularisationPossible = false; + break; + default: + // No static gate -> no modularisation applicable + modularisationPossible = false; + break; + } + + if(modularisationPossible) { + STORM_LOG_TRACE("Recursive CHECK Call"); + bool firstTime = true; + std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> composedModel; + for (auto const ft : dfts) { + STORM_LOG_INFO("Building Model via parallel composition..."); + std::chrono::high_resolution_clock::time_point checkingStart = std::chrono::high_resolution_clock::now(); + + // Find symmetries + std::map<size_t, std::vector<std::vector<size_t>>> emptySymmetry; + storm::storage::DFTIndependentSymmetries symmetries(emptySymmetry); + if(symred) { + auto colouring = ft.colourDFT(); + symmetries = ft.findSymmetries(colouring); + STORM_LOG_INFO("Found " << symmetries.groups.size() << " symmetries."); + STORM_LOG_TRACE("Symmetries: " << std::endl << symmetries); + } + + // Build a single CTMC + STORM_LOG_INFO("Building Model..."); + storm::builder::ExplicitDFTModelBuilderApprox<ValueType> builder(ft, symmetries, enableDC); + typename storm::builder::ExplicitDFTModelBuilderApprox<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula + builder.buildModel(labeloptions, 0, 0.0); + std::shared_ptr<storm::models::sparse::Model<ValueType>> model = builder.getModel(); + //model->printModelInformationToStream(std::cout); + STORM_LOG_INFO("No. states (Explored): " << model->getNumberOfStates()); + STORM_LOG_INFO("No. transitions (Explored): " << model->getNumberOfTransitions()); + explorationTime += std::chrono::high_resolution_clock::now() - checkingStart; + + STORM_LOG_THROW(model->isOfType(storm::models::ModelType::Ctmc), storm::exceptions::NotSupportedException, "Parallel composition only applicable for CTMCs"); + std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> ctmc = model->template as<storm::models::sparse::Ctmc<ValueType>>(); + + ctmc = storm::performDeterministicSparseBisimulationMinimization<storm::models::sparse::Ctmc<ValueType>>(ctmc, {formula}, storm::storage::BisimulationType::Weak)->template as<storm::models::sparse::Ctmc<ValueType>>(); + + if (firstTime) { + composedModel = ctmc; + firstTime = false; + } else { + composedModel = storm::builder::ParallelCompositionBuilder<ValueType>::compose(composedModel, ctmc, isAnd); + } + + // Apply bisimulation + std::chrono::high_resolution_clock::time_point bisimulationStart = std::chrono::high_resolution_clock::now(); + composedModel = storm::performDeterministicSparseBisimulationMinimization<storm::models::sparse::Ctmc<ValueType>>(composedModel, {formula}, storm::storage::BisimulationType::Weak)->template as<storm::models::sparse::Ctmc<ValueType>>(); + std::chrono::high_resolution_clock::time_point bisimulationEnd = std::chrono::high_resolution_clock::now(); + bisimulationTime += bisimulationEnd - bisimulationStart; + + STORM_LOG_INFO("No. states (Composed): " << composedModel->getNumberOfStates()); + STORM_LOG_INFO("No. transitions (Composed): " << composedModel->getNumberOfTransitions()); + if (composedModel->getNumberOfStates() <= 15) { + STORM_LOG_TRACE("Transition matrix: " << std::endl << composedModel->getTransitionMatrix()); + } else { + STORM_LOG_TRACE("Transition matrix: too big to print"); + } + + } + return composedModel; + } + } + + // If we are here, no composition was possible + STORM_LOG_ASSERT(!modularisationPossible, "Modularisation should not be possible."); + std::chrono::high_resolution_clock::time_point checkingStart = std::chrono::high_resolution_clock::now(); + // Find symmetries + std::map<size_t, std::vector<std::vector<size_t>>> emptySymmetry; + storm::storage::DFTIndependentSymmetries symmetries(emptySymmetry); + if(symred) { + auto colouring = dft.colourDFT(); + symmetries = dft.findSymmetries(colouring); + STORM_LOG_INFO("Found " << symmetries.groups.size() << " symmetries."); + STORM_LOG_TRACE("Symmetries: " << std::endl << symmetries); + } + // Build a single CTMC + STORM_LOG_INFO("Building Model..."); + + + storm::builder::ExplicitDFTModelBuilderApprox<ValueType> builder(dft, symmetries, enableDC); + typename storm::builder::ExplicitDFTModelBuilderApprox<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula + builder.buildModel(labeloptions, 0, 0.0); + std::shared_ptr<storm::models::sparse::Model<ValueType>> model = builder.getModel(); + //model->printModelInformationToStream(std::cout); + STORM_LOG_INFO("No. states (Explored): " << model->getNumberOfStates()); + STORM_LOG_INFO("No. transitions (Explored): " << model->getNumberOfTransitions()); + explorationTime += std::chrono::high_resolution_clock::now() - checkingStart; + STORM_LOG_THROW(model->isOfType(storm::models::ModelType::Ctmc), storm::exceptions::NotSupportedException, "Parallel composition only applicable for CTMCs"); + + return model->template as<storm::models::sparse::Ctmc<ValueType>>(); + } + + template<typename ValueType> + typename DFTModelChecker<ValueType>::dft_result DFTModelChecker<ValueType>::checkDFT(storm::storage::DFT<ValueType> const& dft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool enableDC, double approximationError) { + std::chrono::high_resolution_clock::time_point checkingStart = std::chrono::high_resolution_clock::now(); + + // Find symmetries + std::map<size_t, std::vector<std::vector<size_t>>> emptySymmetry; + storm::storage::DFTIndependentSymmetries symmetries(emptySymmetry); + if(symred) { + auto colouring = dft.colourDFT(); + symmetries = dft.findSymmetries(colouring); + STORM_LOG_INFO("Found " << symmetries.groups.size() << " symmetries."); + STORM_LOG_TRACE("Symmetries: " << std::endl << symmetries); + } + + if (approximationError > 0.0) { + // Comparator for checking the error of the approximation + storm::utility::ConstantsComparator<ValueType> comparator; + // Build approximate Markov Automata for lower and upper bound + approximation_result approxResult = std::make_pair(storm::utility::zero<ValueType>(), storm::utility::zero<ValueType>()); + std::shared_ptr<storm::models::sparse::Model<ValueType>> model; + storm::builder::ExplicitDFTModelBuilderApprox<ValueType> builder(dft, symmetries, enableDC); + typename storm::builder::ExplicitDFTModelBuilderApprox<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula + + bool probabilityFormula = formula->isProbabilityOperatorFormula(); + STORM_LOG_ASSERT((formula->isTimeOperatorFormula() && !probabilityFormula) || (!formula->isTimeOperatorFormula() && probabilityFormula), "Probability formula not initialized correctly"); + size_t iteration = 0; + do { + // Iteratively build finer models + std::chrono::high_resolution_clock::time_point explorationStart = std::chrono::high_resolution_clock::now(); + STORM_LOG_INFO("Building model..."); + // TODO Matthias refine model using existing model and MC results + builder.buildModel(labeloptions, iteration, approximationError); + std::chrono::high_resolution_clock::time_point explorationEnd = std::chrono::high_resolution_clock::now(); + explorationTime += explorationEnd - explorationStart; + + // TODO Matthias: possible to do bisimulation on approximated model and not on concrete one? + + // Build model for lower bound + STORM_LOG_INFO("Getting model for lower bound..."); + model = builder.getModelApproximation(probabilityFormula ? false : true); + // We only output the info from the lower bound as the info for the upper bound is the same + STORM_LOG_INFO("No. states: " << model->getNumberOfStates()); + STORM_LOG_INFO("No. transitions: " << model->getNumberOfTransitions()); + buildingTime += std::chrono::high_resolution_clock::now() - explorationEnd; + + // Check lower bound + std::unique_ptr<storm::modelchecker::CheckResult> result = checkModel(model, formula); + result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); + ValueType newResult = result->asExplicitQuantitativeCheckResult<ValueType>().getValueMap().begin()->second; + STORM_LOG_ASSERT(iteration == 0 || !comparator.isLess(newResult, approxResult.first), "New under-approximation " << newResult << " is smaller than old result " << approxResult.first); + approxResult.first = newResult; + + // Build model for upper bound + STORM_LOG_INFO("Getting model for upper bound..."); + explorationEnd = std::chrono::high_resolution_clock::now(); + model = builder.getModelApproximation(probabilityFormula ? true : false); + buildingTime += std::chrono::high_resolution_clock::now() - explorationEnd; + // Check upper bound + result = checkModel(model, formula); + result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); + newResult = result->asExplicitQuantitativeCheckResult<ValueType>().getValueMap().begin()->second; + STORM_LOG_ASSERT(iteration == 0 || !comparator.isLess(approxResult.second, newResult), "New over-approximation " << newResult << " is greater than old result " << approxResult.second); + approxResult.second = newResult; + + ++iteration; + STORM_LOG_INFO("Result after iteration " << iteration << ": (" << std::setprecision(10) << approxResult.first << ", " << approxResult.second << ")"); + totalTime = std::chrono::high_resolution_clock::now() - totalStart; + printTimings(); + STORM_LOG_THROW(!storm::utility::isInfinity<ValueType>(approxResult.first) && !storm::utility::isInfinity<ValueType>(approxResult.second), storm::exceptions::NotSupportedException, "Approximation does not work if result might be infinity."); + } while (!isApproximationSufficient(approxResult.first, approxResult.second, approximationError, probabilityFormula)); + + STORM_LOG_INFO("Finished approximation after " << iteration << " iteration" << (iteration > 1 ? "s." : ".")); + return approxResult; + } else { + // Build a single Markov Automaton + STORM_LOG_INFO("Building Model..."); + std::shared_ptr<storm::models::sparse::Model<ValueType>> model; + // TODO Matthias: use only one builder if everything works again + if (storm::settings::getModule<storm::settings::modules::DFTSettings>().isApproximationErrorSet()) { + storm::builder::ExplicitDFTModelBuilderApprox<ValueType> builder(dft, symmetries, enableDC); + typename storm::builder::ExplicitDFTModelBuilderApprox<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula + builder.buildModel(labeloptions, 0, 0.0); + model = builder.getModel(); + } else { + storm::builder::ExplicitDFTModelBuilder<ValueType> builder(dft, symmetries, enableDC); + typename storm::builder::ExplicitDFTModelBuilder<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula + model = builder.buildModel(labeloptions); + } + //model->printModelInformationToStream(std::cout); + STORM_LOG_INFO("No. states (Explored): " << model->getNumberOfStates()); + STORM_LOG_INFO("No. transitions (Explored): " << model->getNumberOfTransitions()); + explorationTime += std::chrono::high_resolution_clock::now() - checkingStart; + + // Model checking + std::unique_ptr<storm::modelchecker::CheckResult> result = checkModel(model, formula); + result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); + return result->asExplicitQuantitativeCheckResult<ValueType>().getValueMap().begin()->second; + } + } + + template<typename ValueType> + std::unique_ptr<storm::modelchecker::CheckResult> DFTModelChecker<ValueType>::checkModel(std::shared_ptr<storm::models::sparse::Model<ValueType>>& model, std::shared_ptr<const storm::logic::Formula> const& formula) { + // Bisimulation + std::chrono::high_resolution_clock::time_point bisimulationStart = std::chrono::high_resolution_clock::now(); + if (model->isOfType(storm::models::ModelType::Ctmc) && storm::settings::getModule<storm::settings::modules::GeneralSettings>().isBisimulationSet()) { + STORM_LOG_INFO("Bisimulation..."); + model = storm::performDeterministicSparseBisimulationMinimization<storm::models::sparse::Ctmc<ValueType>>(model->template as<storm::models::sparse::Ctmc<ValueType>>(), {formula}, storm::storage::BisimulationType::Weak)->template as<storm::models::sparse::Ctmc<ValueType>>(); + STORM_LOG_INFO("No. states (Bisimulation): " << model->getNumberOfStates()); + STORM_LOG_INFO("No. transitions (Bisimulation): " << model->getNumberOfTransitions()); + } + std::chrono::high_resolution_clock::time_point bisimulationEnd = std::chrono::high_resolution_clock::now(); + bisimulationTime += bisimulationEnd - bisimulationStart; + + // Check the model + STORM_LOG_INFO("Model checking..."); + std::unique_ptr<storm::modelchecker::CheckResult> result(storm::verifySparseModel(model, formula)); + STORM_LOG_INFO("Model checking done."); + STORM_LOG_ASSERT(result, "Result does not exist."); + modelCheckingTime += std::chrono::high_resolution_clock::now() - bisimulationEnd; + return result; + } + + template<typename ValueType> + bool DFTModelChecker<ValueType>::isApproximationSufficient(ValueType lowerBound, ValueType upperBound, double approximationError, bool relative) { + STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Approximation works only for double."); + } + + template<> + bool DFTModelChecker<double>::isApproximationSufficient(double lowerBound, double upperBound, double approximationError, bool relative) { + STORM_LOG_THROW(!std::isnan(lowerBound) && !std::isnan(upperBound), storm::exceptions::NotSupportedException, "Approximation does not work if result is NaN."); + if (relative) { + return upperBound - lowerBound <= approximationError; + } else { + return upperBound - lowerBound <= approximationError * (lowerBound + upperBound) / 2; + } + } + + template<typename ValueType> + void DFTModelChecker<ValueType>::printTimings(std::ostream& os) { + os << "Times:" << std::endl; + os << "Exploration:\t" << explorationTime.count() << std::endl; + os << "Building:\t" << buildingTime.count() << std::endl; + os << "Bisimulation:\t" << bisimulationTime.count() << std::endl; + os << "Modelchecking:\t" << modelCheckingTime.count() << std::endl; + os << "Total:\t\t" << totalTime.count() << std::endl; + } + + template<typename ValueType> + void DFTModelChecker<ValueType>::printResult(std::ostream& os) { + os << "Result: ["; + if (this->approximationError > 0.0) { + approximation_result result = boost::get<approximation_result>(checkResult); + os << "(" << result.first << ", " << result.second << ")"; + } else { + os << boost::get<ValueType>(checkResult); + } + os << "]" << std::endl; + } + + + template class DFTModelChecker<double>; + +#ifdef STORM_HAVE_CARL + template class DFTModelChecker<storm::RationalFunction>; +#endif + } +} diff --git a/src/modelchecker/dft/DFTModelChecker.h b/src/modelchecker/dft/DFTModelChecker.h new file mode 100644 index 000000000..67deb2af0 --- /dev/null +++ b/src/modelchecker/dft/DFTModelChecker.h @@ -0,0 +1,140 @@ +#ifndef STORM_MODELCHECKER_DFT_DFTMODELCHECKER_H_ + +#include "src/logic/Formula.h" +#include "src/modelchecker/results/CheckResult.h" +#include "src/storage/dft/DFT.h" +#include "src/utility/storm.h" + +#include <chrono> + + +namespace storm { + namespace modelchecker { + + /*! + * Analyser for DFTs. + */ + template<typename ValueType> + class DFTModelChecker { + + typedef std::pair<ValueType, ValueType> approximation_result; + typedef boost::variant<ValueType, approximation_result> dft_result; + + public: + + /*! + * Constructor. + */ + DFTModelChecker(); + + /*! + * Main method for checking DFTs. + * + * @param origDft Original DFT + * @param formula Formula to check for + * @param symred Flag indicating if symmetry reduction should be used + * @param allowModularisation Flag indication if modularisation is allowed + * @param enableDC Flag indicating if dont care propagation should be used + * @param approximationError Error allowed for approximation. Value 0 indicates no approximation + */ + void check(storm::storage::DFT<ValueType> const& origDft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred = true, bool allowModularisation = true, bool enableDC = true, double approximationError = 0.0); + + /*! + * Print timings of all operations to stream. + * + * @param os Output stream to write to. + */ + void printTimings(std::ostream& os = std::cout); + + /*! + * Print result to stream. + * + * @param os Output stream to write to. + */ + void printResult(std::ostream& os = std::cout); + + private: + + // Timing values + std::chrono::duration<double> buildingTime = std::chrono::duration<double>::zero(); + std::chrono::duration<double> explorationTime = std::chrono::duration<double>::zero(); + std::chrono::duration<double> bisimulationTime = std::chrono::duration<double>::zero(); + std::chrono::duration<double> modelCheckingTime = std::chrono::duration<double>::zero(); + std::chrono::duration<double> totalTime = std::chrono::duration<double>::zero(); + std::chrono::high_resolution_clock::time_point totalStart; + + // Model checking result + dft_result checkResult; + + // Allowed error bound for approximation + double approximationError; + + /*! + * Internal helper for model checking a DFT. + * + * @param dft DFT + * @param formula Formula to check for + * @param symred Flag indicating if symmetry reduction should be used + * @param allowModularisation Flag indication if modularisation is allowed + * @param enableDC Flag indicating if dont care propagation should be used + * @param approximationError Error allowed for approximation. Value 0 indicates no approximation + * + * @return Model checking result (or in case of approximation two results for lower and upper bound) + */ + dft_result checkHelper(storm::storage::DFT<ValueType> const& dft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool allowModularisation, bool enableDC, double approximationError); + + /*! + * Internal helper for building a CTMC from a DFT via parallel composition. + * + * @param dft DFT + * @param formula Formula to check for + * @param symred Flag indicating if symmetry reduction should be used + * @param allowModularisation Flag indication if modularisation is allowed + * @param enableDC Flag indicating if dont care propagation should be used + * @param approximationError Error allowed for approximation. Value 0 indicates no approximation + * + * @return CTMC representing the DFT + */ + std::shared_ptr<storm::models::sparse::Ctmc<ValueType>> buildModelComposition(storm::storage::DFT<ValueType> const& dft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool allowModularisation, bool enableDC); + + /*! + * Check model generated from DFT. + * + * @param dft The DFT + * @param formula Formula to check for + * @param symred Flag indicating if symmetry reduction should be used + * @param enableDC Flag indicating if dont care propagation should be used + * @param approximationError Error allowed for approximation. Value 0 indicates no approximation + * + * @return Model checking result + */ + dft_result checkDFT(storm::storage::DFT<ValueType> const& dft, std::shared_ptr<const storm::logic::Formula> const& formula, bool symred, bool enableDC, double approximationError = 0.0); + + /*! + * Check the given markov model for the given property. + * + * @param model Model to check + * @param formula Formula to check for + * + * @return Model checking result + */ + std::unique_ptr<storm::modelchecker::CheckResult> checkModel(std::shared_ptr<storm::models::sparse::Model<ValueType>>& model, std::shared_ptr<const storm::logic::Formula> const& formula); + + /*! + * Checks if the computed approximation is sufficient, i.e. + * upperBound - lowerBound <= approximationError * mean(lowerBound, upperBound). + * + * @param lowerBound The lower bound on the result. + * @param upperBound The upper bound on the result. + * @param approximationError The allowed error for approximating. + * @param relative Flag indicating if the error should be relative to 1 or + to the mean of lower and upper bound. + * + * @return True, if the approximation is sufficient. + */ + bool isApproximationSufficient(ValueType lowerBound, ValueType upperBound, double approximationError, bool relative); + + }; + } +} +#endif /* STORM_MODELCHECKER_DFT_DFTMODELCHECKER_H_ */ diff --git a/src/models/sparse/MarkovAutomaton.cpp b/src/models/sparse/MarkovAutomaton.cpp index e50d10bef..32525ad5e 100644 --- a/src/models/sparse/MarkovAutomaton.cpp +++ b/src/models/sparse/MarkovAutomaton.cpp @@ -69,8 +69,8 @@ namespace storm { std::unordered_map<std::string, RewardModelType>&& rewardModels, boost::optional<std::vector<LabelSet>>&& optionalChoiceLabeling) : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, std::move(transitionMatrix), std::move(stateLabeling), std::move(rewardModels), std::move(optionalChoiceLabeling)), markovianStates(markovianStates), exitRates(std::move(exitRates)), closed(this->checkIsClosed()) { - assert(probabilities); - assert(this->getTransitionMatrix().isProbabilistic()); + STORM_LOG_ASSERT(probabilities, "Matrix must be probabilistic."); + STORM_LOG_ASSERT(this->getTransitionMatrix().isProbabilistic(), "Matrix must be probabilistic."); } template <typename ValueType, typename RewardModelType> @@ -279,10 +279,10 @@ namespace storm { // Get number of choices in current state uint_fast64_t numberChoices = this->getTransitionMatrix().getRowGroupIndices()[state + 1] - this->getTransitionMatrix().getRowGroupIndices()[state]; if (isMarkovianState(state)) { - assert(numberChoices == 1); + STORM_LOG_ASSERT(numberChoices == 1, "Wrong number of choices for markovian state."); } if (numberChoices > 1) { - assert(isProbabilisticState(state)); + STORM_LOG_ASSERT(isProbabilisticState(state), "State is not probabilistic."); return false; } } @@ -300,7 +300,7 @@ namespace storm { } template <typename ValueType, typename RewardModelType> - std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> MarkovAutomaton<ValueType, RewardModelType>::convertToCTMC() { + std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> MarkovAutomaton<ValueType, RewardModelType>::convertToCTMC() const { STORM_LOG_TRACE("MA matrix:" << std::endl << this->getTransitionMatrix()); STORM_LOG_TRACE("Markovian states: " << getMarkovianStates()); @@ -311,7 +311,7 @@ namespace storm { storm::solver::stateelimination::StateEliminator<ValueType> stateEliminator(flexibleMatrix, flexibleBackwardTransitions); for (uint_fast64_t state = 0; state < this->getNumberOfStates(); ++state) { - assert(!this->isHybridState(state)); + STORM_LOG_ASSERT(!this->isHybridState(state), "State is hybrid."); if (this->isProbabilisticState(state)) { // Eliminate this probabilistic state stateEliminator.eliminateState(state, true); @@ -328,7 +328,7 @@ namespace storm { // State is eliminated and can be discarded keepStates.set(state, false); } else { - assert(this->isMarkovianState(state)); + STORM_LOG_ASSERT(this->isMarkovianState(state), "State is not markovian."); // Copy transitions for (uint_fast64_t row = flexibleMatrix.getRowGroupIndices()[state]; row < flexibleMatrix.getRowGroupIndices()[state + 1]; ++row) { for (auto const& entry : flexibleMatrix.getRow(row)) { diff --git a/src/models/sparse/MarkovAutomaton.h b/src/models/sparse/MarkovAutomaton.h index 71d35b894..a11d3f320 100644 --- a/src/models/sparse/MarkovAutomaton.h +++ b/src/models/sparse/MarkovAutomaton.h @@ -186,7 +186,12 @@ namespace storm { bool hasOnlyTrivialNondeterminism() const; - std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> convertToCTMC(); + /*! + * Convert the MA into a MA by eliminating all states with probabilistic choices. + * + * @return Ctmc. + */ + std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> convertToCTMC() const; virtual void writeDotToStream(std::ostream& outStream, bool includeLabeling = true, storm::storage::BitVector const* subsystem = nullptr, std::vector<ValueType> const* firstValue = nullptr, std::vector<ValueType> const* secondValue = nullptr, std::vector<uint_fast64_t> const* stateColoring = nullptr, std::vector<std::string> const* colors = nullptr, std::vector<uint_fast64_t>* scheduler = nullptr, bool finalizeOutput = true) const override; diff --git a/src/models/sparse/Model.cpp b/src/models/sparse/Model.cpp index 14f9fa8a3..34ec5d2cc 100644 --- a/src/models/sparse/Model.cpp +++ b/src/models/sparse/Model.cpp @@ -91,7 +91,7 @@ namespace storm { template<typename ValueType, typename RewardModelType> RewardModelType& Model<ValueType, RewardModelType>::rewardModel(std::string const& rewardModelName) { - assert(this->hasRewardModel(rewardModelName)); + STORM_LOG_ASSERT(this->hasRewardModel(rewardModelName), "Model has no reward model."); return this->rewardModels.find(rewardModelName)->second; } @@ -117,7 +117,7 @@ namespace storm { if(this->hasRewardModel(rewardModelName)) { STORM_LOG_THROW(!(this->hasRewardModel(rewardModelName)), storm::exceptions::IllegalArgumentException, "A reward model with the given name '" << rewardModelName << "' already exists."); } - assert(newRewardModel.isCompatible(this->getNumberOfStates(), this->getTransitionMatrix().getRowCount())); + STORM_LOG_ASSERT(newRewardModel.isCompatible(this->getNumberOfStates(), this->getTransitionMatrix().getRowCount()), "New reward model is not compatible."); this->rewardModels.emplace(rewardModelName, newRewardModel); } diff --git a/src/models/sparse/StandardRewardModel.cpp b/src/models/sparse/StandardRewardModel.cpp index a539cfd15..5f72f97ea 100644 --- a/src/models/sparse/StandardRewardModel.cpp +++ b/src/models/sparse/StandardRewardModel.cpp @@ -37,13 +37,13 @@ namespace storm { template<typename ValueType> std::vector<ValueType> const& StandardRewardModel<ValueType>::getStateRewardVector() const { - assert(this->hasStateRewards()); + STORM_LOG_ASSERT(this->hasStateRewards(), "No state rewards available."); return this->optionalStateRewardVector.get(); } template<typename ValueType> std::vector<ValueType>& StandardRewardModel<ValueType>::getStateRewardVector() { - assert(this->hasStateRewards()); + STORM_LOG_ASSERT(this->hasStateRewards(), "No state rewards available."); return this->optionalStateRewardVector.get(); } @@ -54,16 +54,16 @@ namespace storm { template<typename ValueType> ValueType const& StandardRewardModel<ValueType>::getStateReward(uint_fast64_t state) const { - assert(this->hasStateRewards()); - assert(state < this->optionalStateRewardVector.get().size()); + STORM_LOG_ASSERT(this->hasStateRewards(), "No state rewards available."); + STORM_LOG_ASSERT(state < this->optionalStateRewardVector.get().size(), "Invalid state."); return this->optionalStateRewardVector.get()[state]; } template<typename ValueType> template<typename T> void StandardRewardModel<ValueType>::setStateReward(uint_fast64_t state, T const & newReward) { - assert(this->hasStateRewards()); - assert(state < this->optionalStateRewardVector.get().size()); + STORM_LOG_ASSERT(this->hasStateRewards(), "No state rewards available."); + STORM_LOG_ASSERT(state < this->optionalStateRewardVector.get().size(), "Invalid state."); this->optionalStateRewardVector.get()[state] = newReward; } @@ -74,28 +74,28 @@ namespace storm { template<typename ValueType> std::vector<ValueType> const& StandardRewardModel<ValueType>::getStateActionRewardVector() const { - assert(this->hasStateActionRewards()); + STORM_LOG_ASSERT(this->hasStateActionRewards(), "No state action rewards available."); return this->optionalStateActionRewardVector.get(); } template<typename ValueType> std::vector<ValueType>& StandardRewardModel<ValueType>::getStateActionRewardVector() { - assert(this->hasStateActionRewards()); + STORM_LOG_ASSERT(this->hasStateActionRewards(), "No state action rewards available."); return this->optionalStateActionRewardVector.get(); } template<typename ValueType> ValueType const& StandardRewardModel<ValueType>::getStateActionReward(uint_fast64_t choiceIndex) const { - assert(this->hasStateActionRewards()); - assert(choiceIndex < this->optionalStateActionRewardVector.get().size()); + STORM_LOG_ASSERT(this->hasStateActionRewards(), "No state action rewards available."); + STORM_LOG_ASSERT(choiceIndex < this->optionalStateActionRewardVector.get().size(), "Invalid choiceIndex."); return this->optionalStateActionRewardVector.get()[choiceIndex]; } template<typename ValueType> template<typename T> void StandardRewardModel<ValueType>::setStateActionReward(uint_fast64_t choiceIndex, T const &newValue) { - assert(this->hasStateActionRewards()); - assert(choiceIndex < this->optionalStateActionRewardVector.get().size()); + STORM_LOG_ASSERT(this->hasStateActionRewards(), "No state action rewards available."); + STORM_LOG_ASSERT(choiceIndex < this->optionalStateActionRewardVector.get().size(), "Invalid choiceIndex."); this->optionalStateActionRewardVector.get()[choiceIndex] = newValue; } diff --git a/src/models/sparse/StateLabeling.cpp b/src/models/sparse/StateLabeling.cpp index 7957e6a99..70561c8cb 100644 --- a/src/models/sparse/StateLabeling.cpp +++ b/src/models/sparse/StateLabeling.cpp @@ -29,7 +29,7 @@ namespace storm { return true; } - StateLabeling StateLabeling::getSubLabeling(storm::storage::BitVector const& states) { + StateLabeling StateLabeling::getSubLabeling(storm::storage::BitVector const& states) const { StateLabeling result(states.getNumberOfSetBits()); for (auto const& labelIndexPair : nameToLabelingIndexMap) { result.addLabel(labelIndexPair.first, labelings[labelIndexPair.second] % states); diff --git a/src/models/sparse/StateLabeling.h b/src/models/sparse/StateLabeling.h index c7b3411a7..ac8ac730a 100644 --- a/src/models/sparse/StateLabeling.h +++ b/src/models/sparse/StateLabeling.h @@ -49,7 +49,7 @@ namespace storm { * * @param states The selected set of states. */ - StateLabeling getSubLabeling(storm::storage::BitVector const& states); + StateLabeling getSubLabeling(storm::storage::BitVector const& states) const; /*! * Adds a new label to the labelings. Initially, no state is labeled with this label. diff --git a/src/permissivesched/MILPPermissiveSchedulers.h b/src/permissivesched/MILPPermissiveSchedulers.h index b0f6b4529..fe5a48e11 100644 --- a/src/permissivesched/MILPPermissiveSchedulers.h +++ b/src/permissivesched/MILPPermissiveSchedulers.h @@ -51,13 +51,13 @@ class MilpPermissiveSchedulerComputation : public PermissiveSchedulerComputation bool foundSolution() const override { - assert(mCalledOptimizer); + STORM_LOG_ASSERT(mCalledOptimizer, "Optimizer not called."); return !solver.isInfeasible(); } SubMDPPermissiveScheduler<RM> getScheduler() const override { - assert(mCalledOptimizer); - assert(foundSolution()); + STORM_LOG_ASSERT(mCalledOptimizer, "Optimizer not called."); + STORM_LOG_ASSERT(foundSolution(), "Solution not found."); SubMDPPermissiveScheduler<RM> result(this->mdp, true); @@ -103,7 +103,7 @@ class MilpPermissiveSchedulerComputation : public PermissiveSchedulerComputation */ void createVariables(PermissiveSchedulerPenalties const& penalties, storm::storage::BitVector const& relevantStates) { // We need the unique initial state later, so we get that one before looping. - assert(this->mdp.getInitialStates().getNumberOfSetBits() == 1); + STORM_LOG_ASSERT(this->mdp.getInitialStates().getNumberOfSetBits() == 1, "No unique initial state."); uint_fast64_t initialStateIndex = this->mdp.getInitialStates().getNextSetIndex(0); storm::expressions::Variable var; @@ -151,9 +151,9 @@ class MilpPermissiveSchedulerComputation : public PermissiveSchedulerComputation // (5) and (7) are omitted on purpose (-- we currenty do not support controllability of actions -- ) // (1) - assert(this->mdp.getInitialStates().getNumberOfSetBits() == 1); + STORM_LOG_ASSERT(this->mdp.getInitialStates().getNumberOfSetBits() == 1, "No unique initial state."); uint_fast64_t initialStateIndex = this->mdp.getInitialStates().getNextSetIndex(0); - assert(relevantStates[initialStateIndex]); + STORM_LOG_ASSERT(relevantStates[initialStateIndex], "Initial state not relevant."); if(lowerBound) { solver.addConstraint("c1", mProbVariables[initialStateIndex] >= solver.getConstant(boundary)); } else { @@ -206,9 +206,9 @@ class MilpPermissiveSchedulerComputation : public PermissiveSchedulerComputation std::string satstring = to_string(sat); // (8) if(relevantStates[entry.getColumn()]) { - assert(mGammaVariables.count(entry.getColumn()) > 0); - assert(mGammaVariables.count(s) > 0); - assert(mBetaVariables.count(sat) > 0); + STORM_LOG_ASSERT(mGammaVariables.count(entry.getColumn()) > 0, "Entry not found."); + STORM_LOG_ASSERT(mGammaVariables.count(s) > 0, "Entry not found."); + STORM_LOG_ASSERT(mBetaVariables.count(sat) > 0, "Entry not found."); solver.addConstraint("c8-" + satstring, mGammaVariables[entry.getColumn()] < mGammaVariables[s] + (solver.getConstant(1) - mBetaVariables[sat])); } } diff --git a/src/permissivesched/PermissiveSchedulerPenalty.h b/src/permissivesched/PermissiveSchedulerPenalty.h index bb680defa..f25f61ac1 100644 --- a/src/permissivesched/PermissiveSchedulerPenalty.h +++ b/src/permissivesched/PermissiveSchedulerPenalty.h @@ -30,7 +30,7 @@ namespace storm { } void set(uint_fast64_t state, uint_fast64_t action, double penalty) { - assert(penalty >= 1.0); + STORM_LOG_ASSERT(penalty >= 1.0, "Penalty too low."); if(penalty == 1.0) { auto it = mPenalties.find(std::make_pair(state, action)); if(it != mPenalties.end()) { diff --git a/src/permissivesched/PermissiveSchedulers.cpp b/src/permissivesched/PermissiveSchedulers.cpp index 816c92d55..004230296 100644 --- a/src/permissivesched/PermissiveSchedulers.cpp +++ b/src/permissivesched/PermissiveSchedulers.cpp @@ -17,7 +17,7 @@ namespace storm { template<typename RM> boost::optional<SubMDPPermissiveScheduler<RM>> computePermissiveSchedulerViaMILP(storm::models::sparse::Mdp<double, RM> const& mdp, storm::logic::ProbabilityOperatorFormula const& safeProp) { storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Mdp<double, RM>> propMC(mdp); - assert(safeProp.getSubformula().isEventuallyFormula()); + STORM_LOG_ASSERT(safeProp.getSubformula().isEventuallyFormula(), "No eventually formula."); auto backwardTransitions = mdp.getBackwardTransitions(); storm::storage::BitVector goalstates = propMC.check(safeProp.getSubformula().asEventuallyFormula().getSubformula())->asExplicitQualitativeCheckResult().getTruthValuesVector(); goalstates = storm::utility::graph::performProb1A(mdp, backwardTransitions, storm::storage::BitVector(goalstates.size(), true), goalstates); @@ -44,7 +44,7 @@ namespace storm { template<typename RM> boost::optional<SubMDPPermissiveScheduler<RM>> computePermissiveSchedulerViaSMT(storm::models::sparse::Mdp<double, RM> const& mdp, storm::logic::ProbabilityOperatorFormula const& safeProp) { storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Mdp<double, RM>> propMC(mdp); - assert(safeProp.getSubformula().isEventuallyFormula()); + STORM_LOG_ASSERT(safeProp.getSubformula().isEventuallyFormula(), "No eventually formula."); auto backwardTransitions = mdp.getBackwardTransitions(); storm::storage::BitVector goalstates = propMC.check(safeProp.getSubformula().asEventuallyFormula().getSubformula())->asExplicitQualitativeCheckResult().getTruthValuesVector(); goalstates = storm::utility::graph::performProb1A(mdp, backwardTransitions, storm::storage::BitVector(goalstates.size(), true), goalstates); diff --git a/src/permissivesched/PermissiveSchedulers.h b/src/permissivesched/PermissiveSchedulers.h index 7d0bdb3f7..eb09d4b24 100644 --- a/src/permissivesched/PermissiveSchedulers.h +++ b/src/permissivesched/PermissiveSchedulers.h @@ -32,7 +32,7 @@ namespace storm { } void disable(uint_fast64_t choiceIndex) { - assert(choiceIndex < enabledChoices.size()); + STORM_LOG_ASSERT(choiceIndex < enabledChoices.size(), "Invalid choiceIndex."); enabledChoices.set(choiceIndex, false); } diff --git a/src/permissivesched/SmtBasedPermissiveSchedulers.h b/src/permissivesched/SmtBasedPermissiveSchedulers.h index bfc1530b8..ed183866c 100644 --- a/src/permissivesched/SmtBasedPermissiveSchedulers.h +++ b/src/permissivesched/SmtBasedPermissiveSchedulers.h @@ -34,12 +34,12 @@ namespace storm { } bool foundSolution() const override { - assert(mPerformedSmtLoop); + STORM_LOG_ASSERT(mPerformedSmtLoop, "SMT loop not performed."); return mFoundSolution; } SubMDPPermissiveScheduler<RM> getScheduler() const override { - assert(foundSolution()); + STORM_LOG_ASSERT(foundSolution(), "Solution not found."); SubMDPPermissiveScheduler<RM> result(this->mdp, true); for(auto const& entry : multistrategyVariables) { if(!multistrategyVariablesToTakenMap.at(entry.second)) { @@ -98,9 +98,9 @@ namespace storm { // (4) and (7) are omitted on purpose (-- we currenty do not support controllability of actions -- ) // (1) - assert(this->mdp.getInitialStates().getNumberOfSetBits() == 1); + STORM_LOG_ASSERT(this->mdp.getInitialStates().getNumberOfSetBits() == 1, "No unique initial state."); uint_fast64_t initialStateIndex = this->mdp.getInitialStates().getNextSetIndex(0); - assert(relevantStates[initialStateIndex]); + STORM_LOG_ASSERT(relevantStates[initialStateIndex], "Initial state not relevant."); if(lowerBound) { solver.add(mProbVariables[initialStateIndex] >= manager.rational(boundary)); } else { @@ -162,9 +162,9 @@ namespace storm { // std::string satstring = to_string(sat); // // (8) // if(relevantStates[entry.getColumn()]) { -// assert(mGammaVariables.count(entry.getColumn()) > 0); -// assert(mGammaVariables.count(s) > 0); -// assert(mBetaVariables.count(sat) > 0); +// STORM_LOG_ASSERT(mGammaVariables.count(entry.getColumn()) > 0, "Entry not found."); +// STORM_LOG_ASSERT(mGammaVariables.count(s) > 0, "Entry not found."); +// STORM_LOG_ASSERT(mBetaVariables.count(sat) > 0, "Entry not found."); // solver.addConstraint("c8-" + satstring, mGammaVariables[entry.getColumn()] < mGammaVariables[s] + (solver.getConstant(1) - mBetaVariables[sat]) + mProbVariables[s]); // With rewards, we have to change this. // } // } diff --git a/src/settings/modules/DFTSettings.cpp b/src/settings/modules/DFTSettings.cpp index 365ba0271..ab707b5a4 100644 --- a/src/settings/modules/DFTSettings.cpp +++ b/src/settings/modules/DFTSettings.cpp @@ -9,7 +9,6 @@ #include "src/exceptions/InvalidSettingsException.h" - namespace storm { namespace settings { namespace modules { @@ -21,24 +20,35 @@ namespace storm { const std::string DFTSettings::symmetryReductionOptionShortName = "symred"; const std::string DFTSettings::modularisationOptionName = "modularisation"; const std::string DFTSettings::disableDCOptionName = "disabledc"; + const std::string DFTSettings::approximationErrorOptionName = "approximation"; + const std::string DFTSettings::approximationErrorOptionShortName = "approx"; + const std::string DFTSettings::approximationHeuristicOptionName = "approximationheuristic"; const std::string DFTSettings::propExpectedTimeOptionName = "expectedtime"; const std::string DFTSettings::propExpectedTimeOptionShortName = "mttf"; const std::string DFTSettings::propProbabilityOptionName = "probability"; const std::string DFTSettings::propTimeBoundOptionName = "timebound"; const std::string DFTSettings::minValueOptionName = "min"; const std::string DFTSettings::maxValueOptionName = "max"; - +#ifdef STORM_HAVE_Z3 + const std::string DFTSettings::solveWithSmtOptionName = "smt"; +#endif + DFTSettings::DFTSettings() : ModuleSettings(moduleName) { this->addOption(storm::settings::OptionBuilder(moduleName, dftFileOptionName, false, "Parses the model given in the Galileo format.").setShortName(dftFileOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The name of the file from which to read the DFT model.").addValidationFunctionString(storm::settings::ArgumentValidators::existingReadableFileValidator()).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, symmetryReductionOptionName, false, "Exploit symmetric structure of model.").setShortName(symmetryReductionOptionShortName).build()); this->addOption(storm::settings::OptionBuilder(moduleName, modularisationOptionName, false, "Use modularisation (not applicable for expected time).").build()); this->addOption(storm::settings::OptionBuilder(moduleName, disableDCOptionName, false, "Disable Dont Care propagation.").build()); + this->addOption(storm::settings::OptionBuilder(moduleName, approximationErrorOptionName, false, "Approximation error allowed.").setShortName(approximationErrorOptionShortName).addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("error", "The relative approximation error to use.").addValidationFunctionDouble(storm::settings::ArgumentValidators::doubleGreaterValidatorIncluding(0.0)).build()).build()); + this->addOption(storm::settings::OptionBuilder(moduleName, approximationHeuristicOptionName, false, "Set the heuristic used for approximation.").addArgument(storm::settings::ArgumentBuilder::createStringArgument("heuristic", "Sets which heuristic is used for approximation. Must be in {depth, probability}. Default is").setDefaultValueString("depth").addValidationFunctionString(storm::settings::ArgumentValidators::stringInListValidator({"depth", "rateratio"})).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, propExpectedTimeOptionName, false, "Compute expected time of system failure.").setShortName(propExpectedTimeOptionShortName).build()); this->addOption(storm::settings::OptionBuilder(moduleName, propProbabilityOptionName, false, "Compute probability of system failure.").build()); this->addOption(storm::settings::OptionBuilder(moduleName, propTimeBoundOptionName, false, "Compute probability of system failure up to given timebound.").addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("time", "The timebound to use.").addValidationFunctionDouble(storm::settings::ArgumentValidators::doubleGreaterValidatorExcluding(0.0)).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, minValueOptionName, false, "Compute minimal value in case of non-determinism.").build()); this->addOption(storm::settings::OptionBuilder(moduleName, maxValueOptionName, false, "Compute maximal value in case of non-determinism.").build()); +#ifdef STORM_HAVE_Z3 + this->addOption(storm::settings::OptionBuilder(moduleName, solveWithSmtOptionName, true, "Solve the DFT with SMT.").build()); +#endif } bool DFTSettings::isDftFileSet() const { @@ -60,7 +70,29 @@ namespace storm { bool DFTSettings::isDisableDC() const { return this->getOption(disableDCOptionName).getHasOptionBeenSet(); } - + + bool DFTSettings::isApproximationErrorSet() const { + return this->getOption(approximationErrorOptionName).getHasOptionBeenSet(); + } + + double DFTSettings::getApproximationError() const { + return this->getOption(approximationErrorOptionName).getArgumentByName("error").getValueAsDouble(); + } + + storm::builder::ApproximationHeuristic DFTSettings::getApproximationHeuristic() const { + if (!isApproximationErrorSet() || getApproximationError() == 0.0) { + // No approximation is done + return storm::builder::ApproximationHeuristic::NONE; + } + std::string heuristicAsString = this->getOption(approximationHeuristicOptionName).getArgumentByName("heuristic").getValueAsString(); + if (heuristicAsString == "depth") { + return storm::builder::ApproximationHeuristic::DEPTH; + } else if (heuristicAsString == "probability") { + return storm::builder::ApproximationHeuristic::PROBABILITY; + } + STORM_LOG_THROW(false, storm::exceptions::IllegalArgumentValueException, "Illegal value '" << heuristicAsString << "' set as heuristic for approximation."); + } + bool DFTSettings::usePropExpectedTime() const { return this->getOption(propExpectedTimeOptionName).getHasOptionBeenSet(); } @@ -85,8 +117,13 @@ namespace storm { return this->getOption(maxValueOptionName).getHasOptionBeenSet(); } +#ifdef STORM_HAVE_Z3 + bool DFTSettings::solveWithSMT() const { + return this->getOption(solveWithSmtOptionName).getHasOptionBeenSet(); + } +#endif + void DFTSettings::finalize() { - } bool DFTSettings::check() const { @@ -104,4 +141,4 @@ namespace storm { } // namespace modules } // namespace settings -} // namespace storm \ No newline at end of file +} // namespace storm diff --git a/src/settings/modules/DFTSettings.h b/src/settings/modules/DFTSettings.h index 844953c2c..e580a5d05 100644 --- a/src/settings/modules/DFTSettings.h +++ b/src/settings/modules/DFTSettings.h @@ -3,6 +3,7 @@ #include "storm-config.h" #include "src/settings/modules/ModuleSettings.h" +#include "src/builder/DftExplorationHeuristic.h" namespace storm { namespace settings { @@ -32,10 +33,7 @@ namespace storm { * @return The name of the file that contains the dft specification. */ std::string getDftFilename() const; - - //expectedtime, probability, timebound, prob - //min, max - + /*! * Retrieves whether the option to use symmetry reduction is set. * @@ -57,6 +55,27 @@ namespace storm { */ bool isDisableDC() const; + /*! + * Retrieves whether the option to compute an approximation is set. + * + * @return True iff the option was set. + */ + bool isApproximationErrorSet() const; + + /*! + * Retrieves the relative error allowed for approximating the model checking result. + * + * @return The allowed errorbound. + */ + double getApproximationError() const; + + /*! + * Retrieves the heuristic used for approximation. + * + * @return The heuristic to use. + */ + storm::builder::ApproximationHeuristic getApproximationHeuristic() const; + /*! * Retrieves whether the property expected time should be used. * @@ -99,6 +118,15 @@ namespace storm { */ double getPropTimebound() const; +#ifdef STORM_HAVE_Z3 + /*! + * Retrieves whether the DFT should be checked via SMT. + * + * @return True iff the option was set. + */ + bool solveWithSMT() const; +#endif + bool check() const override; void finalize() override; @@ -113,12 +141,18 @@ namespace storm { static const std::string symmetryReductionOptionShortName; static const std::string modularisationOptionName; static const std::string disableDCOptionName; + static const std::string approximationErrorOptionName; + static const std::string approximationErrorOptionShortName; + static const std::string approximationHeuristicOptionName; static const std::string propExpectedTimeOptionName; static const std::string propExpectedTimeOptionShortName; static const std::string propProbabilityOptionName; static const std::string propTimeBoundOptionName; static const std::string minValueOptionName; static const std::string maxValueOptionName; +#ifdef STORM_HAVE_Z3 + static const std::string solveWithSmtOptionName; +#endif }; diff --git a/src/solver/AbstractEquationSolver.h b/src/solver/AbstractEquationSolver.h index 6cc03d56d..491570637 100644 --- a/src/solver/AbstractEquationSolver.h +++ b/src/solver/AbstractEquationSolver.h @@ -1,7 +1,6 @@ #ifndef STORM_SOLVER_ABSTRACTEQUATIONSOLVER_H_ #define STORM_SOLVER_ABSTRACTEQUATIONSOLVER_H_ -#include <memory> #include "src/solver/TerminationCondition.h" #include <memory> diff --git a/src/solver/MinMaxLinearEquationSolver.h b/src/solver/MinMaxLinearEquationSolver.h index 14457e778..6be77c93b 100644 --- a/src/solver/MinMaxLinearEquationSolver.h +++ b/src/solver/MinMaxLinearEquationSolver.h @@ -14,6 +14,7 @@ #include "src/solver/OptimizationDirection.h" #include "src/exceptions/InvalidSettingsException.h" +#include "src/utility/macros.h" namespace storm { namespace storage { diff --git a/src/solver/OptimizationDirection.cpp b/src/solver/OptimizationDirection.cpp index cd88657b9..71d0760ed 100644 --- a/src/solver/OptimizationDirection.cpp +++ b/src/solver/OptimizationDirection.cpp @@ -1,6 +1,6 @@ #include "OptimizationDirection.h" #include <iostream> -#include <cassert> +#include "src/utility/macros.h" namespace storm { namespace solver { @@ -18,7 +18,7 @@ namespace storm { } OptimizationDirection convert(OptimizationDirectionSetting s) { - assert(isSet(s)); + STORM_LOG_ASSERT(isSet(s), "Setting is not set."); return static_cast<OptimizationDirection>(s); } diff --git a/src/storage/BucketPriorityQueue.cpp b/src/storage/BucketPriorityQueue.cpp new file mode 100644 index 000000000..7870888f0 --- /dev/null +++ b/src/storage/BucketPriorityQueue.cpp @@ -0,0 +1,212 @@ +#include "src/storage/BucketPriorityQueue.h" +#include "src/utility/macros.h" +#include "src/adapters/CarlAdapter.h" + +#include <cmath> + +namespace storm { + namespace storage { + + template<typename ValueType> + BucketPriorityQueue<ValueType>::BucketPriorityQueue(size_t nrBuckets, double lowerValue, double ratio) : lowerValue(lowerValue), logBase(std::log(ratio)), nrBuckets(nrBuckets), nrUnsortedItems(0), buckets(nrBuckets), currentBucket(nrBuckets) { + compare = ([this](HeuristicPointer a, HeuristicPointer b) { + return *a < *b; + }); + } + + template<typename ValueType> + void BucketPriorityQueue<ValueType>::fix() { + if (currentBucket < nrBuckets && nrUnsortedItems > buckets[currentBucket].size() / 10) { + // Fix current bucket + std::make_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare); + nrUnsortedItems = 0; + } + } + + template<typename ValueType> + bool BucketPriorityQueue<ValueType>::empty() const { + return currentBucket == nrBuckets && immediateBucket.empty(); + } + + template<typename ValueType> + std::size_t BucketPriorityQueue<ValueType>::size() const { + size_t size = immediateBucket.size(); + for (size_t i = currentBucket; currentBucket < nrBuckets; ++i) { + size += buckets[i].size(); + } + return size; + } + + 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"); + return buckets[currentBucket].front(); + } + + 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; + nrUnsortedItems = 0; + } + buckets[bucket].push_back(item); + if (bucket == currentBucket) { + // Insert in first bucket + if (AUTOSORT) { + std::push_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare); + } else { + ++nrUnsortedItems; + } + } + } + + 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); + + if (oldBucket == newBucket) { + if (currentBucket < newBucket) { + // No change as the bucket is not sorted yet + } else { + if (AUTOSORT) { + // Sort first bucket + fix(); + } else { + ++nrUnsortedItems; + } + } + } else { + // Move to new bucket + STORM_LOG_ASSERT(newBucket < oldBucket, "Will update to higher bucket"); + if (newBucket < currentBucket) { + currentBucket = newBucket; + nrUnsortedItems = 0; + } + // Remove old entry by swap-and-pop + if (buckets[oldBucket].size() >= 2) { + // Find old index by linear search + // Notice: using a map to rememeber index was not efficient + size_t oldIndex = 0; + for ( ; oldIndex < buckets[oldBucket].size(); ++oldIndex) { + if (buckets[oldBucket][oldIndex]->getId() == item->getId()) { + break; + } + } + STORM_LOG_ASSERT(oldIndex < buckets[oldBucket].size(), "Id " << item->getId() << " not found"); + std::iter_swap(buckets[oldBucket].begin() + oldIndex, buckets[oldBucket].end() - 1); + } + buckets[oldBucket].pop_back(); + // Insert new element + buckets[newBucket].push_back(item); + if (newBucket == currentBucket) { + if (AUTOSORT) { + // Sort in first bucket + std::push_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare); + } else { + ++nrUnsortedItems; + } + } + } + } + + + template<typename ValueType> + void BucketPriorityQueue<ValueType>::pop() { + if (!immediateBucket.empty()) { + immediateBucket.pop_back(); + return; + } + STORM_LOG_ASSERT(!empty(), "BucketPriorityQueue is empty"); + std::pop_heap(buckets[currentBucket].begin(), buckets[currentBucket].end(), compare); + buckets[currentBucket].pop_back(); + if (buckets[currentBucket].empty()) { + // Find next bucket with elements + for ( ; currentBucket < nrBuckets; ++currentBucket) { + if (!buckets[currentBucket].empty()) { + nrUnsortedItems = buckets[currentBucket].size(); + if (AUTOSORT) { + fix(); + } + break; + } + } + } + } + + template<typename ValueType> + typename BucketPriorityQueue<ValueType>::HeuristicPointer BucketPriorityQueue<ValueType>::popTop() { + HeuristicPointer item = top(); + pop(); + return item; + } + + template<typename ValueType> + size_t BucketPriorityQueue<ValueType>::getBucket(double priority) const { + STORM_LOG_ASSERT(priority >= lowerValue, "Priority " << priority << " is too low"); + + // For possible values greater 1 + double tmpBucket = std::log(priority - lowerValue) / logBase; + tmpBucket += buckets.size() / 2; + if (tmpBucket < 0) { + tmpBucket = 0; + } + size_t newBucket = tmpBucket; + // For values ensured to be lower 1 + //size_t newBucket = std::log(priority - lowerValue) / logBase; + if (newBucket >= nrBuckets) { + newBucket = nrBuckets - 1; + } + if (!HIGHER) { + newBucket = nrBuckets-1 - newBucket; + } + //std::cout << "get Bucket: " << priority << ", " << newBucket << std::endl; + STORM_LOG_ASSERT(newBucket < nrBuckets, "Priority " << priority << " is too high"); + return newBucket; + } + + template<typename ValueType> + void BucketPriorityQueue<ValueType>::print(std::ostream& out) const { + out << "Bucket priority queue with size " << buckets.size() << ", lower value: " << lowerValue << " and logBase: " << logBase << 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()) { + out << "Bucket " << bucket << ":" << std::endl; + for (HeuristicPointer heuristic : buckets[bucket]) { + out << "\t" << heuristic->getId() << ": " << heuristic->getPriority() << std::endl; + } + } + } + } + + template<typename ValueType> + void BucketPriorityQueue<ValueType>::printSizes(std::ostream& out) const { + out << "Bucket sizes: " << immediateBucket.size() << " | "; + for (size_t bucket = 0; bucket < buckets.size(); ++bucket) { + out << buckets[bucket].size() << " "; + } + std::cout << std::endl; + } + + // Template instantiations + template class BucketPriorityQueue<double>; + +#ifdef STORM_HAVE_CARL + template class BucketPriorityQueue<storm::RationalFunction>; +#endif + } +} diff --git a/src/storage/BucketPriorityQueue.h b/src/storage/BucketPriorityQueue.h new file mode 100644 index 000000000..75f065bdc --- /dev/null +++ b/src/storage/BucketPriorityQueue.h @@ -0,0 +1,73 @@ +#ifndef STORM_STORAGE_BUCKETPRIORITYQUEUE_H_ +#define STORM_STORAGE_BUCKETPRIORITYQUEUE_H_ + +#include "src/builder/DftExplorationHeuristic.h" +#include <algorithm> +#include <functional> +#include <unordered_map> +#include <vector> + +namespace storm { + namespace storage { + + template<typename ValueType> + class BucketPriorityQueue { + + using HeuristicPointer = std::shared_ptr<storm::builder::DFTExplorationHeuristicBoundDifference<ValueType>>; + + public: + explicit BucketPriorityQueue(size_t nrBuckets, double lowerValue, double ratio); + + void fix(); + + bool empty() const; + + std::size_t size() const; + + HeuristicPointer const& top() const; + + void push(HeuristicPointer const& item); + + void update(HeuristicPointer const& item, double oldPriority); + + void pop(); + + HeuristicPointer popTop(); + + void print(std::ostream& out) const; + + void printSizes(std::ostream& out) const; + + private: + + size_t getBucket(double priority) const; + + const double lowerValue; + + const bool HIGHER = true; + + const bool AUTOSORT = false; + + const double logBase; + + const size_t nrBuckets; + + size_t nrUnsortedItems; + + // 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; + + std::function<bool(HeuristicPointer, HeuristicPointer)> compare; + + }; + + } +} + +#endif // STORM_STORAGE_BUCKETPRIORITYQUEUE_H_ diff --git a/src/storage/DynamicPriorityQueue.h b/src/storage/DynamicPriorityQueue.h new file mode 100644 index 000000000..9df114fbf --- /dev/null +++ b/src/storage/DynamicPriorityQueue.h @@ -0,0 +1,71 @@ +#ifndef STORM_STORAGE_DYNAMICPRIORITYQUEUE_H_ +#define STORM_STORAGE_DYNAMICPRIORITYQUEUE_H_ + +#include <algorithm> +#include <vector> + +namespace storm { + namespace storage { + + template<typename T, typename Container = std::vector<T>, typename Compare = std::less<T>> + class DynamicPriorityQueue { + + private: + Container container; + Compare compare; + + public: + explicit DynamicPriorityQueue(Compare const& compare) : container(), compare(compare) { + // Intentionally left empty + } + + explicit DynamicPriorityQueue(Container&& container, Compare const& compare) : container(std::move(container)), compare(compare) { + std::make_heap(container.begin(), container.end(), compare); + } + + void fix() { + std::make_heap(container.begin(), container.end(), compare); + } + + bool empty() const { + return container.empty(); + } + + std::size_t size() const { + return container.size(); + } + + const T& top() const { + return container.front(); + } + + void push(const T& item) { + container.push_back(item); + std::push_heap(container.begin(), container.end(), compare); + } + + void push(T&& item) { + container.push_back(std::move(item)); + std::push_heap(container.begin(), container.end(), compare); + } + + void pop() { + std::pop_heap(container.begin(), container.end(), compare); + container.pop_back(); + } + + T popTop() { + T item = top(); + pop(); + return item; + } + + Container getContainer() const { + return container; + } + + }; + } +} + +#endif // STORM_STORAGE_DYNAMICPRIORITYQUEUE_H_ diff --git a/src/storage/FlexibleSparseMatrix.cpp b/src/storage/FlexibleSparseMatrix.cpp index 102d20304..563f10378 100644 --- a/src/storage/FlexibleSparseMatrix.cpp +++ b/src/storage/FlexibleSparseMatrix.cpp @@ -75,15 +75,15 @@ namespace storm { template<typename ValueType> typename FlexibleSparseMatrix<ValueType>::row_type& FlexibleSparseMatrix<ValueType>::getRow(index_type rowGroup, index_type offset) { - assert(rowGroup < this->getRowGroupCount()); - assert(offset < this->getRowGroupSize(rowGroup)); + STORM_LOG_ASSERT(rowGroup < this->getRowGroupCount(), "Invalid rowGroup."); + STORM_LOG_ASSERT(offset < this->getRowGroupSize(rowGroup), "Invalid offset."); return getRow(rowGroupIndices[rowGroup] + offset); } template<typename ValueType> typename FlexibleSparseMatrix<ValueType>::row_type const& FlexibleSparseMatrix<ValueType>::getRow(index_type rowGroup, index_type offset) const { - assert(rowGroup < this->getRowGroupCount()); - assert(offset < this->getRowGroupSize(rowGroup)); + STORM_LOG_ASSERT(rowGroup < this->getRowGroupCount(), "Invalid rowGroup."); + STORM_LOG_ASSERT(offset < this->getRowGroupSize(rowGroup), "Invalid offset."); return getRow(rowGroupIndices[rowGroup] + offset); } @@ -136,7 +136,7 @@ namespace storm { this->columnCount = 0; for (auto const& row : this->data) { for (auto const& element : row) { - assert(!storm::utility::isZero(element.getValue())); + STORM_LOG_ASSERT(!storm::utility::isZero(element.getValue()), "Entry is 0."); ++this->nonzeroEntryCount; this->columnCount = std::max(element.getColumn() + 1, this->columnCount); } diff --git a/src/storage/SparseMatrix.cpp b/src/storage/SparseMatrix.cpp index 401fc4082..0bb4ebe10 100644 --- a/src/storage/SparseMatrix.cpp +++ b/src/storage/SparseMatrix.cpp @@ -22,6 +22,8 @@ #include "src/utility/macros.h" +#include <iterator> + namespace storm { namespace storage { @@ -267,66 +269,37 @@ namespace storm { } template<typename ValueType> - bool SparseMatrixBuilder<ValueType>::replaceColumns(std::vector<index_type> const& replacements, index_type offset) { - bool changed = false; + void SparseMatrixBuilder<ValueType>::replaceColumns(std::vector<index_type> const& replacements, index_type offset) { index_type maxColumn = 0; - for (auto& elem : columnsAndValues) { - if (elem.getColumn() >= offset) { - elem.setColumn(replacements[elem.getColumn() - offset]); - changed = true; - } - maxColumn = std::max(maxColumn, elem.getColumn()); - } - assert(changed || highestColumn == maxColumn); - highestColumn = maxColumn; - assert(changed || lastColumn == columnsAndValues[columnsAndValues.size() - 1].getColumn()); - lastColumn = columnsAndValues[columnsAndValues.size() - 1].getColumn(); - - if (changed) { - fixColumns(); - } - return changed; - } - template<typename ValueType> - void SparseMatrixBuilder<ValueType>::fixColumns() { - // Sort columns per row - typename SparseMatrix<ValueType>::index_type endGroups; - typename SparseMatrix<ValueType>::index_type endRows; - - if (hasCustomRowGrouping) { - for (index_type group = 0; group < rowGroupIndices.get().size(); ++group) { - endGroups = group < rowGroupIndices.get().size()-1 ? rowGroupIndices.get()[group+1] : rowIndications.size(); - for (index_type i = rowGroupIndices.get()[group]; i < endGroups; ++i) { - endRows = i < rowIndications.size()-1 ? rowIndications[i+1] : columnsAndValues.size(); - // Sort the row - std::sort(columnsAndValues.begin() + rowIndications[i], columnsAndValues.begin() + endRows, - [](MatrixEntry<index_type, value_type> const& a, MatrixEntry<index_type, value_type> const& b) { - return a.getColumn() < b.getColumn(); - }); - // Assert no equal elements - assert(std::is_sorted(columnsAndValues.begin() + rowIndications[i], columnsAndValues.begin() + endRows, - [](MatrixEntry<index_type, value_type> const& a, MatrixEntry<index_type, value_type> const& b) { - return a.getColumn() <= b.getColumn(); - })); + for (index_type row = 0; row < rowIndications.size(); ++row) { + bool changed = false; + auto startRow = std::next(columnsAndValues.begin(), rowIndications[row]); + auto endRow = row < rowIndications.size()-1 ? std::next(columnsAndValues.begin(), rowIndications[row+1]) : columnsAndValues.end(); + for (auto entry = startRow; entry != endRow; ++entry) { + if (entry->getColumn() >= offset) { + // Change column + entry->setColumn(replacements[entry->getColumn() - offset]); + changed = true; } + maxColumn = std::max(maxColumn, entry->getColumn()); } - } else { - for (index_type i = 0; i < rowIndications.size(); ++i) { - endRows = i < rowIndications.size()-1 ? rowIndications[i+1] : columnsAndValues.size(); - // Sort the row - std::sort(columnsAndValues.begin() + rowIndications[i], columnsAndValues.begin() + endRows, + if (changed) { + // Sort columns in row + std::sort(startRow, endRow, [](MatrixEntry<index_type, value_type> const& a, MatrixEntry<index_type, value_type> const& b) { return a.getColumn() < b.getColumn(); }); // Assert no equal elements - assert(std::is_sorted(columnsAndValues.begin() + rowIndications[i], columnsAndValues.begin() + endRows, - [](MatrixEntry<index_type, value_type> const& a, MatrixEntry<index_type, value_type> const& b) { - return a.getColumn() <= b.getColumn(); - })); + STORM_LOG_ASSERT(std::is_sorted(startRow, endRow, + [](MatrixEntry<index_type, value_type> const& a, MatrixEntry<index_type, value_type> const& b) { + return a.getColumn() < b.getColumn(); + }), "Columns not sorted."); } - } + + highestColumn = maxColumn; + lastColumn = columnsAndValues[columnsAndValues.size() - 1].getColumn(); } template<typename ValueType> @@ -688,6 +661,9 @@ namespace storm { *writeIt = std::move(intermediateRowEntry); ++writeIt; } + } else { + // skip the intermediate rows + writeIt = getRow(smallerRow).begin(); } // Write the larger row to its new position. @@ -720,6 +696,9 @@ namespace storm { *writeIt = std::move(*intermediateRowEntryIt); --writeIt; } + } else { + // skip the intermediate rows + writeIt = getRow(smallerRow).end() - 1; } // Write the larger row to its new position. @@ -1139,12 +1118,12 @@ namespace storm { #ifdef STORM_HAVE_CARL template<> - typename std::pair<storm::storage::SparseMatrix<RationalFunction>, std::vector<RationalFunction>> SparseMatrix<RationalFunction>::getJacobiDecomposition() const { + typename std::pair<storm::storage::SparseMatrix<Interval>, std::vector<Interval>> SparseMatrix<Interval>::getJacobiDecomposition() const { STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "This operation is not supported."); } template<> - typename std::pair<storm::storage::SparseMatrix<Interval>, std::vector<Interval>> SparseMatrix<Interval>::getJacobiDecomposition() const { + typename std::pair<storm::storage::SparseMatrix<RationalFunction>, std::vector<RationalFunction>> SparseMatrix<RationalFunction>::getJacobiDecomposition() const { STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "This operation is not supported."); } #endif @@ -1538,7 +1517,7 @@ namespace storm { void SparseMatrix<ValueType>::printAsMatlabMatrix(std::ostream& out) const { // Iterate over all row groups. for (typename SparseMatrix<ValueType>::index_type group = 0; group < this->getRowGroupCount(); ++group) { - assert(this->getRowGroupSize(group) == 1); + STORM_LOG_ASSERT(this->getRowGroupSize(group) == 1, "Incorrect row group size."); for (typename SparseMatrix<ValueType>::index_type i = this->getRowGroupIndices()[group]; i < this->getRowGroupIndices()[group + 1]; ++i) { typename SparseMatrix<ValueType>::index_type nextIndex = this->rowIndications[i]; diff --git a/src/storage/SparseMatrix.h b/src/storage/SparseMatrix.h index 5389e6e85..b9cf3831c 100644 --- a/src/storage/SparseMatrix.h +++ b/src/storage/SparseMatrix.h @@ -11,6 +11,7 @@ #include <boost/optional.hpp> #include "src/utility/OsDetection.h" +#include "src/utility/macros.h" #include "src/adapters/CarlAdapter.h" // Forward declaration for adapter classes. @@ -229,9 +230,8 @@ namespace storm { * * @param replacements Mapping indicating the replacements from offset+i -> value of i. * @param offset Offset to add to each id in vector index. - * @return True if replacement took place, False if nothing changed. */ - bool replaceColumns(std::vector<index_type> const& replacements, index_type offset); + void replaceColumns(std::vector<index_type> const& replacements, index_type offset); private: // A flag indicating whether a row count was set upon construction. @@ -294,11 +294,6 @@ namespace storm { // Stores the currently active row group. This is used for correctly constructing the row grouping of the // matrix. index_type currentRowGroup; - - /*! - * Fixes the matrix by sorting the columns to gain increasing order again. - */ - void fixColumns(); }; /*! diff --git a/src/storage/bisimulation/BisimulationDecomposition.h b/src/storage/bisimulation/BisimulationDecomposition.h index e3322444c..436d42660 100644 --- a/src/storage/bisimulation/BisimulationDecomposition.h +++ b/src/storage/bisimulation/BisimulationDecomposition.h @@ -118,7 +118,7 @@ namespace storm { } OptimizationDirection getOptimizationDirection() const { - assert(optimalityType); + STORM_LOG_ASSERT(optimalityType, "Optimality type not set."); return optimalityType.get(); } // A flag that indicates whether a measure driven initial partition is to be used. If this flag is set @@ -217,10 +217,10 @@ namespace storm { /*! * Refines the partition by considering the given splitter. All blocks that become potential splitters - * because of this refinement, are marked as splitters and inserted into the splitter queue. + * because of this refinement, are marked as splitters and inserted into the splitter vector. * * @param splitter The splitter to use. - * @param splitterQueue The queue into which to insert the newly discovered potential splitters. + * @param splitterVector The vector into which to insert the newly discovered potential splitters. */ virtual void refinePartitionBasedOnSplitter(bisimulation::Block<BlockDataType>& splitter, std::vector<bisimulation::Block<BlockDataType>*>& splitterQueue) = 0; diff --git a/src/storage/bisimulation/DeterministicModelBisimulationDecomposition.cpp b/src/storage/bisimulation/DeterministicModelBisimulationDecomposition.cpp index b030b7cee..8693aa9aa 100644 --- a/src/storage/bisimulation/DeterministicModelBisimulationDecomposition.cpp +++ b/src/storage/bisimulation/DeterministicModelBisimulationDecomposition.cpp @@ -181,7 +181,7 @@ namespace storm { }); - // If the predecessor block was split, we need to insert it into the splitter queue if it is not already + // If the predecessor block was split, we need to insert it into the splitter vector if it is not already // marked as a splitter. if (split && !blockToRefineProbabilistically->data().splitter()) { splitterQueue.emplace_back(blockToRefineProbabilistically); diff --git a/src/storage/dft/DFT.cpp b/src/storage/dft/DFT.cpp index 5c70c475d..88be1df73 100644 --- a/src/storage/dft/DFT.cpp +++ b/src/storage/dft/DFT.cpp @@ -14,7 +14,7 @@ namespace storm { template<typename ValueType> DFT<ValueType>::DFT(DFTElementVector const& elements, DFTElementPointer const& tle) : mElements(elements), mNrOfBEs(0), mNrOfSpares(0), mTopLevelIndex(tle->id()), mMaxSpareChildCount(0) { - assert(elementIndicesCorrect()); + STORM_LOG_ASSERT(elementIndicesCorrect(), "Ids incorrect."); size_t nrRepresentatives = 0; for (auto& elem : mElements) { @@ -100,9 +100,9 @@ namespace storm { } else { // Generate information according to symmetries for (size_t symmetryIndex : symmetries.sortedSymmetries) { - assert(!visited[symmetryIndex]); + STORM_LOG_ASSERT(!visited[symmetryIndex], "Element already considered for symmetry."); auto const& symmetryGroup = symmetries.groups.at(symmetryIndex); - assert(!symmetryGroup.empty()); + STORM_LOG_ASSERT(!symmetryGroup.empty(), "No symmetry available."); // Insert all elements of first subtree of each symmetry size_t groupIndex = stateIndex; @@ -117,14 +117,14 @@ namespace storm { // Mirror symmetries size_t noSymmetricElements = symmetryGroup.front().size(); - assert(noSymmetricElements > 1); + STORM_LOG_ASSERT(noSymmetricElements > 1, "No symmetry available."); for (std::vector<size_t> symmetricElements : symmetryGroup) { - assert(symmetricElements.size() == noSymmetricElements); + STORM_LOG_ASSERT(symmetricElements.size() == noSymmetricElements, "No. of symmetric elements do not coincide."); if (visited[symmetricElements[1]]) { // Elements already mirrored for (size_t index : symmetricElements) { - assert(visited[index]); + STORM_LOG_ASSERT(visited[index], "Element not mirrored."); } continue; } @@ -143,13 +143,13 @@ namespace storm { generationInfo.addStateIndex(symmetricElement, index + offset * i); stateIndex += 2; - assert((activationIndex > 0) == isRepresentative(symmetricElement)); + STORM_LOG_ASSERT((activationIndex > 0) == isRepresentative(symmetricElement), "Bits for representative incorrect."); if (activationIndex > 0) { generationInfo.addSpareActivationIndex(symmetricElement, activationIndex + offset * i); ++stateIndex; } - assert((usageIndex > 0) == mElements[symmetricElement]->isSpareGate()); + STORM_LOG_ASSERT((usageIndex > 0) == mElements[symmetricElement]->isSpareGate(), "Bits for usage incorrect."); if (usageIndex > 0) { generationInfo.addSpareUsageIndex(symmetricElement, usageIndex + offset * i); stateIndex += generationInfo.usageInfoBits(); @@ -187,15 +187,15 @@ namespace storm { generationInfo.generateSymmetries(symmetries); STORM_LOG_TRACE(generationInfo); - assert(stateIndex == mStateVectorSize); - assert(visited.full()); + STORM_LOG_ASSERT(stateIndex == mStateVectorSize, "Id incorrect."); + STORM_LOG_ASSERT(visited.full(), "Not all elements considered."); return generationInfo; } template<typename ValueType> size_t DFT<ValueType>::generateStateInfo(DFTStateGenerationInfo& generationInfo, size_t id, storm::storage::BitVector& visited, size_t stateIndex) const { - assert(!visited[id]); + STORM_LOG_ASSERT(!visited[id], "Element already visited."); visited.set(id); // Reserve bits for element @@ -238,19 +238,19 @@ namespace storm { template<typename ValueType> std::vector<DFT<ValueType>> DFT<ValueType>::topModularisation() const { - assert(isGate(mTopLevelIndex)); + STORM_LOG_ASSERT(isGate(mTopLevelIndex), "Top level element is no gate."); auto const& children = getGate(mTopLevelIndex)->children(); std::map<size_t, std::vector<size_t>> subdfts; for(auto const& child : children) { std::vector<size_t> isubdft; - if(child->nrParents() > 1 || child->hasOutgoingDependencies()) { - STORM_LOG_TRACE("child " << child->name() << "does not allow modularisation."); + if(child->nrParents() > 1 || child->hasOutgoingDependencies() || child->hasRestrictions()) { + STORM_LOG_TRACE("child " << child->name() << " does not allow modularisation."); return {*this}; } if (isGate(child->id())) { isubdft = getGate(child->id())->independentSubDft(false); } else { - assert(isBasicElement(child->id())); + STORM_LOG_ASSERT(isBasicElement(child->id()), "Child is no BE."); if(getBasicElement(child->id())->hasIngoingDependencies()) { STORM_LOG_TRACE("child " << child->name() << "does not allow modularisation."); return {*this}; @@ -309,16 +309,16 @@ namespace storm { // Add rewritten elements for (std::vector<size_t> rewrites : rewriteIds) { - assert(rewrites.size() > 1); - assert(mElements[rewrites[1]]->hasParents()); - assert(mElements[rewrites[1]]->parents().front()->isGate()); + STORM_LOG_ASSERT(rewrites.size() > 1, "No rewritten elements."); + STORM_LOG_ASSERT(mElements[rewrites[1]]->hasParents(), "Rewritten elements has no parents."); + STORM_LOG_ASSERT(mElements[rewrites[1]]->parents().front()->isGate(), "Rewritten element has no parent gate."); DFTGatePointer originalParent = std::static_pointer_cast<DFTGate<ValueType>>(mElements[rewrites[1]]->parents().front()); std::string newParentName = builder.getUniqueName(originalParent->name()); // Accumulate children names std::vector<std::string> childrenNames; for (size_t i = 1; i < rewrites.size(); ++i) { - assert(mElements[rewrites[i]]->parents().front()->id() == originalParent->id()); // Children have the same father + STORM_LOG_ASSERT(mElements[rewrites[i]]->parents().front()->id() == originalParent->id(), "Children have the same father"); childrenNames.push_back(mElements[rewrites[i]]->name()); } @@ -390,7 +390,7 @@ namespace storm { stream << "}" << std::endl; return stream.str(); } - assert(it != mTopModule.end()); + STORM_LOG_ASSERT(it != mTopModule.end(), "Element not found."); stream << mElements[(*it)]->name(); ++it; while(it != mTopModule.end()) { @@ -403,7 +403,7 @@ namespace storm { stream << "[" << mElements[spareModule.first]->name() << "] = {"; if (!spareModule.second.empty()) { std::vector<size_t>::const_iterator it = spareModule.second.begin(); - assert(it != spareModule.second.end()); + STORM_LOG_ASSERT(it != spareModule.second.end(), "Element not found."); stream << mElements[(*it)]->name(); ++it; while(it != spareModule.second.end()) { @@ -481,8 +481,7 @@ namespace storm { } else { useId = getChild(elem->id(), nrUsedChild); } - bool isActive = status[stateGenerationInfo.getSpareActivationIndex(useId)]; - if(useId == elem->id() || isActive) { + if(useId == elem->id() || status[stateGenerationInfo.getSpareActivationIndex(useId)]) { stream << "actively "; } stream << "using " << useId << "]"; @@ -494,13 +493,13 @@ namespace storm { template<typename ValueType> size_t DFT<ValueType>::getChild(size_t spareId, size_t nrUsedChild) const { - assert(mElements[spareId]->isSpareGate()); + STORM_LOG_ASSERT(mElements[spareId]->isSpareGate(), "Element is no spare."); return getGate(spareId)->children()[nrUsedChild]->id(); } template<typename ValueType> size_t DFT<ValueType>::getNrChild(size_t spareId, size_t childId) const { - assert(mElements[spareId]->isSpareGate()); + STORM_LOG_ASSERT(mElements[spareId]->isSpareGate(), "Element is no spare."); DFTElementVector children = getGate(spareId)->children(); for (size_t nrChild = 0; nrChild < children.size(); ++nrChild) { if (children[nrChild]->id() == childId) { @@ -528,6 +527,11 @@ namespace storm { } } + template<typename ValueType> + bool DFT<ValueType>::canHaveNondeterminism() const { + return !getDependencies().empty(); + } + template<typename ValueType> DFTColouring<ValueType> DFT<ValueType>::colourDFT() const { return DFTColouring<ValueType>(*this); @@ -551,8 +555,8 @@ namespace storm { } } - assert(isGate(index1)); - assert(isGate(index2)); + STORM_LOG_ASSERT(isGate(index1), "Element is no gate."); + STORM_LOG_ASSERT(isGate(index2), "Element is no gate."); std::vector<size_t> isubdft1 = getGate(index1)->independentSubDft(false); std::vector<size_t> isubdft2 = getGate(index2)->independentSubDft(false); if(isubdft1.empty() || isubdft2.empty() || isubdft1.size() != isubdft2.size()) { @@ -591,7 +595,7 @@ namespace storm { size_t childLeftId = spareLeft->children().at(i)->id(); size_t childRightId = spareRight->children().at(i)->id(); - assert(bijection.count(childLeftId) == 0); + STORM_LOG_ASSERT(bijection.count(childLeftId) == 0, "Child already part of bijection."); if (childLeftId == childRightId) { // Ignore shared child continue; @@ -692,6 +696,8 @@ namespace storm { // suitable parent gate! - Lets check the independent submodules of the children auto const& children = std::static_pointer_cast<DFTGate<ValueType>>(e)->children(); for(auto const& child : children) { + + auto ISD = std::static_pointer_cast<DFTGate<ValueType>>(child)->independentSubDft(true); // In the ISD, check for other children: diff --git a/src/storage/dft/DFT.h b/src/storage/dft/DFT.h index 6fdfd08ef..7620a059a 100644 --- a/src/storage/dft/DFT.h +++ b/src/storage/dft/DFT.h @@ -118,7 +118,7 @@ namespace storm { if(representativeId == mTopLevelIndex) { return mTopModule; } else { - assert(mSpareModules.count(representativeId)>0); + STORM_LOG_ASSERT(mSpareModules.count(representativeId) > 0, "Representative not found."); return mSpareModules.find(representativeId)->second; } } @@ -145,7 +145,7 @@ namespace storm { * @param index The id of the element */ DFTElementCPointer getElement(size_t index) const { - assert(index < nrElements()); + STORM_LOG_ASSERT(index < nrElements(), "Index invalid."); return mElements[index]; } @@ -166,7 +166,7 @@ namespace storm { } std::shared_ptr<DFTBE<ValueType> const> getBasicElement(size_t index) const { - assert(isBasicElement(index)); + STORM_LOG_ASSERT(isBasicElement(index), "Element is no BE."); return std::static_pointer_cast<DFTBE<ValueType> const>(mElements[index]); } @@ -175,17 +175,17 @@ namespace storm { } std::shared_ptr<DFTGate<ValueType> const> getGate(size_t index) const { - assert(isGate(index)); + STORM_LOG_ASSERT(isGate(index), "Element is no gate."); return std::static_pointer_cast<DFTGate<ValueType> const>(mElements[index]); } std::shared_ptr<DFTDependency<ValueType> const> getDependency(size_t index) const { - assert(isDependency(index)); + STORM_LOG_ASSERT(isDependency(index), "Element is no dependency."); return std::static_pointer_cast<DFTDependency<ValueType> const>(mElements[index]); } std::shared_ptr<DFTRestriction<ValueType> const> getRestriction(size_t index) const { - assert(isRestriction(index)); + STORM_LOG_ASSERT(isRestriction(index), "Element is no restriction."); return std::static_pointer_cast<DFTRestriction<ValueType> const>(mElements[index]); } @@ -198,6 +198,8 @@ namespace storm { } return elements; } + + bool canHaveNondeterminism() const; std::vector<DFT<ValueType>> topModularisation() const; @@ -214,9 +216,9 @@ namespace storm { return mRepresentants.find(id) != mRepresentants.end(); } - DFTElementCPointer getRepresentant(size_t id) const { - assert(hasRepresentant(id)); - return getElement(mRepresentants.find(id)->second); + size_t getRepresentant(size_t id) const { + STORM_LOG_ASSERT(hasRepresentant(id), "Element has no representant."); + return mRepresentants.find(id)->second; } bool hasFailed(DFTStatePointer const& state) const { diff --git a/src/storage/dft/DFTBuilder.cpp b/src/storage/dft/DFTBuilder.cpp index 019ebe6e1..256bd95db 100644 --- a/src/storage/dft/DFTBuilder.cpp +++ b/src/storage/dft/DFTBuilder.cpp @@ -32,7 +32,7 @@ namespace storm { // Child not found -> find first dependent event to assure that child is dependency // TODO: Not sure whether this is the intended behaviour? auto itFind = mElements.find(child + "_1"); - STORM_LOG_ASSERT(itFind != mElements.end(), "Child not found."); + STORM_LOG_ASSERT(itFind != mElements.end(), "Child '" << child << "' for gate '" << gate->name() << "' not found."); STORM_LOG_ASSERT(itFind->second->isDependency(), "Child is no dependency."); STORM_LOG_TRACE("Ignore functional dependency " << child << " in gate " << gate->name()); } diff --git a/src/storage/dft/DFTElementState.h b/src/storage/dft/DFTElementState.h index d191cc2d9..548eb6c88 100644 --- a/src/storage/dft/DFTElementState.h +++ b/src/storage/dft/DFTElementState.h @@ -2,7 +2,7 @@ #ifndef DFTELEMENTSTATE_H #define DFTELEMENTSTATE_H -#include <cassert> +#include "src/utility/macros.h" namespace storm { namespace storage { @@ -18,8 +18,10 @@ namespace storm { return os << "Failsafe"; case DFTElementState::DontCare: return os << "Don't Care"; + default: + STORM_LOG_ASSERT(false, "Element state not known."); + return os; } - assert(false); } inline char toChar(DFTElementState st) { @@ -32,8 +34,10 @@ namespace storm { return 'S'; case DFTElementState::DontCare: return '-'; + default: + STORM_LOG_ASSERT(false, "Element state not known."); + return ' '; } - assert(false); } enum class DFTDependencyState {Passive = 0, Unsuccessful = 1, Successful = 2, DontCare = 3}; @@ -48,8 +52,10 @@ namespace storm { return os << "Unsuccessful"; case DFTDependencyState::DontCare: return os << "Don't Care"; + default: + STORM_LOG_ASSERT(false, "Element state not known."); + return os; } - assert(false); } inline char toChar(DFTDependencyState st) { @@ -62,8 +68,10 @@ namespace storm { return 'U'; case DFTDependencyState::DontCare: return '-'; + default: + STORM_LOG_ASSERT(false, "Element state not known."); + return ' '; } - assert(false); } } diff --git a/src/storage/dft/DFTIsomorphism.h b/src/storage/dft/DFTIsomorphism.h index d94640d9c..1f66e5449 100644 --- a/src/storage/dft/DFTIsomorphism.h +++ b/src/storage/dft/DFTIsomorphism.h @@ -1,6 +1,5 @@ #pragma once -#include <cassert> #include <vector> #include <unordered_map> #include <utility> @@ -197,7 +196,7 @@ namespace storage { } else if(dft.isDependency(id)) { colourize(dft.getDependency(id)); } else { - assert(dft.isRestriction(id)); + STORM_LOG_ASSERT(dft.isRestriction(id), "Element is no restriction."); colourize(dft.getRestriction(id)); } } @@ -233,7 +232,7 @@ namespace storage { res.pdepCandidates[depColour.at(index)] = std::vector<size_t>({index}); } } else { - assert(dft.isRestriction(index)); + STORM_LOG_ASSERT(dft.isRestriction(index), "Element is no restriction."); auto it = res.restrictionCandidates.find(restrictionColour.at(index)); if(it != res.restrictionCandidates.end()) { it->second.push_back(index); @@ -342,14 +341,14 @@ namespace storage { * Construct the initial bijection. */ void constructInitialBijection() { - assert(candidatesCompatible); + STORM_LOG_ASSERT(candidatesCompatible, "Candidates are not compatible."); // We first construct the currentPermutations, which helps to determine the current state of the check. initializePermutationsAndTreatTrivialGroups(bleft.beCandidates, bright.beCandidates, currentPermutations.beCandidates); initializePermutationsAndTreatTrivialGroups(bleft.gateCandidates, bright.gateCandidates, currentPermutations.gateCandidates); initializePermutationsAndTreatTrivialGroups(bleft.pdepCandidates, bright.pdepCandidates, currentPermutations.pdepCandidates); initializePermutationsAndTreatTrivialGroups(bleft.restrictionCandidates, bright.restrictionCandidates, currentPermutations.restrictionCandidates); STORM_LOG_TRACE(bijection.size() << " vs. " << bleft.size() << " vs. " << bright.size()); - assert(bijection.size() == bleft.size()); + STORM_LOG_ASSERT(bijection.size() == bleft.size(), "No. of bijection elements do not match."); } @@ -358,7 +357,7 @@ namespace storage { * @return true if a next bijection exists. */ bool findNextBijection() { - assert(candidatesCompatible); + STORM_LOG_ASSERT(candidatesCompatible, "Candidates are not compatible."); bool foundNext = false; if(!currentPermutations.beCandidates.empty()) { auto it = currentPermutations.beCandidates.begin(); @@ -395,28 +394,28 @@ namespace storage { if(foundNext) { for(auto const& colour : bleft.beCandidates) { if (colour.second.size() > 1) { - assert(currentPermutations.beCandidates.find(colour.first) != currentPermutations.beCandidates.end()); + STORM_LOG_ASSERT(currentPermutations.beCandidates.find(colour.first) != currentPermutations.beCandidates.end(), "Colour not found."); zipVectorsIntoMap(colour.second, currentPermutations.beCandidates.find(colour.first)->second, bijection); } } for(auto const& colour : bleft.gateCandidates) { if (colour.second.size() > 1) { - assert(currentPermutations.gateCandidates.find(colour.first) != currentPermutations.gateCandidates.end()); + STORM_LOG_ASSERT(currentPermutations.gateCandidates.find(colour.first) != currentPermutations.gateCandidates.end(), "Colour not found."); zipVectorsIntoMap(colour.second, currentPermutations.gateCandidates.find(colour.first)->second, bijection); } } for(auto const& colour : bleft.pdepCandidates) { if (colour.second.size() > 1) { - assert(currentPermutations.pdepCandidates.find(colour.first) != currentPermutations.pdepCandidates.end()); + STORM_LOG_ASSERT(currentPermutations.pdepCandidates.find(colour.first) != currentPermutations.pdepCandidates.end(), "Colour not found."); zipVectorsIntoMap(colour.second, currentPermutations.pdepCandidates.find(colour.first)->second, bijection); } } for(auto const& colour : bleft.restrictionCandidates) { if (colour.second.size() > 1) { - assert(currentPermutations.restrictionCandidates.find(colour.first) != currentPermutations.restrictionCandidates.end()); + STORM_LOG_ASSERT(currentPermutations.restrictionCandidates.find(colour.first) != currentPermutations.restrictionCandidates.end(), "Colour not found."); zipVectorsIntoMap(colour.second, currentPermutations.restrictionCandidates.find(colour.first)->second, bijection); } } @@ -430,13 +429,13 @@ namespace storage { * */ bool check() const { - assert(bijection.size() == bleft.size()); + STORM_LOG_ASSERT(bijection.size() == bleft.size(), "No. of bijection elements do not match."); // We can skip BEs, as they are identified by they're homomorphic if they are in the same class for(auto const& indexpair : bijection) { // Check type first. Colouring takes care of a lot, but not necesarily everything (e.g. voting thresholds) equalType(*dft.getElement(indexpair.first), *dft.getElement(indexpair.second)); if(dft.isGate(indexpair.first)) { - assert(dft.isGate(indexpair.second)); + STORM_LOG_ASSERT(dft.isGate(indexpair.second), "Element is no gate."); auto const& lGate = dft.getGate(indexpair.first); auto const& rGate = dft.getGate(indexpair.second); if(lGate->isDynamicGate()) { @@ -483,7 +482,7 @@ namespace storage { } else if(dft.isDependency(indexpair.first)) { - assert(dft.isDependency(indexpair.second)); + STORM_LOG_ASSERT(dft.isDependency(indexpair.second), "Element is no dependency."); auto const& lDep = dft.getDependency(indexpair.first); auto const& rDep = dft.getDependency(indexpair.second); if(bijection.at(lDep->triggerEvent()->id()) != rDep->triggerEvent()->id()) { @@ -493,7 +492,7 @@ namespace storage { return false; } } else if(dft.isRestriction(indexpair.first)) { - assert(dft.isRestriction(indexpair.second)); + STORM_LOG_ASSERT(dft.isRestriction(indexpair.second), "Element is no restriction."); auto const& lRestr = dft.getRestriction(indexpair.first); std::vector<size_t> childrenLeftMapped; for(auto const& child : lRestr->children() ) { @@ -521,8 +520,8 @@ namespace storage { } } else { - assert(dft.isBasicElement(indexpair.first)); - assert(dft.isBasicElement(indexpair.second)); + STORM_LOG_ASSERT(dft.isBasicElement(indexpair.first), "Element is no BE."); + STORM_LOG_ASSERT(dft.isBasicElement(indexpair.second), "Element is no BE."); // No operations required. } } @@ -590,12 +589,12 @@ namespace storage { for(auto const& colour : right) { if(colour.second.size()>1) { auto it = permutations.insert(colour); - assert(it.second); + STORM_LOG_ASSERT(it.second, "Element already contained."); std::sort(it.first->second.begin(), it.first->second.end()); zipVectorsIntoMap(left.at(colour.first), it.first->second, bijection); } else { - assert(colour.second.size() == 1); - assert(bijection.count(left.at(colour.first).front()) == 0); + STORM_LOG_ASSERT(colour.second.size() == 1, "No elements for colour."); + STORM_LOG_ASSERT(bijection.count(left.at(colour.first).front()) == 0, "Element already contained."); bijection[left.at(colour.first).front()] = colour.second.front(); } } @@ -606,7 +605,7 @@ namespace storage { */ void zipVectorsIntoMap(std::vector<size_t> const& a, std::vector<size_t> const& b, std::map<size_t, size_t>& map) const { // Assert should pass due to compatibility check - assert(a.size() == b.size()); + STORM_LOG_ASSERT(a.size() == b.size(), "Sizes do not match."); auto it = b.cbegin(); for(size_t lIndex : a) { map[lIndex] = *it; @@ -638,4 +637,4 @@ namespace std { return hasher(p.first) ^ hasher(p.second); } }; -} \ No newline at end of file +} diff --git a/src/storage/dft/DFTState.cpp b/src/storage/dft/DFTState.cpp index 25edf1d38..cc0eaabec 100644 --- a/src/storage/dft/DFTState.cpp +++ b/src/storage/dft/DFTState.cpp @@ -6,33 +6,43 @@ namespace storm { namespace storage { template<typename ValueType> - DFTState<ValueType>::DFTState(DFT<ValueType> const& dft, DFTStateGenerationInfo const& stateGenerationInfo, size_t id) : mStatus(dft.stateVectorSize()), mId(id), mDft(dft), mStateGenerationInfo(stateGenerationInfo) { + DFTState<ValueType>::DFTState(DFT<ValueType> const& dft, DFTStateGenerationInfo const& stateGenerationInfo, size_t id) : mStatus(dft.stateVectorSize()), mId(id), mPseudoState(false), mDft(dft), mStateGenerationInfo(stateGenerationInfo) { + // TODO Matthias: use construct() // Initialize uses - for(size_t id : mDft.getSpareIndices()) { - std::shared_ptr<DFTGate<ValueType> const> elem = mDft.getGate(id); - assert(elem->isSpareGate()); - assert(elem->nrChildren() > 0); - this->setUses(id, elem->children()[0]->id()); + for(size_t spareId : mDft.getSpareIndices()) { + std::shared_ptr<DFTGate<ValueType> const> elem = mDft.getGate(spareId); + STORM_LOG_ASSERT(elem->isSpareGate(), "Element is no spare gate."); + STORM_LOG_ASSERT(elem->nrChildren() > 0, "Element has no child."); + this->setUses(spareId, elem->children()[0]->id()); } - + // Initialize activation propagateActivation(mDft.getTopLevelIndex()); - std::vector<size_t> alwaysActiveBEs = dft.nonColdBEs(); - mIsCurrentlyFailableBE.insert(mIsCurrentlyFailableBE.end(), alwaysActiveBEs.begin(), alwaysActiveBEs.end()); - sortFailableBEs(); + std::vector<size_t> alwaysActiveBEs = mDft.nonColdBEs(); + mCurrentlyFailableBE.insert(mCurrentlyFailableBE.end(), alwaysActiveBEs.begin(), alwaysActiveBEs.end()); + } + + template<typename ValueType> + DFTState<ValueType>::DFTState(storm::storage::BitVector const& status, DFT<ValueType> const& dft, DFTStateGenerationInfo const& stateGenerationInfo, size_t id) : mStatus(status), mId(id), mPseudoState(true), mDft(dft), mStateGenerationInfo(stateGenerationInfo) { + // Intentionally left empty } template<typename ValueType> - DFTState<ValueType>::DFTState(storm::storage::BitVector const& status, DFT<ValueType> const& dft, DFTStateGenerationInfo const& stateGenerationInfo, size_t id) : mStatus(status), mId(id), mDft(dft), mStateGenerationInfo(stateGenerationInfo) { - + void DFTState<ValueType>::construct() { + STORM_LOG_TRACE("Construct concrete state from pseudo state " << mDft.getStateString(mStatus, mStateGenerationInfo, mId)); + // Clear information from pseudo state + mCurrentlyFailableBE.clear(); + mFailableDependencies.clear(); + mUsedRepresentants.clear(); + STORM_LOG_ASSERT(mPseudoState, "Only pseudo states can be constructed."); for(size_t index = 0; index < mDft.nrElements(); ++index) { // Initialize currently failable BE if (mDft.isBasicElement(index) && isOperational(index)) { std::shared_ptr<const DFTBE<ValueType>> be = mDft.getBasicElement(index); - if ((!be->isColdBasicElement() && be->canFail()) || !mDft.hasRepresentant(index) || isActive(mDft.getRepresentant(index)->id())) { - mIsCurrentlyFailableBE.push_back(index); + if ((!be->isColdBasicElement() && be->canFail()) || !mDft.hasRepresentant(index) || isActive(mDft.getRepresentant(index))) { + mCurrentlyFailableBE.push_back(index); STORM_LOG_TRACE("Currently failable: " << mDft.getBasicElement(index)->toString()); } } else if (mDft.getElement(index)->isSpareGate()) { @@ -42,17 +52,22 @@ namespace storm { STORM_LOG_TRACE("Spare " << index << " uses " << useId); } } - sortFailableBEs(); - + // Initialize failable dependencies for (size_t dependencyId : mDft.getDependencies()) { std::shared_ptr<DFTDependency<ValueType> const> dependency = mDft.getDependency(dependencyId); - assert(dependencyId == dependency->id()); + STORM_LOG_ASSERT(dependencyId == dependency->id(), "Ids do not match."); if (hasFailed(dependency->triggerEvent()->id()) && getElementState(dependency->dependentEvent()->id()) == DFTElementState::Operational) { mFailableDependencies.push_back(dependencyId); STORM_LOG_TRACE("New dependency failure: " << dependency->toString()); } } + mPseudoState = false; + } + + template<typename ValueType> + std::shared_ptr<DFTState<ValueType>> DFTState<ValueType>::copy() const { + return std::make_shared<storm::storage::DFTState<ValueType>>(*this); } template<typename ValueType> @@ -169,9 +184,9 @@ namespace storm { template<typename ValueType> void DFTState<ValueType>::beNoLongerFailable(size_t id) { - auto it = std::find(mIsCurrentlyFailableBE.begin(), mIsCurrentlyFailableBE.end(), id); - if(it != mIsCurrentlyFailableBE.end()) { - mIsCurrentlyFailableBE.erase(it); + auto it = std::find(mCurrentlyFailableBE.begin(), mCurrentlyFailableBE.end(), id); + if (it != mCurrentlyFailableBE.end()) { + mCurrentlyFailableBE.erase(it); } } @@ -182,9 +197,9 @@ namespace storm { } for (auto dependency : mDft.getElement(id)->outgoingDependencies()) { - assert(dependency->triggerEvent()->id() == id); + STORM_LOG_ASSERT(dependency->triggerEvent()->id() == id, "Ids do not match."); if (getElementState(dependency->dependentEvent()->id()) == DFTElementState::Operational) { - assert(!isFailsafe(dependency->dependentEvent()->id())); + STORM_LOG_ASSERT(!isFailsafe(dependency->dependentEvent()->id()), "Dependent event is failsafe."); mFailableDependencies.push_back(dependency->id()); STORM_LOG_TRACE("New dependency failure: " << dependency->toString()); } @@ -194,22 +209,40 @@ namespace storm { template<typename ValueType> void DFTState<ValueType>::updateDontCareDependencies(size_t id) { - assert(mDft.isBasicElement(id)); - assert(hasFailed(id)); + STORM_LOG_ASSERT(mDft.isBasicElement(id), "Element is no BE."); + STORM_LOG_ASSERT(hasFailed(id), "Element has not failed."); for (auto dependency : mDft.getBasicElement(id)->ingoingDependencies()) { - assert(dependency->dependentEvent()->id() == id); + STORM_LOG_ASSERT(dependency->dependentEvent()->id() == id, "Ids do not match."); setDependencyDontCare(dependency->id()); } } template<typename ValueType> - std::pair<std::shared_ptr<DFTBE<ValueType> const>, bool> DFTState<ValueType>::letNextBEFail(size_t index) - { + ValueType DFTState<ValueType>::getBERate(size_t id) const { + STORM_LOG_ASSERT(mDft.isBasicElement(id), "Element is no BE."); + + if (mDft.hasRepresentant(id) && !isActive(mDft.getRepresentant(id))) { + // Return passive failure rate + return mDft.getBasicElement(id)->passiveFailureRate(); + } else { + // Return active failure rate + return mDft.getBasicElement(id)->activeFailureRate(); + } + } + + template<typename ValueType> + ValueType DFTState<ValueType>::getFailableBERate(size_t index) const { + STORM_LOG_ASSERT(index < nrFailableBEs(), "Index invalid."); + return getBERate(mCurrentlyFailableBE[index]); + } + + template<typename ValueType> + std::pair<std::shared_ptr<DFTBE<ValueType> const>, bool> DFTState<ValueType>::letNextBEFail(size_t index) { STORM_LOG_TRACE("currently failable: " << getCurrentlyFailableString()); if (nrFailableDependencies() > 0) { // Consider failure due to dependency - assert(index < nrFailableDependencies()); + STORM_LOG_ASSERT(index < nrFailableDependencies(), "Index invalid."); std::shared_ptr<DFTDependency<ValueType> const> dependency = mDft.getDependency(mFailableDependencies[index]); std::pair<std::shared_ptr<DFTBE<ValueType> const>,bool> res(mDft.getBasicElement(dependency->dependentEvent()->id()), true); mFailableDependencies.erase(mFailableDependencies.begin() + index); @@ -218,10 +251,10 @@ namespace storm { return res; } else { // Consider "normal" failure - assert(index < nrFailableBEs()); - std::pair<std::shared_ptr<DFTBE<ValueType> const>,bool> res(mDft.getBasicElement(mIsCurrentlyFailableBE[index]), false); - assert(res.first->canFail()); - mIsCurrentlyFailableBE.erase(mIsCurrentlyFailableBE.begin() + index); + STORM_LOG_ASSERT(index < nrFailableBEs(), "Index invalid."); + std::pair<std::shared_ptr<DFTBE<ValueType> const>,bool> res(mDft.getBasicElement(mCurrentlyFailableBE[index]), false); + STORM_LOG_ASSERT(res.first->canFail(), "Element cannot fail."); + mCurrentlyFailableBE.erase(mCurrentlyFailableBE.begin() + index); setFailed(res.first->id()); return res; } @@ -229,7 +262,7 @@ namespace storm { template<typename ValueType> void DFTState<ValueType>::letDependencyBeUnsuccessful(size_t index) { - assert(nrFailableDependencies() > 0 && index < nrFailableDependencies()); + STORM_LOG_ASSERT(nrFailableDependencies() > 0 && index < nrFailableDependencies(), "Index invalid."); std::shared_ptr<DFTDependency<ValueType> const> dependency = mDft.getDependency(getDependencyId(index)); mFailableDependencies.erase(mFailableDependencies.begin() + index); setDependencyUnsuccessful(dependency->id()); @@ -243,7 +276,7 @@ namespace storm { template<typename ValueType> bool DFTState<ValueType>::isActive(size_t id) const { - assert(mDft.isRepresentative(id)); + STORM_LOG_ASSERT(mDft.isRepresentative(id), "Element is no representative."); return mStatus[mStateGenerationInfo.getSpareActivationIndex(id)]; } @@ -256,13 +289,13 @@ namespace storm { if(mDft.isBasicElement(elem) && isOperational(elem)) { std::shared_ptr<const DFTBE<ValueType>> be = mDft.getBasicElement(elem); if (be->isColdBasicElement() && be->canFail()) { - mIsCurrentlyFailableBE.push_back(elem); + // Add to failable BEs + mCurrentlyFailableBE.push_back(elem); } } else if (mDft.getElement(elem)->isSpareGate() && !isActive(uses(elem))) { propagateActivation(uses(elem)); } } - sortFailableBEs(); } template<typename ValueType> @@ -277,7 +310,7 @@ namespace storm { template<typename ValueType> uint_fast64_t DFTState<ValueType>::extractUses(size_t from) const { - assert(mStateGenerationInfo.usageInfoBits() < 64); + STORM_LOG_ASSERT(mStateGenerationInfo.usageInfoBits() < 64, "UsageInfoBit size too large."); return mStatus.getAsInt(from, mStateGenerationInfo.usageInfoBits()); } @@ -294,14 +327,14 @@ namespace storm { template<typename ValueType> void DFTState<ValueType>::finalizeUses(size_t spareId) { - assert(hasFailed(spareId)); + STORM_LOG_ASSERT(hasFailed(spareId), "Spare has not failed."); mStatus.setFromInt(mStateGenerationInfo.getSpareUsageIndex(spareId), mStateGenerationInfo.usageInfoBits(), mDft.getMaxSpareChildCount()); } template<typename ValueType> bool DFTState<ValueType>::hasOperationalPostSeqElements(size_t id) const { - assert(!mDft.isDependency(id)); - assert(!mDft.isRestriction(id)); + STORM_LOG_ASSERT(!mDft.isDependency(id), "Element is dependency."); + STORM_LOG_ASSERT(!mDft.isRestriction(id), "Element is restriction."); auto const& postIds = mStateGenerationInfo.seqRestrictionPostElements(id); for(size_t id : postIds) { if(isOperational(id)) { @@ -315,7 +348,7 @@ namespace storm { bool DFTState<ValueType>::claimNew(size_t spareId, size_t currentlyUses, std::vector<std::shared_ptr<DFTElement<ValueType>>> const& children) { auto it = children.begin(); while ((*it)->id() != currentlyUses) { - assert(it != children.end()); + STORM_LOG_ASSERT(it != children.end(), "Currently used element not found."); ++it; } ++it; @@ -355,18 +388,11 @@ namespace storm { n = tmp; } while (n > 0); } + if (changed) { + mPseudoState = true; + } return changed; } - - template<typename ValueType> - void DFTState<ValueType>::sortFailableBEs() { - std::sort( mIsCurrentlyFailableBE.begin( ), mIsCurrentlyFailableBE.end( ), [&](size_t const& lhs, size_t const& rhs) { - ValueType leftRate = mDft.getBasicElement(lhs)->activeFailureRate(); - ValueType rightRate = mDft.getBasicElement(rhs)->activeFailureRate(); - // Sort decreasing - return rightRate < leftRate; - }); - } // Explicitly instantiate the class. @@ -377,4 +403,4 @@ namespace storm { #endif } -} \ No newline at end of file +} diff --git a/src/storage/dft/DFTState.h b/src/storage/dft/DFTState.h index f647aae17..b077b58b0 100644 --- a/src/storage/dft/DFTState.h +++ b/src/storage/dft/DFTState.h @@ -1,12 +1,12 @@ #ifndef DFTSTATE_H #define DFTSTATE_H -#include "../BitVector.h" -#include "DFTElementState.h" +#include "src/storage/dft/DFTElementState.h" +#include "src/storage/BitVector.h" +#include "src/builder/DftExplorationHeuristic.h" #include <sstream> #include <memory> -#include <cassert> namespace storm { namespace storage { @@ -26,21 +26,41 @@ namespace storm { // Status is bitvector where each element has two bits with the meaning according to DFTElementState storm::storage::BitVector mStatus; size_t mId; - std::vector<size_t> mIsCurrentlyFailableBE; + std::vector<size_t> mCurrentlyFailableBE; std::vector<size_t> mFailableDependencies; std::vector<size_t> mUsedRepresentants; + bool mPseudoState; bool mValid = true; const DFT<ValueType>& mDft; const DFTStateGenerationInfo& mStateGenerationInfo; - + public: + /** + * Construct the initial state. + * + * @param dft DFT + * @param stateGenerationInfo General information for state generation + * @param id State id + */ DFTState(DFT<ValueType> const& dft, DFTStateGenerationInfo const& stateGenerationInfo, size_t id); - + /** - * Construct state from underlying bitvector. + * Construct temporary pseudo state. The actual state is constructed later. + * + * @param status BitVector representing the status of the state. + * @param dft DFT + * @param stateGenerationInfo General information for state generation + * @param id Pseudo state id */ DFTState(storm::storage::BitVector const& status, DFT<ValueType> const& dft, DFTStateGenerationInfo const& stateGenerationInfo, size_t id); - + + /** + * Construct concerete state from pseudo state by using the underlying bitvector. + */ + void construct(); + + std::shared_ptr<DFTState<ValueType>> copy() const; + DFTElementState getElementState(size_t id) const; DFTDependencyState getDependencyState(size_t id) const; @@ -96,7 +116,11 @@ namespace storm { bool isInvalid() const { return !mValid; } - + + bool isPseudoState() const { + return mPseudoState; + } + storm::storage::BitVector const& status() const { return mStatus; } @@ -135,17 +159,49 @@ namespace storm { * @param spareId Id of the spare */ void finalizeUses(size_t spareId); - + + /** + * Claim a new spare child for the given spare gate. + * + * @param spareId Id of the spare gate. + * @param currentlyUses Id of the currently used spare child. + * @param children List of children of this spare. + * + * @return True, if claiming was successful. + */ bool claimNew(size_t spareId, size_t currentlyUses, std::vector<std::shared_ptr<DFTElement<ValueType>>> const& children); - bool hasOutgoingEdges() const { - return !mIsCurrentlyFailableBE.empty(); - } - + /** + * Get number of currently failable BEs. + * + * @return Number of failable BEs. + */ size_t nrFailableBEs() const { - return mIsCurrentlyFailableBE.size(); + return mCurrentlyFailableBE.size(); } + /** + * Get the failure rate of the currently failable BE on the given index. + * + * @param index Index of BE in list of currently failable BEs. + * + * @return Failure rate of the BE. + */ + ValueType getFailableBERate(size_t index) const; + + /** + * Get the current failure rate of the given BE. + * + * @param id Id of BE. + * + * @return Failure rate of the BE. + */ + ValueType getBERate(size_t id) const; + + /** Get number of currently failable dependencies. + * + * @return Number of failable dependencies. + */ size_t nrFailableDependencies() const { return mFailableDependencies.size(); } @@ -156,7 +212,7 @@ namespace storm { * @return Id of the dependency */ size_t getDependencyId(size_t index) const { - assert(index < nrFailableDependencies()); + STORM_LOG_ASSERT(index < nrFailableDependencies(), "Index invalid."); return mFailableDependencies[index]; } @@ -214,13 +270,13 @@ namespace storm { } stream << "} "; } else { - auto it = mIsCurrentlyFailableBE.begin(); + auto it = mCurrentlyFailableBE.begin(); stream << "{"; - if(it != mIsCurrentlyFailableBE.end()) { + if(it != mCurrentlyFailableBE.end()) { stream << *it; } ++it; - while(it != mIsCurrentlyFailableBE.end()) { + while(it != mCurrentlyFailableBE.end()) { stream << ", " << *it; ++it; } @@ -235,11 +291,6 @@ namespace storm { private: void propagateActivation(size_t representativeId); - - /*! - * Sort failable BEs in decreasing order of their active failure rate. - */ - void sortFailableBEs(); }; diff --git a/src/storage/dft/DFTStateSpaceGenerationQueues.h b/src/storage/dft/DFTStateSpaceGenerationQueues.h index 6ded7cb38..bd71b5fbf 100644 --- a/src/storage/dft/DFTStateSpaceGenerationQueues.h +++ b/src/storage/dft/DFTStateSpaceGenerationQueues.h @@ -55,7 +55,7 @@ namespace storm { } DFTRestrictionPointer nextRestrictionCheck() { - assert(!restrictionChecksDone()); + STORM_LOG_ASSERT(!restrictionChecksDone(), "All restriction checks done already."); DFTRestrictionPointer next = restrictionChecks.back(); restrictionChecks.pop_back(); return next; diff --git a/src/storage/expressions/Variable.cpp b/src/storage/expressions/Variable.cpp index aeb432f2a..5c9660b83 100644 --- a/src/storage/expressions/Variable.cpp +++ b/src/storage/expressions/Variable.cpp @@ -1,6 +1,5 @@ #include "src/storage/expressions/Variable.h" #include "src/storage/expressions/ExpressionManager.h" -#include <cassert> namespace storm { namespace expressions { @@ -41,7 +40,7 @@ namespace storm { } ExpressionManager const& Variable::getManager() const { - assert(manager != nullptr); + STORM_LOG_ASSERT(manager != nullptr, "Manager is null."); return *manager; } @@ -65,4 +64,4 @@ namespace storm { return this->getType().isNumericalType(); } } -} \ No newline at end of file +} diff --git a/src/storage/expressions/Variable.h b/src/storage/expressions/Variable.h index 3ade039ed..79b430f65 100644 --- a/src/storage/expressions/Variable.h +++ b/src/storage/expressions/Variable.h @@ -6,6 +6,7 @@ #include <functional> #include "src/utility/OsDetection.h" +#include "src/utility/macros.h" namespace storm { namespace expressions { @@ -155,4 +156,4 @@ namespace std { }; } -#endif /* STORM_STORAGE_EXPRESSIONS_VARIABLE_H_ */ \ No newline at end of file +#endif /* STORM_STORAGE_EXPRESSIONS_VARIABLE_H_ */ diff --git a/src/storage/prism/Assignment.cpp b/src/storage/prism/Assignment.cpp index 40b37d4a9..e0ccf7a92 100644 --- a/src/storage/prism/Assignment.cpp +++ b/src/storage/prism/Assignment.cpp @@ -1,5 +1,4 @@ #include "Assignment.h" -#include <cassert> namespace storm { namespace prism { @@ -25,7 +24,7 @@ namespace storm { bool Assignment::isIdentity() const { if(this->expression.isVariable()) { - assert(this->expression.getVariables().size() == 1); + STORM_LOG_ASSERT(this->expression.getVariables().size() == 1, "Invalid number of variables."); return variable == *(this->expression.getVariables().begin()); } return false; diff --git a/src/storage/prism/Command.cpp b/src/storage/prism/Command.cpp index 5cad634e7..a5a59522d 100644 --- a/src/storage/prism/Command.cpp +++ b/src/storage/prism/Command.cpp @@ -1,5 +1,4 @@ #include "Command.h" -#include <cassert> namespace storm { namespace prism { @@ -32,7 +31,7 @@ namespace storm { } storm::prism::Update const& Command::getUpdate(uint_fast64_t index) const { - assert(index < getNumberOfUpdates()); + STORM_LOG_ASSERT(index < getNumberOfUpdates(), "Invalid index."); return this->updates[index]; } diff --git a/src/storage/prism/Program.cpp b/src/storage/prism/Program.cpp index f187ecc08..fb624742b 100644 --- a/src/storage/prism/Program.cpp +++ b/src/storage/prism/Program.cpp @@ -1597,10 +1597,6 @@ namespace storm { return Command(newCommandIndex, false, actionIndex, actionName, newGuard, newUpdates, this->getFilename(), 0); } - uint_fast64_t Program::getNumberOfActions() const { - return this->actions.size(); - } - storm::jani::Model Program::toJani(bool allVariablesGlobal) const { ToJaniConverter converter; return converter.convert(*this, allVariablesGlobal); diff --git a/src/storage/prism/Program.h b/src/storage/prism/Program.h index 9e4617913..b1c2a4dc5 100644 --- a/src/storage/prism/Program.h +++ b/src/storage/prism/Program.h @@ -574,8 +574,6 @@ namespace storm { std::unordered_map<uint_fast64_t, std::string> buildActionIndexToActionNameMap() const; - uint_fast64_t getNumberOfActions() const; - /*! * Converts the PRISM model into an equivalent JANI model. */ diff --git a/src/storage/sparse/StateStorage.cpp b/src/storage/sparse/StateStorage.cpp index f86bbcd48..4ea26edb8 100644 --- a/src/storage/sparse/StateStorage.cpp +++ b/src/storage/sparse/StateStorage.cpp @@ -15,6 +15,7 @@ namespace storm { } template class StateStorage<uint32_t>; + template class StateStorage<uint_fast64_t>; } } } \ No newline at end of file diff --git a/src/storm-dyftee.cpp b/src/storm-dyftee.cpp index 1bad1c994..58f464f07 100644 --- a/src/storm-dyftee.cpp +++ b/src/storm-dyftee.cpp @@ -1,11 +1,14 @@ -#include "logic/Formula.h" -#include "utility/initialize.h" -#include "utility/storm.h" -#include "modelchecker/DFTAnalyser.h" +#include "src/logic/Formula.h" +#include "src/utility/initialize.h" +#include "src/utility/storm.h" +#include "src/parser/DFTGalileoParser.h" +#include "src/modelchecker/dft/DFTModelChecker.h" + +#include "src/modelchecker/dft/DFTASFChecker.h" #include "src/cli/cli.h" #include "src/exceptions/BaseException.h" #include "src/utility/macros.h" -#include <boost/lexical_cast.hpp> +#include "src/builder/DftSmtBuilder.h" #include "src/settings/modules/GeneralSettings.h" #include "src/settings/modules/DFTSettings.h" @@ -15,6 +18,7 @@ //#include "src/settings/modules/CuddSettings.h" //#include "src/settings/modules/SylvanSettings.h" #include "src/settings/modules/GmmxxEquationSolverSettings.h" +#include "src/settings/modules/MinMaxEquationSolverSettings.h" #include "src/settings/modules/NativeEquationSolverSettings.h" //#include "src/settings/modules/BisimulationSettings.h" //#include "src/settings/modules/GlpkSettings.h" @@ -23,6 +27,8 @@ //#include "src/settings/modules/ParametricSettings.h" #include "src/settings/modules/EliminationSettings.h" +#include <boost/lexical_cast.hpp> + /*! * Load DFT from filename, build corresponding Model and check against given property. * @@ -30,18 +36,36 @@ * @param property PCTC formula capturing the property to check. */ template <typename ValueType> -void analyzeDFT(std::string filename, std::string property, bool symred = false, bool allowModularisation = false, bool enableDC = true) { +void analyzeDFT(std::string filename, std::string property, bool symred, bool allowModularisation, bool enableDC, double approximationError) { std::cout << "Running DFT analysis on file " << filename << " with property " << property << std::endl; storm::parser::DFTGalileoParser<ValueType> parser; storm::storage::DFT<ValueType> dft = parser.parseDFT(filename); std::vector<std::shared_ptr<storm::logic::Formula const>> formulas = storm::parseFormulasForExplicit(property); STORM_LOG_ASSERT(formulas.size() == 1, "Wrong number of formulas."); + + storm::modelchecker::DFTModelChecker<ValueType> modelChecker; + modelChecker.check(dft, formulas[0], symred, allowModularisation, enableDC, approximationError); + modelChecker.printTimings(); + modelChecker.printResult(); +} + +/*! + * Analyze the DFT with use of SMT solving. + * + * @param filename Path to DFT file in Galileo format. + */ +template<typename ValueType> +void analyzeWithSMT(std::string filename) { + std::cout << "Running DFT analysis on file " << filename << " with use of SMT" << std::endl; - DFTAnalyser<ValueType> analyser; - analyser.check(dft, formulas[0], symred, allowModularisation, enableDC); - analyser.printTimings(); - analyser.printResult(); + storm::parser::DFTGalileoParser<ValueType> parser; + storm::storage::DFT<ValueType> dft = parser.parseDFT(filename); + storm::modelchecker::DFTASFChecker asfChecker(dft); + asfChecker.convert(); + asfChecker.toFile("test.smt2"); + //bool sat = dftSmtBuilder.check(); + //std::cout << "SMT result: " << sat << std::endl; } /*! @@ -59,6 +83,7 @@ void initializeSettings() { //storm::settings::addModule<storm::settings::modules::CuddSettings>(); //storm::settings::addModule<storm::settings::modules::SylvanSettings>(); storm::settings::addModule<storm::settings::modules::GmmxxEquationSolverSettings>(); + storm::settings::addModule<storm::settings::modules::MinMaxEquationSolverSettings>(); storm::settings::addModule<storm::settings::modules::NativeEquationSolverSettings>(); //storm::settings::addModule<storm::settings::modules::BisimulationSettings>(); //storm::settings::addModule<storm::settings::modules::GlpkSettings>(); @@ -92,6 +117,24 @@ int main(const int argc, const char** argv) { STORM_LOG_THROW(false, storm::exceptions::InvalidSettingsException, "No input model."); } + bool parametric = false; +#ifdef STORM_HAVE_CARL + parametric = generalSettings.isParametricSet(); +#endif + +#ifdef STORM_HAVE_Z3 + if (dftSettings.solveWithSMT()) { + // Solve with SMT + if (parametric) { + // analyzeWithSMT<storm::RationalFunction>(dftSettings.getDftFilename()); + } else { + analyzeWithSMT<double>(dftSettings.getDftFilename()); + } + storm::utility::cleanUp(); + return 0; + } +#endif + // Set min or max bool minimal = true; if (dftSettings.isComputeMaximalValue()) { @@ -101,7 +144,6 @@ int main(const int argc, const char** argv) { // Construct pctlFormula std::string pctlFormula = ""; - bool allowModular = true; std::string operatorType = ""; std::string targetFormula = ""; @@ -112,7 +154,6 @@ int main(const int argc, const char** argv) { STORM_LOG_THROW(!dftSettings.usePropProbability() && !dftSettings.usePropTimebound(), storm::exceptions::InvalidSettingsException, "More than one property given."); operatorType = "T"; targetFormula = "F \"failed\""; - allowModular = false; } else if (dftSettings.usePropProbability()) { STORM_LOG_THROW(!dftSettings.usePropTimebound(), storm::exceptions::InvalidSettingsException, "More than one property given."); operatorType = "P";; @@ -132,20 +173,20 @@ int main(const int argc, const char** argv) { STORM_LOG_ASSERT(!pctlFormula.empty(), "Pctl formula empty."); - bool parametric = false; -#ifdef STORM_HAVE_CARL - parametric = generalSettings.isParametricSet(); -#endif - + double approximationError = 0.0; + if (dftSettings.isApproximationErrorSet()) { + approximationError = dftSettings.getApproximationError(); + } + // From this point on we are ready to carry out the actual computations. if (parametric) { #ifdef STORM_HAVE_CARL - analyzeDFT<storm::RationalFunction>(dftSettings.getDftFilename(), pctlFormula, dftSettings.useSymmetryReduction(), allowModular && dftSettings.useModularisation(), !dftSettings.isDisableDC() ); + analyzeDFT<storm::RationalFunction>(dftSettings.getDftFilename(), pctlFormula, dftSettings.useSymmetryReduction(), dftSettings.useModularisation(), !dftSettings.isDisableDC(), approximationError); #else STORM_LOG_THROW(false, storm::exceptions::NotSupportedException, "Parameters are not supported in this build."); #endif } else { - analyzeDFT<double>(dftSettings.getDftFilename(), pctlFormula, dftSettings.useSymmetryReduction(), allowModular && dftSettings.useModularisation(), !dftSettings.isDisableDC()); + analyzeDFT<double>(dftSettings.getDftFilename(), pctlFormula, dftSettings.useSymmetryReduction(), dftSettings.useModularisation(), !dftSettings.isDisableDC(), approximationError); } // All operations have now been performed, so we clean up everything and terminate. diff --git a/src/utility/bitoperations.h b/src/utility/bitoperations.h index bec6a6414..18b42a1cd 100644 --- a/src/utility/bitoperations.h +++ b/src/utility/bitoperations.h @@ -10,6 +10,6 @@ inline size_t smallestIntWithNBitsSet(size_t n) { inline size_t nextBitPermutation(size_t v) { if (v==0) return static_cast<size_t>(0); // From https://graphics.stanford.edu/~seander/bithacks.html#NextBitPermutation - unsigned int t = (v | (v - 1)) + 1; + size_t t = (v | (v - 1)) + 1; return t | ((((t & -t) / (v & -v)) >> 1) - 1); } diff --git a/src/utility/constants.cpp b/src/utility/constants.cpp index 2739bfa1f..99a4d527b 100644 --- a/src/utility/constants.cpp +++ b/src/utility/constants.cpp @@ -42,13 +42,18 @@ namespace storm { bool isConstant(ValueType const& a) { return true; } - + template<typename ValueType> bool isInteger(ValueType const& number) { ValueType iPart; ValueType result = std::modf(number, &iPart); return result = zero<ValueType>(); } + + template<typename ValueType> + bool isInfinity(ValueType const& a) { + return a == infinity<ValueType>(); + } template<> bool isInteger(int const& number) { @@ -129,7 +134,13 @@ namespace storm { // FIXME: this should be treated more properly. return storm::RationalFunction(-1.0); } - + + template<> + bool isInfinity(storm::RationalFunction const& a) { + // FIXME: this should be treated more properly. + return a == infinity<storm::RationalFunction>(); + } + template<> storm::RationalNumber infinity() { // FIXME: this should be treated more properly. @@ -347,10 +358,11 @@ namespace storm { template bool isOne(double const& value); template bool isZero(double const& value); template bool isConstant(double const& value); - - template double one(); - template double zero(); - template double infinity(); + template bool isInfinity(double const& value); + + template double one(); + template double zero(); + template double infinity(); template double pow(double const& value, uint_fast64_t exponent); @@ -367,10 +379,11 @@ namespace storm { template bool isOne(float const& value); template bool isZero(float const& value); template bool isConstant(float const& value); - - template float one(); - template float zero(); - template float infinity(); + template bool isInfinity(float const& value); + + template float one(); + template float zero(); + template float infinity(); template float pow(float const& value, uint_fast64_t exponent); template bool isInteger(float const& number); @@ -389,6 +402,7 @@ namespace storm { template bool isOne(int const& value); template bool isZero(int const& value); template bool isConstant(int const& value); + template bool isInfinity(int const& value); template int one(); template int zero(); @@ -407,6 +421,7 @@ namespace storm { template bool isOne(storm::storage::sparse::state_type const& value); template bool isZero(storm::storage::sparse::state_type const& value); template bool isConstant(storm::storage::sparse::state_type const& value); + template bool isInfinity(storm::storage::sparse::state_type const& value); template uint32_t one(); template uint32_t zero(); @@ -471,6 +486,7 @@ namespace storm { template bool isOne(RationalFunction const& value); template bool isZero(RationalFunction const& value); template bool isConstant(RationalFunction const& value); + template bool isInfinity(RationalFunction const& value); template RationalFunction one(); template RationalFunction zero(); @@ -487,6 +503,7 @@ namespace storm { template bool isOne(Interval const& value); template bool isZero(Interval const& value); template bool isConstant(Interval const& value); + template bool isInfinity(Interval const& value); template Interval one(); template Interval zero(); diff --git a/src/utility/constants.h b/src/utility/constants.h index c384a4312..75cac92e2 100644 --- a/src/utility/constants.h +++ b/src/utility/constants.h @@ -43,6 +43,9 @@ namespace storm { template<typename ValueType> bool isConstant(ValueType const& a); + template<typename ValueType> + bool isInfinity(ValueType const& a); + template<typename ValueType> ValueType pow(ValueType const& value, uint_fast64_t exponent); @@ -82,4 +85,4 @@ namespace storm { } } -#endif /* STORM_UTILITY_CONSTANTS_H_ */ \ No newline at end of file +#endif /* STORM_UTILITY_CONSTANTS_H_ */ diff --git a/src/utility/logging.h b/src/utility/logging.h index 2fc58a008..951935e46 100644 --- a/src/utility/logging.h +++ b/src/utility/logging.h @@ -1,6 +1,16 @@ #ifndef STORM_UTILITY_LOGGING_H_ #define STORM_UTILITY_LOGGING_H_ +// Include config to know whether CARL is available or not. +#include "storm-config.h" +#ifdef STORM_HAVE_CARL +// Load streaming operator from CARL +#include <carl/io/streamingOperators.h> +namespace l3pp { + using carl::operator<<; +} +#endif + #include <l3pp.h> #if !defined(STORM_LOG_DISABLE_DEBUG) && !defined(STORM_LOG_DISABLE_TRACE) diff --git a/src/utility/vector.cpp b/src/utility/vector.cpp deleted file mode 100644 index 20af3c628..000000000 --- a/src/utility/vector.cpp +++ /dev/null @@ -1,17 +0,0 @@ -#include "src/utility/vector.h" - -//template<typename ValueType> -//std::ostream& operator<<(std::ostream& out, std::vector<ValueType> const& vector) { -std::ostream& operator<<(std::ostream& out, std::vector<double> const& vector) { - out << "vector (" << vector.size() << ") [ "; - for (uint_fast64_t i = 0; i < vector.size() - 1; ++i) { - out << vector[i] << ", "; - } - out << vector.back(); - out << " ]"; - return out; -} - -// Explicitly instantiate functions. -//template std::ostream& operator<<(std::ostream& out, std::vector<double> const& vector); -//template std::ostream& operator<<(std::ostream& out, std::vector<uint_fast64_t> const& vector); \ No newline at end of file diff --git a/src/utility/vector.h b/src/utility/vector.h index c973b4250..c7759cfa3 100644 --- a/src/utility/vector.h +++ b/src/utility/vector.h @@ -15,11 +15,6 @@ #include "src/utility/macros.h" #include "src/solver/OptimizationDirection.h" -// Template was causing problems as Carl has the same function -//template<typename ValueType> -//std::ostream& operator<<(std::ostream& out, std::vector<ValueType> const& vector); -std::ostream& operator<<(std::ostream& out, std::vector<double> const& vector); - namespace storm { namespace utility { namespace vector {