logos-blockchain/research
1
0
Fork 0
mirror of https://github.com/logos-blockchain/research.git synced 2026-08-07 03:33:33 +00:00
Code Issues Packages Projects Releases Wiki Activity

Importing tsi-sim v2

Browse Source
This commit is contained in:
Marcin Pawlowski 2026-07-30 18:52:01 +02:00
parent 97a4e8cc30
commit e86bb0cb6c
No known key found for this signature in database
44 changed files with 2755 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

18
tools/simulators/tsi/tsi-sim-mc/.gitignore vendored Normal file
View File

@ -0,0 +1,18 @@
# Generated artifacts
runs/
results/
figures/
!results/.gitkeep
!figures/.gitkeep
# Python
__pycache__/
*.py[cod]
*.egg-info/
.venv/
venv/
.pytest_cache/
.mypy_cache/
.ruff_cache/
*.parquet
*.csv

57
tools/simulators/tsi/tsi-sim-mc/Makefile Normal file
View File

@ -0,0 +1,57 @@
.PHONY: install smoke sweep sweep-fullscale figures verify test test-slow test-all lint clean
VENV ?= .venv
PY := $(VENV)/bin/python
PIP := $(VENV)/bin/pip
# Pin BLAS/OpenMP to a single thread so joblib's N worker processes don't oversubscribe
# the cores (numpy in each worker would otherwise each spawn a full BLAS thread pool).
export OMP_NUM_THREADS := 1
export OPENBLAS_NUM_THREADS := 1
export MKL_NUM_THREADS := 1
export NUMEXPR_NUM_THREADS := 1
$(VENV):
python3 -m venv $(VENV)
$(PIP) install --upgrade pip
install: $(VENV)
$(PIP) install -e ".[dev]"
# Every sweep writes results + figures into a fresh dated folder under runs/ so new
# runs never overwrite old ones (runs/<YYYY-MM-DD_HHMMSS>_<label>/).
# Fast end-to-end check (tiny scaled-k grid); uses all cores across configs.
smoke: install
$(PY) scripts/run_sweep.py --config configs/smoke.yaml --label smoke
# Full scaled-k parameter sweep (multicore across configs; auto-generates figures).
sweep: install
$(PY) scripts/run_sweep.py --config configs/default.yaml --label default
# Full-scale (true k) grid; small + heavy, so load-balance one config at a time.
sweep-fullscale: install
$(PY) scripts/run_sweep.py --config configs/fullscale.yaml --label fullscale --batch-size 1
# Re-render figures from a run's parquet into a fresh dated figures/ folder:
# make figures RESULTS=runs/<dir>/results.parquet
figures: install
$(PY) scripts/make_figures.py --results $(RESULTS)
# Analytic sanity checks (simulator vs closed-form theory; replicates run across cores).
verify: install
$(PY) scripts/verify.py
test: install # fast subset (addopts already excludes slow)
$(PY) -m pytest
test-slow: install
$(PY) -m pytest -m slow
test-all: install
$(PY) -m pytest -m ''
lint: install
$(VENV)/bin/ruff check src scripts tests
clean:
rm -rf runs/* results/*.parquet figures/* .pytest_cache .ruff_cache .mypy_cache
find . -name __pycache__ -type d -prune -exec rm -rf {} +

102
tools/simulators/tsi/tsi-sim-mc/README.md Normal file
View File

@ -0,0 +1,102 @@
# tsi-sim-mc — Cryptarchia TSI simulator (multicore build)
> **This is the multicore-optimised, reviewed copy of `../tsi-sim/`.** The original is left
> untouched. Versus the original it adds: a **sparse lottery sampler** (~30–100× faster per
> epoch), a **hardened multicore sweep** (loky + single-thread BLAS), an **opt-in parallel
> chunked lottery**, an optional **spec fixed-point mode**, corrected figures, config
> validation, and a larger test suite. See "Performance & reproducibility" below.
Monte-Carlo simulation framework used to choose safe values for the **uncle-reference**
parameters of Cryptarchia's Total Stake Inference (TSI):
- `U` — max uncles referenced per block (`MAX_UNCLES`); `U=0` is the no-uncle baseline.
- `W` — uncle reference window in slots (spec default 300).
- swept against network size `N`, stake distribution `S` (uniform / Pareto), and
network **latency** `L` (in slots — deliberately *not* `D`, which denotes the stake estimate).
It measures how well the inferred total active stake `D` tracks the true total stake, and
whether uncle references recover the active slots that network latency loses to forks.
> This lives under a `raw/` docs path, so it is invisible to the repository's
> markdown-lint CI. It is a standalone Python package with its own tooling.
## Model
Reduced **canonical-chain-with-orphans** model: we simulate the global winning-slot
sequence (stake-weighted φ lottery), build a real block tree with latency- and
multi-winner-induced forks, resolve the canonical chain (honest longest-chain), let
canonical blocks reference uncles per the spec's selection rules, and count TSI density
`m = honest-chain blocks + deduplicated referenced uncles` in the measurement window.
All honest nodes converge to the same deep chain (k-finality), so a single per-epoch
`D` is faithful. A full per-node model is the planned next phase (`per_node_dest` flag
scaffolds it).
See the sibling spec `../` and `../../cryptarchia-total-stake-inference.md` for the math.
## Quick start
```bash
make install # create .venv and install (editable) with dev deps
make test # unit tests + fast theory checks
make verify # simulator vs closed-form analytic checks
make smoke # tiny scaled-k sweep + figures (end-to-end smoke test)
make sweep figures # full scaled-k parameter sweep + academic figures
```
Outputs: `results/*.parquet` (one row per config×epoch) and `figures/*.{pdf,png}`
(both git-ignored).
## Scale
True constants (`k=2160`, `f=1/30`) give 648,000-slot epochs. Mean accuracy is provably
`k`-invariant (only variance scales `~1/T`), so sweeps use a **scaled `k`**
(`configs/default.yaml`); the final accuracy/variance figures re-run at true `k`
(`configs/fullscale.yaml`, `make sweep-fullscale`). `configs/smoke.yaml` is a tiny dev grid.
## Performance & reproducibility
- **Sparse lottery (the main win).** The number of slots a node wins is `Binomial(n_slots,
p_i)` and the won slots are a uniform distinct subset — distributionally identical to an
independent Bernoulli per slot, but without the dense `(n_nodes, n_slots)` array that was
~95% of runtime. Full-scale (`k=2160`) epoch: **~3.3 s → ~0.1 s (~30×)**; at `k=256`,
**~0.39 s → ~0.009 s (~40×)**.
- **Multicore across configs (the main lever).** `run_trajectory` is a pure function of a
hash-seeded, immutable `SimConfig`, so the sweep is order-independent and embarrassingly
parallel. `run_sweep` uses joblib's process-based **loky** backend with
`inner_max_num_threads=1` (and the Makefile pins `*_NUM_THREADS=1`) to use all cores
without BLAS oversubscription. `--n-jobs -1` (default) uses every core; `--batch-size 1`
suits the small heavy full-scale grid. Measured on a 14-core box, the full scaled-`k`
sweep (`configs/default.yaml`, 3888 configs) runs in **~58 s parallel vs ~561 s serial
(9.7×)**.
- **Opt-in within-config parallelism.** `lottery_chunks > 1` splits the per-slot lottery
across slot-chunks with independent `SeedSequence.spawn` children. After the sparse fix the
lottery is a small fraction of an epoch, so this rarely helps — it exists for a single
isolated config with an enormous `n_slots`. **Reproducibility caveat:** the chunked stream
differs from the serial stream and *changes with `n_chunks`*, so `lottery_chunks` must be a
pinned, recorded parameter, never derived from the core count.
- **RNG reproducibility.** Every draw is a deterministic spawn off `SeedSequence(hash(config))`
— child 0 draws stake, child `e+1` drives epoch `e`, which spawns lottery/aux sub-streams.
Results are identical regardless of parallel scheduling order. The sparse sampler consumes
the RNG differently from the original dense one, so committed baselines here were
regenerated against the sparse sampler.
## Faithfulness notes (from review)
- `update_D(..., fixed_point=True)` mirrors the spec's integer `f`-truncation
(`int(f·1000)/1000 = 0.033`), reproducing the on-chain estimator's ~1% systematic
overestimate. Default `False` keeps the exact-`f`, analysis-faithful behaviour.
- TSI counts *blocks* (wins), so at full uncle recovery the estimate equilibrates at the
**block-count ceiling** `-ln(1-f)/f ≈ 1.017`, not 1.0 — a deterministic overshoot, not
noise. Figures overlay this ceiling; `theory.block_count_ceiling` computes it.
## Layout
```
src/tsi_sim/ constants config rng stake lottery latency blocktree uncles tsi epoch
engine metrics theory sweep verify plotting/{style,figures,make_figures}
scripts/ run_sweep.py make_figures.py verify.py (thin shims; installed as
tsi-sweep / tsi-figures / tsi-verify)
configs/ smoke.yaml default.yaml fullscale.yaml
tests/ test_{lottery,uncles,tsi_counting,blocktree,config,rng,stake,theory,
latency,theory_convergence}.py
```

16
tools/simulators/tsi/tsi-sim-mc/configs/default.yaml Normal file
View File

@ -0,0 +1,16 @@
# Full scaled-k parameter sweep. Mean accuracy / q_eff / convergence are k-invariant,
# so k is scaled down for tractability; re-run headline + variance figures at full scale
# with fullscale.yaml.
# Sweep axes (cartesian product x replicates); every value is a list.
n_nodes: [1000, 2000, 4000] # number of nodes / stake holders
stake_dist: [uniform, pareto] # stake distribution (uniform = equal, pareto = heavy-tailed)
latency: [0, 1, 2, 4, 6, 8, 10, 15, 20, 30] # L: network latency in slots; block visible at t+L
max_uncles: [0, 1, 2, 3, 4, 8, 10] # U: max uncle references per block (0 = baseline)
uncle_strategy: [oldest, random] # uncle selection: oldest-first fill vs random coin-flip
replicates: 12 # independent RNG replicates per grid cell
base: # per-run settings shared by every cell (not swept)
k: 256 # security parameter (scaled; T = 6*floor(256/f) = 46080 slots)
epochs: 75 # epochs simulated per trajectory
f: 0.03333333333333333 # slot activation coefficient (default 1/30); configurable
genesis_d_factor: 0.5 # genesis D_est = factor x true total stake
pareto_shape: 1.16 # Pareto (Lomax) tail index for the pareto stake distribution

16
tools/simulators/tsi/tsi-sim-mc/configs/fullscale.yaml Normal file
View File

@ -0,0 +1,16 @@
# Full-scale (true k=2160) confirmation runs. Same latency/uncle grid as default.yaml but
# at the true security parameter (k not scaled down); heavier per config, so N and reps are
# kept modest. Note: T = 6*floor(2160/f) = 388,800 slots per epoch.
# Sweep axes (cartesian product x replicates); every value is a list.
n_nodes: [1000] # number of nodes / stake holders
stake_dist: [pareto] # heavy-tailed (realistic) stake distribution
latency: [0, 1, 2, 4, 6, 8, 10, 15, 20, 30] # L: network latency in slots; block visible at t+L
max_uncles: [0, 1, 2, 3, 4, 8, 10] # U: max uncle references per block (0 = baseline)
uncle_strategy: [oldest] # uncle selection: oldest-first fill
replicates: 10 # independent RNG replicates per grid cell
base: # per-run settings shared by every cell (not swept)
k: 2160 # true security parameter (not scaled)
epochs: 60 # epochs simulated per trajectory
f: 0.03333333333333333 # slot activation coefficient (default 1/30); configurable
genesis_d_factor: 0.001 # genesis D_est = factor x true total stake
pareto_shape: 1.16 # Pareto (Lomax) tail index for the pareto stake distribution

14
tools/simulators/tsi/tsi-sim-mc/configs/smoke.yaml Normal file
View File

@ -0,0 +1,14 @@
# Tiny end-to-end smoke grid (fast; for dev + CI-style checks).
# Sweep axes (cartesian product x replicates); every value is a list.
n_nodes: [1000] # number of nodes / stake holders
stake_dist: [uniform, pareto] # stake distribution (uniform = equal, pareto = heavy-tailed)
latency: [0, 4, 8] # L: network latency in slots; block visible to others at t+L
max_uncles: [0, 1, 4] # U: max uncle references per block (0 = baseline, no uncles)
uncle_strategy: [oldest] # uncle selection: oldest-first fill
replicates: 4 # independent RNG replicates per grid cell
base: # per-run settings shared by every cell (not swept)
k: 16 # security parameter (scaled; T = 6*floor(16/f) = 2880 slots)
epochs: 25 # epochs simulated per trajectory
f: 0.03333333333333333 # slot activation coefficient (default 1/30); configurable
genesis_d_factor: 0.5 # genesis D_est = factor x true total stake
pareto_shape: 1.16 # Pareto (Lomax) tail index for the pareto stake distribution

47
tools/simulators/tsi/tsi-sim-mc/pyproject.toml Normal file
View File

@ -0,0 +1,47 @@
[build-system]
requires = ["hatchling"]
build-backend = "hatchling.build"
[project]
name = "tsi-sim"
version = "0.1.0"
description = "Monte-Carlo simulation of Cryptarchia Total Stake Inference with uncle references"
readme = "README.md"
requires-python = ">=3.11"
license = { text = "MIT" }
authors = [{ name = "Logos" }]
dependencies = [
"numpy>=1.26",
"pandas>=2.1",
"pyarrow>=14",
"matplotlib>=3.8",
"scipy>=1.11",
"pyyaml>=6.0",
"tqdm>=4.66",
"joblib>=1.3",
]
[project.optional-dependencies]
dev = ["pytest>=8.0", "pytest-xdist>=3.5", "ruff>=0.5", "mypy>=1.8"]
accel = ["numba>=0.60"] # optional; only used if a single-config tree build is a bottleneck
[project.scripts]
tsi-sweep = "tsi_sim.sweep:main"
tsi-verify = "tsi_sim.verify:main"
tsi-figures = "tsi_sim.plotting.make_figures:main"
[tool.hatch.build.targets.wheel]
packages = ["src/tsi_sim"]
[tool.ruff]
line-length = 100
target-version = "py311"
[tool.ruff.lint]
select = ["E", "F", "I", "UP", "B", "NPY"]
ignore = ["E741"] # allow single-char names like L (latency), q, m — domain symbols
[tool.pytest.ini_options]
testpaths = ["tests"]
markers = ["slow: long-running (full-scale k) checks"]
addopts = "-q -m 'not slow'" # fast subset by default; `make test-all` runs everything

1
tools/simulators/tsi/tsi-sim-mc/requirements-dev.txt Normal file
View File

@ -0,0 +1 @@
-e .[dev]

1
tools/simulators/tsi/tsi-sim-mc/requirements.txt Normal file
View File

@ -0,0 +1 @@
-e .

14
tools/simulators/tsi/tsi-sim-mc/scripts/make_figures.py Normal file
View File

@ -0,0 +1,14 @@
#!/usr/bin/env python
"""Thin shim so `python scripts/make_figures.py` works without installing."""
from __future__ import annotations
import sys
from pathlib import Path
sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "src"))
from tsi_sim.plotting.make_figures import main # noqa: E402
if __name__ == "__main__":
main()

14
tools/simulators/tsi/tsi-sim-mc/scripts/run_sweep.py Normal file
View File

@ -0,0 +1,14 @@
#!/usr/bin/env python
"""Run a TSI parameter sweep from a YAML config -> results parquet."""
from __future__ import annotations
import sys
from pathlib import Path
sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "src"))
from tsi_sim.sweep import main # noqa: E402
if __name__ == "__main__":
main()

14
tools/simulators/tsi/tsi-sim-mc/scripts/verify.py Normal file
View File

@ -0,0 +1,14 @@
#!/usr/bin/env python
"""Thin shim so `python scripts/verify.py` works without installing; see tsi_sim.verify."""
from __future__ import annotations
import sys
from pathlib import Path
sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "src"))
from tsi_sim.verify import main # noqa: E402
if __name__ == "__main__":
raise SystemExit(main())

3
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/__init__.py Normal file
View File

@ -0,0 +1,3 @@
"""Cryptarchia Total Stake Inference simulator (uncle references)."""
__version__ = "0.1.0"

127
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/blocktree.py Normal file
View File

@ -0,0 +1,127 @@
"""Block tree, latency-driven forks, and honest longest-chain fork choice.
Blocks are stored in parallel arrays (id == index). A virtual genesis is block 0 at
slot -1, height 0. Every real block is produced at an active slot by one winning node
and points at the best tip *visible to that node at production time*, which is what makes
network latency (and same-slot multi-winners) produce forks.
Fork choice is honest longest-chain with a first-seen tie-break (prefer higher height,
then earlier slot, then lower id) — no adversary is modelled, so the spec's density /
deep-fork rules never engage.
"""
from __future__ import annotations
import heapq
from dataclasses import dataclass
import numpy as np
from .latency import LatencyModel
GENESIS = 0
@dataclass
class BlockTree:
slot: np.ndarray # int64, slot of each block (genesis = -1)
parent: np.ndarray # int64, parent id (genesis = -1)
height: np.ndarray # int64, chain height (genesis = 0)
leader: np.ndarray # int64, producing node id (genesis = -1)
uncles: list[tuple[int, ...]] # referenced uncle ids per block (filled later)
@property
def n_blocks(self) -> int:
return self.slot.shape[0]
def ancestors(self, block_id: int) -> list[int]:
"""Ancestor chain of ``block_id`` from itself down to (excluding) genesis."""
out: list[int] = []
b = block_id
while b > GENESIS:
out.append(b)
b = int(self.parent[b])
return out
def canonical_chain(self) -> list[int]:
"""Honest longest-chain: ancestors of the best tip over the whole tree.
Returns real block ids (genesis excluded), tip-first.
"""
tip = self._best_over_all()
return self.ancestors(tip)
def _rank(self, bid: int) -> tuple[int, int, int]:
# Preference order for "better tip": higher height, earlier slot, lower id.
return (int(self.height[bid]), -int(self.slot[bid]), -bid)
def _best_over_all(self) -> int:
best = GENESIS
best_rank = self._rank(GENESIS)
for bid in range(1, self.n_blocks):
r = self._rank(bid)
if r > best_rank:
best_rank, best = r, bid
return best
def build_tree(
active_slots: np.ndarray,
winners_per_slot: list[np.ndarray],
latency: LatencyModel,
rng: np.random.Generator,
) -> BlockTree:
"""Construct the block tree from grouped lottery winners under a latency model."""
# Preallocate with genesis in slot 0.
slot = [-1]
parent = [-1]
height = [0]
leader = [-1]
# global_best = best publicly-visible tip so far, as (height, slot, id).
def better(a: tuple[int, int, int], b: tuple[int, int, int]) -> tuple[int, int, int]:
# higher height, then earlier slot, then lower id
ah, as_, ai = a
bh, bs, bi = b
if ah != bh:
return a if ah > bh else b
if as_ != bs:
return a if as_ < bs else b
return a if ai < bi else b
global_best = (0, -1, GENESIS)
own_best: dict[int, tuple[int, int, int]] = {}
# min-heap of (visible_at, block_id) awaiting public visibility
pending: list[tuple[int, int]] = []
next_id = 1
for si in range(active_slots.shape[0]):
t = int(active_slots[si])
# advance visibility frontier to slot t
while pending and pending[0][0] <= t:
_, bid = heapq.heappop(pending)
cand = (height[bid], slot[bid], bid)
global_best = better(global_best, cand)
for v in winners_per_slot[si].tolist():
gb = global_best
ob = own_best.get(v, (0, -1, GENESIS))
chosen = better(gb, ob)
p_id = chosen[2]
h = chosen[0] + 1
bid = next_id
next_id += 1
slot.append(t)
parent.append(p_id)
height.append(h)
leader.append(v)
own_best[v] = (h, t, bid)
va = latency.visible_at(t, rng)
heapq.heappush(pending, (va, bid))
return BlockTree(
slot=np.asarray(slot, np.int64),
parent=np.asarray(parent, np.int64),
height=np.asarray(height, np.int64),
leader=np.asarray(leader, np.int64),
uncles=[() for _ in range(next_id)],
)

62
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/concurrency.py Normal file
View File

@ -0,0 +1,62 @@
"""Concurrent block-proposal analysis.
Every lottery win is a block *proposal*. Two proposals produced within ``L`` slots of each
other cannot see one another (a block becomes visible only after the network latency ``L``),
so they are mutually concurrent — competing forks. Bucketing the timeline into
non-overlapping windows of ``L`` slots and counting proposals per bucket gives a direct
view of how many proposals are concurrent, and the busiest bucket is the peak number of
concurrent proposals. This is the quantity that bounds how many uncles can appear, so it
informs the ``MAX_UNCLES`` choice.
"""
from __future__ import annotations
from dataclasses import replace
import numpy as np
from .config import SimConfig
from .lottery import sample_wins, win_probs
from .rng import seedseq_for
from .stake import make_stake
def window_counts(winner_slots: np.ndarray, epoch_len: int, bucket: int) -> np.ndarray:
"""Proposals per non-overlapping ``bucket``-slot window over ``[0, epoch_len)``."""
if epoch_len <= 0:
return np.empty(0, np.int64)
bucket = max(int(bucket), 1)
n_windows = (epoch_len + bucket - 1) // bucket
if winner_slots.size == 0:
return np.zeros(n_windows, np.int64)
return np.bincount(winner_slots // bucket, minlength=n_windows)
def proposal_slots(config: SimConfig, replicate: int = 0) -> np.ndarray:
"""Simulate one epoch of block proposals at the *true* lottery difficulty (D=D_true).
Returns the sorted slots at which proposals (all lottery winners, including forks) occur.
Independent of the TSI trajectory — the proposal process depends only on stake, ``f``,
and ``epoch_len`` — so this is a cheap, self-contained re-simulation for the plots.
"""
cfg = replace(config, replicate=replicate)
root = seedseq_for(cfg)
children = root.spawn(2)
stake = make_stake(cfg, np.random.default_rng(children[0]))
p_win = win_probs(stake, float(stake.sum()), cfg.f)
winner_slots, _ = sample_wins(p_win, cfg.epoch_len, np.random.default_rng(children[1]))
return winner_slots
def concurrency_stats(config: SimConfig, replicate: int = 0) -> dict:
"""Per-bucket proposal-count stats for one simulated epoch (bucket = ``max(L, 1)``)."""
ws = proposal_slots(config, replicate)
bucket = max(config.latency, 1)
counts = window_counts(ws, config.epoch_len, bucket)
return {
"bucket": bucket,
"counts": counts,
"max": int(counts.max()) if counts.size else 0,
"mean": float(counts.mean()) if counts.size else 0.0,
"p99": float(np.percentile(counts, 99)) if counts.size else 0.0,
}

163
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/config.py Normal file
View File

@ -0,0 +1,163 @@
"""Configuration dataclasses for single runs and parameter sweeps."""
from __future__ import annotations
import itertools
from dataclasses import dataclass, field, replace
from typing import Any, Literal
from . import constants
StakeDist = Literal["uniform", "pareto"]
UncleStrategy = Literal["oldest", "random"]
@dataclass(frozen=True)
class SimConfig:
"""A single fully-specified simulation run (one grid cell, one replicate)."""
# --- network / stake ---
n_nodes: int = 1000
stake_dist: StakeDist = "uniform"
pareto_shape: float = 1.16 # Pareto (Lomax) tail index; ~80/20 by default
uniform_random: bool = False # if True, draw i.i.d. uniform stakes; else equal
total_stake: float = 1.0e9 # FIXED across distributions for comparability
# --- network latency (slots) ---
latency: int = 0 # L: block visible to others at t + L
latency_stochastic: bool = False # if True, L is the mean of a stochastic model
# --- uncle references ---
uncle_window: int = constants.W_DEFAULT # W
max_uncles: int = 0 # U (0 = baseline, no uncles)
uncle_strategy: UncleStrategy = "oldest"
# Coin-flip inclusion prob for the "random" strategy. Only 0.5 reproduces the spec's
# unbiased coin (cryptarchia-v1-protocol.md); other values are a deliberate, non-spec
# sensitivity knob, not protocol behaviour.
uncle_random_p: float = 0.5
# --- consensus / TSI ---
f: float = constants.F
beta: float = constants.BETA_DEFAULT
k: int = 64 # scaled by default; full scale = 2160
genesis_d_factor: float = 0.5 # genesis D = factor * true total stake
epochs: int = 40
# If True, mirror the spec's integer fixed-point f-truncation (f_p = int(f*1000)/1000),
# which the on-chain estimator uses; this reproduces its ~1% systematic overestimate.
# Default False keeps the analysis-faithful exact-f behaviour.
fixed_point: bool = False
per_node_dest: bool = False # Phase-2 hook: per-node D_est (unused in reduced model)
# --- performance ---
# >1 parallelises the per-slot lottery across slot-chunks (opt-in; must be pinned and
# recorded because it changes the RNG stream — see lottery.sample_wins_chunked).
lottery_chunks: int = 1
# --- bookkeeping ---
replicate: int = 0
root_seed: int = 12345
def __post_init__(self) -> None:
# frozen dataclass: validation only (no attribute assignment)
if self.stake_dist not in ("uniform", "pareto"):
raise ValueError(f"stake_dist must be uniform|pareto, got {self.stake_dist!r}")
if self.uncle_strategy not in ("oldest", "random"):
raise ValueError(f"uncle_strategy must be oldest|random, got {self.uncle_strategy!r}")
checks = {
"n_nodes": self.n_nodes >= 1,
"k": self.k >= 1,
"epochs": self.epochs >= 1,
"latency": self.latency >= 0,
"max_uncles": self.max_uncles >= 0,
"uncle_window": self.uncle_window >= 1,
"lottery_chunks": self.lottery_chunks >= 1,
"uncle_random_p": 0.0 <= self.uncle_random_p <= 1.0,
"f": 0.0 < self.f < 1.0,
"beta": self.beta > 0.0,
"genesis_d_factor": self.genesis_d_factor > 0.0,
"pareto_shape": self.pareto_shape > 0.0,
"total_stake": self.total_stake > 0.0,
}
bad = [name for name, ok in checks.items() if not ok]
if bad:
raise ValueError(f"invalid SimConfig field(s): {bad}")
# derived geometry -------------------------------------------------------
@property
def epoch_len(self) -> int:
return constants.epoch_len(self.k, self.f)
@property
def period_T(self) -> int:
return constants.period_T(self.k, self.f)
def key(self) -> tuple:
"""Hashable identity used to seed the RNG deterministically.
Must include EVERY field that affects the run (guarded by test_rng), otherwise two
distinct configs would share an RNG stream.
"""
return (
self.n_nodes, self.stake_dist, self.pareto_shape, self.uniform_random,
self.total_stake, self.latency, self.latency_stochastic, self.uncle_window,
self.max_uncles, self.uncle_strategy, self.uncle_random_p, self.f, self.beta,
self.k, self.genesis_d_factor, self.epochs, self.fixed_point, self.per_node_dest,
self.lottery_chunks, self.replicate,
)
@dataclass
class SweepConfig:
"""A cartesian grid of runs plus replicates, all sharing ``base`` settings."""
n_nodes: list[int] = field(default_factory=lambda: [1000])
stake_dist: list[StakeDist] = field(default_factory=lambda: ["uniform", "pareto"])
latency: list[int] = field(default_factory=lambda: [0, 1, 2, 4, 8])
max_uncles: list[int] = field(default_factory=lambda: [0, 1, 2, 3, 4])
uncle_strategy: list[UncleStrategy] = field(default_factory=lambda: ["oldest", "random"])
f: list[float] = field(default_factory=lambda: [constants.F])
replicates: int = 8
base: dict[str, Any] = field(default_factory=dict)
def expand(self) -> list[SimConfig]:
"""Materialise every ``SimConfig`` in the grid × replicates."""
base = SimConfig(**self.base)
cells: list[SimConfig] = []
axes = itertools.product(
self.n_nodes, self.stake_dist, self.latency, self.max_uncles,
self.uncle_strategy, self.f,
)
for n, dist, lat, u, strat, fval in axes:
# U=0 is strategy-independent; keep only one strategy to avoid duplicate work.
if u == 0 and strat != self.uncle_strategy[0]:
continue
for rep in range(self.replicates):
cells.append(
replace(
base,
n_nodes=n,
stake_dist=dist,
latency=lat,
max_uncles=u,
uncle_strategy=strat,
f=fval,
replicate=rep,
)
)
return cells
@classmethod
def from_dict(cls, d: dict[str, Any]) -> SweepConfig:
d = dict(d)
base = d.pop("base", {})
known = {
"n_nodes", "stake_dist", "latency", "max_uncles",
"uncle_strategy", "f", "replicates",
}
unknown = set(d) - known
if unknown:
raise ValueError(
f"unknown sweep keys: {sorted(unknown)} (did you mean one of {sorted(known)}? "
"per-run settings belong under 'base:')"
)
return cls(base=base, **d)

39
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/constants.py Normal file
View File

@ -0,0 +1,39 @@
"""Protocol constants and epoch/window geometry.
All slot geometry derives from the pair ``(k, f)`` so a scaled-down ``k`` (used for
parameter sweeps) automatically shrinks the epoch and measurement window. See
``cryptarchia-v1-protocol.md`` and ``cryptarchia-total-stake-inference.md``.
"""
from __future__ import annotations
# --- True protocol values (full scale) -------------------------------------
K_TRUE = 2160 # security parameter (blocks)
F = 1.0 / 30.0 # slot activation coefficient (fixed; never tuned)
W_DEFAULT = 300 # uncle reference window w_u (slots)
BETA_DEFAULT = 1.0 # TSI learning rate
SLOT_SECONDS = 1 # slot length
def floor_k_over_f(k: int, f: float = F) -> int:
"""``floor(k / f)`` — the base quantum of the epoch schedule."""
return int(k / f)
def epoch_len(k: int, f: float = F) -> int:
"""Epoch length in slots: ``10 * floor(k/f)``."""
return 10 * floor_k_over_f(k, f)
def period_T(k: int, f: float = F) -> int:
"""TSI measurement window length ``T`` in slots: ``6 * floor(k/f)``.
This is the first ``6*floor(k/f)`` slots of the (previous) epoch over which the
block density is measured.
"""
return 6 * floor_k_over_f(k, f)
def expected_blocks_in_window(k: int, f: float = F) -> float:
"""Expected honest-chain block count in the measurement window at equilibrium."""
return period_T(k, f) * f

39
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/engine.py Normal file
View File

@ -0,0 +1,39 @@
"""Multi-epoch trajectory driver for a single config."""
from __future__ import annotations
from typing import Any
import numpy as np
from . import tsi
from .config import SimConfig
from .epoch import simulate_epoch
from .metrics import metric_row
from .rng import seedseq_for
from .stake import make_stake
def run_trajectory(config: SimConfig) -> list[dict[str, Any]]:
"""Run ``config.epochs`` epochs of TSI, returning one metric row per epoch.
Stake is drawn once (it is fixed; only ``D_est`` evolves). ``D_est`` starts at the
hardcoded genesis value ``genesis_d_factor * D_true`` and is updated each epoch from
the measured density. The RNG is a spawn hierarchy off the config's root SeedSequence:
child 0 draws the stake, child ``e+1`` drives epoch ``e`` — so every draw is a
deterministic, order-independent function of the config identity.
"""
root = seedseq_for(config)
children = root.spawn(config.epochs + 1)
stake = make_stake(config, np.random.default_rng(children[0]))
d_true = float(stake.sum())
d_est = config.genesis_d_factor * d_true
T = config.period_T
rows: list[dict[str, Any]] = []
for epoch in range(config.epochs):
er = simulate_epoch(config, stake, d_est, children[epoch + 1])
d_next = tsi.update_D(d_est, er.m, T, config.f, config.beta, config.fixed_point)
rows.append(metric_row(config, epoch, d_est, d_next, d_true, er))
d_est = d_next
return rows

84
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/epoch.py Normal file
View File

@ -0,0 +1,84 @@
"""Single-epoch simulation: lottery -> block tree -> uncles -> density counting."""
from __future__ import annotations
from dataclasses import dataclass
import numpy as np
from . import concurrency, lottery, tsi
from .blocktree import build_tree
from .config import SimConfig
from .latency import make_latency
@dataclass
class EpochResult:
m: int # TSI block count in window
q: float # honest active-slot fraction (window)
q_eff: float # uncle-recovered active-slot fraction (window)
n_active: int # active slots in window
n_honest: int # honest slots in window
n_recovered: int # orphan slots recovered by uncles in window
total_winners_window: int # total lottery wins in window (incl. multi-winner)
n_blocks: int # real blocks produced this epoch
n_canonical: int # canonical chain length
n_orphans: int # orphaned blocks
max_concurrent: int # most block proposals in any latency-sized (max(L,1)) bucket
mean_concurrent: float # mean proposals per latency-sized bucket
def simulate_epoch(
config: SimConfig, stake: np.ndarray, d_est: float, epoch_ss: np.random.SeedSequence
) -> EpochResult:
f = config.f
T = config.period_T
# independent sub-streams: one for the lottery, one for the tree/uncle auxiliary draws
lottery_ss, aux_ss = epoch_ss.spawn(2)
aux_rng = np.random.default_rng(aux_ss)
p_win = lottery.win_probs(stake, d_est, f)
if config.lottery_chunks > 1:
winner_slots, winner_nodes = lottery.sample_wins_chunked(
p_win, config.epoch_len, lottery_ss, config.lottery_chunks
)
else:
winner_slots, winner_nodes = lottery.sample_wins(
p_win, config.epoch_len, np.random.default_rng(lottery_ss)
)
active_slots, groups = lottery.group_by_slot(winner_slots, winner_nodes)
latency = make_latency(config)
tree = build_tree(active_slots, groups, latency, aux_rng)
canonical = tree.canonical_chain()
from .uncles import annotate_uncles
annotate_uncles(tree, canonical, config, aux_rng)
m = tsi.density_m(tree, canonical, T)
ref = tsi.referenced_uncle_ids(tree, canonical)
ss = tsi.slot_stats(tree, canonical, ref, active_slots, T)
total_winners_window = int((winner_slots < T).sum())
n_real = tree.n_blocks - 1
# Concurrent proposals: bucket the whole epoch into latency-sized windows and count
# block proposals (every winner is a proposal) per bucket. The max bucket is the peak
# number of mutually-concurrent proposals (they cannot see each other within L slots).
bucket = max(config.latency, 1)
counts = concurrency.window_counts(winner_slots, config.epoch_len, bucket)
max_concurrent = int(counts.max()) if counts.size else 0
mean_concurrent = float(counts.mean()) if counts.size else 0.0
return EpochResult(
m=m,
q=ss.q,
q_eff=ss.q_eff,
n_active=ss.n_active,
n_honest=ss.n_honest,
n_recovered=ss.n_recovered,
total_winners_window=total_winners_window,
n_blocks=n_real,
n_canonical=len(canonical),
n_orphans=n_real - len(canonical),
max_concurrent=max_concurrent,
mean_concurrent=mean_concurrent,
)

54
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/latency.py Normal file
View File

@ -0,0 +1,54 @@
"""Network latency models.
Latency ``L`` is the number of slots between a block being produced and it becoming
visible to the rest of the network. ``L`` is deliberately named to avoid clashing with
``D`` (the stake estimate). A leader at slot ``t`` can only build on blocks whose
``visible_at <= t`` (its own block is visible to itself immediately), which is what
produces latency-induced forks.
"""
from __future__ import annotations
from typing import Protocol
import numpy as np
class LatencyModel(Protocol):
def visible_at(self, produced_slot: int, rng: np.random.Generator) -> int:
"""Slot at which a block produced at ``produced_slot`` becomes visible to others."""
...
class FixedSlotLatency:
"""Deterministic integer-slot latency: visible to all others at ``t + L``."""
def __init__(self, latency: int) -> None:
self.latency = int(latency)
def visible_at(self, produced_slot: int, rng: np.random.Generator) -> int:
return produced_slot + self.latency
class RealisticLatency:
"""Stochastic latency with mean ``L`` slots (optional sensitivity model).
Rounds an exponential draw (mean ``L``) up to whole slots. A stand-in for the
reference notebook's blend/broadcast delay model; not used by the primary sweep.
"""
def __init__(self, mean_latency: float) -> None:
self.mean_latency = float(mean_latency)
def visible_at(self, produced_slot: int, rng: np.random.Generator) -> int:
if self.mean_latency <= 0:
return produced_slot
draw = rng.exponential(self.mean_latency)
return produced_slot + int(np.ceil(draw))
def make_latency(config) -> LatencyModel: # noqa: ANN001 - avoid import cycle with config
"""Build the latency model for a config."""
if config.latency_stochastic:
return RealisticLatency(config.latency)
return FixedSlotLatency(config.latency)

114
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/lottery.py Normal file
View File

@ -0,0 +1,114 @@
"""Stake-weighted slot lottery (sparse sampler).
Per node ``i`` and slot, an independent Bernoulli win with probability
``phi_f(alpha_i) = 1 - (1 - f)^alpha_i`` where ``alpha_i = w_i / D_est``. Multiple winners
in a slot are possible (a guaranteed fork).
The number of slots a node wins is exactly ``Binomial(n_slots, p_i)``, and the won slots
are a uniformly-random distinct subset — this is *distributionally identical* to drawing an
independent Bernoulli(p_i) in every slot, but avoids materialising the dense
``(n_nodes, n_slots)`` array (which was ~95% of the whole simulator's runtime). Winners are
returned as sparse ``(winner_slots, winner_nodes)`` coordinates, sorted by slot.
"""
from __future__ import annotations
import numpy as np
from joblib import Parallel, delayed
def phi(f: float, alpha: np.ndarray | float) -> np.ndarray | float:
"""Leader-lottery win probability ``1 - (1 - f)^alpha``."""
return 1.0 - (1.0 - f) ** alpha
def win_probs(stake: np.ndarray, d_est: float, f: float) -> np.ndarray:
"""Per-node win probability ``phi_f(w_i / D_est)``."""
return phi(f, stake / d_est)
def _winners_from_counts(
counts: np.ndarray, offset: int, span: int, rng: np.random.Generator
) -> tuple[list[np.ndarray], list[np.ndarray]]:
"""For each winning node, sample ``counts[i]`` distinct slots in ``[offset, offset+span)``."""
nz = np.nonzero(counts)[0]
slot_parts: list[np.ndarray] = []
node_parts: list[np.ndarray] = []
for i in nz:
c = int(counts[i])
slots_i = rng.choice(span, size=c, replace=False).astype(np.int64) + offset
slot_parts.append(slots_i)
node_parts.append(np.full(c, i, np.int64))
return slot_parts, node_parts
def _finalize(
slot_parts: list[np.ndarray], node_parts: list[np.ndarray]
) -> tuple[np.ndarray, np.ndarray]:
if not slot_parts:
return np.empty(0, np.int64), np.empty(0, np.int64)
winner_slots = np.concatenate(slot_parts)
winner_nodes = np.concatenate(node_parts)
order = np.argsort(winner_slots, kind="stable")
return winner_slots[order], winner_nodes[order]
def sample_wins(
p_win: np.ndarray, n_slots: int, rng: np.random.Generator, chunk: int = 8192
) -> tuple[np.ndarray, np.ndarray]:
"""Sample lottery wins over ``n_slots`` slots (sparse; ``chunk`` kept for API compat)."""
counts = rng.binomial(n_slots, p_win)
slot_parts, node_parts = _winners_from_counts(counts, 0, n_slots, rng)
return _finalize(slot_parts, node_parts)
def sample_wins_chunked(
p_win: np.ndarray,
n_slots: int,
seedseq: np.random.SeedSequence,
n_chunks: int,
n_jobs: int = -1,
) -> tuple[np.ndarray, np.ndarray]:
"""Parallel sparse lottery: partition slots into ``n_chunks`` independent ranges.
Correct because ``Binomial(n_slots, p) = sum_c Binomial(L_c, p)`` and per-chunk distinct
subsets are independent. Deterministic given ``(seedseq, n_chunks)`` — but the exact
winner identities differ from the serial sampler and *change with* ``n_chunks``, so
``n_chunks`` must be a pinned, recorded config parameter, never derived from core count.
Note: after the sparse rewrite the lottery is a small fraction of an epoch, so this has
little ROI versus across-config parallelism; it exists for the rare isolated config with
an enormous ``n_slots`` and no across-config work to fill cores.
"""
if n_chunks <= 1:
return sample_wins(p_win, n_slots, np.random.default_rng(seedseq))
bounds = np.linspace(0, n_slots, n_chunks + 1).astype(np.int64)
children = seedseq.spawn(n_chunks)
def one_chunk(c: int) -> tuple[list[np.ndarray], list[np.ndarray]]:
lo, hi = int(bounds[c]), int(bounds[c + 1])
span = hi - lo
rng = np.random.default_rng(children[c])
counts = rng.binomial(span, p_win)
return _winners_from_counts(counts, lo, span, rng)
results = Parallel(n_jobs=n_jobs, prefer="threads")(
delayed(one_chunk)(c) for c in range(n_chunks)
)
slot_parts: list[np.ndarray] = []
node_parts: list[np.ndarray] = []
for sp, npar in results:
slot_parts.extend(sp)
node_parts.extend(npar)
return _finalize(slot_parts, node_parts)
def group_by_slot(
winner_slots: np.ndarray, winner_nodes: np.ndarray
) -> tuple[np.ndarray, list[np.ndarray]]:
"""Group sorted winner coordinates into ``(active_slots, winners_per_active_slot)``."""
if winner_slots.size == 0:
return np.empty(0, np.int64), []
active_slots, starts = np.unique(winner_slots, return_index=True)
groups = np.split(winner_nodes, starts[1:])
return active_slots, groups

76
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/metrics.py Normal file
View File

@ -0,0 +1,76 @@
"""Per-epoch metric rows and equilibrium summaries."""
from __future__ import annotations
from typing import Any
import numpy as np
from .config import SimConfig
from .epoch import EpochResult
# Config fields recorded on every row for grouping/plotting.
_CONFIG_FIELDS = (
"n_nodes", "stake_dist", "pareto_shape", "latency", "uncle_window",
"max_uncles", "uncle_strategy", "uncle_random_p", "f", "beta", "k",
"genesis_d_factor", "epochs", "fixed_point", "replicate",
)
def metric_row(
config: SimConfig,
epoch: int,
d_in: float,
d_out: float,
d_true: float,
er: EpochResult,
) -> dict[str, Any]:
row: dict[str, Any] = {field: getattr(config, field) for field in _CONFIG_FIELDS}
row.update(
epoch=epoch,
d_in=d_in,
d_out=d_out,
d_true=d_true,
ratio=d_out / d_true,
m=er.m,
measured_density=er.m / config.period_T,
q=er.q,
q_eff=er.q_eff,
n_active=er.n_active,
n_honest=er.n_honest,
n_recovered=er.n_recovered,
total_winners_window=er.total_winners_window,
n_blocks=er.n_blocks,
n_canonical=er.n_canonical,
n_orphans=er.n_orphans,
orphan_rate=er.n_orphans / er.n_blocks if er.n_blocks else float("nan"),
max_concurrent=er.max_concurrent,
mean_concurrent=er.mean_concurrent,
)
return row
def equilibrium_stats(ratios: np.ndarray, burn_in: int) -> dict[str, float]:
"""Mean/variance of the stake ratio after ``burn_in`` epochs."""
ratios = np.asarray(ratios, dtype=float)
if ratios.size == 0:
return {"mean_ratio": float("nan"), "var_ratio": float("nan"), "std_ratio": float("nan")}
tail = ratios[burn_in:]
if tail.size == 0:
tail = ratios[-1:]
return {
"mean_ratio": float(np.mean(tail)),
"var_ratio": float(np.var(tail)),
"std_ratio": float(np.std(tail)),
}
def epochs_to_within(ratios: np.ndarray, target: float, eps: float) -> int:
"""First epoch index after which ``|ratio - target| <= eps`` holds for the rest."""
within = np.abs(np.asarray(ratios, dtype=float) - target) <= eps
n = within.size
if n == 0:
return 0
# first index i such that within[i:] are all True == 1 + last index that is False
false_idx = np.flatnonzero(~within)
return int(false_idx[-1] + 1) if false_idx.size else 0

1
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/plotting/__init__.py Normal file
View File

@ -0,0 +1 @@
"""Academic-quality figure generation."""

339
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/plotting/figures.py Normal file
View File

@ -0,0 +1,339 @@
"""Figure builders. Each takes the results DataFrame and returns a Matplotlib Figure.
Rows are per (config, epoch). We summarise each config to its *equilibrium* by averaging
the tail epochs (after a burn-in), keeping replicates so we can draw percentile bands.
"""
from __future__ import annotations
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from .. import concurrency, theory
from ..config import SimConfig
from . import style
CONFIG_COLS = ["n_nodes", "stake_dist", "latency", "max_uncles", "uncle_strategy", "k"]
def _ceiling_line(ax, f: float) -> None:
"""Overlay the intrinsic block-count overshoot ceiling ``-ln(1-f)/f`` (~1.017)."""
c = float(theory.block_count_ceiling(f))
ax.axhline(c, color="0.55", lw=0.9, ls=":", zorder=0)
ax.text(0.99, c, f" full-recovery ceiling {c:.3f}", transform=ax.get_yaxis_transform(),
ha="right", va="bottom", fontsize=7, color="0.4")
def equilibrium(df: pd.DataFrame, burn_frac: float = 0.5) -> pd.DataFrame:
"""Per-(config, replicate) equilibrium means over the tail epochs."""
cutoff = df["epochs"] * burn_frac
tail = df[df["epoch"] >= cutoff]
keys = [*CONFIG_COLS, "replicate"]
return (
tail.groupby(keys, as_index=False)
.agg(
ratio=("ratio", "mean"),
q=("q", "mean"),
q_eff=("q_eff", "mean"),
var_ratio=("ratio", "var"),
orphan_rate=("orphan_rate", "mean"),
n_active=("n_active", "mean"),
)
)
def _series_over_replicates(eq: pd.DataFrame, xcol: str, ycol: str, xvals) -> np.ndarray:
"""Build a ``(n_replicates, len(xvals))`` matrix of ``ycol`` for band plots."""
reps = sorted(eq["replicate"].unique())
mat = np.full((len(reps), len(xvals)), np.nan)
for i, r in enumerate(reps):
sub = eq[eq["replicate"] == r].set_index(xcol)[ycol]
for j, x in enumerate(xvals):
if x in sub.index:
mat[i, j] = sub.loc[x]
return mat
def _provenance(df: pd.DataFrame) -> str:
k = [int(x) for x in sorted(df["k"].unique())]
n = [int(x) for x in sorted(df["n_nodes"].unique())]
reps = int(df["replicate"].nunique())
return f"tsi-sim | k={k} N={n} reps={reps} f={float(df['f'].iloc[0]):.4g}"
# --- Figure 1: accuracy vs U per latency ------------------------------------
def accuracy_vs_u(df: pd.DataFrame, stake_dist: str, strategy: str = "oldest") -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["stake_dist"] == stake_dist) & (eq["n_nodes"] == df["n_nodes"].max())]
eq = eq[(eq["uncle_strategy"] == strategy) | (eq["max_uncles"] == 0)]
latencies = sorted(eq["latency"].unique())
uvals = sorted(eq["max_uncles"].unique())
fig, ax = plt.subplots()
for i, lat in enumerate(latencies):
sub = eq[eq["latency"] == lat]
mat = _series_over_replicates(sub, "max_uncles", "ratio", uvals)
style.band_plot(ax, uvals, mat, color=style.color_for(i), label=f"L={lat}")
ax.axhline(1.0, color="0.4", lw=1.0, ls="--", zorder=0)
_ceiling_line(ax, float(df["f"].iloc[0]))
ax.set_xlabel("max uncles per block $U$")
ax.set_ylabel(r"inferred / true stake $\langle \hat D / D_{\mathrm{true}} \rangle$")
ax.set_title(f"TSI accuracy vs uncle cap ({stake_dist} stake)")
ax.set_xticks(uvals)
ax.legend(title="latency (slots)", ncol=2)
return fig
# --- Figure 2: q_eff vs U per latency ---------------------------------------
def qeff_vs_u(df: pd.DataFrame, stake_dist: str, strategy: str = "oldest") -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["stake_dist"] == stake_dist) & (eq["n_nodes"] == df["n_nodes"].max())]
eq = eq[(eq["uncle_strategy"] == strategy) | (eq["max_uncles"] == 0)]
latencies = sorted(eq["latency"].unique())
uvals = sorted(eq["max_uncles"].unique())
fig, ax = plt.subplots()
for i, lat in enumerate(latencies):
sub = eq[eq["latency"] == lat]
mat = _series_over_replicates(sub, "max_uncles", "q_eff", uvals)
style.band_plot(ax, uvals, mat, color=style.color_for(i), label=f"L={lat}")
base_q = sub[sub["max_uncles"] == 0]["q"].mean()
if not np.isnan(base_q):
ax.axhline(base_q, color=style.color_for(i), lw=0.8, ls=":", alpha=0.6)
ax.axhline(1.0, color="0.4", lw=1.0, ls="--", zorder=0)
ax.set_xlabel("max uncles per block $U$")
ax.set_ylabel(r"effective utilisation $q_{\mathrm{eff}}$")
ax.set_title(f"Active-slot recovery vs uncle cap ({stake_dist} stake)")
ax.set_xticks(uvals)
ax.legend(title="latency (slots)", ncol=2)
return fig
# --- Figure 3: convergence (ratio vs epoch) per U ---------------------------
def convergence(df: pd.DataFrame, stake_dist: str, latency: int,
strategy: str = "oldest") -> plt.Figure:
style.apply_style()
sub = df[(df["stake_dist"] == stake_dist) & (df["latency"] == latency)
& (df["n_nodes"] == df["n_nodes"].max())]
sub = sub[(sub["uncle_strategy"] == strategy) | (sub["max_uncles"] == 0)]
uvals = sorted(sub["max_uncles"].unique())
epochs = sorted(sub["epoch"].unique())
fig, ax = plt.subplots()
for i, u in enumerate(uvals):
su = sub[sub["max_uncles"] == u]
mat = np.full((su["replicate"].nunique(), len(epochs)), np.nan)
for ri, r in enumerate(sorted(su["replicate"].unique())):
s = su[su["replicate"] == r].set_index("epoch")["ratio"]
for j, e in enumerate(epochs):
if e in s.index:
mat[ri, j] = s.loc[e]
style.band_plot(ax, epochs, mat, color=style.color_for(i), label=f"U={u}")
ax.axhline(1.0, color="0.4", lw=1.0, ls="--", zorder=0)
if any(u > 0 for u in uvals):
_ceiling_line(ax, float(sub["f"].iloc[0]))
ax.set_xlabel("epoch")
ax.set_ylabel(r"$\hat D / D_{\mathrm{true}}$")
ax.set_title(f"Convergence ({stake_dist}, latency L={latency})")
ax.legend(title="uncle cap", ncol=2)
return fig
# --- Figure 4: accuracy heatmap over latency x U ----------------------------
def heatmap_accuracy(df: pd.DataFrame, stake_dist: str, n_nodes: int,
strategy: str = "oldest") -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["stake_dist"] == stake_dist) & (eq["n_nodes"] == n_nodes)]
eq = eq[(eq["uncle_strategy"] == strategy) | (eq["max_uncles"] == 0)]
piv = eq.groupby(["latency", "max_uncles"])["ratio"].mean().unstack("max_uncles")
latencies = piv.index.to_numpy()
uvals = piv.columns.to_numpy()
data = piv.to_numpy()
fig, ax = plt.subplots()
vmax = np.nanmax(np.abs(data - 1.0))
im = ax.imshow(data, origin="lower", aspect="auto", cmap=style.DIVERGING_CMAP,
vmin=1 - vmax, vmax=1 + vmax)
ax.set_xticks(range(len(uvals)), uvals)
ax.set_yticks(range(len(latencies)), latencies)
ax.set_xlabel("max uncles per block $U$")
ax.set_ylabel("network latency $L$ (slots)")
ax.set_title(f"Accuracy $\\hat D/D_{{\\mathrm{{true}}}}$ ({stake_dist}, N={n_nodes})")
for yi in range(len(latencies)):
for xi in range(len(uvals)):
v = data[yi, xi]
if not np.isnan(v):
safe = 0.98 <= v <= 1.02
ax.text(xi, yi, f"{v:.2f}", ha="center", va="center", fontsize=7,
color="black", fontweight="bold" if safe else "normal")
fig.colorbar(im, ax=ax, label=r"$\hat D / D_{\mathrm{true}}$")
ax.grid(False)
return fig
# --- Figure 5: distribution comparison --------------------------------------
def dist_comparison(df: pd.DataFrame, latency: int, strategy: str = "oldest") -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["latency"] == latency) & (eq["n_nodes"] == df["n_nodes"].max())]
eq = eq[(eq["uncle_strategy"] == strategy) | (eq["max_uncles"] == 0)]
uvals = sorted(eq["max_uncles"].unique())
fig, ax = plt.subplots()
for i, dist in enumerate(sorted(eq["stake_dist"].unique())):
sub = eq[eq["stake_dist"] == dist]
mat = _series_over_replicates(sub, "max_uncles", "ratio", uvals)
style.band_plot(ax, uvals, mat, color=style.color_for(i), label=dist)
ax.axhline(1.0, color="0.4", lw=1.0, ls="--", zorder=0)
_ceiling_line(ax, float(df["f"].iloc[0]))
ax.set_xlabel("max uncles per block $U$")
ax.set_ylabel(r"$\hat D / D_{\mathrm{true}}$")
ax.set_title(f"Stake-distribution comparison (latency L={latency})")
ax.set_xticks(uvals)
ax.legend(title="distribution")
return fig
# --- Figure 6: variance vs U ------------------------------------------------
def variance_vs_u(df: pd.DataFrame, stake_dist: str, latency: int,
strategy: str = "oldest") -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["stake_dist"] == stake_dist) & (eq["latency"] == latency)
& (eq["n_nodes"] == df["n_nodes"].max())]
eq = eq[(eq["uncle_strategy"] == strategy) | (eq["max_uncles"] == 0)]
uvals = sorted(eq["max_uncles"].unique())
# per-epoch tail variance (mean of each replicate's within-tail variance), matching the
# marginal Var[D/D_true] that theory.variance_ratio denotes -- NOT the variance across
# replicate means (which is ~T smaller).
var_by_u = eq.groupby("max_uncles")["var_ratio"].mean()
qeff_by_u = eq.groupby("max_uncles")["q_eff"].mean()
f = float(df["f"].iloc[0])
T = int(round(6 * float(df["k"].iloc[0]) / f)) # measurement window length
fig, ax = plt.subplots()
ax.plot(uvals, [var_by_u.get(u, np.nan) for u in uvals], "o-",
color=style.color_for(0), label="empirical per-epoch tail variance")
theo = [float(theory.variance_ratio(f, min(qeff_by_u.get(u, np.nan), 1.0), T)) for u in uvals]
ax.plot(uvals, theo, "s--", color=style.color_for(1), label=r"theory $(q_{\mathrm{eff}})$")
ax.set_xlabel("max uncles per block $U$")
ax.set_ylabel(r"$\mathrm{Var}[\hat D / D_{\mathrm{true}}]$")
ax.set_title(f"Estimate variance vs uncle cap ({stake_dist}, L={latency})")
ax.set_xticks(uvals)
ax.legend()
return fig
# --- Figure 7: strategy comparison ------------------------------------------
def strategy_comparison(df: pd.DataFrame, stake_dist: str, latency: int) -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["stake_dist"] == stake_dist) & (eq["latency"] == latency)
& (eq["n_nodes"] == df["n_nodes"].max()) & (eq["max_uncles"] > 0)]
uvals = sorted(eq["max_uncles"].unique())
fig, ax = plt.subplots()
for i, strat in enumerate(sorted(eq["uncle_strategy"].unique())):
sub = eq[eq["uncle_strategy"] == strat]
mat = _series_over_replicates(sub, "max_uncles", "ratio", uvals)
style.band_plot(ax, uvals, mat, color=style.color_for(i), label=strat)
ax.axhline(1.0, color="0.4", lw=1.0, ls="--", zorder=0)
ax.set_xlabel("max uncles per block $U$")
ax.set_ylabel(r"$\hat D / D_{\mathrm{true}}$")
ax.set_title(f"Uncle-selection strategy ({stake_dist}, L={latency})")
ax.set_xticks(uvals)
ax.legend(title="strategy")
return fig
# --- Figure 8: diagnostics (orphan rate vs latency) -------------------------
def orphan_diagnostics(df: pd.DataFrame, stake_dist: str) -> plt.Figure:
style.apply_style()
eq = equilibrium(df)
eq = eq[(eq["stake_dist"] == stake_dist) & (eq["max_uncles"] == 0)
& (eq["n_nodes"] == df["n_nodes"].max())]
g = eq.groupby("latency")["orphan_rate"].agg(["mean", "std"])
fig, ax = plt.subplots()
ax.errorbar(g.index, g["mean"], yerr=g["std"], marker="o", color=style.color_for(1),
capsize=3)
ax.set_xlabel("network latency $L$ (slots)")
ax.set_ylabel("orphan rate (orphans / blocks)")
ax.set_title(f"Fork/orphan rate vs latency ({stake_dist} stake)")
return fig
def _config_from_df(df: pd.DataFrame, stake_dist: str, n_nodes: int, latency: int) -> SimConfig:
"""Reconstruct a single-epoch SimConfig from the swept parameters for re-simulation."""
r = df.iloc[0]
return SimConfig(
n_nodes=int(n_nodes), stake_dist=stake_dist, latency=int(latency),
k=int(r["k"]), pareto_shape=float(r["pareto_shape"]), f=float(r["f"]), epochs=1,
)
# --- Figure 9: peak concurrent block proposals vs latency -------------------
def concurrency_vs_latency(df: pd.DataFrame, stake_dist: str, n_nodes: int,
reps: int = 6) -> plt.Figure:
"""Proposals per latency-sized bucket vs L: the max bucket is the peak concurrency."""
style.apply_style()
latencies = sorted(int(x) for x in df["latency"].unique())
f = float(df["f"].iloc[0])
maxc, p99c, meanc = [], [], []
for lat in latencies:
cfg = _config_from_df(df, stake_dist, n_nodes, lat)
stats = [concurrency.concurrency_stats(cfg, replicate=r) for r in range(reps)]
maxc.append(max(s["max"] for s in stats))
p99c.append(float(np.mean([s["p99"] for s in stats])))
meanc.append(float(np.mean([s["mean"] for s in stats])))
fig, ax = plt.subplots()
ax.plot(latencies, maxc, "o-", color=style.color_for(1), label="max (peak concurrency)")
ax.plot(latencies, p99c, "s--", color=style.color_for(0), label="99th percentile bucket")
ax.plot(latencies, meanc, "^:", color=style.color_for(2), label="mean per bucket")
# expected proposals per bucket ~ L * (-ln(1-f)); reference for the mean
ref = [max(lat, 1) * (-np.log(1 - f)) for lat in latencies]
ax.plot(latencies, ref, color="0.6", lw=0.9, ls="-", zorder=0,
label=r"expected mean $=L\,(-\ln(1-f))$")
ax.set_xlabel("network latency $L$ (slots) = bucket size")
ax.set_ylabel("block proposals per $L$-slot bucket")
ax.set_title(f"Concurrent block proposals per latency-window ({stake_dist}, N={n_nodes})")
ax.legend()
return fig
# --- Figure 10: proposals-per-bucket across time (small multiples) ----------
def concurrency_timeseries(df: pd.DataFrame, stake_dist: str, n_nodes: int,
window_slots: int = 3000) -> plt.Figure:
"""Block proposals per latency-sized bucket across time, one panel per latency.
Each panel steps through the proposals-per-bucket over the first ``window_slots``; the
dashed line marks the peak concurrency observed over the *whole* epoch for that latency.
"""
style.apply_style()
latencies = sorted(int(x) for x in df["latency"].unique() if x > 0)
if len(latencies) > 4:
idx = np.linspace(0, len(latencies) - 1, 4).round().astype(int)
latencies = [latencies[i] for i in idx]
fig, axes = plt.subplots(len(latencies), 1, sharex=True,
figsize=(6.4, 1.5 * len(latencies) + 0.5))
axes = np.atleast_1d(axes)
for i, (lat, ax) in enumerate(zip(latencies, axes, strict=True)):
cfg = _config_from_df(df, stake_dist, n_nodes, lat)
ws = concurrency.proposal_slots(cfg, replicate=0)
epoch_max = int(concurrency.window_counts(ws, cfg.epoch_len, lat).max())
span = min(window_slots, cfg.epoch_len)
counts = concurrency.window_counts(ws[ws < span], span, lat)
centers = (np.arange(counts.size) + 0.5) * lat
ax.step(centers, counts, where="mid", color=style.color_for(i), lw=1.2)
ax.axhline(epoch_max, color="0.45", ls="--", lw=0.9)
ax.text(0.995, 0.92, f"L={lat} · epoch peak = {epoch_max}", transform=ax.transAxes,
ha="right", va="top", fontsize=8)
ax.set_ylim(0, epoch_max + 1)
ax.set_ylabel("proposals")
axes[-1].set_xlabel("slot (time)")
axes[0].set_title(
f"Block proposals across time, bucketed by latency ({stake_dist}, N={n_nodes})"
)
return fig

76
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/plotting/make_figures.py Normal file
View File

@ -0,0 +1,76 @@
"""Render academic figures from a results frame (importable + CLI ``tsi-figures``)."""
from __future__ import annotations
import argparse
import datetime as _dt
from pathlib import Path
import pandas as pd
from . import figures as F
from . import style
def timestamped_figdir(outdir: str | Path, label: str) -> Path:
"""A fresh ``<outdir>/<YYYY-MM-DD_HHMMSS>_<label>`` folder (never overwrites)."""
ts = _dt.datetime.now().strftime("%Y-%m-%d_%H%M%S")
base = Path(outdir) / f"{ts}_{label}"
d, i = base, 2
while d.exists():
d = base.with_name(f"{base.name}_{i}")
i += 1
return d
def render(df: pd.DataFrame, out: str | Path) -> int:
"""Render every applicable figure from ``df`` into directory ``out``; return the count."""
out = Path(out)
prov = F._provenance(df)
dists = sorted(df["stake_dist"].unique())
latencies = sorted(df["latency"].unique())
lat_focus = latencies[len(latencies) // 2] if latencies else 0
written: list[Path] = []
for dist in dists:
written += style.save(F.accuracy_vs_u(df, dist), out / f"01_accuracy_vs_u_{dist}", prov)
written += style.save(F.qeff_vs_u(df, dist), out / f"02_qeff_vs_u_{dist}", prov)
written += style.save(F.convergence(df, dist, lat_focus),
out / f"03_convergence_{dist}_L{lat_focus}", prov)
written += style.save(F.orphan_diagnostics(df, dist), out / f"08_orphans_{dist}", prov)
for n in sorted(df["n_nodes"].unique()):
n = int(n)
written += style.save(F.heatmap_accuracy(df, dist, n),
out / f"04_heatmap_{dist}_N{n}", prov)
written += style.save(F.concurrency_vs_latency(df, dist, n),
out / f"09_concurrency_vs_latency_{dist}_N{n}", prov)
written += style.save(F.concurrency_timeseries(df, dist, n),
out / f"10_concurrency_timeseries_{dist}_N{n}", prov)
written += style.save(F.variance_vs_u(df, dist, lat_focus),
out / f"06_variance_vs_u_{dist}_L{lat_focus}", prov)
if df[(df["stake_dist"] == dist) & (df["max_uncles"] > 0)]["uncle_strategy"].nunique() > 1:
written += style.save(F.strategy_comparison(df, dist, lat_focus),
out / f"07_strategy_{dist}_L{lat_focus}", prov)
written += style.save(
F.dist_comparison(df, lat_focus), out / f"05_dist_comparison_L{lat_focus}", prov
)
return len(written)
def main(argv: list[str] | None = None) -> None:
ap = argparse.ArgumentParser(description="Render TSI figures")
ap.add_argument("--results", required=True)
ap.add_argument("--outdir", default="figures",
help="parent dir; a dated sub-folder is created (never overwrites)")
ap.add_argument("--label", default=None, help="run label (default: results file stem)")
ap.add_argument("--out", default=None, help="explicit output dir (skips the dated folder)")
args = ap.parse_args(argv)
label = args.label or Path(args.results).stem
out = Path(args.out) if args.out else timestamped_figdir(args.outdir, label)
n = render(pd.read_parquet(args.results), out)
print(f"wrote {n} figures -> {out}")
if __name__ == "__main__": # pragma: no cover
main()

92
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/plotting/style.py Normal file
View File

@ -0,0 +1,92 @@
"""Shared academic matplotlib theme, palette, and helpers.
Palette: Okabe-Ito — the standard colorblind-safe qualitative set — for categorical
series (the ``U`` lines); ``cividis`` (perceptually uniform, CVD-safe) for heatmaps.
Figures are saved as vector PDF (primary) plus 300-dpi PNG preview.
"""
from __future__ import annotations
from collections.abc import Sequence
from pathlib import Path
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
# Okabe-Ito colorblind-safe qualitative palette
OKABE_ITO = [
"#0072B2", # blue
"#D55E00", # vermillion
"#009E73", # bluish green
"#CC79A7", # reddish purple
"#E69F00", # orange
"#56B4E9", # sky blue
"#F0E442", # yellow
"#000000", # black
]
SEQUENTIAL_CMAP = "cividis"
# Diverging map centered at ratio = 1 (accuracy heatmaps)
DIVERGING_CMAP = "RdBu_r"
def apply_style() -> None:
"""Install the shared rcParams theme (idempotent)."""
mpl.rcParams.update(
{
"figure.dpi": 150,
"savefig.dpi": 300,
"savefig.bbox": "tight",
"figure.figsize": (6.4, 4.0),
"font.size": 10,
"font.family": "sans-serif",
"axes.titlesize": 11,
"axes.labelsize": 10,
"axes.spines.top": False,
"axes.spines.right": False,
"axes.grid": True,
"grid.alpha": 0.25,
"grid.linewidth": 0.6,
"legend.frameon": False,
"legend.fontsize": 8.5,
"lines.linewidth": 1.8,
"axes.prop_cycle": mpl.cycler(color=OKABE_ITO),
"mathtext.default": "regular",
}
)
def color_for(index: int) -> str:
return OKABE_ITO[index % len(OKABE_ITO)]
def band_plot(
ax: plt.Axes,
x: Sequence[float],
series: np.ndarray,
*,
color: str,
label: str | None = None,
percentiles: tuple[float, float] = (10, 90),
) -> None:
"""Plot the mean of ``series`` (shape ``(n_replicates, len(x))``) with a percentile band."""
x = np.asarray(x, dtype=float)
mean = np.nanmean(series, axis=0)
lo = np.nanpercentile(series, percentiles[0], axis=0)
hi = np.nanpercentile(series, percentiles[1], axis=0)
ax.plot(x, mean, color=color, label=label)
ax.fill_between(x, lo, hi, color=color, alpha=0.18, linewidth=0)
def save(fig: plt.Figure, out_stem: str | Path, provenance: str | None = None) -> list[Path]:
"""Save ``fig`` as a 300-dpi PNG. ``out_stem`` has no suffix. Returns written paths."""
out_stem = Path(out_stem)
out_stem.parent.mkdir(parents=True, exist_ok=True)
if provenance:
# leave headroom so the footnote never collides with the x-axis label
fig.subplots_adjust(bottom=0.16)
fig.text(0.99, 0.005, provenance, fontsize=6, alpha=0.5, va="bottom", ha="right")
p = out_stem.with_suffix(".png")
fig.savefig(p)
plt.close(fig)
return [p]

32
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/rng.py Normal file
View File

@ -0,0 +1,32 @@
"""Deterministic, order-independent RNG derivation.
Each ``SimConfig`` maps to an independent root ``SeedSequence`` seeded from a hash of its
identity plus the global root seed. This guarantees a config yields the same random stream
regardless of the order in which parallel workers run it. The engine ``spawn``\\s children
of this root — one per epoch, and independent sub-streams within an epoch — so every draw
(including the optional parallel chunked lottery) is a deterministic function of the root.
"""
from __future__ import annotations
import hashlib
import numpy as np
from .config import SimConfig
def _entropy(config: SimConfig) -> int:
payload = repr((config.root_seed, config.key())).encode()
digest = hashlib.blake2b(payload, digest_size=16).digest()
return int.from_bytes(digest, "big")
def seedseq_for(config: SimConfig) -> np.random.SeedSequence:
"""Return the reproducible root ``SeedSequence`` for this exact config+replicate."""
return np.random.SeedSequence(_entropy(config))
def rng_for(config: SimConfig) -> np.random.Generator:
"""Return the reproducible ``Generator`` for this exact config+replicate."""
return np.random.default_rng(seedseq_for(config))

32
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/stake.py Normal file
View File

@ -0,0 +1,32 @@
"""Stake distribution generation.
Total stake is held FIXED across distributions (via renormalisation) so that accuracy
comparisons between ``uniform`` and ``pareto`` isolate the *shape* effect on the lottery
(winner multiplicity / forking), not a difference in aggregate stake.
"""
from __future__ import annotations
import numpy as np
from .config import SimConfig
def make_stake(config: SimConfig, rng: np.random.Generator) -> np.ndarray:
"""Return an ``(n_nodes,)`` non-negative stake vector summing to ``total_stake``."""
n = config.n_nodes
if config.stake_dist == "uniform":
if config.uniform_random:
w = rng.random(n)
else:
w = np.ones(n)
elif config.stake_dist == "pareto":
# numpy.pareto draws Lomax = Pareto(shape) - 1, heavy-tailed for small shape.
w = rng.pareto(config.pareto_shape, n) + 1.0
else: # pragma: no cover - guarded by Literal typing
raise ValueError(f"unknown stake_dist: {config.stake_dist}")
total = w.sum()
if total <= 0: # pragma: no cover - degenerate
raise ValueError("stake vector summed to zero")
return w * (config.total_stake / total)

115
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/sweep.py Normal file
View File

@ -0,0 +1,115 @@
"""Parameter-sweep expansion, parallel execution, and result persistence.
Across-config parallelism is the main multicore lever: ``run_trajectory`` is a pure
function of an immutable, hash-seeded ``SimConfig``, so results are bitwise
order-independent and the grid is embarrassingly parallel. We use joblib's process-based
**loky** backend (this workload is CPU-bound Python that holds the GIL, so threads would
serialise), and pin each worker to a single BLAS thread to avoid oversubscription.
"""
from __future__ import annotations
import argparse
import datetime as _dt
import os
from pathlib import Path
import pandas as pd
import yaml
from joblib import Parallel, delayed
from tqdm import tqdm
from .config import SweepConfig
from .engine import run_trajectory
# Keep numpy/BLAS single-threaded inside each worker process (belt-and-braces alongside
# joblib's inner_max_num_threads); prevents N_workers x N_blas_threads oversubscription.
for _var in ("OMP_NUM_THREADS", "OPENBLAS_NUM_THREADS", "MKL_NUM_THREADS", "NUMEXPR_NUM_THREADS"):
os.environ.setdefault(_var, "1")
def run_sweep(
sweep: SweepConfig, n_jobs: int = -1, progress: bool = True, batch_size: str | int = "auto"
) -> pd.DataFrame:
"""Expand the grid, run every config across cores (loky), return one big frame.
``n_jobs=-1`` uses all logical cores. ``batch_size`` defaults to joblib "auto" (good for
many tiny scaled-k tasks); pass ``1`` for a small grid of heavy full-scale configs so
they load-balance rather than pre-batch.
"""
configs = sweep.expand()
runner = Parallel(
n_jobs=n_jobs,
backend="loky",
inner_max_num_threads=1,
batch_size=batch_size,
return_as="generator",
)(delayed(run_trajectory)(c) for c in configs)
if progress:
runner = tqdm(runner, total=len(configs), desc="configs")
rows: list[dict] = []
for traj in runner:
rows.extend(traj)
return pd.DataFrame(rows)
def persist(df: pd.DataFrame, path: str | Path) -> None:
path = Path(path)
path.parent.mkdir(parents=True, exist_ok=True)
if path.suffix == ".parquet":
df.to_parquet(path, index=False)
else:
df.to_csv(path, index=False)
def load_sweep_yaml(path: str | Path) -> SweepConfig:
with open(path) as fh:
data = yaml.safe_load(fh)
return SweepConfig.from_dict(data)
def new_run_dir(outdir: str | Path, label: str) -> Path:
"""Create a fresh timestamped run directory so runs never overwrite each other."""
ts = _dt.datetime.now().strftime("%Y-%m-%d_%H%M%S")
base = Path(outdir) / f"{ts}_{label}"
run_dir, i = base, 2
while run_dir.exists(): # avoid same-second collisions
run_dir = base.with_name(f"{base.name}_{i}")
i += 1
run_dir.mkdir(parents=True)
return run_dir
def main(argv: list[str] | None = None) -> None:
parser = argparse.ArgumentParser(description="Run a TSI parameter sweep + figures")
parser.add_argument("--config", required=True, help="sweep YAML")
parser.add_argument("--outdir", default="runs",
help="parent dir; a dated sub-folder is created per run")
parser.add_argument("--label", default=None, help="run label (default: config name)")
parser.add_argument("--n-jobs", type=int, default=-1, help="-1 = all logical cores")
parser.add_argument("--batch-size", default="auto",
help="'auto' (default) or an int; use 1 for a small heavy grid")
parser.add_argument("--no-figures", action="store_true",
help="skip auto figure generation")
args = parser.parse_args(argv)
batch_size = int(args.batch_size) if args.batch_size != "auto" else "auto"
label = args.label or Path(args.config).stem
run_dir = new_run_dir(args.outdir, label)
sweep = load_sweep_yaml(args.config)
df = run_sweep(sweep, n_jobs=args.n_jobs, batch_size=batch_size)
results_path = run_dir / "results.parquet"
persist(df, results_path)
n_cfg = df[["n_nodes", "stake_dist", "latency", "max_uncles", "uncle_strategy", "replicate"]]
print(f"wrote {len(df)} rows ({len(n_cfg.drop_duplicates())} configs) -> {results_path}")
if not args.no_figures:
from .plotting.make_figures import render
n_fig = render(df, run_dir / "figures")
print(f"wrote {n_fig} figures -> {run_dir / 'figures'}")
print(f"run dir: {run_dir}")
if __name__ == "__main__": # pragma: no cover
main()

60
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/theory.py Normal file
View File

@ -0,0 +1,60 @@
"""Closed-form TSI results from ``analysis-total-stake-inference.md``.
Used both as figure overlays and as ground truth for the verification checks. ``q`` is
the honest active-slot utilisation; with uncle references, substitute the effective
``q_eff`` to predict the improved accuracy.
"""
from __future__ import annotations
import numpy as np
ArrayLike = np.ndarray | float
def expected_ratio(f: float, q: ArrayLike) -> ArrayLike:
"""Equilibrium ``E[D_inf] / D_true = log(1-f) / log(1 - f/q)`` for ``q in (f, 1]``."""
q = np.asarray(q, dtype=float)
return np.log(1.0 - f) / np.log(1.0 - f / q)
def block_count_ceiling(f: float) -> float:
"""Equilibrium ratio at *full* uncle recovery for the equal-stake limit.
TSI counts blocks (all lottery wins, rate ``-ln(1-f)`` per slot in the small-stake
limit), whereas ``f`` is the *active-slot* rate. So even with every orphan recovered the
estimate equilibrates at ``-ln(1-f)/f`` (~1.017 for f=1/30), not 1.0. This is a
deterministic overshoot floor, not noise; concentrated (Pareto) stake gives a smaller
value because of the concavity of ``phi``.
"""
return -np.log(1.0 - f) / f
def fixed_point_bias(f: float, precision: int = 1000) -> float:
"""Extra multiplicative bias from the spec's integer f-truncation ``int(f*P)/P``."""
f_p = int(f * precision) / precision
return f / f_p
def variance_ratio(f: float, q: ArrayLike, T: int, beta: float = 1.0) -> ArrayLike:
"""Equilibrium ``Var[D_inf / D_true]``."""
q = np.asarray(q, dtype=float)
er = expected_ratio(f, q)
return (beta / f) ** 2 * (q / T) * er**2 * (1.0 - f) * f
def variance_bound(f: float, T: int, beta: float = 1.0) -> float:
"""Upper bound on ``Var[D_inf / D_true]`` at ``q = 1`` (perfect network)."""
return (beta / f) ** 2 / T * (1.0 - f) * f
def beta_stability_bound(f: float, q: ArrayLike) -> ArrayLike:
"""Stability threshold: convergence requires ``beta < 2f/((q-f) log(1/(1-f/q)))``."""
q = np.asarray(q, dtype=float)
return 2.0 * f / ((q - f) * np.log(1.0 / (1.0 - f / q)))
def optimal_beta(f: float, q: ArrayLike) -> ArrayLike:
"""Convergence-optimal learning rate ``f/((q-f) log(1/(1-f/q)))``."""
q = np.asarray(q, dtype=float)
return f / ((q - f) * np.log(1.0 / (1.0 - f / q)))

91
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/tsi.py Normal file
View File

@ -0,0 +1,91 @@
"""Total Stake Inference: density counting and the per-epoch estimate update.
The estimate update counts *blocks* exactly as the spec's ``density_over_slots`` does:
``m = honest-chain blocks in window + deduplicated referenced uncle blocks (by their own
slot) in window``. We additionally report slot-based ``q`` / ``q_eff`` (honest and
uncle-recovered active-slot fractions) for comparison against the closed-form theory,
which is written in terms of active-slot utilisation. The two differ only at rare
multi-winner slots.
"""
from __future__ import annotations
from dataclasses import dataclass
import numpy as np
from .blocktree import BlockTree
def referenced_uncle_ids(tree: BlockTree, canonical_ids: list[int]) -> set[int]:
"""Deduplicated set of uncle ids referenced by the canonical chain."""
ref: set[int] = set()
for b in canonical_ids:
ref.update(tree.uncles[b])
return ref
def _in_window(slot: int, T: int) -> bool:
return 0 <= slot < T
def density_m(tree: BlockTree, canonical_ids: list[int], T: int) -> int:
"""Block count ``m`` for the TSI update (honest blocks + deduped uncles, in window)."""
s = tree.slot[canonical_ids]
honest = int(((s >= 0) & (s < T)).sum())
ref = referenced_uncle_ids(tree, canonical_ids)
uncle = sum(1 for u in ref if _in_window(int(tree.slot[u]), T))
return honest + uncle
PRECISION = 1000 # spec on-chain fixed-point scale (cryptarchia-total-stake-inference.md)
def update_D(
d_prev: float, m: int, T: int, f: float, beta: float, fixed_point: bool = False
) -> float:
"""Spec TSI recursion: ``max(1, D_prev * (1 - beta*(f_eff - m/T)/f_eff))``.
With ``fixed_point=True`` the target rate ``f`` is quantised exactly as the on-chain
algorithm does (``f_p = int(f*PRECISION)/PRECISION`` = 0.033 for f=1/30), which drives
the estimate to a measured density of 0.033 instead of 1/30 — a ~1% systematic
overestimate the deployed estimator exhibits but the exact-``f`` float model omits.
The remaining integer divisions in the spec (``tse/PRECISION``) are negligible at
realistic stakes and are not reproduced.
"""
f_eff = (int(f * PRECISION) / PRECISION) if fixed_point else f
measured_density = m / T
d_new = d_prev * (1.0 - beta * (f_eff - measured_density) / f_eff)
return max(d_new, 1.0)
@dataclass
class SlotStats:
n_active: int # active slots (>=1 winner) in window
n_honest: int # honest-chain-occupied slots in window
n_recovered: int # orphan-only slots recovered via referenced uncles
q: float # n_honest / n_active
q_eff: float # (n_honest + n_recovered) / n_active
def slot_stats(
tree: BlockTree,
canonical_ids: list[int],
ref_uncle_ids: set[int],
active_slots: np.ndarray,
T: int,
) -> SlotStats:
"""Slot-based utilisation stats used for theory overlays."""
active_in = active_slots[(active_slots >= 0) & (active_slots < T)]
n_active = int(active_in.size)
honest_slots = {int(tree.slot[b]) for b in canonical_ids if _in_window(int(tree.slot[b]), T)}
recovered: set[int] = set()
for u in ref_uncle_ids:
su = int(tree.slot[u])
if _in_window(su, T) and su not in honest_slots:
recovered.add(su)
n_honest = len(honest_slots)
n_rec = len(recovered)
q = n_honest / n_active if n_active else float("nan")
q_eff = (n_honest + n_rec) / n_active if n_active else float("nan")
return SlotStats(n_active=n_active, n_honest=n_honest, n_recovered=n_rec, q=q, q_eff=q_eff)

82
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/uncles.py Normal file
View File

@ -0,0 +1,82 @@
"""Proposer-local uncle selection.
For each canonical block ``B`` (processed oldest-first so ancestors' references are
known), candidates are orphan (non-canonical) blocks ``U`` with
``0 < slot_B - slot_U <= W`` that have not already been referenced by an ancestor of
``B``. Two strategies match the spec: deterministic oldest-first, and random (oldest-first
order, a coin of probability ``uncle_random_p`` per candidate, capped at ``U``). The spec's
coin is unbiased (``uncle_random_p = 0.5``, the default); other values are a non-spec
sensitivity knob. Dedup across ancestors is enforced by threading a ``referenced`` set down
the canonical chain.
"""
from __future__ import annotations
import numpy as np
from .blocktree import BlockTree
from .config import SimConfig
def _orphans_sorted(tree: BlockTree, canonical_ids: list[int]) -> tuple[np.ndarray, np.ndarray]:
"""Return orphan block ids sorted by (slot, id) and their slots."""
canonical = np.zeros(tree.n_blocks, dtype=bool)
canonical[canonical_ids] = True
all_real = np.arange(1, tree.n_blocks)
orphan_ids = all_real[~canonical[1:]]
orphan_slots = tree.slot[orphan_ids]
order = np.lexsort((orphan_ids, orphan_slots)) # by slot, then id
return orphan_ids[order], orphan_slots[order]
def annotate_uncles(
tree: BlockTree, canonical_ids: list[int], config: SimConfig, rng: np.random.Generator
) -> None:
"""Fill ``tree.uncles[B]`` for every canonical block ``B`` per the selection rule."""
u_max = config.max_uncles
if u_max <= 0:
return
w = config.uncle_window
orphan_ids, orphan_slots = _orphans_sorted(tree, canonical_ids)
if orphan_ids.size == 0:
return
referenced: set[int] = set()
# oldest canonical block first
for b in reversed(canonical_ids):
sb = int(tree.slot[b])
lo = int(np.searchsorted(orphan_slots, sb - w, side="left")) # slot_U >= sb - W
hi = int(np.searchsorted(orphan_slots, sb, side="left")) # slot_U < sb
if hi <= lo:
continue
window_ids = orphan_ids[lo:hi] # already oldest-first
selected = _select(window_ids, referenced, config, rng)
if selected:
tree.uncles[b] = tuple(selected)
referenced.update(selected)
def _select(
window_ids: np.ndarray, referenced: set[int], config: SimConfig, rng: np.random.Generator
) -> list[int]:
u_max = config.max_uncles
out: list[int] = []
if config.uncle_strategy == "oldest":
for bid in window_ids.tolist():
if bid in referenced:
continue
out.append(bid)
if len(out) >= u_max:
break
elif config.uncle_strategy == "random":
p = config.uncle_random_p
for bid in window_ids.tolist():
if bid in referenced:
continue
if rng.random() < p:
out.append(bid)
if len(out) >= u_max:
break
else: # pragma: no cover - guarded by Literal
raise ValueError(config.uncle_strategy)
return out

116
tools/simulators/tsi/tsi-sim-mc/src/tsi_sim/verify.py Normal file
View File

@ -0,0 +1,116 @@
"""Analytic sanity checks: simulator vs closed-form theory (``tsi-verify``).
Replicate runs are evaluated across CPU cores. Exits non-zero if any check fails.
"""
from __future__ import annotations
from dataclasses import replace
import numpy as np
import pandas as pd
from joblib import Parallel, delayed
from .config import SimConfig
from .engine import run_trajectory
from .epoch import simulate_epoch
from .rng import seedseq_for
from .stake import make_stake
from .theory import block_count_ceiling, expected_ratio, fixed_point_bias
F = 1.0 / 30.0
K = 128 # scaled: T = 6*floor(128/f) = 23040 slots
EPOCHS = 45
REPS = 12
BURN = 25
def tail_mean(cfg: SimConfig, col: str, reps: int = REPS, n_jobs: int = -1) -> tuple[float, float]:
"""Mean over ``reps`` replicates of each trajectory's post-burn-in tail mean of ``col``."""
def one(r: int) -> float:
df = pd.DataFrame(run_trajectory(replace(cfg, replicate=r)))
return float(df[col].iloc[BURN:].mean())
vals = np.array(Parallel(n_jobs=n_jobs, backend="loky", inner_max_num_threads=1)(
delayed(one)(r) for r in range(reps)
))
return float(vals.mean()), float(vals.std() / np.sqrt(reps))
def check(name: str, ok: bool, detail: str) -> bool:
print(f"[{'PASS' if ok else 'FAIL'}] {name}: {detail}")
return ok
def main() -> int:
results = []
# 1. Active-slot rate ~= f when D_est = D_true, L=0, U=0.
cfg = SimConfig(n_nodes=2000, stake_dist="uniform", latency=0, max_uncles=0,
k=K, epochs=1, genesis_d_factor=1.0)
ss = seedseq_for(cfg)
children = ss.spawn(2)
stake = make_stake(cfg, np.random.default_rng(children[0]))
er = simulate_epoch(cfg, stake, float(stake.sum()), children[1])
active_rate = er.n_active / cfg.period_T
results.append(check("active-slot rate ~ f (L=0,U=0,D=D_true)",
abs(active_rate - F) / F < 0.05,
f"active_rate={active_rate:.5f} f={F:.5f}"))
# 2. U=0 equilibrium ratio ~= expected_ratio(f, measured q).
cfg = SimConfig(n_nodes=1000, stake_dist="uniform", latency=4, max_uncles=0,
k=K, epochs=EPOCHS, genesis_d_factor=0.5)
ratio, se = tail_mean(cfg, "ratio")
q, _ = tail_mean(cfg, "q")
pred = float(expected_ratio(F, q))
results.append(check("U=0 ratio ~ theory(q)",
abs(ratio - pred) < 0.02 + 2 * se,
f"sim={ratio:.4f}±{se:.4f} theory(q={q:.3f})={pred:.4f}"))
# 3. Underestimate at higher latency (q < 1 => ratio < 1), U=0.
cfg = SimConfig(n_nodes=1000, stake_dist="uniform", latency=10, max_uncles=0,
k=K, epochs=EPOCHS, genesis_d_factor=0.5)
ratio, se = tail_mean(cfg, "ratio")
q, _ = tail_mean(cfg, "q")
results.append(check("U=0 underestimates true stake at latency",
ratio < 0.98 and q < 0.98,
f"ratio={ratio:.4f} q={q:.3f}"))
# 4. q_eff -> 1 and ratio -> block-count ceiling as U grows (uncles recover forks).
# The residual |ratio-1| ~ 0.017 is the intrinsic winners-vs-active-slots floor
# (density_m counts blocks), NOT a convergence failure -- see check 5.
base = dict(n_nodes=1000, stake_dist="uniform", latency=8, uncle_strategy="oldest",
k=K, epochs=EPOCHS, genesis_d_factor=0.5)
r0, _ = tail_mean(SimConfig(max_uncles=0, **base), "ratio")
r4, _ = tail_mean(SimConfig(max_uncles=4, **base), "ratio")
qe4, _ = tail_mean(SimConfig(max_uncles=4, **base), "q_eff")
results.append(check("uncles recover accuracy (q_eff->1, |ratio-1| shrinks)",
qe4 > 0.99 and abs(r4 - 1) < abs(r0 - 1) and abs(r4 - 1) < 0.03,
f"ratio U0={r0:.4f} -> U4={r4:.4f}; q_eff(U4)={qe4:.4f}"))
# 5. Full recovery equilibrates at the block-count ceiling -ln(1-f)/f (~1.017), not 1.
ceiling = float(block_count_ceiling(F))
r4_val, se4 = tail_mean(SimConfig(max_uncles=4, **base), "ratio")
results.append(check("full-recovery ratio ~ block-count ceiling",
abs(r4_val - ceiling) < 0.015 + 2 * se4,
f"ratio(U4)={r4_val:.4f} ceiling=-ln(1-f)/f={ceiling:.4f}"))
# 6. Fixed-point mode adds the spec's ~1% f-truncation overestimate.
fp_base = dict(n_nodes=1000, stake_dist="uniform", latency=0, max_uncles=0,
k=K, epochs=EPOCHS, genesis_d_factor=1.0)
r_float, _ = tail_mean(SimConfig(fixed_point=False, **fp_base), "ratio")
r_fixed, _ = tail_mean(SimConfig(fixed_point=True, **fp_base), "ratio")
bias = float(fixed_point_bias(F))
results.append(check("fixed_point mode adds ~1% f-truncation bias",
r_fixed > r_float and abs(r_fixed / r_float - bias) < 0.01,
f"float={r_float:.4f} fixed={r_fixed:.4f} ratio={r_fixed/r_float:.4f} "
f"expected~{bias:.4f}"))
print()
n_pass = sum(results)
print(f"{n_pass}/{len(results)} checks passed")
return 0 if n_pass == len(results) else 1
if __name__ == "__main__": # pragma: no cover
raise SystemExit(main())

54
tools/simulators/tsi/tsi-sim-mc/tests/test_blocktree.py Normal file
View File

@ -0,0 +1,54 @@
import numpy as np
from tsi_sim.blocktree import build_tree
from tsi_sim.latency import FixedSlotLatency
def _winners(*groups):
return [np.array(g, dtype=np.int64) for g in groups]
def test_latency_induces_fork_and_longest_chain():
# L=2. slot0: node0; slot1: node1 (can't see block1 yet -> forks on genesis);
# slot3: node2 (sees both, builds on the earlier-slot tip -> block1).
active = np.array([0, 1, 3], dtype=np.int64)
winners = _winners([0], [1], [2])
tree = build_tree(active, winners, FixedSlotLatency(2), np.random.default_rng(0))
assert tree.n_blocks == 4 # genesis + 3
assert tree.height.tolist() == [0, 1, 1, 2]
# block 3 built on block 1 (earlier slot wins the height-1 tie)
assert tree.parent[3] == 1
assert tree.canonical_chain() == [3, 1] # tip-first
# block 2 is the orphan
canon = set(tree.canonical_chain())
orphans = [b for b in range(1, tree.n_blocks) if b not in canon]
assert orphans == [2]
def test_same_slot_multiwinner_forks_even_at_zero_latency():
active = np.array([0, 1], dtype=np.int64)
winners = _winners([0, 1], [2]) # two winners in slot 0 -> guaranteed fork
tree = build_tree(active, winners, FixedSlotLatency(0), np.random.default_rng(0))
# blocks 1 and 2 are siblings at height 1 on genesis
assert tree.parent[1] == 0 and tree.parent[2] == 0
assert tree.height[1] == 1 and tree.height[2] == 1
# block 3 at slot 1 extends one of them (height 2)
assert tree.height[3] == 2
assert len(tree.canonical_chain()) == 2
def test_self_extension_within_latency():
# A single node winning consecutive slots builds on its own block despite latency.
active = np.array([0, 1], dtype=np.int64)
winners = _winners([5], [5])
tree = build_tree(active, winners, FixedSlotLatency(10), np.random.default_rng(0))
assert tree.parent[2] == 1 # node 5 self-extends
assert tree.height.tolist() == [0, 1, 2]
def test_ancestors():
active = np.array([0, 1, 2], dtype=np.int64)
winners = _winners([0], [0], [0]) # one node, clean chain
tree = build_tree(active, winners, FixedSlotLatency(0), np.random.default_rng(0))
assert tree.ancestors(3) == [3, 2, 1]

40
tools/simulators/tsi/tsi-sim-mc/tests/test_concurrency.py Normal file
View File

@ -0,0 +1,40 @@
import numpy as np
from tsi_sim import concurrency
from tsi_sim.config import SimConfig
def test_window_counts_partitions_all_proposals():
ws = np.array([0, 1, 2, 5, 5, 9], dtype=np.int64)
counts = concurrency.window_counts(ws, epoch_len=10, bucket=5)
assert counts.tolist() == [3, 3] # slots 0-4 -> 3, slots 5-9 -> 3
assert counts.sum() == ws.size
def test_window_counts_bucket_one_is_per_slot():
ws = np.array([0, 0, 3], dtype=np.int64)
counts = concurrency.window_counts(ws, epoch_len=4, bucket=1)
assert counts.tolist() == [2, 0, 0, 1] # slot 0 has 2 concurrent, slot 3 has 1
def test_window_counts_empty():
assert concurrency.window_counts(np.empty(0, np.int64), 10, 2).tolist() == [0, 0, 0, 0, 0]
def test_concurrency_stats_scale_with_latency():
# Bigger latency bucket => more proposals per bucket (max/mean grow roughly with L).
small = concurrency.concurrency_stats(SimConfig(n_nodes=1000, latency=2, k=32, epochs=1))
large = concurrency.concurrency_stats(SimConfig(n_nodes=1000, latency=20, k=32, epochs=1))
assert large["mean"] > small["mean"]
assert large["max"] >= small["max"]
# mean per bucket ~ bucket * (-ln(1-f))
expected = large["bucket"] * (-np.log(1 - 1 / 30))
assert abs(large["mean"] - expected) < 0.25 * expected
def test_proposal_slots_reproducible_and_sorted():
cfg = SimConfig(n_nodes=500, latency=4, k=16, epochs=1)
a = concurrency.proposal_slots(cfg, replicate=0)
b = concurrency.proposal_slots(cfg, replicate=0)
np.testing.assert_array_equal(a, b)
assert np.all(np.diff(a) >= 0)

81
tools/simulators/tsi/tsi-sim-mc/tests/test_config.py Normal file
View File

@ -0,0 +1,81 @@
import dataclasses
import pytest
from tsi_sim.config import SimConfig, SweepConfig
def test_expand_cardinality_and_u0_strategy_dedup():
sweep = SweepConfig(
n_nodes=[1000, 2000],
stake_dist=["uniform"],
latency=[0, 4],
max_uncles=[0, 1, 2],
uncle_strategy=["oldest", "random"],
replicates=3,
base={"k": 16, "epochs": 5},
)
configs = sweep.expand()
# U=0 keeps only the first strategy; U>0 keeps both.
# per (n,dist,lat): U0 x1 strat + U1 x2 + U2 x2 = 5 strat-U combos, x3 reps = 15
# x 2 n_nodes x 1 dist x 2 lat = 60
assert len(configs) == 60
u0 = [c for c in configs if c.max_uncles == 0]
assert all(c.uncle_strategy == "oldest" for c in u0)
assert {c.replicate for c in configs} == {0, 1, 2}
def test_base_propagation():
sweep = SweepConfig(n_nodes=[500], stake_dist=["pareto"], latency=[2], max_uncles=[0],
uncle_strategy=["oldest"], replicates=1,
base={"k": 32, "epochs": 7, "fixed_point": True})
(c,) = sweep.expand()
assert c.k == 32 and c.epochs == 7 and c.fixed_point is True and c.stake_dist == "pareto"
@pytest.mark.parametrize("kwargs", [
{"k": 0}, {"epochs": 0}, {"n_nodes": 0}, {"latency": -1}, {"max_uncles": -1},
{"uncle_window": 0}, {"lottery_chunks": 0}, {"uncle_random_p": 1.5}, {"f": 0.0},
{"f": 1.0}, {"beta": 0.0}, {"genesis_d_factor": 0.0}, {"pareto_shape": 0.0},
{"stake_dist": "zipf"}, {"uncle_strategy": "newest"},
])
def test_validation_rejects_bad_fields(kwargs):
with pytest.raises(ValueError):
SimConfig(**kwargs)
def test_unknown_sweep_key_rejected():
with pytest.raises(ValueError, match="unknown sweep keys"):
SweepConfig.from_dict({"latencies": [0, 1], "base": {}}) # typo: latencies vs latency
def test_from_dict_roundtrip_ok():
sw = SweepConfig.from_dict({"latency": [0, 3], "max_uncles": [0, 2], "replicates": 2,
"base": {"k": 8, "epochs": 4}})
assert sw.latency == [0, 3] and sw.replicates == 2 and sw.base["k"] == 8
def test_key_covers_every_field():
# Guard against the silent shared-RNG bug: key() must reflect all run-affecting fields.
ignored = {"root_seed"} # root_seed enters _entropy separately, not via key()
names = {f.name for f in dataclasses.fields(SimConfig)} - ignored
a = SimConfig()
for name in names:
cur = getattr(a, name)
alt = _perturb(cur)
b = dataclasses.replace(a, **{name: alt})
assert a.key() != b.key(), f"key() does not distinguish field {name!r}"
def _perturb(v):
if isinstance(v, bool):
return not v
if isinstance(v, int):
return v + 1
if isinstance(v, float):
return v + 0.001
if v == "uniform":
return "pareto"
if v == "oldest":
return "random"
return v

35
tools/simulators/tsi/tsi-sim-mc/tests/test_latency.py Normal file
View File

@ -0,0 +1,35 @@
import numpy as np
from tsi_sim.config import SimConfig
from tsi_sim.engine import run_trajectory
from tsi_sim.latency import FixedSlotLatency, RealisticLatency, make_latency
def test_fixed_slot_latency():
lat = FixedSlotLatency(3)
assert lat.visible_at(10, np.random.default_rng(0)) == 13
def test_realistic_latency_zero_mean_immediate():
lat = RealisticLatency(0.0)
assert lat.visible_at(5, np.random.default_rng(0)) == 5
def test_realistic_latency_positive_delays():
lat = RealisticLatency(4.0)
rng = np.random.default_rng(0)
vals = [lat.visible_at(100, rng) for _ in range(200)]
assert all(v >= 101 for v in vals) # strictly after production
assert np.mean([v - 100 for v in vals]) > 1.0
def test_make_latency_dispatch():
assert isinstance(make_latency(SimConfig(latency=2)), FixedSlotLatency)
assert isinstance(make_latency(SimConfig(latency=2, latency_stochastic=True)), RealisticLatency)
def test_stochastic_latency_integration_runs():
cfg = SimConfig(n_nodes=500, k=8, epochs=5, latency=6, latency_stochastic=True, max_uncles=2)
rows = run_trajectory(cfg)
assert len(rows) == 5
assert all(np.isfinite(r["ratio"]) for r in rows)

77
tools/simulators/tsi/tsi-sim-mc/tests/test_lottery.py Normal file
View File

@ -0,0 +1,77 @@
import numpy as np
from tsi_sim import lottery
def test_phi_bounds():
f = 1 / 30
assert lottery.phi(f, 0.0) == 0.0
# a single all-stake node (alpha=1) wins at exactly rate f
assert abs(lottery.phi(f, 1.0) - f) < 1e-12
def test_win_probs_monotone_in_stake():
stake = np.array([1.0, 2.0, 3.0])
p = lottery.win_probs(stake, d_est=6.0, f=1 / 30)
assert np.all(np.diff(p) > 0)
def test_sample_wins_sorted_and_rate():
rng = np.random.default_rng(0)
n, slots = 500, 4000
p = np.full(n, 0.001)
ws, wn = lottery.sample_wins(p, slots, rng, chunk=512)
assert np.all(np.diff(ws) >= 0) # sorted by slot
assert ws.shape == wn.shape
assert np.all((wn >= 0) & (wn < n))
# expected wins ~ n * slots * p
assert abs(ws.size - n * slots * 0.001) < 4 * np.sqrt(n * slots * 0.001)
def test_group_by_slot():
ws = np.array([0, 0, 2, 5, 5, 5])
wn = np.array([3, 7, 1, 2, 4, 9])
active, groups = lottery.group_by_slot(ws, wn)
assert list(active) == [0, 2, 5]
assert [g.tolist() for g in groups] == [[3, 7], [1], [2, 4, 9]]
def test_group_by_slot_empty():
active, groups = lottery.group_by_slot(np.empty(0, int), np.empty(0, int))
assert active.size == 0 and groups == []
def test_sparse_per_node_distinct_slots():
# Each node wins any slot at most once (independent Bernoulli-per-slot invariant).
rng = np.random.default_rng(1)
p = np.full(300, 0.05)
ws, wn = lottery.sample_wins(p, 2000, rng)
for node in np.unique(wn):
slots = ws[wn == node]
assert slots.size == np.unique(slots).size # no duplicate (node, slot)
def test_sparse_preserves_multiwinner_slots():
# With high p, some slots have >1 distinct winner (guaranteed forks) — must be possible.
rng = np.random.default_rng(2)
p = np.full(50, 0.3)
ws, _ = lottery.sample_wins(p, 500, rng)
_, counts = np.unique(ws, return_counts=True)
assert counts.max() >= 2
def test_chunked_matches_serial_distribution():
# Chunked sampler is deterministic given (seedseq, n_chunks) and statistically matches
# serial. NOTE: SeedSequence.spawn is stateful, so each call needs a FRESH SeedSequence.
p = np.full(400, 0.02)
a_s, a_n = lottery.sample_wins_chunked(p, 100000, np.random.SeedSequence(123), 4, n_jobs=1)
b_s, b_n = lottery.sample_wins_chunked(p, 100000, np.random.SeedSequence(123), 4, n_jobs=1)
np.testing.assert_array_equal(a_s, b_s) # deterministic
np.testing.assert_array_equal(a_n, b_n)
assert np.all(np.diff(a_s) >= 0) # sorted
for node in np.unique(a_n): # per-node distinct slots preserved
s = a_s[a_n == node]
assert s.size == np.unique(s).size
serial_n = lottery.sample_wins(p, 100000, np.random.default_rng(7))[0].size
assert abs(a_s.size - serial_n) < 6 * np.sqrt(serial_n) # same rate

31
tools/simulators/tsi/tsi-sim-mc/tests/test_rng.py Normal file
View File

@ -0,0 +1,31 @@
import numpy as np
from tsi_sim.config import SimConfig
from tsi_sim.engine import run_trajectory
from tsi_sim.rng import rng_for, seedseq_for
def test_seedseq_and_rng_deterministic():
cfg = SimConfig(k=8, epochs=3)
a = np.random.default_rng(seedseq_for(cfg)).random(5)
b = rng_for(cfg).random(5)
np.testing.assert_array_equal(a, b)
def test_distinct_configs_get_distinct_streams():
c0 = SimConfig(k=8, epochs=3, latency=0)
c1 = SimConfig(k=8, epochs=3, latency=1)
assert not np.array_equal(rng_for(c0).random(4), rng_for(c1).random(4))
def test_trajectory_is_order_independent_and_reproducible():
cfg = SimConfig(n_nodes=500, k=8, epochs=6, latency=4, max_uncles=2)
r1 = run_trajectory(cfg)
r2 = run_trajectory(cfg)
assert [row["ratio"] for row in r1] == [row["ratio"] for row in r2]
def test_replicates_differ():
a = run_trajectory(SimConfig(n_nodes=500, k=8, epochs=6, latency=4, replicate=0))
b = run_trajectory(SimConfig(n_nodes=500, k=8, epochs=6, latency=4, replicate=1))
assert a[-1]["ratio"] != b[-1]["ratio"]

40
tools/simulators/tsi/tsi-sim-mc/tests/test_stake.py Normal file
View File

@ -0,0 +1,40 @@
import numpy as np
from tsi_sim.config import SimConfig
from tsi_sim.stake import make_stake
def _rng():
return np.random.default_rng(0)
def test_uniform_equal_and_sum():
cfg = SimConfig(n_nodes=100, stake_dist="uniform", total_stake=1e9)
w = make_stake(cfg, _rng())
assert w.shape == (100,)
assert np.allclose(w, w[0]) # equal weights
assert abs(w.sum() - 1e9) < 1e-3
def test_uniform_random_varies_but_sums():
cfg = SimConfig(n_nodes=200, stake_dist="uniform", uniform_random=True, total_stake=5e8)
w = make_stake(cfg, _rng())
assert w.std() > 0
assert abs(w.sum() - 5e8) < 1e-2
assert np.all(w >= 0)
def test_pareto_sum_fixed_and_heavier_tailed():
n = 5000
uni = make_stake(SimConfig(n_nodes=n, stake_dist="uniform", total_stake=1e9), _rng())
par = make_stake(SimConfig(n_nodes=n, stake_dist="pareto", pareto_shape=1.16,
total_stake=1e9), _rng())
assert abs(par.sum() - 1e9) < 1.0
assert par.max() > uni.max() * 5 # heavy tail: richest holds far more
assert np.all(par >= 0)
def test_total_stake_fixed_across_distributions():
uni = make_stake(SimConfig(n_nodes=1000, stake_dist="uniform", total_stake=7e8), _rng())
par = make_stake(SimConfig(n_nodes=1000, stake_dist="pareto", total_stake=7e8), _rng())
assert abs(uni.sum() - par.sum()) < 1.0 # comparability guarantee

40
tools/simulators/tsi/tsi-sim-mc/tests/test_theory.py Normal file
View File

@ -0,0 +1,40 @@
import numpy as np
from tsi_sim import theory
F = 1 / 30
T = 10000
def test_expected_ratio_unbiased_at_q1():
assert abs(float(theory.expected_ratio(F, 1.0)) - 1.0) < 1e-12
def test_expected_ratio_monotone_in_q():
qs = np.linspace(0.5, 1.0, 20)
er = theory.expected_ratio(F, qs)
assert np.all(np.diff(er) > 0) # accuracy improves as q -> 1
assert np.all(er <= 1.0 + 1e-12) # always an underestimate
def test_variance_bound_matches_at_q1():
v = float(theory.variance_ratio(F, 1.0, T))
assert abs(v - theory.variance_bound(F, T)) < 1e-15
def test_optimal_beta_is_half_stability_bound():
for q in (0.7, 0.85, 0.95):
opt = float(theory.optimal_beta(F, q))
bound = float(theory.beta_stability_bound(F, q))
assert abs(opt - bound / 2) < 1e-12
def test_block_count_ceiling_above_one():
c = theory.block_count_ceiling(F)
assert 1.015 < c < 1.02 # -ln(1-1/30)/(1/30) ~ 1.01705
def test_fixed_point_bias_about_one_percent():
b = theory.fixed_point_bias(F)
assert abs(b - (F / (33 / 1000))) < 1e-12
assert 1.005 < b < 1.02

49
tools/simulators/tsi/tsi-sim-mc/tests/test_theory_convergence.py Normal file
View File

@ -0,0 +1,49 @@
"""End-to-end statistical checks against closed-form theory (scaled k)."""
import numpy as np
import pandas as pd
import pytest
from tsi_sim.config import SimConfig
from tsi_sim.engine import run_trajectory
from tsi_sim.theory import expected_ratio
F = 1 / 30
def _tail_mean(col, **cfg):
reps = cfg.pop("reps", 6)
burn = cfg["epochs"] // 2
vals = []
for r in range(reps):
df = pd.DataFrame(run_trajectory(SimConfig(replicate=r, **cfg)))
vals.append(df[col].iloc[burn:].mean())
return float(np.mean(vals))
def test_baseline_exact_without_forks():
# No latency, no uncles: active-slot rate == f, so the estimate is unbiased.
ratio = _tail_mean("ratio", n_nodes=1000, stake_dist="uniform", latency=0,
max_uncles=0, k=64, epochs=35, genesis_d_factor=0.5, reps=6)
assert abs(ratio - 1.0) < 0.02
@pytest.mark.slow
def test_u0_matches_expected_ratio():
cfg = dict(n_nodes=1000, stake_dist="uniform", latency=6, max_uncles=0,
k=96, epochs=45, genesis_d_factor=0.5)
ratio = _tail_mean("ratio", reps=8, **cfg)
q = _tail_mean("q", reps=8, **cfg)
assert abs(ratio - float(expected_ratio(F, q))) < 0.02
@pytest.mark.slow
def test_uncles_recover_accuracy():
common = dict(n_nodes=1000, stake_dist="uniform", latency=8, uncle_strategy="oldest",
k=96, epochs=45, genesis_d_factor=0.5)
r0 = _tail_mean("ratio", max_uncles=0, reps=8, **common)
r4 = _tail_mean("ratio", max_uncles=4, reps=8, **common)
qe4 = _tail_mean("q_eff", max_uncles=4, reps=8, **common)
assert qe4 > 0.99
assert abs(r4 - 1) < abs(r0 - 1) # uncles reduce the error
assert abs(r4 - 1) < 0.03 # residual is the small block-count overshoot

87
tools/simulators/tsi/tsi-sim-mc/tests/test_tsi_counting.py Normal file
View File

@ -0,0 +1,87 @@
import numpy as np
from tsi_sim.blocktree import BlockTree
from tsi_sim.tsi import density_m, referenced_uncle_ids, slot_stats, update_D
def make_tree(slots, parents, heights, uncles):
n = len(slots)
return BlockTree(
slot=np.array(slots, np.int64),
parent=np.array(parents, np.int64),
height=np.array(heights, np.int64),
leader=np.zeros(n, np.int64),
uncles=uncles,
)
def test_density_counts_honest_plus_deduped_uncles_in_window():
# canonical 1(slot0),4(slot3); orphans 2(slot1),3(slot2). block4 refs uncles 2 and 3.
tree = make_tree(
slots=[-1, 0, 1, 2, 3],
parents=[-1, 0, 0, 0, 1],
heights=[0, 1, 1, 1, 2],
uncles=[(), (), (), (), (2, 3)],
)
canonical = [4, 1] # tip-first
# window T=10 includes all slots
assert density_m(tree, canonical, T=10) == 4 # 2 honest (slots 0,3) + 2 uncles (1,2)
# window T=2 excludes slots 2,3 -> honest {slot0}=1, uncle slot1=1 (slot2 excluded)
assert density_m(tree, canonical, T=2) == 2
def test_uncle_counted_by_own_slot_and_deduped():
tree = make_tree(
slots=[-1, 0, 5, 1],
parents=[-1, 0, 1, 0],
heights=[0, 1, 2, 1],
uncles=[(), (), (3,), ()], # block2 (slot5) references orphan 3 (slot1)
)
canonical = [2, 1]
assert referenced_uncle_ids(tree, canonical) == {3}
# uncle counted by its OWN slot (1), so window T=2 includes it
assert density_m(tree, canonical, T=2) == 2 # honest slot0 + uncle slot1
# window that excludes the uncle's own slot
assert density_m(tree, canonical, T=1) == 1 # only honest slot0
def test_slot_stats_q_and_qeff():
# active slots 0,1,2 in window; honest occupies 0,2; orphan at slot1 recovered by uncle
tree = make_tree(
slots=[-1, 0, 1, 2],
parents=[-1, 0, 0, 1],
heights=[0, 1, 1, 2],
uncles=[(), (), (), (2,)], # block3 refs orphan 2 (slot1)
)
canonical = [3, 1]
active = np.array([0, 1, 2], np.int64)
ref = referenced_uncle_ids(tree, canonical)
ss = slot_stats(tree, canonical, ref, active, T=10)
assert ss.n_active == 3
assert ss.n_honest == 2 # slots 0 and 2
assert ss.n_recovered == 1 # slot 1 recovered via uncle
assert abs(ss.q - 2 / 3) < 1e-9
assert abs(ss.q_eff - 1.0) < 1e-9
def test_update_D_fixed_point():
f, T = 1 / 30, 3000
m = int(round(T * f)) # measured density == f -> D unchanged
assert abs(update_D(1000.0, m, T, f, beta=1.0) - 1000.0) < 1e-6
# measured below f -> estimate drops
assert update_D(1000.0, m - 20, T, f, 1.0) < 1000.0
# clamp at 1
assert update_D(1.0, 0, T, f, 1.0) == 1.0
def test_update_D_fixed_point_mode_targets_truncated_f():
# In fixed-point mode the target rate is int(f*1000)/1000 = 0.033, not 1/30.
f, T = 1 / 30, 10000
m = 340
# For the same measured density, fixed-point (lower target 0.033) raises the estimate
# more than exact-f, i.e. it is systematically higher.
assert (update_D(1000.0, m, T, f, 1.0, fixed_point=True)
> update_D(1000.0, m, T, f, 1.0, fixed_point=False))
# A density of exactly 0.033 is the fixed-point fixed point (estimate unchanged).
m_trunc = int(round(0.033 * T)) # 330
assert abs(update_D(1000.0, m_trunc, T, f, 1.0, fixed_point=True) - 1000.0) < 1e-6

110
tools/simulators/tsi/tsi-sim-mc/tests/test_uncles.py Normal file
View File

@ -0,0 +1,110 @@
import numpy as np
from tsi_sim.blocktree import BlockTree
from tsi_sim.config import SimConfig
from tsi_sim.uncles import annotate_uncles
def make_tree(slots, parents, heights, leaders):
n = len(slots)
return BlockTree(
slot=np.array(slots, np.int64),
parent=np.array(parents, np.int64),
height=np.array(heights, np.int64),
leader=np.array(leaders, np.int64),
uncles=[() for _ in range(n)],
)
def _canonical_and_orphan_tree():
# genesis(0); canonical chain 1(slot0)->3(slot3)->4(slot5); orphan 2(slot1)
tree = make_tree(
slots=[-1, 0, 1, 3, 5],
parents=[-1, 0, 0, 1, 3],
heights=[0, 1, 1, 2, 3],
leaders=[-1, 0, 1, 2, 3],
)
canonical = [4, 3, 1] # tip-first
return tree, canonical
def test_oldest_selection_and_window():
tree, canonical = _canonical_and_orphan_tree()
cfg = SimConfig(max_uncles=1, uncle_window=300, uncle_strategy="oldest")
annotate_uncles(tree, canonical, cfg, np.random.default_rng(0))
# orphan 2 (slot1) is within window of block 3 (slot3) -> referenced there
referenced = {u for b in canonical for u in tree.uncles[b]}
assert referenced == {2}
def test_no_uncles_when_u_zero():
tree, canonical = _canonical_and_orphan_tree()
cfg = SimConfig(max_uncles=0)
annotate_uncles(tree, canonical, cfg, np.random.default_rng(0))
assert all(tree.uncles[b] == () for b in canonical)
def test_window_excludes_out_of_range_orphan():
tree, canonical = _canonical_and_orphan_tree()
cfg = SimConfig(max_uncles=1, uncle_window=1, uncle_strategy="oldest")
annotate_uncles(tree, canonical, cfg, np.random.default_rng(0))
# orphan 2 at slot1; nearest canonical after it is block3 at slot3 -> gap 2 > W=1
referenced = {u for b in canonical for u in tree.uncles[b]}
assert referenced == set()
def _wide_orphan_tree():
# canonical 1(0)->6(6); orphans 2,3,4,5 at slots 1,2,3,4 (all within window of block6)
tree = make_tree(
slots=[-1, 0, 1, 2, 3, 4, 6],
parents=[-1, 0, 0, 0, 0, 0, 1],
heights=[0, 1, 1, 1, 1, 1, 2],
leaders=[-1, 0, 1, 2, 3, 4, 0],
)
return tree, [6, 1] # tip-first
def test_random_strategy_deterministic_and_capped():
import numpy as np
tree_a, canon = _wide_orphan_tree()
tree_b, _ = _wide_orphan_tree()
cfg = SimConfig(max_uncles=2, uncle_window=300, uncle_strategy="random", uncle_random_p=0.5)
annotate_uncles(tree_a, canon, cfg, np.random.default_rng(7))
annotate_uncles(tree_b, canon, cfg, np.random.default_rng(7))
assert tree_a.uncles == tree_b.uncles # same seed -> identical
total = sum(len(tree_a.uncles[b]) for b in canon)
assert total <= cfg.max_uncles # capped
def test_random_p_one_matches_oldest():
import numpy as np
tree_r, canon = _wide_orphan_tree()
tree_o, _ = _wide_orphan_tree()
annotate_uncles(tree_r, canon, SimConfig(max_uncles=2, uncle_strategy="random",
uncle_random_p=1.0), np.random.default_rng(1))
annotate_uncles(tree_o, canon, SimConfig(max_uncles=2, uncle_strategy="oldest"),
np.random.default_rng(1))
assert tree_r.uncles == tree_o.uncles # p=1 deterministically takes oldest-first
def test_random_p_zero_selects_nothing():
import numpy as np
tree, canon = _wide_orphan_tree()
annotate_uncles(tree, canon, SimConfig(max_uncles=4, uncle_strategy="random",
uncle_random_p=0.0), np.random.default_rng(1))
assert all(tree.uncles[b] == () for b in canon)
def test_dedup_across_ancestors():
# Two canonical blocks both within window of the single orphan: only one references it.
tree = make_tree(
slots=[-1, 0, 1, 2, 3],
parents=[-1, 0, 0, 1, 3],
heights=[0, 1, 1, 2, 3],
leaders=[-1, 0, 9, 2, 3],
)
canonical = [4, 3, 1]
cfg = SimConfig(max_uncles=4, uncle_window=300, uncle_strategy="oldest")
annotate_uncles(tree, canonical, cfg, np.random.default_rng(0))
counts = sum(len(tree.uncles[b]) for b in canonical)
assert counts == 1 # orphan 2 referenced exactly once despite two eligible blocks