@ -84,7 +84,7 @@ namespace storm {
if ( comparator . isZero ( lowerBound ) & & comparator . isInfinity ( upperBound ) ) {
return this - > computeUntilProbabilitiesHelper ( this - > getModel ( ) . getTransitionMatrix ( ) , this - > getModel ( ) . getBackwardTransitions ( ) , phiStates , psiStates , qualitative , * this - > linearEquationSolver ) ;
}
// From this point on, we know that we have to solve a more complicated problem [t, t'] with either t != 0
// or t' != inf.
@ -107,7 +107,7 @@ namespace storm {
if ( comparator . isZero ( lowerBound ) ) {
// In this case, the interval is of the form [0, t].
// Note that this excludes [0, inf] since this is untimed reachability and we considered this case earlier.
// Find the maximal rate of all 'maybe' states to take it as the uniformization rate.
ValueType uniformizationRate = 0 ;
for ( auto const & state : statesWithProbabilityGreater0NonPsi ) {
@ -118,12 +118,12 @@ namespace storm {
// Compute the uniformized matrix.
storm : : storage : : SparseMatrix < ValueType > uniformizedMatrix = this - > computeUniformizedMatrix ( this - > getModel ( ) . getTransitionMatrix ( ) , statesWithProbabilityGreater0 , psiStates , uniformizationRate , exitRates ) ;
// Finally compute the transient probabilities.
std : : vector < ValueType > psiStateValues ( statesWithProbabilityGreater0 . getNumberOfSetBits ( ) , storm : : utility : : zero < ValueType > ( ) ) ;
storm : : utility : : vector : : setVectorValues ( psiStateValues , psiStates % statesWithProbabilityGreater0 , storm : : utility : : one < ValueType > ( ) ) ;
std : : vector < ValueType > subresult = this - > computeTransientProbabilities ( uniformizedMatrix , uniformizationRate * upperBound , psiStateValues , * this - > linearEquationSolver ) ;
std : : vector < ValueType > subresult = this - > computeTransientProbabilities ( uniformizedMatrix , upperBound , uniformizationRate , psiStateValues , * this - > linearEquationSolver ) ;
storm : : utility : : vector : : setVectorValues ( result , statesWithProbabilityGreater0 , subresult ) ;
} else if ( comparator . isInfinity ( upperBound ) ) {
// In this case, the interval is of the form [t, inf] with t != 0.
@ -149,7 +149,7 @@ namespace storm {
storm : : storage : : SparseMatrix < ValueType > uniformizedMatrix = this - > computeUniformizedMatrix ( this - > getModel ( ) . getTransitionMatrix ( ) , storm : : storage : : BitVector ( this - > getModel ( ) . getNumberOfStates ( ) , true ) , absorbingStates , uniformizationRate , exitRates ) ;
// Finally compute the transient probabilities.
result = this - > computeTransientProbabilities ( uniformizedMatrix , uniformizationRate * lowerBound , result , * this - > linearEquationSolver ) ;
result = this - > computeTransientProbabilities ( uniformizedMatrix , lowerBound , uniformizationRate , result , * this - > linearEquationSolver ) ;
} else {
// In this case, the interval is of the form [t, t'] with t != 0 and t' != inf.
@ -167,7 +167,7 @@ namespace storm {
// Start by computing the transient probabilities of reaching a psi state in time t' - t.
std : : vector < ValueType > psiStateValues ( statesWithProbabilityGreater0 . getNumberOfSetBits ( ) , storm : : utility : : zero < ValueType > ( ) ) ;
storm : : utility : : vector : : setVectorValues ( psiStateValues , psiStates % statesWithProbabilityGreater0 , storm : : utility : : one < ValueType > ( ) ) ;
std : : vector < ValueType > subresult = this - > computeTransientProbabilities ( uniformizedMatrix , uniformizationRate * ( upperBound - lowerBound ) , psiStateValues , * this - > linearEquationSolver ) ;
std : : vector < ValueType > subresult = this - > computeTransientProbabilities ( uniformizedMatrix , upperBound - lowerBound , uniformizationRate , psiStateValues , * this - > linearEquationSolver ) ;
storm : : utility : : vector : : setVectorValues ( result , statesWithProbabilityGreater0 , subresult ) ;
// Then compute the transient probabilities of being in such a state after t time units. For this,
@ -178,7 +178,7 @@ namespace storm {
uniformizationRate = std : : max ( uniformizationRate , exitRates [ state ] ) ;
}
uniformizationRate * = 1.02 ;
// Set the result to zero for all states that are known to violate phi.
storm : : utility : : vector : : setVectorValues ( result , absorbingStates , storm : : utility : : zero < ValueType > ( ) ) ;
@ -186,7 +186,7 @@ namespace storm {
// Finally, we compute the second set of transient probabilities.
uniformizedMatrix = this - > computeUniformizedMatrix ( this - > getModel ( ) . getTransitionMatrix ( ) , storm : : storage : : BitVector ( this - > getModel ( ) . getNumberOfStates ( ) , true ) , absorbingStates , uniformizationRate , exitRates ) ;
result = this - > computeTransientProbabilities ( uniformizedMatrix , uniformizationRate * lowerBound , result , * this - > linearEquationSolver ) ;
result = this - > computeTransientProbabilities ( uniformizedMatrix , lowerBound , uniformizationRate , result , * this - > linearEquationSolver ) ;
}
}
}
@ -227,7 +227,10 @@ namespace storm {
}
template < class ValueType >
std : : vector < ValueType > SparseCtmcCslModelChecker < ValueType > : : computeTransientProbabilities ( storm : : storage : : SparseMatrix < ValueType > const & uniformizedMatrix , ValueType const & lambda , std : : vector < ValueType > values , storm : : solver : : LinearEquationSolver < ValueType > const & linearEquationSolver ) const {
template < bool computeCumulativeReward >
std : : vector < ValueType > SparseCtmcCslModelChecker < ValueType > : : computeTransientProbabilities ( storm : : storage : : SparseMatrix < ValueType > const & uniformizedMatrix , ValueType timeBound , ValueType uniformizationRate , std : : vector < ValueType > values , storm : : solver : : LinearEquationSolver < ValueType > const & linearEquationSolver ) const {
ValueType lambda = timeBound * uniformizationRate ;
// If no time can pass, the current values are the result.
if ( lambda = = storm : : utility : : zero < ValueType > ( ) ) {
return values ;
@ -241,6 +244,16 @@ namespace storm {
element / = std : : get < 2 > ( foxGlynnResult ) ;
}
// If the cumulative reward is to be computed, we need to adjust the weights.
if ( computeCumulativeReward ) {
ValueType sum = storm : : utility : : zero < ValueType > ( ) ;
for ( auto & element : std : : get < 3 > ( foxGlynnResult ) ) {
sum + = element ;
element = ( 1 - sum ) / uniformizationRate ;
}
}
// Initialize result.
std : : vector < ValueType > result ;
uint_fast64_t startingIteration = std : : get < 0 > ( foxGlynnResult ) ;
@ -249,13 +262,26 @@ namespace storm {
storm : : utility : : vector : : scaleVectorInPlace ( result , std : : get < 3 > ( foxGlynnResult ) [ 0 ] ) ;
+ + startingIteration ;
} else {
result = std : : vector < ValueType > ( values . size ( ) ) ;
if ( computeCumulativeReward ) {
result = std : : vector < ValueType > ( values . size ( ) ) ;
storm : : utility : : vector : : applyPointwiseInPlace ( result , [ & uniformizationRate ] ( ValueType const & a ) { return a / uniformizationRate ; } ) ;
} else {
result = std : : vector < ValueType > ( values . size ( ) ) ;
}
}
std : : vector < ValueType > multiplicationResult ( result . size ( ) ) ;
// Perform the matrix-vector multiplications (without adding).
if ( std : : get < 0 > ( foxGlynnResult ) > 1 ) {
if ( ! computeCumulativeReward & & std : : get < 0 > ( foxGlynnResult ) > 1 ) {
// Perform the matrix-vector multiplications (without adding).
linearEquationSolver . performMatrixVectorMultiplication ( uniformizedMatrix , values , nullptr , std : : get < 0 > ( foxGlynnResult ) - 1 , & multiplicationResult ) ;
} else if ( computeCumulativeReward ) {
std : : function < ValueType ( ValueType const & , ValueType const & ) > addAndScale = [ & uniformizationRate ] ( ValueType const & a , ValueType const & b ) { return a + b / uniformizationRate ; } ;
// For the iterations below the left truncation point, we need to add and scale the result with the uniformization rate.
for ( uint_fast64_t index = 1 ; index < startingIteration ; + + index ) {
linearEquationSolver . performMatrixVectorMultiplication ( uniformizedMatrix , values , nullptr , 1 , & multiplicationResult ) ;
storm : : utility : : vector : : applyPointwiseInPlace ( result , values , addAndScale ) ;
}
}
// For the indices that fall in between the truncation points, we need to perform the matrix-vector
@ -277,6 +303,60 @@ namespace storm {
return SparseDtmcPrctlModelChecker < ValueType > : : computeUntilProbabilitiesHelper ( transitionMatrix , backwardTransitions , phiStates , psiStates , qualitative , linearEquationSolver ) ;
}
template < class ValueType >
std : : unique_ptr < CheckResult > SparseCtmcCslModelChecker < ValueType > : : computeInstantaneousRewards ( storm : : logic : : InstantaneousRewardFormula const & rewardPathFormula , bool qualitative , boost : : optional < storm : : logic : : OptimalityType > const & optimalityType ) {
return std : : unique_ptr < CheckResult > ( new ExplicitQuantitativeCheckResult < ValueType > ( this - > computeInstantaneousRewardsHelper ( rewardPathFormula . getContinuousTimeBound ( ) ) ) ) ;
}
template < class ValueType >
std : : vector < ValueType > SparseCtmcCslModelChecker < ValueType > : : computeInstantaneousRewardsHelper ( double timeBound ) const {
// Only compute the result if the model has a state-based reward this->getModel().
STORM_LOG_THROW ( this - > getModel ( ) . hasStateRewards ( ) , storm : : exceptions : : InvalidPropertyException , " Missing reward model for formula. Skipping formula. " ) ;
// Initialize result to state rewards of the this->getModel().
std : : vector < ValueType > result ( this - > getModel ( ) . getStateRewardVector ( ) ) ;
// If the time-bound is not zero, we need to perform a transient analysis.
if ( timeBound > 0 ) {
ValueType uniformizationRate = 0 ;
for ( auto const & rate : this - > getModel ( ) . getExitRateVector ( ) ) {
uniformizationRate = std : : max ( uniformizationRate , rate ) ;
}
storm : : storage : : SparseMatrix < ValueType > uniformizedMatrix = this - > computeUniformizedMatrix ( this - > getModel ( ) . getTransitionMatrix ( ) , storm : : storage : : BitVector ( this - > getModel ( ) . getNumberOfStates ( ) , true ) , storm : : storage : : BitVector ( this - > getModel ( ) . getNumberOfStates ( ) ) , uniformizationRate , this - > getModel ( ) . getExitRateVector ( ) ) ;
result = this - > computeTransientProbabilities ( uniformizedMatrix , timeBound , uniformizationRate , result , * this - > linearEquationSolver ) ;
}
return result ;
}
template < class ValueType >
std : : unique_ptr < CheckResult > SparseCtmcCslModelChecker < ValueType > : : computeCumulativeRewards ( storm : : logic : : CumulativeRewardFormula const & rewardPathFormula , bool qualitative , boost : : optional < storm : : logic : : OptimalityType > const & optimalityType ) {
return std : : unique_ptr < CheckResult > ( new ExplicitQuantitativeCheckResult < ValueType > ( this - > computeCumulativeRewardsHelper ( rewardPathFormula . getContinuousTimeBound ( ) ) ) ) ;
}
template < class ValueType >
std : : vector < ValueType > SparseCtmcCslModelChecker < ValueType > : : computeCumulativeRewardsHelper ( double timeBound ) const {
// Only compute the result if the model has a state-based reward this->getModel().
STORM_LOG_THROW ( this - > getModel ( ) . hasStateRewards ( ) , storm : : exceptions : : InvalidPropertyException , " Missing reward model for formula. Skipping formula. " ) ;
// Initialize result to state rewards of the this->getModel().
std : : vector < ValueType > result ( this - > getModel ( ) . getStateRewardVector ( ) ) ;
// If the time-bound is not zero, we need to perform a transient analysis.
if ( timeBound > 0 ) {
ValueType uniformizationRate = 0 ;
for ( auto const & rate : this - > getModel ( ) . getExitRateVector ( ) ) {
uniformizationRate = std : : max ( uniformizationRate , rate ) ;
}
storm : : storage : : SparseMatrix < ValueType > uniformizedMatrix = this - > computeUniformizedMatrix ( this - > getModel ( ) . getTransitionMatrix ( ) , storm : : storage : : BitVector ( this - > getModel ( ) . getNumberOfStates ( ) , true ) , storm : : storage : : BitVector ( this - > getModel ( ) . getNumberOfStates ( ) ) , uniformizationRate , this - > getModel ( ) . getExitRateVector ( ) ) ;
result = this - > computeTransientProbabilities < true > ( uniformizedMatrix , timeBound , uniformizationRate , result , * this - > linearEquationSolver ) ;
}
return result ;
}
template < class ValueType >
storm : : storage : : SparseMatrix < ValueType > SparseCtmcCslModelChecker < ValueType > : : computeProbabilityMatrix ( storm : : storage : : SparseMatrix < ValueType > const & rateMatrix , std : : vector < ValueType > const & exitRates ) {
// Turn the rates into probabilities by scaling each row with the exit rate of the state.