Markovian model for Leios efficiency (#559)

* Implemented model
* Added documentation
* Ran example
* Updated logbook

Author: bwbush

.github/workflows/lean_action_ci.yml

@@ -0,0 +1,20 @@
+name: Lean Action CI
+
+on:
+  pull_request:
+    paths:
+      - "analysis/markov/**"
+  push:
+    branches:
+      - main
+    paths:
+      - "analysis/markov/**"
+
+jobs:
+  build:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v4
+      - uses: leanprover/lean-action@v1
+        with:
+          lake-package-directory: analysis/markov

Logbook.md

@@ -1,5 +1,18 @@
 # Leios logbook
 
+## 2025-10-02
+
+### Markovian simulator for Linear Leios
+
+A new [Markovian simulation of Linear Leios](analysis/markov/ReadMe.md) computes the probability of EB certifications as RBs are produced. It probabilistically models four key transitions in Linear Leios and computes the efficiency of RB and EB production under potential adversarial conditions.
+
+1. *Forge RB:* create a new RB.
+2. *Certify:* include a certificate in the RB.
+3. *Forge EB:* create a new EB.
+4. *Vote:* cast votes to reach a quorum.
+
+![Example results](analysis/markov/example-results.png)
+
 ## 2025-09-19
 
 ### Antithesis configuration for Rust simulator

analysis/markov/.gitignore

@@ -0,0 +1 @@
+/.lake

analysis/markov/Linleios.lean

@@ -0,0 +1 @@
+import Linleios.Evolve

analysis/markov/Linleios/Evolve.lean

@@ -0,0 +1,122 @@
+
+import Std.Data.HashMap
+import Batteries.Lean.HashMap
+import Lean.Data.Json.FromToJson
+
+import Linleios.Types
+
+
+open Lean (Json)
+open Lean.ToJson (toJson)
+open Std (HashMap)
+
+
+def forgeRb {env : Environment} (state : State) : Outcomes :=
+  let state' :=
+    {
+      state with
+      clock := state.clock + 1
+    }
+  let p := 1 - env.fAdversary
+  [
+    ({state' with hasRb := true, rbCount := state.rbCount + 1}, p)
+  , ({state' with hasRb := false, canCertify := false}, 1 - p)
+  ]
+
+def certify {env : Environment} (state : State) : Outcomes :=
+  if state.hasRb && state.canCertify
+    then let p := env.pSpacingOkay
+         [
+           ⟨{state with ebCount := state.ebCount + 1}, p⟩
+         , ⟨state, 1 - p⟩
+         ]
+    else [(state, 1)]
+
+def forgeEb {env : Environment} (state : State) : Outcomes :=
+  if state.hasRb
+    then let p := 1 - env.pLate
+         [
+           ({state with canCertify := true}, p)
+         , ({state with canCertify := false}, 1 - p)
+         ]
+    else [(state, 1)]
+
+def vote {env : Environment} (state : State) : Outcomes :=
+  if state.hasRb
+    then let p := env.pQuorum
+         [
+           (state, p)
+         , ({state with canCertify := false}, 1 - p)
+         ]
+    else [(state, 1)]
+
+def step {env : Environment} : List (State → Outcomes) :=
+  [@forgeRb env, @certify env, @forgeEb env, @vote env]
+
+
+def prune (ε : Float) : Probabilities → Probabilities :=
+  HashMap.filter (fun _ p => p > ε)
+
+def evolve (transition : State → Outcomes) : Probabilities → Probabilities :=
+  HashMap.fold
+    (
+      fun acc state p =>
+        HashMap.mergeWith (fun _ => Add.add) acc
+          ∘ HashMap.map (fun _ => Mul.mul p ∘ List.sum ∘ List.map Prod.snd)
+          ∘ List.groupByKey Prod.fst
+          $ transition state
+    )
+    ∅
+
+def chain (transitions : List (State → Outcomes)) : Probabilities → Probabilities :=
+  match transitions with
+  | [] => id
+  | (t :: ts) => chain ts ∘ evolve t
+
+def simulate (env : Environment) (start : Probabilities) (ε : Float) : Nat → Probabilities
+| 0     => start
+| n + 1 => let state' := prune ε $ chain (@step env) start
+           simulate env state' ε n
+
+def totalProbability (states : Probabilities) : Probability :=
+  states.values.sum
+
+def ebDistribution : Probabilities → HashMap Nat Probability :=
+  HashMap.fold
+    (
+      fun acc state p =>
+        HashMap.mergeWith (fun _ => Add.add) acc
+          $ singleton ⟨ state.ebCount , p ⟩
+    )
+    ∅
+
+def ebDistributionJson : Probabilities → Json :=
+  Json.mkObj
+    ∘ List.map (fun ⟨k, v⟩ => ⟨toString k, toJson v⟩)
+    ∘ (fun z => z.mergeSort (by exact fun x y => x.fst ≤ y.fst))
+    ∘ HashMap.toList
+    ∘ ebDistribution
+
+def rbEfficiency (states : Probabilities) : Float :=
+  let clock := states.keys.head!.clock
+  let rbCount :=
+    HashMap.fold
+      (fun acc state p =>acc + state.rbCount.toFloat * p)
+      0
+      states
+  rbCount / clock.toFloat
+
+def ebEfficiency (states : Probabilities) : Float :=
+  let clock := states.keys.head!.clock
+  let ebCount :=
+    HashMap.fold
+      (fun acc state p =>acc + state.ebCount.toFloat * p)
+      0
+      states
+  ebCount / (clock.toFloat - 1)
+
+def efficiency (states : Probabilities) : Float :=
+  let rb := rbEfficiency states
+  let eb := ebEfficiency states
+  let ρ := 12.5 / 0.9
+  (rb + ρ * eb) / (1 + ρ)

analysis/markov/Linleios/Probability.lean

@@ -0,0 +1,27 @@
+
+import Linleios.Util
+
+
+def nPools : Nat := 2500
+
+def stakeDistribution (nPools : Nat) : List Float :=
+  (List.range nPools).map (fun k => ((k + 1).toFloat / nPools.toFloat)^10 - (k.toFloat / nPools.toFloat)^10)
+
+private def calibrateCommittee(committeeSize : Float) : Float :=
+  let stakes : List Float := stakeDistribution nPools
+  let f (m : Float) : Float :=
+    let ps : List Float := stakes.map (fun s => 1 - Float.exp (- m * s))
+    ps.sum - committeeSize
+  bisectionSearch f committeeSize nPools.toFloat 0.5 10
+
+private def committeeDistribution (pSuccessfulVote committeeSize : Float) : Float × Float :=
+  let stakes : List Float := stakeDistribution nPools
+  let m := calibrateCommittee committeeSize
+  let ps : List Float := stakes.map (fun s => pSuccessfulVote * (1 - Float.exp (- m * s)))
+  let μ := ps.sum
+  let σ := (ps.map (fun p => p * (1 - p))).sum.sqrt
+  ⟨ μ , σ ⟩
+
+def pQuorum (pSuccessfulVote committeeSize τ : Float) : Float :=
+  let ⟨ μ , σ ⟩ := committeeDistribution pSuccessfulVote committeeSize
+  1 - cdfGaussian (τ * committeeSize) μ σ

analysis/markov/Linleios/Types.lean

@@ -0,0 +1,60 @@
+import Std.Data.HashMap
+import Batteries.Lean.HashMap
+
+import Linleios.Probability
+
+
+open Std (HashMap)
+
+
+abbrev Probability := Float
+
+
+structure Environment where
+  activeSlotCoefficient : Probability
+  Lheader : Nat
+  Lvote : Nat
+  Ldiff : Nat
+  pSpacingOkay : Probability
+  pQuorum : Probability
+  pLate : Probability
+  fAdversary : Float
+deriving Repr
+
+def makeEnvironment (Lheader Lvote Ldiff : Nat) (activeSlotCoefficient committeeSize τ pRbHeaderArrives pEbValidates pEbUnvalidated fAdversary : Float) : Environment :=
+  {
+    activeSlotCoefficient := activeSlotCoefficient
+    Lheader := Lheader
+    Lvote := Lvote
+    Ldiff := Ldiff
+    pSpacingOkay := (1 - activeSlotCoefficient).pow (3 * Lheader + Lvote + Ldiff - 1).toFloat
+    pQuorum := pQuorum (pRbHeaderArrives * pEbValidates * (1 - fAdversary)) committeeSize τ
+    pLate := pEbUnvalidated * (1 - committeeSize / nPools.toFloat)
+    fAdversary := fAdversary
+  }
+
+
+structure State where
+  clock : Nat
+  rbCount : Nat
+  ebCount : Nat
+  hasRb : Bool
+  canCertify : Bool
+deriving Repr, BEq, Hashable, Inhabited
+
+theorem genesis : (default : State).clock = 0 ∧ (default : State).rbCount = 0 ∧ (default : State).ebCount = 0 := by
+  constructor
+  rfl
+  constructor
+  rfl
+  rfl
+
+
+def Probabilities := HashMap State Probability
+deriving Repr, EmptyCollection
+
+instance : Inhabited Probabilities where
+  default := (∅ : Probabilities).insert Inhabited.default 1
+
+
+abbrev Outcomes := List (State × Probability)

analysis/markov/Linleios/Util.lean

@@ -0,0 +1,30 @@
+
+partial def erf (x : Float) : Float :=
+  if x < 0
+    then - erf (- x)
+    else
+      let p := 0.3275911
+      let a₁ := 0.254829592
+      let a₂ := -0.284496736
+      let a₃ := 1.421413741
+      let a₄ := -1.453152027
+      let a₅ := 1.061405429
+      let t := 1 / (1 + p * x)
+      1 - (a₁ * t + a₂ * t^2 + a₃ * t^3 + a₄ * t^4 + a₅ * t^5) * Float.exp (- x^2)
+
+def cdfGaussian (x μ σ : Float) : Float :=
+  (1 + erf ((x - μ) / σ / Float.sqrt 2)) / 2
+
+def bisectionSearch (f : Float → Float) (low high : Float) (ε : Float) (maxIter : Nat) : Float :=
+  match maxIter with
+  | 0 => (low + high) / 2
+  | maxIter' + 1 =>
+    let mid := (low + high) / 2
+    let fmid := f mid
+    if high - low < ε || Float.abs fmid < ε then
+      mid
+    else if f low * fmid < 0 then
+      bisectionSearch f low mid ε maxIter'
+    else
+      bisectionSearch f mid high ε maxIter'
+termination_by maxIter

analysis/markov/Main.lean

@@ -0,0 +1,77 @@
+
+import Cli
+import Linleios
+import Parser.Char.Numeric
+import Parser.Stream
+
+open Cli
+open Lean (Json)
+
+
+instance : ParseableType Float where
+  name     := "Float"
+  parse? x := match (Parser.Char.ASCII.parseFloat : SimpleParser Substring Char Float).run x with
+              | .ok _ y => some y
+              | _       => none
+
+
+def runMarkovCmd (p : Parsed) : IO UInt32 :=
+  do
+    let activeSlotCoefficient : Float := p.flag! "active-slot-coefficient" |>.as! Float
+    let Lheader : Nat := p.flag! "l-header" |>.as! Nat
+    let Lvote : Nat := p.flag! "l-vote" |>.as! Nat
+    let Ldiff : Nat := p.flag! "l-diff" |>.as! Nat
+    let committeeSize : Nat := p.flag! "committee-size" |>.as! Nat
+    let τ : Float := p.flag! "quorum-fraction" |>.as! Float
+    let pRbHeaderArrives : Float := p.flag! "p-rb-header-arrives" |>.as! Float
+    let pEbValidates : Float := p.flag! "p-eb-validates" |>.as! Float
+    let pLateDiffusion : Float := p.flag! "p-late-diffusion" |>.as! Float
+    let fAdversary : Float := p.flag! "adversary-fraction" |>.as! Float
+    let ε : Float := p.flag! "tolerance" |>.as! Float
+    let rbCount : Nat := p.flag! "rb-count" |>.as! Nat
+    let env := makeEnvironment Lheader Lvote Ldiff activeSlotCoefficient committeeSize.toFloat τ pRbHeaderArrives pEbValidates pLateDiffusion fAdversary
+    let sn := simulate env default ε rbCount
+    if p.hasFlag "output-file"
+      then IO.FS.writeFile (p.flag! "output-file" |>.as! String) (Json.pretty $ ebDistributionJson sn)
+    IO.println s!"RB efficiency: {(reprPrec (rbEfficiency sn) 0).pretty}"
+    IO.println s!"EB efficiency: {(reprPrec (ebEfficiency sn) 0).pretty}"
+    IO.println s!"Overall efficiency: {(reprPrec (efficiency sn) 0).pretty}"
+    IO.eprintln s!"Missing probability: {1 - totalProbability sn}"
+    pure 0
+
+def markov : Cmd := `[Cli|
+  markov VIA runMarkovCmd; ["0.0.1"]
+  "Run a Markov simulation of Linear Leios."
+  FLAGS:
+    "active-slot-coefficient" : Float  ; "The active slot coefficient for RBs, in probability per slot."
+    "l-header"                : Nat    ; "L_header protocol parameter, in slots."
+    "l-vote"                  : Nat    ; "L_vote protocol parameter, in slots."
+    "l-diff"                  : Nat    ; "L_diff protocol parameter, in slots."
+    "committee-size"          : Nat    ; "Expected number of voters in the committee, out of 2500 stakepools total."
+    "quorum-fraction"         : Float  ; "τ protocol parameter, in %/100."
+    "p-rb-header-arrives"     : Float  ; "Probability that the RB header arrives before L_header."
+    "p-eb-validates"          : Float  ; "Probability that the EB is fully validated before 3 L_header + L_vote."
+    "p-late-diffusion"        : Float  ; "Probability than a RB producer hasn't validated the previous EB."
+    "adversary-fraction"      : Float  ; "Fraction of stake that is adversarial: the adversary never produces RBs or EBs and never votes."
+    "tolerance"               : Float  ; "Discard states with less than this probability."
+    "rb-count"                : Nat    ; "Number of potential RBs to simulate."
+    "output-file"             : String ; "Path to the JSON output file for the EB distribution."
+  EXTENSIONS:
+    defaultValues! #[
+      ("active-slot-coefficient", "0.05")
+    , ("l-header"               , "1"   )
+    , ("l-vote"                 , "4"   )
+    , ("l-diff"                 , "7"   )
+    , ("committee-size"         , "600" )
+    , ("quorum-fraction"        , "0.75")
+    , ("p-rb-header-arrives"    , "0.95")
+    , ("p-eb-validates"         , "0.90")
+    , ("p-late-diffusion"       , "0.05")
+    , ("adversary-fraction"     , "0"   )
+    , ("tolerance"              , "1e-8")
+    , ("rb-count"               , "100" )
+    ]
+]
+
+def main (args : List String) : IO UInt32 :=
+  markov.validate args