Merge pull request #87 from Zinoex/fm/product_bellman

Bellman operator for ProductProcess

Author: Zinoex

README.md

@@ -16,9 +16,9 @@ with great efficiency. Furthermore, it provides methods for accelerating the com
 certificate using CUDA hardware. 
 
 ## Features
-- O-maximization and value iteration
+- Value iteration (Bellman operator via O-maximization)
 - Dense and sparse matrix support
-- Parametric probability types for customizable precision
+- Parametric probability and value types for customizable precision including rationals and floating-point numbers
 - Multithreaded CPU and CUDA-accelerated value iteration
 - Data loading and writing in formats by various tools (PRISM, bmdp-tool, IMDP.jl)
 
@@ -72,7 +72,7 @@ V, k, residual = value_iteration(problem)
 See [Usage](https://www.baymler.com/IntervalMDP.jl/dev/usage/) for more information about different specifications, using sparse matrices, and CUDA.
 
 ## Ecosystem
-Building upon IntervalMDP.jl, we are designing an ecosystem of tools, which current consists of:
+Building upon IntervalMDP.jl, we are designing an ecosystem of tools, which currently consists of:
 
 - [IntervalMDPAbstractions.jl](https://github.com/Zinoex/IntervalMDPAbstractions.jl) - constructing abstractions of stochastic dynamical systems to verify properties.
 

docs/src/reference/specifications.md

@@ -8,21 +8,30 @@ strategy
 Specification
 system_property
 Property
+BasicProperty
+ProductProperty
 satisfaction_mode
 SatisfactionMode
 strategy_mode
 StrategyMode
 ```
 
-## Temporal logic
+## DFA Reachability
 
 ```@docs
-LTLFormula
-isfinitetime(prop::LTLFormula)
-LTLfFormula
-isfinitetime(prop::LTLfFormula)
-time_horizon(prop::LTLfFormula)
-PCTLFormula
+AbstractDFAReachability
+
+FiniteTimeDFAReachability
+isfinitetime(prop::FiniteTimeDFAReachability)
+terminal_states(prop::FiniteTimeDFAReachability)
+reach(prop::FiniteTimeDFAReachability)
+time_horizon(prop::FiniteTimeDFAReachability)
+
+InfiniteTimeDFAReachability
+isfinitetime(prop::InfiniteTimeDFAReachability)
+terminal_states(prop::InfiniteTimeDFAReachability)
+reach(prop::InfiniteTimeDFAReachability)
+convergence_eps(prop::InfiniteTimeDFAReachability)
 ```
 
 ## Reachability

docs/src/reference/systems.md

@@ -16,16 +16,16 @@ MixtureIntervalMarkovDecisionProcess
 MixtureIntervalMarkovChain
 stateptr(mdp::MixtureIntervalMarkovDecisionProcess)
 DFA
-letters2alphabet(letters::AbstractVector{String})
-alphabet2index(alphabet::AbstractVector{String})
+num_states(dfa::DFA)
+num_labels(dfa::DFA)
 transition(dfa::DFA)
-alphabetptr(dfa::DFA)
+labelmap(dfa::DFA)
 initial_state(dfa::DFA)
 accepting_states(dfa::DFA)
-ProductIntervalMarkovDecisionProcessDFA
-imdp(md::ProductIntervalMarkovDecisionProcessDFA)
-automaton(md::ProductIntervalMarkovDecisionProcessDFA)
-labelling_function(md::ProductIntervalMarkovDecisionProcessDFA)
+ProductProcess
+markov_process(proc::ProductProcess)
+automaton(proc::ProductProcess)
+labelling_function(proc::ProductProcess)
 ```
 
 ## Probability representation
@@ -75,12 +75,14 @@ weighting_probs
 ### Labelling of IMDP states to Automaton alphabet
 ```@docs
 LabellingFunction
-count_mapping(map::AbstractArray)
 mapping(labelling_func::LabellingFunction)
+num_labels(labelling_func::LabellingFunction)
 ```
 
 ### Transition function for DFA
 ```@docs
 TransitionFunction
 transition(transition_func::TransitionFunction)
+num_states(tf::TransitionFunction)
+num_labels(tf::TransitionFunction)
 ```

ext/IntervalMDPCudaExt.jl

@@ -117,14 +117,19 @@ function IntervalMDP.checkdevice(b::AbstractGPUArray, A::AbstractMatrix)
     @assert false "The reward vector is on the GPU ($(typeof(b))) and the transition matrix is a CPU matrix ($(typeof(A)))."
 end
 
+IntervalMDP.arrayfactory(
+    ::MR,
+    T,
+    num_states,
+) where {R, MR <: Union{CuSparseMatrixCSC{R}, CuMatrix{R}}} = CUDA.zeros(T, num_states)
+
 include("cuda/utils.jl")
 include("cuda/array.jl")
 include("cuda/sorting.jl")
 include("cuda/workspace.jl")
 include("cuda/strategy.jl")
 include("cuda/bellman/dense.jl")
 include("cuda/bellman/sparse.jl")
-include("cuda/value_iteration.jl")
 include("cuda/interval_probabilities.jl")
 include("cuda/specification.jl")
 

ext/cuda/strategy.jl

@@ -1,10 +1,3 @@
-function IntervalMDP.construct_action_cache(
-    ::IntervalProbabilities{R, VR},
-    dims,
-) where {R <: Real, VR <: AbstractGPUVector{R}}
-    return CUDA.zeros(Int32, dims)
-end
-
 abstract type ActiveCache end
 abstract type OptimizingActiveCache <: ActiveCache end
 

ext/cuda/value_iteration.jl

@@ -19,8 +19,9 @@ export StationaryStrategy, TimeVaryingStrategy
 export construct_strategy_cache, time_length
 
 include("specification.jl")
-export Property, LTLFormula, LTLfFormula, PCTLFormula
+export Property, BasicProperty, ProductProperty
 
+export AbstractDFAReachability, FiniteTimeDFAReachability, InfiniteTimeDFAReachability
 export AbstractReachability,
     FiniteTimeReachability, InfiniteTimeReachability, ExactTimeReachability
 export AbstractReachAvoid, FiniteTimeReachAvoid, InfiniteTimeReachAvoid, ExactTimeReachAvoid

src/IntervalMDP.jl

@@ -1,5 +1,5 @@
 """
-    bellman(V, model; upper_bound = false)
+    bellman(V, model; upper_bound = false, maximize = true)
 
 Compute robust Bellman update with the value function `V` and the model `model`, e.g. [`IntervalMarkovDecisionProcess`](@ref),
 that upper or lower bounds the expectation of the value function `V` via O-maximization [1].
@@ -56,7 +56,7 @@ Vcur = bellman(Vprev, model; upper_bound = false)
 
 !!! note
     This function will construct a workspace object and an output vector.
-    For a hot-loop, it is more efficient to use `bellman!` and pass in pre-allocated objects.
+    For a hot-loop, it is more efficient to use [`bellman!`](@ref) and pass in pre-allocated objects.
 
 [1] M. Lahijanian, S. B. Andersson and C. Belta, "Formal Verification and Synthesis for Discrete-Time Stochastic Systems," in IEEE Transactions on Automatic Control, vol. 60, no. 8, pp. 2031-2045, Aug. 2015, doi: 10.1109/TAC.2015.2398883.
 
@@ -67,7 +67,7 @@ function bellman(V, model; upper_bound = false, maximize = true)
 end
 
 """
-    bellman!(workspace, strategy_cache, Vres, V, model, stateptr; upper_bound = false, maximize = true)
+    bellman!(workspace, strategy_cache, Vres, V, model; upper_bound = false, maximize = true)
 
 Compute in-place robust Bellman update with the value function `V` and the model `model`, 
 e.g. [`IntervalMarkovDecisionProcess`](@ref), that upper or lower bounds the expectation of the value function `V` via O-maximization [1].
@@ -155,7 +155,7 @@ function bellman!(
     strategy_cache,
     Vres,
     V,
-    model;
+    model::IntervalMarkovProcess;
     upper_bound = false,
     maximize = true,
 )
@@ -171,6 +171,72 @@ function bellman!(
     )
 end
 
+function bellman!(
+    workspace::ProductWorkspace,
+    strategy_cache,
+    Vres,
+    V,
+    model::ProductProcess;
+    upper_bound = false,
+    maximize = true,
+)
+    mp = markov_process(model)
+    lf = labelling_function(model)
+    dfa = automaton(model)
+
+    W = workspace.intermediate_values
+
+    @inbounds for state in dfa
+        local_strategy_cache = localize_strategy_cache(strategy_cache, state)
+
+        # Select the value function for the current DFA state
+        # according to the appropriate DFA transition function
+        map!(W, CartesianIndices(product_num_states(mp))) do idx
+            return V[idx, dfa[state, lf[idx]]]
+        end
+
+        # For each state in the product process, compute the Bellman operator
+        # for the corresponding Markov process
+        bellman!(
+            workspace.underlying_workspace,
+            local_strategy_cache,
+            selectdim(Vres, ndims(Vres), state),
+            W,
+            mp;
+            upper_bound = upper_bound,
+            maximize = maximize,
+        )
+    end
+
+    return Vres
+end
+
+function localize_strategy_cache(strategy_cache::NoStrategyCache, dfa_state)
+    return strategy_cache
+end
+
+function localize_strategy_cache(strategy_cache::TimeVaryingStrategyCache, dfa_state)
+    return TimeVaryingStrategyCache(
+        selectdim(
+            strategy_cache.cur_strategy,
+            ndims(strategy_cache.cur_strategy),
+            dfa_state,
+        ),
+    )
+end
+
+function localize_strategy_cache(strategy_cache::StationaryStrategyCache, dfa_state)
+    return StationaryStrategyCache(
+        selectdim(strategy_cache.strategy, ndims(strategy_cache.strategy), dfa_state),
+    )
+end
+
+function localize_strategy_cache(strategy_cache::ActiveGivenStrategyCache, dfa_state)
+    return ActiveGivenStrategyCache(
+        selectdim(strategy_cache.strategy, ndims(strategy_cache.strategy), dfa_state),
+    )
+end
+
 ######################################################
 # Bellman operator for IntervalMarkovDecisionProcess #
 ######################################################