Design: warm-session daemon

design/daemon-measurements/README.md

@@ -0,0 +1,19 @@
+# Cold-start measurements behind `../daemon.md`
+
+Three standalone scripts, each run as a **fresh process** from the repo's
+venv on the group's GPU box (2× Quadro GV100, warm Julia depot, package at
+`d3b7a22`). `raw_stack.py` and `user_level.py` were run twice each;
+phase-timing spread between runs was under 5%.
+
+- `raw_stack.py` — the fixed stack cost, no model: Julia boot, each `using`,
+  first CUDA operation. ≈ 60 s to a live CUDA context.
+- `user_level.py` — user-visible phases for Lukšan–Vlček N=100 000 on
+  `backend="cuda"`: cold 244 s to first solution, warm re-solve 0.43 s.
+  The legs after the warm re-solve exercise `set_parameters` and rebuilds;
+  the `set_parameters` leg currently dies with `ModelError: Scalar indexing
+  is disallowed` — the device bug described in the design doc.
+- `opf.py` — AC OPF realism (needs pglib, `PGLIB_DIR` or `~/git/pglib-opf`):
+  `case1354_pegase` cold 308 s / warm re-solve 0.71 s, then
+  `case2869_pegase` in the same process: 8.5 s total despite being a
+  different fingerprint — JIT reuse is by expression type, which is what a
+  warm daemon inherits.

design/daemon-measurements/opf.py

@@ -0,0 +1,51 @@
+"""P0 scenario C: a realistic model — AC OPF from pglib, backend="cuda".
+Fresh process. case1354_pegase cold, then case2869_pegase in the SAME process:
+the second case has a different fingerprint (different index sets/data) but
+identical expression types, so it measures what a warm daemon pays for a new
+case of a known model family — the reuse that matters for real workflows."""
+import pathlib
+import sys
+import time
+
+T0 = time.perf_counter()
+LAST = [T0]
+
+
+def mark(label):
+    now = time.perf_counter()
+    print(f"{label:44s} {now - LAST[0]:8.2f}s  (cum {now - T0:7.2f}s)", flush=True)
+    LAST[0] = now
+
+
+sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[2] / "examples"))
+import matpower  # noqa: E402
+
+import examodels as exa  # noqa: E402,F401
+from ac_opf import ac_opf  # noqa: E402
+
+mark("imports (pure python)")
+
+PGLIB = pathlib.Path.home() / "git/pglib-opf"
+
+data = matpower.read(str(PGLIB / "pglib_opf_case1354_pegase.m"))
+mark("read case1354_pegase")
+
+core, _var = ac_opf(data, backend="cuda")
+mark("trace + Core (Julia boot lands here)")
+
+m = exa.Model(core)
+mark("Model(core)")
+
+s1 = m.solve()
+mark(f"solve #1 cold ({s1.status})")
+
+s2 = m.solve()
+mark(f"solve #2 warm ({s2.status})")
+
+data2 = matpower.read(str(PGLIB / "pglib_opf_case2869_pegase.m"))
+core2, _ = ac_opf(data2, backend="cuda")
+m2 = exa.Model(core2)
+mark("case2869: read + trace + Model")
+
+s3 = m2.solve()
+mark(f"case2869 solve, warm process ({s3.status})")

design/daemon-measurements/raw_stack.py

@@ -0,0 +1,32 @@
+"""P0 scenario A: raw cold-stack decomposition, independent of examodels' surface.
+Fresh process, warm Julia depot, timing each phase once — boot/load/JIT phases
+are one-shot by nature; run-to-run spread is checked by running the script twice."""
+import time
+
+T0 = time.perf_counter()
+LAST = [T0]
+
+
+def mark(label):
+    now = time.perf_counter()
+    print(f"{label:38s} {now - LAST[0]:8.2f}s  (cum {now - T0:7.2f}s)", flush=True)
+    LAST[0] = now
+
+
+from juliacall import Main  # noqa: E402
+
+mark("import juliacall (Julia boot)")
+Main.seval("using ExaModels")
+mark("using ExaModels")
+Main.seval("using NLPModelsIpopt")
+mark("using NLPModelsIpopt")
+Main.seval("using MadNLP")
+mark("using MadNLP")
+Main.seval("using CUDA")
+mark("using CUDA")
+Main.seval("using CUDSS, MadNLPGPU")
+mark("using CUDSS, MadNLPGPU")
+Main.seval("sum(CUDA.zeros(Float64, 8))")
+mark("first CUDA op (context init)")
+Main.seval("sum(CUDA.zeros(Float64, 8))")
+mark("second CUDA op (warm)")

design/daemon-measurements/user_level.py

@@ -0,0 +1,79 @@
+"""P0 scenario B: user-visible cold breakdown + what a warm daemon could serve.
+Fresh process. Model: Luksan-Vlcek at N=100_000 on backend="cuda" (the size
+docs/tutorial/gpu.md quotes 0.75s warm for), objective carrying a parameter
+block so the set_parameters leg is measurable.
+
+Warm-process legs at the end:
+  - re-solve            = live-instance hit, same data
+  - set_parameters+solve = live-instance hit, new parameter values
+  - full rebuild+solve  = replay-per-solve cost, same expression types
+  - new-structure model = marginal cost of an unseen fingerprint in a warm
+                          daemon (iterations capped: JIT overhead is the datum,
+                          not solve quality)
+"""
+import time
+
+T0 = time.perf_counter()
+LAST = [T0]
+
+
+def mark(label):
+    now = time.perf_counter()
+    print(f"{label:38s} {now - LAST[0]:8.2f}s  (cum {now - T0:7.2f}s)", flush=True)
+    LAST[0] = now
+
+
+import examodels as exa  # noqa: E402
+
+mark("import examodels")
+
+N = 100_000
+start = [-1.2 if i % 2 == 0 else 1.0 for i in range(N)]
+
+core = exa.Core(backend="cuda")
+mark("Core(backend='cuda')")
+
+x = core.add_var(N, start=start)
+p = core.add_par([1.0] * N)
+core.add_obj(lambda i: 100 * (x[i - 1] ** 2 - x[i]) ** 2 + (x[i - 1] - p[i]) ** 2,
+             over=range(1, N))
+core.add_con(lambda i: 3 * x[i + 1] ** 3 + 2 * x[i + 2] - 5
+             + exa.sin(x[i + 1] - x[i + 2]) * exa.sin(x[i + 1] + x[i + 2])
+             + 4 * x[i + 1] - x[i] * exa.exp(x[i] - x[i + 1]) - 3,
+             over=range(0, N - 2))
+mark("trace expressions")
+
+m = exa.Model(core)
+mark("Model(core)")
+
+s1 = m.solve()
+mark(f"solve #1 cold ({s1.status})")
+
+s2 = m.solve()
+mark(f"solve #2 warm ({s2.status})")
+
+m.set_parameters(p, [1.001] * N)
+s3 = m.solve()
+mark(f"set_parameters + solve ({s3.status})")
+
+core2 = exa.Core(backend="cuda")
+x2 = core2.add_var(N, start=start)
+p2 = core2.add_par([1.0] * N)
+core2.add_obj(lambda i: 100 * (x2[i - 1] ** 2 - x2[i]) ** 2 + (x2[i - 1] - p2[i]) ** 2,
+              over=range(1, N))
+core2.add_con(lambda i: 3 * x2[i + 1] ** 3 + 2 * x2[i + 2] - 5
+              + exa.sin(x2[i + 1] - x2[i + 2]) * exa.sin(x2[i + 1] + x2[i + 2])
+              + 4 * x2[i + 1] - x2[i] * exa.exp(x2[i] - x2[i + 1]) - 3,
+              over=range(0, N - 2))
+m2 = exa.Model(core2)
+s4 = m2.solve()
+mark(f"rebuild same structure + solve ({s4.status})")
+
+core3 = exa.Core(backend="cuda")
+x3 = core3.add_var(N, start=1.0)
+core3.add_obj(lambda i: (x3[i] - 2) ** 4 + exa.cos(x3[i]), over=range(N))
+core3.add_con(lambda i: exa.tanh(x3[i]) + x3[i + 1] ** 2, over=range(N - 1),
+              lcon=-10.0, ucon=10.0)
+m3 = exa.Model(core3)
+s5 = m3.solve(max_iter=3)
+mark(f"NEW structure build + solve, max_iter=3 ({s5.status})")