Gradient descent fails on bounded reachability properties

[jianlin-herman-li]

Hi, Storm developers,

I am trying to apply gradient descent to the dice model from [1], following the instructions in [2].

storm-pars works for the unbounded reachability property P >= 3/20 [ F d = 2 ], but fails for the bounded variant P >= 3/20 [ F <= 10 d = 2 ].

I am wondering whether the storm-pars --mode feasibility --feasibility:method gd supports bounded reachability properties, or if I am using it incorrectly.

Thanks for helping.

[1] https://arxiv.org/abs/1903.07993

[2] https://www.stormchecker.org/documentation/usage/running-storm-on-parametric-models.html

kydice.prism:

dtmc

const double x;
const double y;

module die

    // local state
    s : [0..7] init 0;
    // value of the dice
    d : [0..6] init 0;
    
    [] s=0 -> x : (s'=1) + 1-x : (s'=2);

    [] s=1 -> y : (s'=3) + 1-y : (s'=4);
    [] s=2 -> y : (s'=5) + 1-y : (s'=6);
    
    [] s=3 -> x : (s'=1)          + 1-x : (s'=7) & (d'=1);
    [] s=4 -> x : (s'=7) & (d'=3) + 1-x : (s'=7) & (d'=2);
    [] s=5 -> x : (s'=2)          + 1-x : (s'=7) & (d'=4);
    [] s=6 -> x : (s'=7) & (d'=6) + 1-x : (s'=7) & (d'=5);
    [] s=7 -> (s'=7);
    
endmodule

label "done" = (s=7);

Output from storm-pars:

> storm-pars --mode feasibility --feasibility:method gd --prism kydice.prism --prop "P>=3/20 [F d = 2]"
Storm-pars 1.11.1

Date: Tue Jan 27 19:37:36 2026
Command line arguments: --mode feasibility '--feasibility:method' gd --prism kydice.prism --prop 'P>=3/20 [F d = 2]'
Current working directory: /home/j2655li/pal-dev/benchmarks/kydice

Time for model input parsing: 0.000s.

Time for model construction: 0.010s.

-------------------------------------------------------------- 
Model type:     DTMC (sparse)
States:         13
Transitions:    20
Reward Models:  none
State Labels:   3 labels
   * deadlock -> 0 item(s)
   * init -> 1 item(s)
   * (d = 2) -> 1 item(s)
Choice Labels:  none
-------------------------------------------------------------- 

Time for model simplification: 0.000s.

-------------------------------------------------------------- 
Model type:     DTMC (sparse)
States:         6
Transitions:    10
Reward Models:  none
State Labels:   4 labels
   * (d = 2) -> 1 item(s)
   * deadlock -> 0 item(s)
   * target -> 1 item(s)
   * init -> 1 item(s)
Choice Labels:  none
-------------------------------------------------------------- 

Time for model preprocessing: 0.000s.

-------------------------------------------------------------- 
Model type:     DTMC (sparse)
States:         6
Transitions:    10
Reward Models:  none
State Labels:   4 labels
   * (d = 2) -> 1 item(s)
   * deadlock -> 0 item(s)
   * target -> 1 item(s)
   * init -> 1 item(s)
Choice Labels:  none
-------------------------------------------------------------- 
Find feasible solution for  P=? [F (d = 2)]
Finding an extremum using Gradient Descent
Trying out a new starting point
Trying initial guess (p->0.5 for every parameter p or set start point)
Starting at {x : 4503608634569751/9007199254740992, y : 4503608634569751/9007199254740992}
Aborting because the bound is satisfied
Result at initial state: 0.1666665556 ( approx. 0.1666665556) at [x=4503608634569751/9007199254740992, y=4503608634569751/9007199254740992].
Time for model checking: 0.001s.
> storm-pars --mode feasibility --feasibility:method gd --prism kydice.prism --prop "P>=3/20 [F<=10 d = 2]"
Storm-pars 1.11.1

Date: Tue Jan 27 19:37:42 2026
Command line arguments: --mode feasibility '--feasibility:method' gd --prism kydice.prism --prop 'P>=3/20 [F<=10 d = 2]'
Current working directory: /home/j2655li/pal-dev/benchmarks/kydice

Time for model input parsing: 0.000s.

Time for model construction: 0.010s.

-------------------------------------------------------------- 
Model type:     DTMC (sparse)
States:         13
Transitions:    20
Reward Models:  none
State Labels:   3 labels
   * deadlock -> 0 item(s)
   * init -> 1 item(s)
   * (d = 2) -> 1 item(s)
Choice Labels:  none
-------------------------------------------------------------- 

Time for model simplification: 0.000s.

-------------------------------------------------------------- 
Model type:     DTMC (sparse)
States:         6
Transitions:    10
Reward Models:  none
State Labels:   4 labels
   * init -> 1 item(s)
   * target -> 1 item(s)
   * deadlock -> 0 item(s)
   * (d = 2) -> 1 item(s)
Choice Labels:  none
-------------------------------------------------------------- 

Time for model preprocessing: 0.000s.

-------------------------------------------------------------- 
Model type:     DTMC (sparse)
States:         6
Transitions:    10
Reward Models:  none
State Labels:   4 labels
   * init -> 1 item(s)
   * target -> 1 item(s)
   * deadlock -> 0 item(s)
   * (d = 2) -> 1 item(s)
Choice Labels:  none
-------------------------------------------------------------- 
Find feasible solution for  P=? [true U<=10 (d = 2)]
Finding an extremum using Gradient Descent
ERROR (storm-pars.cpp:559): An unexpected exception occurred and caused Storm-pars to terminate. The message of this exception is: std::bad_cast

[volkm]

Thanks for the bug report.

The problem occurs around this line because we do not explicitly check for bounded untils.

@linusheck Do you know whether bounded untils should be supported or not? Depending on this we need to throw an exception.

[linusheck]

Hi! Nobody has yet asked for GD with finite-horizon properties, so the code assumes infinite-horizon. I think that you should be able to fix this by changing a couple lines of code (as long as the theory works out: is the bounded property on the derived system of equations the derived bounded probability? Probably yes:)). I’m on vacation right now though, until mid-February.

[linusheck]

Hi all,

I think that computing step-bounded derivatives for Markov chains requires a bit of additional correctness arguments. Storm computes derivatives using the systems of equations described in this paper: https://arxiv.org/pdf/2111.04407. One would need to argue that a step-bounded reachability solution of the derived system of equations is actually the same as the derivative of a step-bounded reachability solution. (Upon further thought, this is probably wrong for the equation system stated in the paper, as one "step" moves from the SOC into the derived SOC.)

Another way you could compute this value: You could see your target value as an “autodifferentiation” of a program computing the step-bounded probability. So, to compute a step-bounded probability, I just need to apply the Bellman operator n times:

set x to be 1 at target states
for i in range(n): x = B(x)
return x

Autodifferentiating this program would give you the value you want.