[examples/docs] GPU two-stage reliable at tol=1e-4 (deterministic Schur assembly)

The deterministic GPU Schur assembly (MadNLP) makes scopf_twostage_model converge reliably at
the default tol=1e-4 (100%, matching CPU obj 5352.103) — no retry, no loosened tolerance.

- examples/scopf.jl: drop the single-as-oracle reference retry from the two-stage and compare
  modes (the GPU :single path is itself unreliable now and is no longer a valid oracle); the
  two-stage solve is direct. --retry stays as a safety net (default 1). Recommended:
    --mode twostage --gpu --inertia based
- CLAUDE.md / SCOPF-GPU-TWOSTAGE-FINDINGS.md: record the resolution (atomic-assembly
  nondeterminism, the fix, 47/47) and the open GPU :single caveat (same class, needs the
  analogous cuda_sparse.jl fix).

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LpUgCcm5tg8uHi6QaA4eqP

Author: michel2323

examples/scopf.jl

@@ -17,10 +17,9 @@
 #   julia --project=./examples examples/scopf.jl --case case9 --mode twostage --form rect
 #   julia --project=./examples examples/scopf.jl --gpu --inertia free
 #
-# GPU two-stage (Schur) converges to the true optimum on case9 with this config — it is a
-# per-run roulette (~60%), so retry; see SCOPF-GPU-TWOSTAGE-FINDINGS.md for why:
-#   julia --project=./examples examples/scopf.jl --case case9 --mode twostage --gpu \
-#         --inertia based --tol 5e-3 --richardson 50 --retry 6
+# GPU two-stage (Schur) converges reliably to the true optimum on case9 at the default
+# tol=1e-4 (deterministic Schur assembly; see SCOPF-GPU-TWOSTAGE-FINDINGS.md):
+#   julia --project=./examples examples/scopf.jl --case case9 --mode twostage --gpu --inertia based
 
 using ExaModelsPower
 using ExaModels: solution
@@ -178,19 +177,10 @@ elseif mode == :twostage
     lin = backend === nothing ? MadNLP.MumpsSolver : MadNLPGPU.CUDSSSolver
     kkt_opts = schur_kkt_options(info, backend, cudss_ir)
 
-    # On GPU the Schur solve can converge to a suboptimal local point, so use the reliable
-    # :single path as a reference objective; the retry then accepts only the true optimum.
-    refobj = nothing
-    if backend !== nothing
-        mref, vref, _ = scopf_model(case, contingencies; form = form, backend = backend)
-        rref = madnlp(mref; tol = 1.0e-4, print_level = MadNLP.ERROR,
-                      kkt_system = MadNLP.SparseCondensedKKTSystem, linear_solver = MadNLPGPU.CUDSSSolver,
-                      inertia_correction_method = inertia, max_iter = max_iter)
-        rref.status == MadNLP.SOLVE_SUCCEEDED && (refobj = rref.objective)
-        @info "single-path reference objective (retry target)" refobj
-    end
-
-    result = solve_with_retry(; ref = refobj) do
+    # With the deterministic GPU Schur assembly the two-stage solve converges reliably at
+    # tol=1e-4 (matching CPU); --retry defaults to 1 (a single solve) and is kept only as a
+    # safety net.
+    result = solve_with_retry() do
         madnlp(model;
             callback = MadNLP.SparseCallback,
             kkt_system = MadNLP.SchurComplementKKTSystem,
@@ -222,20 +212,16 @@ elseif mode == :compare
     # Sparse first-stage Schur complement solver: MUMPS on CPU, cuDSS on GPU. The
     # per-scenario sparse blocks use cuDSS/MUMPS internally.
     lin = backend === nothing ? MadNLP.MumpsSolver : MadNLPGPU.CUDSSSolver
-    # Retry the (GPU-flaky) two-stage solve against the :single objective so the comparison
-    # is over the true optimum, not a suboptimal local point.
-    r2 = solve_with_retry(; ref = r1.objective) do
-        madnlp(m2;
-            callback = MadNLP.SparseCallback,
-            kkt_system = MadNLP.SchurComplementKKTSystem,
-            linear_solver = lin,
-            kkt_options = schur_kkt_options(info, backend, cudss_ir),
-            inertia_correction_method = inertia,
-            max_iter = max_iter,
-            tol = tol, print_level = MadNLP.ERROR,
-            rich_kw()...,
-        )
-    end
+    r2 = madnlp(m2;
+        callback = MadNLP.SparseCallback,
+        kkt_system = MadNLP.SchurComplementKKTSystem,
+        linear_solver = lin,
+        kkt_options = schur_kkt_options(info, backend, cudss_ir),
+        inertia_correction_method = inertia,
+        max_iter = max_iter,
+        tol = tol, print_level = MadNLP.ERROR,
+        rich_kw()...,
+    )
     pg0 = Array(solution(r2, v2.pg0))
 
     obj_gap = abs(r1.objective - r2.objective)