sp
/
tempest
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Only keep track of results from the last iteration (instead of all iterations) for 2 of the 3 vectors

tempestpy_adaptions
Matthias Volk 6 years ago
parent
8ae800b130
commit
f9fb90499d
1 changed files with 48 additions and 29 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Unified View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 75
      src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp

75
src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp
View File

@ -37,28 +37,47 @@ namespace storm {
template<typename ValueType> template<typename ValueType>
struct UnifPlusVectors { struct UnifPlusVectors {
std::vector<std::vector<ValueType>> resLower;
std::vector<std::vector<ValueType>> resUpper;
UnifPlusVectors(uint64_t steps, uint64_t noStates) : numberOfStates(noStates), resLowerOld(numberOfStates, -1), resLowerNew(numberOfStates, -1), resUpperOld(numberOfStates, -1), resUpperNew(numberOfStates, -1) {
// Intentionally empty.
wUpper = std::vector<std::vector<ValueType>>(steps, std::vector<ValueType>(numberOfStates, -1));
}
/**
* Prepare new iteration by setting the new result vectors as old result vectors, and initializing the new result vectors with -1 again.
*/
void prepareNewIteration() {
resLowerOld.swap(resLowerNew);
std::fill(resLowerNew.begin(), resLowerNew.end(), -1);
resUpperOld.swap(resUpperNew);
std::fill(resUpperNew.begin(), resUpperNew.end(), -1);
}
uint64_t numberOfStates;
std::vector<ValueType> resLowerOld;
std::vector<ValueType> resLowerNew;
std::vector<ValueType> resUpperOld;
std::vector<ValueType> resUpperNew;
std::vector<std::vector<ValueType>> wUpper; std::vector<std::vector<ValueType>> wUpper;
}; };
template<typename ValueType> template<typename ValueType>
void calculateUnifPlusVector(Environment const& env, uint64_t k, uint64_t state, bool calcLower, ValueType lambda, uint64_t numberOfProbabilisticChoices, std::vector<std::vector<ValueType>> const & relativeReachability, OptimizationDirection dir, UnifPlusVectors<ValueType>& unifVectors, storm::storage::SparseMatrix<ValueType> const& fullTransitionMatrix, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::unique_ptr<storm::solver::MinMaxLinearEquationSolver<ValueType>> const& solver, storm::utility::numerical::FoxGlynnResult<ValueType> const& poisson, bool cycleFree) { void calculateUnifPlusVector(Environment const& env, uint64_t k, uint64_t state, bool calcLower, ValueType lambda, uint64_t numberOfProbabilisticChoices, std::vector<std::vector<ValueType>> const & relativeReachability, OptimizationDirection dir, UnifPlusVectors<ValueType>& unifVectors, storm::storage::SparseMatrix<ValueType> const& fullTransitionMatrix, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::unique_ptr<storm::solver::MinMaxLinearEquationSolver<ValueType>> const& solver, storm::utility::numerical::FoxGlynnResult<ValueType> const& poisson, bool cycleFree) {
std::vector<std::vector<ValueType>>& resVector = calcLower ? unifVectors.resLower : unifVectors.wUpper;
std::vector<ValueType>& resVectorOld = calcLower ? unifVectors.resLowerOld : unifVectors.wUpper[k+1];
std::vector<ValueType>& resVectorNew = calcLower ? unifVectors.resLowerNew : unifVectors.wUpper[k];
if (resVector[k][state] != -1) {
if (resVectorNew[state] != -1) {
// Result already calculated. // Result already calculated.
return; return;
} }
auto numberOfStates = fullTransitionMatrix.getRowGroupCount(); auto numberOfStates = fullTransitionMatrix.getRowGroupCount();
uint64_t N = resVector.size() - 1;
uint64_t N = unifVectors.wUpper.size() - 1;
auto const& rowGroupIndices = fullTransitionMatrix.getRowGroupIndices(); auto const& rowGroupIndices = fullTransitionMatrix.getRowGroupIndices();
ValueType res; ValueType res;
// First case, k==N, independent from kind of state. // First case, k==N, independent from kind of state.
if (k == N) { if (k == N) {
resVector[k][state] = storm::utility::zero<ValueType>();
resVectorNew[state] = storm::utility::zero<ValueType>();
return; return;
} }
@ -72,10 +91,10 @@ namespace storm {
res += poisson.weights[i - poisson.left]; res += poisson.weights[i - poisson.left];
} }
} }
resVector[k][state] = res;
resVectorNew[state] = res;
} else { } else {
// w upper // w upper
resVector[k][state] = storm::utility::one<ValueType>();
resVectorNew[state] = storm::utility::one<ValueType>();
} }
return; return;
} }
@ -86,12 +105,13 @@ namespace storm {
res = storm::utility::zero<ValueType>(); res = storm::utility::zero<ValueType>();
for (auto const& element : fullTransitionMatrix.getRow(rowGroupIndices[state])) { for (auto const& element : fullTransitionMatrix.getRow(rowGroupIndices[state])) {
uint64_t successor = element.getColumn(); uint64_t successor = element.getColumn();
if (resVector[k+1][successor] == -1) {
if (resVectorOld[successor] == -1) {
STORM_LOG_ASSERT(false, "Need to calculate previous result.");
calculateUnifPlusVector(env, k+1, successor, calcLower, lambda, numberOfProbabilisticChoices, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree); calculateUnifPlusVector(env, k+1, successor, calcLower, lambda, numberOfProbabilisticChoices, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree);
} }
res += element.getValue() * resVector[k+1][successor];
res += element.getValue() * resVectorOld[successor];
} }
resVector[k][state]=res;
resVectorNew[state]=res;
return; return;
} }
@ -110,10 +130,10 @@ namespace storm {
continue; continue;
} }
if (resVector[k][successor] == -1) {
if (resVectorNew[successor] == -1) {
calculateUnifPlusVector(env, k, successor, calcLower, lambda, numberOfProbabilisticChoices, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree); calculateUnifPlusVector(env, k, successor, calcLower, lambda, numberOfProbabilisticChoices, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree);
} }
between += element.getValue() * resVector[k][successor];
between += element.getValue() * resVectorNew[successor];
} }
if (maximize(dir)) { if (maximize(dir)) {
res = storm::utility::max(res, between); res = storm::utility::max(res, between);
@ -125,7 +145,7 @@ namespace storm {
} }
} }
} }
resVector[k][state] = res;
resVectorNew[state] = res;
return; return;
} }
@ -150,10 +170,10 @@ namespace storm {
continue; continue;
} }
if (resVector[k][successor] == -1) {
if (resVectorNew[successor] == -1) {
calculateUnifPlusVector(env, k, successor, calcLower, lambda, numberOfProbabilisticStates, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree); calculateUnifPlusVector(env, k, successor, calcLower, lambda, numberOfProbabilisticStates, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree);
} }
res += relativeReachability[j][stateCount] * resVector[k][successor];
res += relativeReachability[j][stateCount] * resVectorNew[successor];
++stateCount; ++stateCount;
} }
@ -166,28 +186,30 @@ namespace storm {
solver->solveEquations(env, dir, x, b); solver->solveEquations(env, dir, x, b);
// Expand the solution for the probabilistic states to all states. // Expand the solution for the probabilistic states to all states.
storm::utility::vector::setVectorValues(resVector[k], ~markovianStates, x);
storm::utility::vector::setVectorValues(resVectorNew, ~markovianStates, x);
} }
template <typename ValueType> template <typename ValueType>
void calculateResUpper(Environment const& env, std::vector<std::vector<ValueType>> const& relativeReachability, OptimizationDirection dir, uint64_t k, uint64_t state, ValueType lambda, uint64_t numberOfProbabilisticStates, UnifPlusVectors<ValueType>& unifVectors, storm::storage::SparseMatrix<ValueType> const& fullTransitionMatrix, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::unique_ptr<storm::solver::MinMaxLinearEquationSolver<ValueType>> const& solver, storm::utility::numerical::FoxGlynnResult<ValueType> const & poisson, bool cycleFree) { void calculateResUpper(Environment const& env, std::vector<std::vector<ValueType>> const& relativeReachability, OptimizationDirection dir, uint64_t k, uint64_t state, ValueType lambda, uint64_t numberOfProbabilisticStates, UnifPlusVectors<ValueType>& unifVectors, storm::storage::SparseMatrix<ValueType> const& fullTransitionMatrix, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::unique_ptr<storm::solver::MinMaxLinearEquationSolver<ValueType>> const& solver, storm::utility::numerical::FoxGlynnResult<ValueType> const & poisson, bool cycleFree) {
// Avoiding multiple computation of the same value. // Avoiding multiple computation of the same value.
if (unifVectors.resUpper[k][state] != -1) {
if (unifVectors.resUpperNew[state] != -1) {
STORM_LOG_ASSERT(false, "Result was already calculated.");
return; return;
} }
uint64_t N = unifVectors.resUpper.size() - 1;
uint64_t N = unifVectors.wUpper.size() - 1;
ValueType res = storm::utility::zero<ValueType>(); ValueType res = storm::utility::zero<ValueType>();
for (uint64_t i = k; i < N; ++i) { for (uint64_t i = k; i < N; ++i) {
if (unifVectors.wUpper[N-1-(i-k)][state] == -1) { if (unifVectors.wUpper[N-1-(i-k)][state] == -1) {
STORM_LOG_ASSERT(false, "Need to calculate previous result.");
calculateUnifPlusVector(env, N-1-(i-k), state, false, lambda, numberOfProbabilisticStates, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree); calculateUnifPlusVector(env, N-1-(i-k), state, false, lambda, numberOfProbabilisticStates, relativeReachability, dir, unifVectors, fullTransitionMatrix, markovianStates, psiStates, solver, poisson, cycleFree);
} }
if (i >= poisson.left && i <= poisson.right) { if (i >= poisson.left && i <= poisson.right) {
res += poisson.weights[i - poisson.left] * unifVectors.wUpper[N-1-(i-k)][state]; res += poisson.weights[i - poisson.left] * unifVectors.wUpper[N-1-(i-k)][state];
} }
} }
unifVectors.resUpper[k][state] = res;
unifVectors.resUpperNew[state] = res;
} }
template <typename ValueType> template <typename ValueType>
@ -237,7 +259,7 @@ namespace storm {
bool cycleFree = sccDecomposition.empty(); bool cycleFree = sccDecomposition.empty();
// Vectors to store computed vectors. // Vectors to store computed vectors.
UnifPlusVectors<ValueType> unifVectors;
UnifPlusVectors<ValueType> unifVectors(0, 0);
// Transitions from goal states will be ignored. However, we mark them as non-probabilistic to make sure // Transitions from goal states will be ignored. However, we mark them as non-probabilistic to make sure
// we do not apply the MDP algorithm to them. // we do not apply the MDP algorithm to them.
@ -313,7 +335,6 @@ namespace storm {
} }
} }
std::vector<ValueType> init(numberOfStates, -1);
ValueType maxNorm = storm::utility::zero<ValueType>(); ValueType maxNorm = storm::utility::zero<ValueType>();
// Maximal step size // Maximal step size
uint64_t N; uint64_t N;
@ -362,10 +383,7 @@ namespace storm {
} }
// (4) Define vectors/matrices. // (4) Define vectors/matrices.
std::vector<std::vector<ValueType>> v = std::vector<std::vector<ValueType>>(N + 1, init);
unifVectors.resLower = v;
unifVectors.resUpper = v;
unifVectors.wUpper = v;
unifVectors = UnifPlusVectors<ValueType>(N+1, numberOfStates);
// (5) Compute vectors and maxNorm. // (5) Compute vectors and maxNorm.
for (int64_t k = N; k >= 0; --k) { for (int64_t k = N; k >= 0; --k) {
@ -374,12 +392,13 @@ namespace storm {
calculateUnifPlusVector(env, k, i, false, lambda, numberOfProbabilisticChoices, relativeReachabilities, dir, unifVectors, fullTransitionMatrix, markovianAndGoalStates, psiStates, solver, foxGlynnResult, cycleFree); calculateUnifPlusVector(env, k, i, false, lambda, numberOfProbabilisticChoices, relativeReachabilities, dir, unifVectors, fullTransitionMatrix, markovianAndGoalStates, psiStates, solver, foxGlynnResult, cycleFree);
calculateResUpper(env, relativeReachabilities, dir, k, i, lambda, numberOfProbabilisticChoices, unifVectors, fullTransitionMatrix, markovianAndGoalStates, psiStates, solver, foxGlynnResult, cycleFree); calculateResUpper(env, relativeReachabilities, dir, k, i, lambda, numberOfProbabilisticChoices, unifVectors, fullTransitionMatrix, markovianAndGoalStates, psiStates, solver, foxGlynnResult, cycleFree);
} }
unifVectors.prepareNewIteration();
} }
// Only iterate over result vector, as the results can only get more precise. // Only iterate over result vector, as the results can only get more precise.
maxNorm = storm::utility::zero<ValueType>(); maxNorm = storm::utility::zero<ValueType>();
for (uint64_t i = 0; i < numberOfStates; i++){ for (uint64_t i = 0; i < numberOfStates; i++){
ValueType diff = storm::utility::abs(unifVectors.resUpper[0][i] - unifVectors.resLower[0][i]);
ValueType diff = storm::utility::abs(unifVectors.resUpperOld[i] - unifVectors.resLowerOld[i]);
maxNorm = std::max(maxNorm, diff); maxNorm = std::max(maxNorm, diff);
} }
@ -389,7 +408,7 @@ namespace storm {
} while (maxNorm > epsilon * (1 - kappa)); } while (maxNorm > epsilon * (1 - kappa));
return unifVectors.resLower[0];
return unifVectors.resLowerOld;
} }
template <typename ValueType> template <typename ValueType>

Write Preview
Loading…
Cancel
Save