diff --git a/.grok/skills/unsga3-oracle/SKILL.md b/.grok/skills/unsga3-oracle/SKILL.md index cade722..a29b05f 100644 --- a/.grok/skills/unsga3-oracle/SKILL.md +++ b/.grok/skills/unsga3-oracle/SKILL.md @@ -68,8 +68,18 @@ Checks sphere residual, niche occupancy, coverage vs true PF. Known historical r 3. If algorithm code changed: `dotnet test Unsga3.slnx -c Release` must stay green. 4. Do **not** commit `tools/oracle/out/*.csv` (gitignored). +### 4. Constrained problems and mating opt-ins + +```powershell +python tools/oracle/run_new_surfaces.py --seeds 15 +``` + +Writes `docs/NEW-SURFACES-RESULTS.md`. Does **not** write `docs/WILCOXON-RESULTS.md`. C# defaults stay `RankNicheDistance` + `IndependentWithReplacement`. C1-DTLZ1 passes pymoo `n_var=7` (k=5). Every cell is measured or `skip:`. + ## Do not - Claim bit-identical fronts (different RNG streams). - Tighten CI IGD bars without multi-seed evidence. - Run full n=15 when the user only asked for a smoke check. +- Invent IGD, HV, or Wilcoxon numbers. A missing front or a seed with no feasible point is `skip:`. +- Rewrite `docs/WILCOXON-RESULTS.md` from `run_new_surfaces.py` or from a mating-opt-in sweep. diff --git a/CHANGELOG.md b/CHANGELOG.md index dff1c8e..bd5e717 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,6 +9,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/). ### Added +- Measured 15-seed tables for OSY, TNK, and C1-DTLZ1 (feasible non-dominated IGD vs pymoo `UNSGA3`) and for mating opt-ins (`TwoShuffledPasses`, `Algorithm2`, and both) in `docs/NEW-SURFACES-RESULTS.md`. Protocol is in `docs/EQUIVALENCE.md`. Package defaults are unchanged. `docs/WILCOXON-RESULTS.md` is not rewritten. - Optional **TypeSafe / Jev System One** pass (`tools/typesafe-pareto`) for Score + Choice judgments over a small Pareto candidate sample. Additive semantic layer — does not replace NSGA-III / U-NSGA-III objectives. Live calls only when `TYPESAFE_API_KEY` is set; CI uses mocks. Metrics append to `metrics/typesafe-runs.jsonl`. ### Changed diff --git a/README.md b/README.md index 25a6bc1..d76baf6 100644 --- a/README.md +++ b/README.md @@ -8,7 +8,7 @@ > **v0.2.0** — opt-in two-shuffle mating, Algorithm 2 tournament, odd-N pairing, feasible-only niching, equal-CV Deb fix, two-layer Das–Dennis, and OSY/TNK/C1-DTLZ1. Defaults stay unconstrained `RankNicheDistance`, independent mating, and single-layer Das–Dennis. The published 15-seed table is not rewritten. > **15-seed IGD vs pymoo `UNSGA3`:** ZDT1 **median 0.053 vs 0.070** (MWU *p*≈0.05) compares the full C# non-dominated front with pymoo `res.F`. DTLZ2 **median 0.0045 vs 0.0028** (~1.6×) compares C# n_var=12 with pymoo's default n_var=10. Neither pair is a same-set, same-problem ranking. Notes: [`docs/ORACLE-RESULTS.md`](docs/ORACLE-RESULTS.md). -> Details: [`docs/WILCOXON-RESULTS.md`](docs/WILCOXON-RESULTS.md) · single-seed notes: [`docs/ORACLE-RESULTS.md`](docs/ORACLE-RESULTS.md) +> Details: [`docs/WILCOXON-RESULTS.md`](docs/WILCOXON-RESULTS.md) · single-seed notes: [`docs/ORACLE-RESULTS.md`](docs/ORACLE-RESULTS.md) · constrained + mating opt-ins: [`docs/NEW-SURFACES-RESULTS.md`](docs/NEW-SURFACES-RESULTS.md) ```text https://github.com/AppSprout-dev/Unsga3 @@ -123,6 +123,7 @@ API docs used: [HTTP](https://docs.typesafe.ai/api.md) · [Python SDK](https://d |-----|----------| | [docs/ORACLE-RESULTS.md](docs/ORACLE-RESULTS.md) | Single-seed C# vs pymoo | | [docs/WILCOXON-RESULTS.md](docs/WILCOXON-RESULTS.md) | Multi-seed Mann–Whitney / Wilcoxon | +| [docs/NEW-SURFACES-RESULTS.md](docs/NEW-SURFACES-RESULTS.md) | 15-seed constrained IGD and mating opt-ins | | [docs/EQUIVALENCE.md](docs/EQUIVALENCE.md) | Protocol & intentional deltas (ZDT2 quality A/B = gens=250, `PymooCompatible`; matches unsga3-bend) | | [docs/RESEARCH-STANDARDS.md](docs/RESEARCH-STANDARDS.md) | Literature + indicator standards | | [docs/NOTICE.md](docs/NOTICE.md) | Attribution (papers + validation tools) | diff --git a/docs/EQUIVALENCE.md b/docs/EQUIVALENCE.md index d36cdae..2cab62a 100644 --- a/docs/EQUIVALENCE.md +++ b/docs/EQUIVALENCE.md @@ -45,13 +45,13 @@ ZDT2 **gens=100** is an early-stress snapshot (collapse on Bend, C#, and pymoo), | Bi | ZDT1–4, ZDT6 | 2 | yes | | Many | DTLZ1–4, DTLZ7 | 3+ | yes | | Hard | WFG1, WFG2, WFG9 | 3–5 | planned | -| Constrained | OSY, TNK, C1-DTLZ1 | 2–3 | yes (self-tests; no IGD table) | +| Constrained | OSY, TNK, C1-DTLZ1 | 2–3 | yes (self-tests; 15-seed IGD in [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md)) | ## Unit checks - Reference association & niche counts - Non-dominated ranks / constraint domination (equal CV is mutual non-domination) -- Feasible niching, CV fill, and OSY / TNK / C1-DTLZ1 (`ConstrainedSurvivalTests`, `ConstrainedProblemTests`). No IGD table. +- Feasible niching, CV fill, and OSY / TNK / C1-DTLZ1 (`ConstrainedSurvivalTests`, `ConstrainedProblemTests`). Multi-seed IGD is in [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md). - **Normalization intercepts / ASF axis extremes** (`NormalizationTests`) - Tournament pressure - `populationSize: null` ⇒ `|refs|` @@ -95,6 +95,39 @@ pymoo 0.6.2 `TournamentSelection` (`pressure=2`) **is** the two-shuffled-pass po Odd N (Das–Dennis M=2, p=12 → 13; M=3, p=12 → 91): `CreateOffspring` pairs the leftover parent with parent 0, so `parents[N-1]` enters SBX, then the cursor restarts at 0. Even N stays consecutive pairs `(0,1), (2,3), …` and wraps to `(0,1)` with no extra random draw. Published ZDT1 (52) and DTLZ2 (92) are even, so their pair sequence is unchanged. +## New surfaces (constrained + mating opt-ins) + +Measured tables live in [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md). They are written by `tools/oracle/run_new_surfaces.py`. [WILCOXON-RESULTS.md](WILCOXON-RESULTS.md) is not a target of that script. Package defaults stay `RankNicheDistance` and `IndependentWithReplacement`. + +Seeds are 1..15, the same indexes as the published ZDT1 and DTLZ2 tables. C# and pymoo do not share an RNG. Pairing is by seed index. Operators: SBX η=30, p_c=1, PM η=20, p_m=1/n, duplicate elimination on, single-layer Das–Dennis. + +### Constrained + +| Problem | M | n_var | Partitions | Pop | Gens | C# settings | +|---------|---|------:|------------|----:|-----:|-------------| +| OSY | 2 | 6 | 12 | 52 | 250 | constructor defaults | +| TNK | 2 | 2 | 12 | 52 | 250 | constructor defaults | +| C1-DTLZ1 | 3 | **7** (k=5) | 12 | 92 | 150 | constructor defaults | + +Pop and partitions follow the ZDT (M=2, pop=52) and DTLZ2 (M=3, pop=92) oracles. OSY and TNK use the ZDT2 quality generation count (250). C1-DTLZ1 uses the DTLZ2 generation count (150). `C1Dtlz1Problem` defaults to k=5. pymoo `get_problem("c1dtlz1")` defaults to n_var=12; the harness passes **n_var=7**. + +IGD is pymoo's mean nearest distance on the **feasible non-dominated subset of the final population**, both solvers, one shared reference front. OSY and TNK use pymoo's shipped front file. C1-DTLZ1 uses the DTLZ1 simplex at the same Das–Dennis partitions (the C1 inequality does not remove that front). `res.F` is a different set and is not the comparison. A seed with no feasible point is `skip:`, not a filled-in IGD. + +pymoo `UNSGA3` still uses its own tournament and two-shuffle pool. This is default Unsga3 against default pymoo, on a shared front. It is not a claim that the mating operators match. A second pymoo run of the same seed is not always bit-identical with SciPy OpenBLAS. The results file keeps the first run in the table and lists the second run when it differs. The two values are not averaged. + +### Mating opt-ins + +Same populations and generation counts as the quality protocols above. ZDT1 gens=100, ZDT2 gens=250, DTLZ2 gens=150. IGD is the C# feasible non-dominated front against `ParetoFronts` (the published ZDT1 C# definition). + +| Column | Tournament | Pool | +|--------|------------|------| +| default | `RankNicheDistance` | `IndependentWithReplacement` | +| two_shuffled | `RankNicheDistance` | `TwoShuffledPasses` | +| algorithm2 | `Algorithm2` | `IndependentWithReplacement` | +| both | `Algorithm2` | `TwoShuffledPasses` | + +The ZDT2 and DTLZ2 default columns in that table are `RankNicheDistance`. The published DTLZ2 Wilcoxon column is `PymooCompatible`. ZDT2 has no published Wilcoxon column. Neither published table is rewritten. + ## DTLZ2 gap history | Stage | C# pymoo-mode IGD | vs pymoo 0.0035 | diff --git a/docs/NEW-SURFACES-RESULTS.md b/docs/NEW-SURFACES-RESULTS.md new file mode 100644 index 0000000..1e434d0 --- /dev/null +++ b/docs/NEW-SURFACES-RESULTS.md @@ -0,0 +1,292 @@ +# New surfaces: constrained problems and mating opt-ins + +Measured 2026-09-30. pymoo 0.6.2. .NET SDK 10.0.401. + +Generated by `tools/oracle/run_new_surfaces.py`. Every numeric cell is that run. A cell that was not measured is `skip:`. This file does not replace [`WILCOXON-RESULTS.md`](WILCOXON-RESULTS.md). + +## Command + +```bash +python tools/oracle/run_new_surfaces.py +``` + +C# side, one config, all seeds in one process: + +```bash +dotnet run --project tools/OracleCompare -c Release -- \ + --problem zdt1 --partitions 12 --pop 52 --gens 100 --seed 1 --seed-count 15 \ + --tournament rank-niche --mating independent +dotnet run --project tools/OracleCompare -c Release -- \ + --problem osy --partitions 12 --pop 52 --gens 250 --seed 1 --seed-count 15 +dotnet run --project tools/OracleCompare -c Release -- \ + --problem c1dtlz1 --partitions 12 --pop 92 --gens 150 --seed 1 --seed-count 15 +``` + +Omitted `--gens` on `osy` and `tnk` is 250. Omitted `--gens` on `c1dtlz1` is 150. Omitted `--tournament` is `RankNicheDistance`. Omitted `--mating` is `IndependentWithReplacement`. + +## Protocol + +Operators on both sides: SBX η=30, p_c=1, PM η=20, p_m=1/n, duplicate elimination on, single-layer Das–Dennis. Seeds are `1 .. n` (published ZDT1 / DTLZ2 indexes). C# and pymoo do not share an RNG stream. Pairing is by seed index. + +### Constrained + +| Problem | M | n_var | Partitions | Refs | Pop | Gens | C# tournament | C# mating | pymoo | +|---------|---|------:|------------|-----:|----:|-----:|---------------|-----------|-------| +| osy | 2 | 6 | 12 | 13 | 52 | 250 | RankNicheDistance | IndependentWithReplacement | UNSGA3 defaults | +| tnk | 2 | 2 | 12 | 13 | 52 | 250 | RankNicheDistance | IndependentWithReplacement | UNSGA3 defaults | +| c1dtlz1 | 3 | 7 | 12 | 91 | 92 | 150 | RankNicheDistance | IndependentWithReplacement | UNSGA3 defaults | + +Scored set for the comparison: feasible non-dominated subset of the **final population**, both solvers, one `pymoo.indicators.igd.IGD` call, one shared reference front. OSY and TNK use pymoo's shipped Pareto-front file (`pareto_front()`). C1-DTLZ1 uses `pareto_front(das-dennis partitions=12)` at **n_var=7** (k=5). That front is the DTLZ1 simplex (0.5 × the directions). The C1 inequality does not remove it. pymoo's default `get_problem("c1dtlz1")` is n_var=12 and is not this protocol. + +C# constructor defaults are not pymoo's mating. pymoo `UNSGA3` uses `comp_by_rank_and_ref_line_dist` and two shuffled tournament passes. The table is default Unsga3 against default pymoo UNSGA3 on a shared front definition. It is not a mating-matched ablation. The mating section is that ablation, on unconstrained problems. Each pymoo seed is run twice. The table is the first run. A second value is listed only when it differs. The two runs are not averaged. + +`pymoo res.F` is reported in its own column. It is the survival result, not the feasible non-dominated population, and it is not an input to the hypothesis tests. + +A seed with no feasible point on either side is `skip: no feasible points` and is left out of the median and the tests. Tests run only when at least five seeds are paired. + +### Mating opt-ins + +| Problem | Partitions | Pop | Gens | Reference front | Seeds | +|---------|------------|----:|-----:|-----------------|-------| +| zdt1 | 12 | 52 | 100 | library sampler | 1..15 | +| zdt2 | 12 | 52 | 250 | library sampler | 1..15 | +| dtlz2 | 12 | 92 | 150 | library sampler | 1..15 | + +IGD is the C# feasible non-dominated front against `ParetoFronts` (ZDT: 500 points; DTLZ2: Das–Dennis sphere at the same partitions). That is the definition of the published ZDT1 C# column. + +| Column | Tournament | Mating pool | +|--------|------------|-------------| +| default | RankNicheDistance | IndependentWithReplacement | +| two_shuffled | RankNicheDistance | TwoShuffledPasses | +| algorithm2 | Algorithm2 | IndependentWithReplacement | +| both | Algorithm2 | TwoShuffledPasses | + +ZDT2 uses gens=250. Its default column is `RankNicheDistance`, not `PymooCompatible`. There is no published ZDT2 Wilcoxon column. DTLZ2's default column is also `RankNicheDistance`. The published DTLZ2 column is `PymooCompatible` and is not rewritten here. + +Paired tests compare each opt-in with the default column (same seed index). Lower IGD is better. p-values are two-sided and not multiplicity-adjusted. + +## Constrained results + +### OSY + +- C#: n=15 median=8.60511977766423 mean=15.00603704731623 std=18.277259721680554 min=3.1977789885300525 max=78.81046575608664 +- pymoo feasible ND: n=15 median=14.695651997216531 mean=16.390897264845076 std=18.16684552002373 min=2.8772635546031213 max=78.3785704501295 +- Unpaired skips: 0 +- Reference front: n=99 min=[-273.84271673317005, 4.474099944179967] max=[-71.50302980207687, 75.81677090600435] + +| Test | Statistic | p-value | Verdict (α=0.05) | +|------|-----------|---------|------------------| +| Mann–Whitney U (U₁) | 100.0 | 0.6186695570822622 | fail to reject H₀ | +| Wilcoxon signed-rank | 50.0 | 0.5894979267172403 | fail to reject H₀ | + +- Median ratio (C# / pymoo), feasible ND: **0.5855554948697823** +- Better median: **Unsga3** +- n paired: **15** + +| Seed | C# feasible ND IGD | pymoo feasible ND IGD | pymoo res.F IGD | C# n | pymoo n | +|------|--------------------|-----------------------|-----------------|-----:|--------:| +| 1 | 7.393868356452636 | 19.873296739502415 | 20.457323760937363 | 52 | 52 | +| 2 | 19.104206102940815 | 14.695651997216531 | 15.631623678589158 | 52 | 52 | +| 3 | 7.2144869425846885 | 4.100397940273672 | 5.265802898388595 | 52 | 52 | +| 4 | 7.277854060433869 | 4.037204777945615 | 5.413195115963867 | 52 | 52 | +| 5 | 13.829327385395237 | 13.600364067513409 | 14.880338946532621 | 52 | 52 | +| 6 | 4.262666692258077 | 16.853821590938367 | 17.764567157386704 | 52 | 52 | +| 7 | 8.476602751515053 | 15.058042326917555 | 15.869660108196387 | 52 | 52 | +| 8 | 12.693839026382387 | 7.639536880805481 | 9.28309052265498 | 52 | 52 | +| 9 | 16.770967865700595 | 6.8701731181178625 | 8.389030370673947 | 52 | 52 | +| 10 | 16.859101274258233 | 18.70594737453411 | 19.84032979491332 | 52 | 52 | +| 11 | 3.1977789885300525 | 7.438037364425758 | 9.285283780356588 | 52 | 52 | +| 12 | 78.81046575608664 | 78.3785704501295 | 78.5261490731193 | 52 | 52 | +| 13 | 7.599973079682027 | 15.649896295864115 | 16.64522732647078 | 52 | 52 | +| 14 | 8.60511977766423 | 20.085254493888623 | 20.807186832192865 | 52 | 52 | +| 15 | 12.994297649858929 | 2.8772635546031213 | 5.480216887119926 | 52 | 52 | + +A second pymoo run of the same seed is not always bit-identical on this host (SciPy OpenBLAS, dynamic threads). The table cell is the first run. The second run differed for: +- seed 11 feasible ND: table 7.438037364425758, repeat 8.076576440054877; res.F: table 9.285283780356588, repeat 9.212851572401725 +The two runs are not averaged. Hypothesis tests use the table column. + +### TNK + +- C#: n=15 median=0.01096188711215301 mean=0.011118198573506589 std=0.001136409134697483 min=0.00885014717358177 max=0.013056494475681189 +- pymoo feasible ND: n=15 median=0.010597315398511284 mean=0.010934792695545665 std=0.0012548513842829893 min=0.009235158580058828 max=0.013635571499417975 +- Unpaired skips: 0 +- Reference front: n=104 min=[0.0421200313409779, 0.042139109050040514] max=[1.0382306834095432, 1.0382573364803416] + +| Test | Statistic | p-value | Verdict (α=0.05) | +|------|-----------|---------|------------------| +| Mann–Whitney U (U₁) | 129.0 | 0.506914933812496 | fail to reject H₀ | +| Wilcoxon signed-rank | 51.0 | 0.6292607875039373 | fail to reject H₀ | + +- Median ratio (C# / pymoo), feasible ND: **1.0344022707574545** +- Better median: **pymoo** +- n paired: **15** + +| Seed | C# feasible ND IGD | pymoo feasible ND IGD | pymoo res.F IGD | C# n | pymoo n | +|------|--------------------|-----------------------|-----------------|-----:|--------:| +| 1 | 0.00885014717358177 | 0.011745010131922605 | 0.03073598150149818 | 51 | 52 | +| 2 | 0.010714260853631669 | 0.00993246903675877 | 0.031473870797655935 | 52 | 52 | +| 3 | 0.010086882018114977 | 0.010346608934930241 | 0.029738849086012136 | 52 | 52 | +| 4 | 0.011849173306602404 | 0.012291260966643193 | 0.032651987970428724 | 52 | 52 | +| 5 | 0.010681320317344064 | 0.010217866086475881 | 0.03355602313854264 | 52 | 50 | +| 6 | 0.01096188711215301 | 0.012178757428066601 | 0.0310280506119408 | 52 | 52 | +| 7 | 0.012808627329843603 | 0.009838092718677301 | 0.031211880509497254 | 52 | 52 | +| 8 | 0.010935333916166332 | 0.009626535399212209 | 0.030421854840813743 | 52 | 52 | +| 9 | 0.011612955208595 | 0.011339085025311144 | 0.03237259904985194 | 52 | 52 | +| 10 | 0.01217629024160554 | 0.013635571499417975 | 0.030011821324008717 | 52 | 52 | +| 11 | 0.011155615353057518 | 0.010597315398511284 | 0.03109513668033898 | 52 | 52 | +| 12 | 0.013056494475681189 | 0.010021185833450245 | 0.030798243938314294 | 52 | 52 | +| 13 | 0.009500655338888067 | 0.010711898966898985 | 0.032204458363750896 | 52 | 51 | +| 14 | 0.01079258160183822 | 0.009235158580058828 | 0.030087970674572504 | 51 | 52 | +| 15 | 0.011590754355495458 | 0.01230507442684969 | 0.05583771456799684 | 52 | 51 | + +A second pymoo run of each seed matched the table cell (feasible ND and res.F). + +### C1DTLZ1 + +- C#: n=15 median=0.05850096220717 mean=0.10236604544269495 std=0.1062343831726953 min=0.018648262198023267 max=0.36360493779311853 +- pymoo feasible ND: n=15 median=0.06849687222996582 mean=0.08767356938876994 std=0.06407846396051092 min=0.030555609123652613 max=0.25580120269113255 +- Unpaired skips: 0 +- Reference front: n=91 min=[0.0, 0.0, 0.0] max=[0.5, 0.5, 0.5] + +| Test | Statistic | p-value | Verdict (α=0.05) | +|------|-----------|---------|------------------| +| Mann–Whitney U (U₁) | 105.0 | 0.7715511878155722 | fail to reject H₀ | +| Wilcoxon signed-rank | 59.0 | 0.9773446461818073 | fail to reject H₀ | + +- Median ratio (C# / pymoo), feasible ND: **0.8540676428372325** +- Better median: **Unsga3** +- n paired: **15** + +| Seed | C# feasible ND IGD | pymoo feasible ND IGD | pymoo res.F IGD | C# n | pymoo n | +|------|--------------------|-----------------------|-----------------|-----:|--------:| +| 1 | 0.046769488306550026 | 0.10884312699518343 | 0.10978632739039286 | 92 | 92 | +| 2 | 0.06631521343953564 | 0.08057941812059227 | 0.08124534293602614 | 92 | 92 | +| 3 | 0.08626629056523412 | 0.25580120269113255 | 0.2558479934707872 | 92 | 92 | +| 4 | 0.16078833256353278 | 0.09501388981900726 | 0.09652886014516196 | 92 | 92 | +| 5 | 0.03524944210378766 | 0.06465991475755106 | 0.06523079941668194 | 92 | 92 | +| 6 | 0.041427433679500754 | 0.08582115315029125 | 0.08623866061236912 | 92 | 92 | +| 7 | 0.36360493779311853 | 0.09440230787423456 | 0.09589951305982343 | 14 | 92 | +| 8 | 0.02996549372296548 | 0.06029455189147133 | 0.06044303051452107 | 92 | 92 | +| 9 | 0.33063366847721587 | 0.0450795236953868 | 0.04759689097769537 | 92 | 92 | +| 10 | 0.11647500054209174 | 0.030555609123652613 | 0.032307343244459136 | 92 | 92 | +| 11 | 0.03822735177855429 | 0.03476111340226754 | 0.035590039194817916 | 92 | 92 | +| 12 | 0.09311394335572692 | 0.06849687222996582 | 0.07032117093816809 | 92 | 92 | +| 13 | 0.04950486090741722 | 0.20898780381298027 | 0.2090523986557476 | 92 | 92 | +| 14 | 0.05850096220717 | 0.0401023681803815 | 0.040324598632333036 | 92 | 92 | +| 15 | 0.018648262198023267 | 0.041704685087450656 | 0.04181883070364761 | 92 | 92 | + +A second pymoo run of the same seed is not always bit-identical on this host (SciPy OpenBLAS, dynamic threads). The table cell is the first run. The second run differed for: +- seed 14 feasible ND: table 0.0401023681803815, repeat 0.040102390539722504; res.F: table 0.040324598632333036, repeat 0.04032461993491003 +The two runs are not averaged. Hypothesis tests use the table column. + +## Mating opt-in results + +### ZDT1 + +| Config | n | median | mean | std | min | max | +|--------|--:|--------|------|-----|-----|-----| +| default | 15 | 0.05317298862666853 | 0.06221177985400482 | 0.022563211247080415 | 0.03775116533704069 | 0.10566779337111303 | +| two_shuffled | 15 | 0.04131806757185226 | 0.052214014939089524 | 0.024819303989865225 | 0.029066657506991816 | 0.10484173675149552 | +| algorithm2 | 15 | 0.0821518523540782 | 0.09080800483928475 | 0.02895823610889682 | 0.05128490163910914 | 0.15608487161055737 | +| both | 15 | 0.06231785035827368 | 0.07663242028848027 | 0.03807979130308405 | 0.04092846717129279 | 0.19111712541912093 | + +| Opt-in vs default | U₁ | MWU p | W | signed-rank p | median ratio (opt-in / default) | better median | +|-------------------|----|-------|---|----------------|---------------------------------|---------------| +| two_shuffled | 72.0 | 0.09709106759647713 | 37.0 | 0.20127970635570083 | 0.7770499390573936 | two_shuffled | +| algorithm2 | 177.0 | 0.007940336264246587 | 16.0 | 0.013487378597561894 | 1.5449921938913835 | default | +| both | 143.0 | 0.21337370709311845 | 30.0 | 0.09383863903411002 | 1.1719832186942116 | default | + +Published G6 check (default column vs docs/WILCOXON-RESULTS.md ZDT1 Unsga3 column): 15/15 match. The published file was not rewritten. + +| Seed | default | two_shuffled | algorithm2 | both | +|------|------|------|------|------| +| 1 | 0.051430749249856716 | 0.04597518639666193 | 0.09773312575437405 | 0.04092846717129279 | +| 2 | 0.060915075620609 | 0.09331932820152633 | 0.15608487161055737 | 0.06041016403017146 | +| 3 | 0.03775116533704069 | 0.04377755386806775 | 0.07478299329094622 | 0.08695099439793498 | +| 4 | 0.064032991910825 | 0.0880420897255572 | 0.1302724497287653 | 0.08968887162152264 | +| 5 | 0.0735681028410708 | 0.046320181668926146 | 0.10195045782325228 | 0.10517063708489827 | +| 6 | 0.10566779337111303 | 0.029066657506991816 | 0.0821518523540782 | 0.19111712541912093 | +| 7 | 0.10090480098935474 | 0.029403749788930678 | 0.06012295441798301 | 0.04627993573127827 | +| 8 | 0.08367968202294207 | 0.038411391600916917 | 0.11288330471205056 | 0.06838610141334237 | +| 9 | 0.08449622088420745 | 0.07325143777689323 | 0.07258320653829962 | 0.0880376933745096 | +| 10 | 0.049751355380425216 | 0.04131806757185226 | 0.07224559481125653 | 0.10292146729419119 | +| 11 | 0.03906442386490016 | 0.040518398743594475 | 0.06687808552167199 | 0.04431527314273339 | +| 12 | 0.05317298862666853 | 0.03609276165561767 | 0.06732278311341969 | 0.057090444801660735 | +| 13 | 0.038506968538297465 | 0.03233330571475704 | 0.05128490163910914 | 0.05324698420859879 | +| 14 | 0.0422284046353372 | 0.10484173675149552 | 0.10339438588163138 | 0.0526242942776751 | +| 15 | 0.04800597453742421 | 0.04053837711455391 | 0.11242910539187595 | 0.06231785035827368 | + +### ZDT2 + +| Config | n | median | mean | std | min | max | +|--------|--:|--------|------|-----|-----|-----| +| default | 15 | 0.020011241536364334 | 0.019551964653898097 | 0.0018114813062077356 | 0.015230938620873397 | 0.022052340426814396 | +| two_shuffled | 15 | 0.019565881988858867 | 0.019428639231120187 | 0.0016837246541701134 | 0.017074868100335514 | 0.022785591788352672 | +| algorithm2 | 15 | 0.018683811009861915 | 0.021538811349568917 | 0.007743083122715551 | 0.014221988974628922 | 0.04278986447488385 | +| both | 15 | 0.016708931230469345 | 0.02189179932000593 | 0.01730144669591226 | 0.014234999363259494 | 0.08364571898741646 | + +| Opt-in vs default | U₁ | MWU p | W | signed-rank p | median ratio (opt-in / default) | better median | +|-------------------|----|-------|---|----------------|---------------------------------|---------------| +| two_shuffled | 99.0 | 0.5897385933186072 | 56.0 | 0.8424296456125092 | 0.9777445319074201 | two_shuffled | +| algorithm2 | 96.0 | 0.506914933812496 | 58.0 | 0.932106941814807 | 0.933665758614216 | algorithm2 | +| both | 67.0 | 0.06196961991403876 | 41.0 | 0.2933828921517897 | 0.8349772401730274 | both | + +| Seed | default | two_shuffled | algorithm2 | both | +|------|------|------|------|------| +| 1 | 0.018749140703627317 | 0.022785591788352672 | 0.019230070908833854 | 0.01847638538169906 | +| 2 | 0.020011241536364334 | 0.01748536238959677 | 0.017793818608281095 | 0.020725937466819992 | +| 3 | 0.021445915731687564 | 0.01933175149093086 | 0.018683811009861915 | 0.08364571898741646 | +| 4 | 0.018344177623392074 | 0.017074868100335514 | 0.03370223537989431 | 0.016457346675429668 | +| 5 | 0.02097496219036023 | 0.01760943595027763 | 0.029098969876935516 | 0.015455616319794403 | +| 6 | 0.01898988127499011 | 0.019823722687072557 | 0.021193018598070632 | 0.024331026160292524 | +| 7 | 0.015230938620873397 | 0.019975335735137555 | 0.01662155580886381 | 0.020874698703884125 | +| 8 | 0.019409230894157355 | 0.019641226287095873 | 0.021164351260310018 | 0.01602322673379817 | +| 9 | 0.020327198964233743 | 0.019565881988858867 | 0.017039290245269283 | 0.01787489778724865 | +| 10 | 0.02056490642434661 | 0.021003204778074235 | 0.01842873649990813 | 0.016708931230469345 | +| 11 | 0.021271217388244102 | 0.01810857387463844 | 0.04278986447488385 | 0.015934019953539255 | +| 12 | 0.01762216792112181 | 0.019251634790507018 | 0.016716426727679132 | 0.014234999363259494 | +| 13 | 0.022052340426814396 | 0.02116415280684568 | 0.014221988974628922 | 0.017909653026130328 | +| 14 | 0.02057721205568473 | 0.017308251955835963 | 0.017277309176413735 | 0.014828880289716394 | +| 15 | 0.01770893805257369 | 0.021300593843243145 | 0.019120722693699577 | 0.014895651720591088 | + +### DTLZ2 + +| Config | n | median | mean | std | min | max | +|--------|--:|--------|------|-----|-----|-----| +| default | 15 | 0.006992905261819032 | 0.006780046456980142 | 0.0007644600632042623 | 0.005500880505570052 | 0.007988023600015104 | +| two_shuffled | 15 | 0.006752328386064863 | 0.006750090991582382 | 0.0006690204050476533 | 0.005729100723183952 | 0.008001562511716862 | +| algorithm2 | 15 | 0.004511679044049125 | 0.004480297686163159 | 0.0006584518614684514 | 0.00361016139623501 | 0.005666001374006467 | +| both | 15 | 0.0037876499460551736 | 0.003825095830262369 | 0.0005597422611197194 | 0.0028974416801691273 | 0.00484655941562858 | + +| Opt-in vs default | U₁ | MWU p | W | signed-rank p | median ratio (opt-in / default) | better median | +|-------------------|----|-------|---|----------------|---------------------------------|---------------| +| two_shuffled | 107.0 | 0.8357050269952793 | 59.0 | 0.9773446461818073 | 0.9655970062875429 | two_shuffled | +| algorithm2 | 1.0 | 4.1432200975375465e-06 | 0.0 | 0.0007265138800579529 | 0.6451794890862746 | algorithm2 | +| both | 0.0 | 3.391821390825104e-06 | 0.0 | 0.0007265138800579529 | 0.5416418218527258 | both | + +Published G6 check (algorithm2 column vs docs/WILCOXON-RESULTS.md DTLZ2 Unsga3 column (that column is PymooCompatible, not Algorithm2)): 15/15 match. The published file was not rewritten. + +| Seed | default | two_shuffled | algorithm2 | both | +|------|------|------|------|------| +| 1 | 0.006992905261819032 | 0.006029451564963478 | 0.004031675764658275 | 0.0037327492292463776 | +| 2 | 0.007079455149960263 | 0.006096842854814725 | 0.005666001374006467 | 0.003919299025962759 | +| 3 | 0.006105593844136021 | 0.006981580446672411 | 0.005130162013811255 | 0.0032189452906896342 | +| 4 | 0.005901386776532219 | 0.005729100723183952 | 0.004776945104309114 | 0.004186982214529411 | +| 5 | 0.005907718617953248 | 0.008001562511716862 | 0.0046632117731021825 | 0.0038362306742307788 | +| 6 | 0.0062889516896962525 | 0.00753601190614311 | 0.005459990389872513 | 0.0037876499460551736 | +| 7 | 0.006844947662146848 | 0.007210495562316086 | 0.00361016139623501 | 0.002957892594413461 | +| 8 | 0.005500880505570052 | 0.006539627292832637 | 0.005111158912956062 | 0.004283318228317427 | +| 9 | 0.007725669378909481 | 0.006898398477956585 | 0.00466990040193047 | 0.004715996545119468 | +| 10 | 0.0071601600211100794 | 0.006687295805576186 | 0.003815591821094395 | 0.0036246524202702296 | +| 11 | 0.007988023600015104 | 0.0059226572658979725 | 0.0039744680194661 | 0.0036853864968974534 | +| 12 | 0.00719451996479668 | 0.007518337819038473 | 0.004511679044049125 | 0.00484655941562858 | +| 13 | 0.0074065684134115035 | 0.0071594455226918425 | 0.003984364774023208 | 0.0035685714409307745 | +| 14 | 0.0060556600142339006 | 0.006188228733866554 | 0.0041724452871580876 | 0.004114762251474884 | +| 15 | 0.007548255954411433 | 0.006752328386064863 | 0.0036267092157751075 | 0.0028974416801691273 | + +## Notes + +- Package defaults are unchanged: `RankNicheDistance`, `IndependentWithReplacement`, single-layer Das–Dennis. +- Front files under `tools/oracle/out/` are gitignored. +- Hypothesis tests reuse `mannwhitney_u` and `wilcoxon_signed_rank` from `run_multiseed_wilcoxon.py`. + diff --git a/docs/ORACLE-RESULTS.md b/docs/ORACLE-RESULTS.md index da7ec2f..bd90ddc 100644 --- a/docs/ORACLE-RESULTS.md +++ b/docs/ORACLE-RESULTS.md @@ -131,7 +131,7 @@ ZDT2 quality A/B is gens=250 + `PymooCompatible` (not the loose smoke bar). ZDT1 | IGD mean-distance | **aligned** | | ASF / axis intercepts | **aligned** | | Collapsed nadir (span ≤ 1e-6) | **delta**: nadir = ideal + 1 after the worst-of-pop fallback. pymoo 0.6.2 stops at worst-of-pop. Locked by `Collapsed_span_sets_nadir_to_ideal_plus_one` (`{2, 2+1e-8}` → nadir 3). | -| Infeasible points | **feasible hyperplane**. Ideal, worst, and ASF use feasible objectives when any feasible member exists. An all-infeasible generation does not move the hyperplane. Survival niches the feasible set and fills a shortfall by ascending CV. The old whole-pool fixture (feasible (1, 1) vs infeasible (0, 0) → ideal (0, 0)) is replaced. OSY, TNK, and C1-DTLZ1 have self-tests. No constrained IGD table. | +| Infeasible points | **feasible hyperplane**. Ideal, worst, and ASF use feasible objectives when any feasible member exists. An all-infeasible generation does not move the hyperplane. Survival niches the feasible set and fills a shortfall by ascending CV. The old whole-pool fixture (feasible (1, 1) vs infeasible (0, 0) → ideal (0, 0)) is replaced. OSY, TNK, and C1-DTLZ1 have self-tests. 15-seed IGD is in [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md). | | Mating re-association | **delta**: after survival, `PrepareForSelection` normalizes the survivors again and overwrites niche ids. pymoo keeps the survival ids. Fixture: `PrepareForSelection_overwrites_survival_niche_ids`. | | Persistent ideal + ND extremes | **aligned** | | `TournamentMode.PymooCompatible` | **implemented** | diff --git a/docs/RESEARCH-STANDARDS.md b/docs/RESEARCH-STANDARDS.md index 3d70092..b26e1b1 100644 --- a/docs/RESEARCH-STANDARDS.md +++ b/docs/RESEARCH-STANDARDS.md @@ -19,7 +19,7 @@ Sources consulted (2026-08): | Bi-objective | ZDT1–4, ZDT6 | Convex / non-convex / disconnected / multimodal / biased density | | Many-objective | DTLZ1–4, DTLZ7 | Linear simplex, sphere, multimodal sphere, bias, disconnected | | Hard (later) | WFG1–9 | Non-separable, scaled, deceptive | -| Constrained | OSY, TNK, C1-DTLZ1 | In the library, with survival self-tests. No IGD table until measured runs exist | +| Constrained | OSY, TNK, C1-DTLZ1 | In the library. 15-seed IGD vs pymoo is in [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md) | **Seada selection paper** stresses ZDT4 (many local fronts) and DTLZ1 for M=3,5. @@ -100,4 +100,4 @@ Export path: dump final `F` as CSV from both sides; compute IGD in this library. - Fixed-seed smoke IGD bounds (loose) + unit metric tests against hand-checked values - GitHub Packages publish workflow -WFG + constrained + automated pymoo bridge remain follow-ons. +WFG and a fuller automated pymoo bridge remain follow-ons. Constrained IGD for OSY, TNK, and C1-DTLZ1 is recorded in [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md). diff --git a/docs/ROADMAP.md b/docs/ROADMAP.md index d38075d..57183d9 100644 --- a/docs/ROADMAP.md +++ b/docs/ROADMAP.md @@ -19,7 +19,7 @@ Living plan for Unsga3. Issues track concrete work; this page is the narrative. - [ ] nuget.org publish (in addition to GitHub Packages) - [ ] GD+ indicator (IGD+ has a hand-case test; GD and IGD were already implemented) - [x] Constrained demos (OSY, TNK, C1-DTLZ1) with formulation and survival self-tests -- [ ] Constrained IGD / Wilcoxon vs pymoo (no table until measured runs exist) +- [x] Constrained IGD / Wilcoxon vs pymoo — [NEW-SURFACES-RESULTS.md](NEW-SURFACES-RESULTS.md). The published ZDT/DTLZ Wilcoxon file is unchanged. - [ ] API docs site (DocFX or similar) ## Medium term (0.3+) diff --git a/tools/OracleCompare/Program.cs b/tools/OracleCompare/Program.cs index 6e1a627..a54ee37 100644 --- a/tools/OracleCompare/Program.cs +++ b/tools/OracleCompare/Program.cs @@ -7,22 +7,30 @@ using Unsga3.Problems; using Unsga3.Utilities; -// Fixed-protocol C# side of the pymoo oracle (see tools/oracle/run_pymoo_oracle.py). -// IGD is scored on the full feasible non-dominated front. pymoo's script scores res.F. +// Fixed-protocol C# side of the pymoo oracle (see tools/oracle/run_pymoo_oracle.py +// and tools/oracle/run_new_surfaces.py). +// IGD is scored on the full feasible non-dominated front when this process has a +// reference front. pymoo's script scores res.F unless the new-surfaces harness +// re-scores both populations itself. // -// dotnet run --project tools/OracleCompare -- --problem zdt1 --partitions 12 --pop 52 --gens 100 --seed 1 --pymoo-mode +// dotnet run --project tools/OracleCompare -- --problem zdt1 --partitions 12 --pop 52 --gens 100 --seed 1 // # ZDT2 quality protocol (matches unsga3-bend A/B): gens=250 + --pymoo-mode. // # Omitted --gens on --problem zdt2 is 250. gens=100 is an early-stress snapshot. -// # RankNicheDistance (omit --pymoo-mode) is optional, not the ZDT2 A/B default. // dotnet run --project tools/OracleCompare -- --problem zdt2 --partitions 12 --pop 52 --seed 1 --pymoo-mode +// +// Mating opt-ins (package defaults stay rank-niche + independent): +// --tournament algorithm2 --mating two-shuffled --seed 1 --seed-count 15 string problemName = "zdt1"; int partitions = 12; int? pop = null; int gens = 100; bool gensExplicit = false; -int seed = 1; +int seedStart = 1; +int seedCount = 1; bool pymooMode = false; +string? tournamentArg = null; +string? matingArg = null; string? outDir = null; for (int i = 0; i < args.Length; i++) @@ -36,103 +44,269 @@ gens = int.Parse(args[++i], CultureInfo.InvariantCulture); gensExplicit = true; break; - case "--seed": seed = int.Parse(args[++i], CultureInfo.InvariantCulture); break; + case "--seed": seedStart = int.Parse(args[++i], CultureInfo.InvariantCulture); break; + case "--seed-count": seedCount = int.Parse(args[++i], CultureInfo.InvariantCulture); break; case "--pymoo-mode": pymooMode = true; break; + case "--tournament": tournamentArg = args[++i]; break; + case "--mating": matingArg = args[++i]; break; case "--out-dir": outDir = args[++i]; break; + default: + Console.Error.WriteLine($"Unknown argument: {args[i]}"); + return 1; } } +if (seedCount < 1) +{ + Console.Error.WriteLine("--seed-count must be at least 1."); + return 1; +} + // Quality protocol: ZDT2 A/B default is gens=250 (unsga3-bend honesty). -// ZDT1 stays 100; DTLZ2 callers still pass --gens 150. Explicit --gens always wins. -if (!gensExplicit && problemName.Equals("zdt2", StringComparison.OrdinalIgnoreCase)) - gens = 250; - -IProblem problem; -int m; -double[][] pf; -switch (problemName.ToLowerInvariant()) +// Constrained bi-objective uses that same 250. C1-DTLZ1 uses the DTLZ2 budget (150). +// ZDT1 stays 100. Explicit --gens always wins. +if (!gensExplicit) { - case "zdt1": - problem = new Zdt1Problem(); - m = 2; - pf = ParetoFronts.Zdt1(500); - break; - case "zdt2": - problem = new Zdt2Problem(); - m = 2; - pf = ParetoFronts.Zdt2(500); - break; - case "dtlz2": - problem = new Dtlz2Problem(nObjectives: 3, k: 10); - m = 3; - pf = ParetoFronts.Dtlz2(3, partitions); - break; - default: - Console.Error.WriteLine($"Unknown problem: {problemName}"); - return 1; + if (problemName.Equals("zdt2", StringComparison.OrdinalIgnoreCase) + || problemName.Equals("osy", StringComparison.OrdinalIgnoreCase) + || problemName.Equals("tnk", StringComparison.OrdinalIgnoreCase)) + gens = 250; + else if (problemName.Equals("c1dtlz1", StringComparison.OrdinalIgnoreCase)) + gens = 150; +} + +TournamentMode mode = tournamentArg is not null + ? ParseTournament(tournamentArg) + : pymooMode ? TournamentMode.PymooCompatible : TournamentMode.RankNicheDistance; +MatingPoolMode mating = matingArg is not null + ? ParseMating(matingArg) + : MatingPoolMode.IndependentWithReplacement; + +if (!TryCreateProblem(problemName, partitions, out IProblem problem, out int m, out double[][]? pf)) +{ + Console.Error.WriteLine($"Unknown problem: {problemName}"); + return 1; } var dirs = ReferenceDirections.DasDennis(m, partitions); int popSize = pop ?? dirs.Length; -var mode = pymooMode ? TournamentMode.PymooCompatible : TournamentMode.RankNicheDistance; +string tournamentSlug = TournamentSlug(mode); +string matingSlug = MatingSlug(mating); Console.WriteLine( - $"Unsga3 | problem={problemName} M={m} n_var={problem.NumberOfVariables} refs={dirs.Length} pop={popSize} gens={gens} seed={seed} tournament={mode}"); + $"Unsga3 | problem={problemName} M={m} n_var={problem.NumberOfVariables} refs={dirs.Length} pop={popSize} gens={gens} seed={seedStart} seed_count={seedCount} tournament={mode} mating={mating}"); + +string dir = outDir ?? Path.Combine(AppContext.BaseDirectory, "..", "..", "..", "..", "oracle", "out"); +dir = Path.GetFullPath(dir); +Directory.CreateDirectory(dir); + +for (int s = 0; s < seedCount; s++) +{ + int seed = seedStart + s; + var algo = new Unsga3Algorithm( + dirs, + popSize, + seed: seed, + tournamentMode: mode, + matingPool: mating); + var result = algo.Run(problem, gens); -var algo = new Unsga3Algorithm(dirs, popSize, seed: seed, tournamentMode: mode); -var result = algo.Run(problem, gens); + var population = result.FinalPopulation; + int feasibleCount = 0; + var feasibleNd = new List(); + foreach (var individual in result.NonDominatedSolutions) + { + if (!individual.IsFeasible) + continue; + feasibleNd.Add((double[])individual.Objectives.Clone()); + } -var obtained = result.NonDominatedSolutions - .Where(i => i.IsFeasible) - .Select(i => (double[])i.Objectives.Clone()) - .ToArray(); -if (obtained.Length == 0) - obtained = result.FinalPopulation.Select(i => (double[])i.Objectives.Clone()).ToArray(); + foreach (var individual in population) + { + if (individual.IsFeasible) + feasibleCount++; + } -double igd = PerformanceIndicators.InvertedGenerationalDistance(obtained, pf); -double? hv = m == 2 - ? PerformanceIndicators.Hypervolume2D(obtained, new[] { 1.1, 1.1 }) - : null; + // Score only when a local reference front exists. OSY and TNK have none here; + // tools/oracle/run_new_surfaces.py scores those from the population file. + double? igd = null; + if (pf is not null && feasibleNd.Count > 0) + igd = PerformanceIndicators.InvertedGenerationalDistance(feasibleNd, pf); -Console.WriteLine($"front=non_dominated n={obtained.Length}"); -Console.WriteLine($"front_size={obtained.Length}"); -Console.WriteLine($"IGD={igd.ToString("G6", CultureInfo.InvariantCulture)}"); -if (hv is double h) - Console.WriteLine($"HV2(r=1.1,1.1)={h.ToString("G6", CultureInfo.InvariantCulture)}"); + // The (1.1, 1.1) reference is the locked ZDT point. OSY and TNK are not on that scale. + bool zdtHv = problemName.Equals("zdt1", StringComparison.OrdinalIgnoreCase) + || problemName.Equals("zdt2", StringComparison.OrdinalIgnoreCase); + double? hv = zdtHv && feasibleNd.Count > 0 + ? PerformanceIndicators.Hypervolume2D(feasibleNd, new[] { 1.1, 1.1 }) + : null; -string dir = outDir ?? Path.Combine(AppContext.BaseDirectory, "..", "..", "..", "..", "oracle", "out"); -dir = Path.GetFullPath(dir); -Directory.CreateDirectory(dir); -string stem = $"csharp_{problemName}_p{partitions}_pop{popSize}_g{gens}_s{seed}_{(pymooMode ? "pymoo" : "default")}"; -string fPath = Path.Combine(dir, $"{stem}_F.csv"); -using (var sw = new StreamWriter(fPath)) + string stem = $"csharp_{problemName}_p{partitions}_pop{popSize}_g{gens}_s{seed}_{tournamentSlug}_{matingSlug}"; + string popPath = Path.Combine(dir, $"{stem}_pop.csv"); + using (var sw = new StreamWriter(popPath)) + { + foreach (var individual in population) + { + var cells = new string[individual.Objectives.Length + 1]; + for (int j = 0; j < individual.Objectives.Length; j++) + cells[j] = individual.Objectives[j].ToString("G17", CultureInfo.InvariantCulture); + cells[^1] = individual.ConstraintViolation.ToString("G17", CultureInfo.InvariantCulture); + sw.WriteLine(string.Join(",", cells)); + } + } + + string fPath = Path.Combine(dir, $"{stem}_F.csv"); + using (var sw = new StreamWriter(fPath)) + { + foreach (var row in feasibleNd) + sw.WriteLine(string.Join(",", row.Select(v => v.ToString("G17", CultureInfo.InvariantCulture)))); + } + + var meta = new Dictionary + { + ["source"] = "Unsga3", + ["algorithm"] = "Unsga3Algorithm", + ["tournament"] = mode.ToString(), + ["mating"] = mating.ToString(), + ["problem"] = problemName, + ["n_obj"] = m, + ["n_var"] = problem.NumberOfVariables, + ["partitions"] = partitions, + ["n_ref_dirs"] = dirs.Length, + ["pop_size"] = popSize, + ["n_gen"] = gens, + ["seed"] = seed, + ["n_population"] = population.Count, + ["n_feasible"] = feasibleCount, + ["n_solutions"] = feasibleNd.Count, + ["front_definition"] = "non_dominated_feasible", + ["igd"] = igd, + ["igd_reference"] = pf is null ? "none_local" : "library_pareto_front", + ["hv2"] = hv, + ["F_csv"] = Path.GetFileName(fPath), + ["pop_csv"] = Path.GetFileName(popPath), + }; + string metaPath = Path.Combine(dir, $"{stem}_meta.json"); + File.WriteAllText(metaPath, JsonSerializer.Serialize(meta, new JsonSerializerOptions { WriteIndented = true })); + + string igdText = igd is double measured + ? measured.ToString("G17", CultureInfo.InvariantCulture) + : "skip"; + string igdG6 = igd is double measuredG6 + ? measuredG6.ToString("G6", CultureInfo.InvariantCulture) + : "skip"; + if (seedCount == 1) + { + Console.WriteLine($"front=non_dominated n={feasibleNd.Count}"); + Console.WriteLine($"front_size={feasibleNd.Count}"); + Console.WriteLine($"feasible={feasibleCount}"); + Console.WriteLine(igd is null ? "IGD=skip" : $"IGD={igdG6}"); + if (hv is double h) + Console.WriteLine($"HV2(r=1.1,1.1)={h.ToString("G6", CultureInfo.InvariantCulture)}"); + Console.WriteLine($"wrote {fPath}"); + Console.WriteLine($"wrote {popPath}"); + Console.WriteLine($"wrote {metaPath}"); + } + + Console.WriteLine( + $"RESULT seed={seed.ToString(CultureInfo.InvariantCulture)} igd={igdText} igd_g6={igdG6} n={feasibleNd.Count.ToString(CultureInfo.InvariantCulture)} feasible={feasibleCount.ToString(CultureInfo.InvariantCulture)} pop_csv={Path.GetFileName(popPath)}"); +} + +return 0; + +static bool TryCreateProblem(string problemName, int partitions, out IProblem problem, out int m, out double[][]? pf) { - foreach (var row in obtained) - sw.WriteLine(string.Join(",", row.Select(v => v.ToString("G17", CultureInfo.InvariantCulture)))); + switch (problemName.ToLowerInvariant()) + { + case "zdt1": + problem = new Zdt1Problem(); + m = 2; + pf = ParetoFronts.Zdt1(500); + return true; + case "zdt2": + problem = new Zdt2Problem(); + m = 2; + pf = ParetoFronts.Zdt2(500); + return true; + case "dtlz2": + problem = new Dtlz2Problem(nObjectives: 3, k: 10); + m = 3; + pf = ParetoFronts.Dtlz2(3, partitions); + return true; + case "osy": + problem = new OsyProblem(); + m = 2; + pf = null; + return true; + case "tnk": + problem = new TnkProblem(); + m = 2; + pf = null; + return true; + case "c1dtlz1": + // k=5 matches C1Dtlz1Problem / Dtlz1Problem (n_var = 7 when M=3). + // pymoo get_problem("c1dtlz1") defaults to n_var=12. The new-surfaces + // harness passes n_var=7. The reference front is the DTLZ1 simplex; + // the C1 inequality does not remove that front. + problem = new C1Dtlz1Problem(nObjectives: 3, k: 5); + m = 3; + pf = ParetoFronts.Dtlz1(3, partitions); + return true; + default: + problem = null!; + m = 0; + pf = null; + return false; + } } -var meta = new Dictionary +static TournamentMode ParseTournament(string value) { - ["source"] = "Unsga3", - ["algorithm"] = "Unsga3Algorithm", - ["tournament"] = mode.ToString(), - ["problem"] = problemName, - ["n_obj"] = m, - ["n_var"] = problem.NumberOfVariables, - ["partitions"] = partitions, - ["n_ref_dirs"] = dirs.Length, - ["pop_size"] = popSize, - ["n_gen"] = gens, - ["seed"] = seed, - ["n_solutions"] = obtained.Length, - ["front_definition"] = "non_dominated_feasible", - ["igd"] = igd, - ["hv2"] = hv, - ["F_csv"] = Path.GetFileName(fPath), + switch (value.ToLowerInvariant()) + { + case "default": + case "rank-niche": + case "rankniche": + return TournamentMode.RankNicheDistance; + case "pymoo": + case "pymoo-compatible": + return TournamentMode.PymooCompatible; + case "algorithm2": + case "algorithm-2": + return TournamentMode.Algorithm2; + default: + throw new ArgumentException( + $"Unknown --tournament '{value}'. Expected rank-niche, pymoo, or algorithm2."); + } +} + +static MatingPoolMode ParseMating(string value) +{ + switch (value.ToLowerInvariant()) + { + case "default": + case "independent": + return MatingPoolMode.IndependentWithReplacement; + case "two-shuffled": + case "two-shuffle": + case "twoshuffle": + return MatingPoolMode.TwoShuffledPasses; + default: + throw new ArgumentException( + $"Unknown --mating '{value}'. Expected independent or two-shuffled."); + } +} + +static string TournamentSlug(TournamentMode mode) => mode switch +{ + TournamentMode.RankNicheDistance => "rankniche", + TournamentMode.PymooCompatible => "pymoo", + TournamentMode.Algorithm2 => "algorithm2", + _ => throw new ArgumentOutOfRangeException(nameof(mode), mode, "Unhandled tournament mode."), }; -string metaPath = Path.Combine(dir, $"{stem}_meta.json"); -File.WriteAllText(metaPath, JsonSerializer.Serialize(meta, new JsonSerializerOptions { WriteIndented = true })); -Console.WriteLine($"wrote {fPath}"); -Console.WriteLine($"wrote {metaPath}"); -return 0; +static string MatingSlug(MatingPoolMode mode) => mode switch +{ + MatingPoolMode.IndependentWithReplacement => "independent", + MatingPoolMode.TwoShuffledPasses => "twoshuffle", + _ => throw new ArgumentOutOfRangeException(nameof(mode), mode, "Unhandled mating pool."), +}; diff --git a/tools/oracle/run_new_surfaces.py b/tools/oracle/run_new_surfaces.py new file mode 100644 index 0000000..61a616f --- /dev/null +++ b/tools/oracle/run_new_surfaces.py @@ -0,0 +1,904 @@ +#!/usr/bin/env python3 +""" +Measured A/B for surfaces added in v0.2.0. + +Does not write docs/WILCOXON-RESULTS.md and does not change package defaults. + +1. Constrained problems (OSY, TNK, C1-DTLZ1), 15 seeds by default. + C# uses RankNicheDistance + IndependentWithReplacement (the constructor defaults). + pymoo uses UNSGA3 defaults. Both feasible non-dominated final populations are + scored with pymoo's IGD against one shared reference front. +2. Mating opt-ins on ZDT1 (and ZDT2 / DTLZ2 when requested). Four C# configs, + library IGD, same seeds as the published ZDT1 table. The published table is + only used as a G6 string check of the default ZDT1 column. + +Usage: + python tools/oracle/run_new_surfaces.py + python tools/oracle/run_new_surfaces.py --seeds 1 --mating zdt1 --constrained osy +""" +from __future__ import annotations + +import argparse +import json +import math +import statistics +import subprocess +import sys +from dataclasses import dataclass +from datetime import date +from pathlib import Path + +import numpy as np +import pymoo +from pymoo.algorithms.moo.unsga3 import UNSGA3 +from pymoo.indicators.igd import IGD +from pymoo.optimize import minimize +from pymoo.problems import get_problem +from pymoo.util.nds.non_dominated_sorting import NonDominatedSorting +from pymoo.util.ref_dirs import get_reference_directions + +from run_multiseed_wilcoxon import mannwhitney_u, wilcoxon_signed_rank + +ROOT = Path(__file__).resolve().parents[2] +OUT = Path(__file__).resolve().parent / "out" +ORACLE_COMPARE = ROOT / "tools" / "OracleCompare" + +# Published ZDT1 C# column in docs/WILCOXON-RESULTS.md (G6 console strings). +# Used only to check a remeasure of the default config. Not a new measurement. +PUBLISHED_ZDT1_CSHARP_G6 = { + 1: "0.0514307", + 2: "0.0609151", + 3: "0.0377512", + 4: "0.064033", + 5: "0.0735681", + 6: "0.105668", + 7: "0.100905", + 8: "0.0836797", + 9: "0.0844962", + 10: "0.0497514", + 11: "0.0390644", + 12: "0.053173", + 13: "0.038507", + 14: "0.0422284", + 15: "0.048006", +} + +# Published DTLZ2 Unsga3 column is PymooCompatible, not Algorithm2. +# Compared only as a G6 string check. The published file is not rewritten. +PUBLISHED_DTLZ2_PYMOO_COMPAT_G6 = { + 1: "0.00403168", + 2: "0.005666", + 3: "0.00513016", + 4: "0.00477695", + 5: "0.00466321", + 6: "0.00545999", + 7: "0.00361016", + 8: "0.00511116", + 9: "0.0046699", + 10: "0.00381559", + 11: "0.00397447", + 12: "0.00451168", + 13: "0.00398436", + 14: "0.00417245", + 15: "0.00362671", +} + + +@dataclass(frozen=True) +class ConstrainedProtocol: + name: str + partitions: int + pop: int + gens: int + n_obj: int + n_var: int + + +@dataclass(frozen=True) +class MatingProtocol: + name: str + partitions: int + pop: int + gens: int + + +# Bi-objective constrained: same pop and partitions as ZDT, generations of the +# ZDT2 quality budget. C1-DTLZ1: same pop / partitions / generations as DTLZ2, +# but n_var=7 (k=5), which is C1Dtlz1Problem's default. pymoo's own c1dtlz1 +# default is n_var=12 and is not used. +CONSTRAINED: dict[str, ConstrainedProtocol] = { + "osy": ConstrainedProtocol("osy", partitions=12, pop=52, gens=250, n_obj=2, n_var=6), + "tnk": ConstrainedProtocol("tnk", partitions=12, pop=52, gens=250, n_obj=2, n_var=2), + "c1dtlz1": ConstrainedProtocol("c1dtlz1", partitions=12, pop=92, gens=150, n_obj=3, n_var=7), +} + +MATING: dict[str, MatingProtocol] = { + "zdt1": MatingProtocol("zdt1", partitions=12, pop=52, gens=100), + # gens=250 is the quality budget. The tournament in this A/B is still + # RankNicheDistance on the default arm, not the ZDT2 PymooCompatible column + # (that column was never published). + "zdt2": MatingProtocol("zdt2", partitions=12, pop=52, gens=250), + # Default arm is RankNicheDistance, not the published DTLZ2 PymooCompatible run. + "dtlz2": MatingProtocol("dtlz2", partitions=12, pop=92, gens=150), +} + +# (column key, --tournament, --mating). First column is the package default. +MATING_CONFIGS: tuple[tuple[str, str, str], ...] = ( + ("default", "rank-niche", "independent"), + ("two_shuffled", "rank-niche", "two-shuffled"), + ("algorithm2", "algorithm2", "independent"), + ("both", "algorithm2", "two-shuffled"), +) + +MATING_LABELS = { + "default": "RankNicheDistance + IndependentWithReplacement", + "two_shuffled": "RankNicheDistance + TwoShuffledPasses", + "algorithm2": "Algorithm2 + IndependentWithReplacement", + "both": "Algorithm2 + TwoShuffledPasses", +} + + +def fmt(x: float) -> str: + if isinstance(x, float) and math.isnan(x): + return "skip: undefined" + return repr(float(x)) + + +def assert_formulations() -> None: + """Lock the pymoo problems to the fixtures in ConstrainedProblemTests.""" + osy = get_problem("osy") + f_osy, g_osy = osy.evaluate( + np.array([[5.0, 1.0, 5.0, 0.0, 5.0, 0.0]]), + return_values_of=["F", "G"], + ) + if not np.allclose(f_osy, [[-274.0, 76.0]]) or not np.allclose(g_osy, [[-2.0, 0.0, -3.0, 0.0, 0.0, 0.0]]): + raise RuntimeError(f"OSY formulation mismatch: F={f_osy} G={g_osy}") + + tnk = get_problem("tnk") + f_tnk, g_tnk = tnk.evaluate(np.array([[1.0, 1.0]]), return_values_of=["F", "G"]) + if not np.allclose(f_tnk, [[1.0, 1.0]]) or not np.allclose(g_tnk, [[-0.9, 0.0]]): + raise RuntimeError(f"TNK formulation mismatch: F={f_tnk} G={g_tnk}") + + c1 = get_problem("c1dtlz1", n_obj=3, n_var=7) + if int(c1.n_var) != 7: + raise RuntimeError(f"C1-DTLZ1 n_var={c1.n_var}, expected 7") + f_c1, g_c1 = c1.evaluate(np.full((1, 7), 0.5), return_values_of=["F", "G"]) + if not np.allclose(f_c1, [[0.125, 0.125, 0.25]]) or not np.allclose(g_c1, [[-1.0 / 12.0]]): + raise RuntimeError(f"C1-DTLZ1 formulation mismatch: F={f_c1} G={g_c1}") + + +def reference_front(proto: ConstrainedProtocol) -> np.ndarray: + ref_dirs = get_reference_directions("das-dennis", proto.n_obj, n_partitions=proto.partitions) + if len(ref_dirs) != {2: 13, 3: 91}[proto.n_obj]: + raise RuntimeError(f"{proto.name}: unexpected ref count {len(ref_dirs)}") + if proto.name == "c1dtlz1": + problem = get_problem("c1dtlz1", n_obj=proto.n_obj, n_var=proto.n_var) + pf = np.atleast_2d(problem.pareto_front(ref_dirs)) + else: + problem = get_problem(proto.name) + if int(problem.n_var) != proto.n_var: + raise RuntimeError(f"{proto.name}: pymoo n_var={problem.n_var}, expected {proto.n_var}") + pf = np.atleast_2d(problem.pareto_front()) + if pf.size == 0 or not np.isfinite(pf).all(): + raise RuntimeError(f"{proto.name}: Pareto front missing or non-finite") + return pf + + +def feasible_nd(objectives: np.ndarray, cv: np.ndarray) -> np.ndarray | None: + objectives = np.atleast_2d(np.asarray(objectives, dtype=float)) + cv = np.asarray(cv, dtype=float).reshape(-1) + if objectives.shape[0] != cv.shape[0]: + raise RuntimeError(f"F rows {objectives.shape[0]} != CV {cv.shape[0]}") + mask = cv <= 0.0 + if not np.any(mask): + return None + feasible = np.atleast_2d(objectives[mask]) + if feasible.shape[0] == 1: + return feasible + index = NonDominatedSorting().do(feasible, only_non_dominated_front=True) + return np.atleast_2d(feasible[index]) + + +def igd_of(obtained: np.ndarray, pf: np.ndarray) -> float: + return float(IGD(pf)(np.atleast_2d(obtained))) + + +def load_population(path: Path) -> tuple[np.ndarray, np.ndarray]: + data = np.atleast_2d(np.loadtxt(path, delimiter=",")) + if data.shape[1] < 2: + raise RuntimeError(f"{path} has no constraint column") + return data[:, :-1], data[:, -1] + + +def parse_result_line(line: str) -> dict: + fields: dict[str, str] = {} + for part in line.split()[1:]: + key, value = part.split("=", 1) + fields[key] = value + igd = None if fields["igd"] == "skip" else float(fields["igd"]) + return { + "seed": int(fields["seed"]), + "igd": igd, + "igd_g6": fields["igd_g6"], + "n": int(fields["n"]), + "feasible": int(fields["feasible"]), + "pop_csv": fields["pop_csv"], + } + + +def run_csharp( + problem: str, + partitions: int, + pop: int, + gens: int, + seed_start: int, + seed_count: int, + tournament: str, + mating: str, +) -> list[dict]: + args = [ + "dotnet", + "run", + "--project", + str(ORACLE_COMPARE), + "-c", + "Release", + "--no-build", + "--", + "--problem", + problem, + "--partitions", + str(partitions), + "--pop", + str(pop), + "--gens", + str(gens), + "--seed", + str(seed_start), + "--seed-count", + str(seed_count), + "--tournament", + tournament, + "--mating", + mating, + "--out-dir", + str(OUT), + ] + completed = subprocess.run(args, capture_output=True, text=True, cwd=str(ROOT)) + if completed.returncode != 0: + raise RuntimeError( + f"C# {problem} {tournament}/{mating} failed:\n{completed.stdout}\n{completed.stderr}" + ) + rows = [parse_result_line(line) for line in completed.stdout.splitlines() if line.startswith("RESULT ")] + if len(rows) != seed_count: + raise RuntimeError( + f"C# {problem} returned {len(rows)} RESULT lines, expected {seed_count}:\n{completed.stdout}" + ) + return rows + + +def run_pymoo_constrained(proto: ConstrainedProtocol, seed: int, pf: np.ndarray) -> dict: + ref_dirs = get_reference_directions("das-dennis", proto.n_obj, n_partitions=proto.partitions) + if proto.name == "c1dtlz1": + problem = get_problem("c1dtlz1", n_obj=proto.n_obj, n_var=proto.n_var) + else: + problem = get_problem(proto.name) + algorithm = UNSGA3(ref_dirs, pop_size=proto.pop) + result = minimize( + problem, + algorithm, + ("n_gen", proto.gens), + seed=seed, + verbose=False, + save_history=False, + ) + objectives = np.atleast_2d(result.pop.get("F")) + raw_cv = result.pop.get("CV") + if raw_cv is None: + raise RuntimeError(f"{proto.name} seed={seed}: pymoo population has no CV") + cv = np.asarray(raw_cv, dtype=float).reshape(-1) + nd = feasible_nd(objectives, cv) + res_f = np.atleast_2d(result.F) if result.F is not None else np.empty((0, proto.n_obj)) + res_f_ok = res_f.size > 0 and np.isfinite(res_f).all() + stem = f"pymoo_{proto.name}_p{proto.partitions}_pop{proto.pop}_g{proto.gens}_s{seed}" + np.savetxt(OUT / f"{stem}_pop.csv", np.column_stack([objectives, cv]), delimiter=",") + return { + "seed": seed, + "igd": None if nd is None else igd_of(nd, pf), + "n": 0 if nd is None else int(nd.shape[0]), + "feasible": int(np.sum(cv <= 0.0)), + "resf_igd": igd_of(res_f, pf) if res_f_ok else None, + "resf_n": int(res_f.shape[0]) if res_f_ok else 0, + "n_var": int(problem.n_var), + } + + +def summarize(values: list[float]) -> dict: + ordered = sorted(values) + n = len(ordered) + return { + "n": n, + "mean": statistics.fmean(ordered), + "std": statistics.stdev(ordered) if n > 1 else 0.0, + "median": statistics.median(ordered), + "min": ordered[0], + "max": ordered[-1], + } + + +def hypothesis( + left: list[float], + right: list[float], + left_name: str = "Unsga3", + right_name: str = "pymoo", +) -> dict: + u_stat, p_mw, note_mw = mannwhitney_u(left, right) + diffs = [a - b for a, b in zip(left, right)] + w_stat, p_wx, note_wx = wilcoxon_signed_rank(diffs) + left_med = statistics.median(left) + right_med = statistics.median(right) + if left_med < right_med: + better = left_name + elif right_med < left_med: + better = right_name + else: + better = "tie" + + def verdict(p: float) -> str: + if isinstance(p, float) and math.isnan(p): + return "skip: undefined" + if p < 0.05: + return f"reject H₀ (better median = {better})" + return "fail to reject H₀" + + ratio = left_med / right_med if right_med != 0.0 else float("nan") + return { + "mannwhitney_U": u_stat, + "mannwhitney_p": p_mw, + "mannwhitney_note": note_mw, + "mannwhitney_verdict": verdict(p_mw), + "wilcoxon_W": w_stat, + "wilcoxon_p": p_wx, + "wilcoxon_note": note_wx, + "wilcoxon_verdict": verdict(p_wx), + "median_ratio": ratio, + "better_median": better, + "n": len(left), + } + + +def paired_or_skip(rows: list[dict], key_a: str, key_b: str) -> tuple[list[float], list[float], int]: + left: list[float] = [] + right: list[float] = [] + skipped = 0 + for row in rows: + a = row[key_a] + b = row[key_b] + if a is None or b is None: + skipped += 1 + continue + left.append(float(a)) + right.append(float(b)) + return left, right, skipped + + +def measure_constrained(names: list[str], seeds: list[int]) -> dict: + blocks: dict = {} + for name in names: + proto = CONSTRAINED[name] + print(f"\n=== constrained {name} seeds {seeds[0]}..{seeds[-1]} ===", flush=True) + pf = reference_front(proto) + print(f" C# default ...", flush=True) + csharp_rows = run_csharp( + proto.name, + proto.partitions, + proto.pop, + proto.gens, + seeds[0], + len(seeds), + "rank-niche", + "independent", + ) + per_seed = [] + for row, seed in zip(csharp_rows, seeds): + if row["seed"] != seed: + raise RuntimeError(f"{name}: seed order {row['seed']} != {seed}") + objectives, cv = load_population(OUT / row["pop_csv"]) + nd = feasible_nd(objectives, cv) + fair = None if nd is None else igd_of(nd, pf) + per_seed.append( + { + "seed": seed, + "csharp_igd": fair, + "csharp_n": 0 if nd is None else int(nd.shape[0]), + "csharp_feasible": int(np.sum(cv <= 0.0)), + "csharp_library_igd": row["igd"], + "pymoo_igd": None, + "pymoo_n": None, + "pymoo_feasible": None, + "pymoo_resf_igd": None, + "pymoo_resf_n": None, + } + ) + print( + f" C# seed={seed} fair_igd=" + f"{'skip' if fair is None else fmt(fair)} n={per_seed[-1]['csharp_n']}", + flush=True, + ) + + for entry in per_seed: + print(f" py seed={entry['seed']} ...", end=" ", flush=True) + pymoo_row = run_pymoo_constrained(proto, entry["seed"], pf) + repeat = run_pymoo_constrained(proto, entry["seed"], pf) + if pymoo_row["n_var"] != proto.n_var or repeat["n_var"] != proto.n_var: + raise RuntimeError(f"{name}: pymoo ran n_var={pymoo_row['n_var']}") + entry["pymoo_igd"] = pymoo_row["igd"] + entry["pymoo_n"] = pymoo_row["n"] + entry["pymoo_feasible"] = pymoo_row["feasible"] + entry["pymoo_resf_igd"] = pymoo_row["resf_igd"] + entry["pymoo_resf_n"] = pymoo_row["resf_n"] + entry["pymoo_igd_repeat"] = repeat["igd"] + entry["pymoo_resf_igd_repeat"] = repeat["resf_igd"] + entry["pymoo_repeat_match"] = ( + pymoo_row["igd"] == repeat["igd"] and pymoo_row["resf_igd"] == repeat["resf_igd"] + ) + shown = "skip" if pymoo_row["igd"] is None else fmt(pymoo_row["igd"]) + flag = "" if entry["pymoo_repeat_match"] else " repeat-differs" + print(f"fair_igd={shown} n={pymoo_row['n']}{flag}", flush=True) + + left, right, skipped = paired_or_skip(per_seed, "csharp_igd", "pymoo_igd") + tests = None + if len(left) >= 5: + tests = hypothesis(left, right, left_name="Unsga3", right_name="pymoo") + block = { + "protocol": { + "problem": proto.name, + "partitions": proto.partitions, + "pop": proto.pop, + "gens": proto.gens, + "n_obj": proto.n_obj, + "n_var": proto.n_var, + "n_ref_dirs": int(len(get_reference_directions( + "das-dennis", proto.n_obj, n_partitions=proto.partitions + ))), + "pf_points": int(pf.shape[0]), + "pf_min": [float(v) for v in pf.min(axis=0)], + "pf_max": [float(v) for v in pf.max(axis=0)], + "csharp_tournament": "RankNicheDistance", + "csharp_mating": "IndependentWithReplacement", + "pymoo_algorithm": "UNSGA3", + "scored_set": "feasible_nondominated_final_population", + "igd": "pymoo.indicators.igd.IGD", + }, + "seeds": per_seed, + "n_skipped_pairs": skipped, + "tests": tests, + "csharp_summary": summarize(left) if left else None, + "pymoo_summary": summarize(right) if right else None, + } + blocks[name] = block + return blocks + + +def g6_check(rows: list[dict], expected: dict[int, str], compared_to: str) -> dict: + mismatches = [] + for row in rows: + published = expected[row["seed"]] + if row["igd_g6"] != published: + mismatches.append( + {"seed": row["seed"], "measured_g6": row["igd_g6"], "published_g6": published} + ) + return { + "compared_to": compared_to, + "n_match": len(rows) - len(mismatches), + "n_seeds": len(rows), + "mismatches": mismatches, + } + + +def measure_mating(names: list[str], seeds: list[int]) -> dict: + blocks: dict = {} + for name in names: + proto = MATING[name] + print(f"\n=== mating {name} seeds {seeds[0]}..{seeds[-1]} ===", flush=True) + columns: dict[str, list[dict]] = {} + for key, tournament, mating in MATING_CONFIGS: + print(f" C# {key} ({tournament} / {mating}) ...", flush=True) + rows = run_csharp( + proto.name, + proto.partitions, + proto.pop, + proto.gens, + seeds[0], + len(seeds), + tournament, + mating, + ) + if any(row["igd"] is None for row in rows): + raise RuntimeError(f"{name} {key} returned a skip IGD; unconstrained fronts should score") + if [row["seed"] for row in rows] != seeds: + raise RuntimeError(f"{name} {key} seed mismatch") + columns[key] = rows + for row in rows: + print(f" seed={row['seed']} igd={fmt(row['igd'])} g6={row['igd_g6']}", flush=True) + + default_rows = columns["default"] + comparisons = {} + for key, _, _ in MATING_CONFIGS: + if key == "default": + continue + left = [float(row["igd"]) for row in columns[key]] + right = [float(row["igd"]) for row in default_rows] + comparisons[key] = hypothesis(left, right, left_name=key, right_name="default") + + published_checks = [] + if seeds == list(range(1, 16)) and name == "zdt1": + published_checks.append( + g6_check( + default_rows, + PUBLISHED_ZDT1_CSHARP_G6, + "default column vs docs/WILCOXON-RESULTS.md ZDT1 Unsga3 column", + ) + ) + if seeds == list(range(1, 16)) and name == "dtlz2": + published_checks.append( + g6_check( + columns["algorithm2"], + PUBLISHED_DTLZ2_PYMOO_COMPAT_G6, + "algorithm2 column vs docs/WILCOXON-RESULTS.md DTLZ2 Unsga3 column " + "(that column is PymooCompatible, not Algorithm2)", + ) + ) + + per_seed = [] + for index, seed in enumerate(seeds): + entry = {"seed": seed} + for key, _, _ in MATING_CONFIGS: + entry[key] = columns[key][index]["igd"] + entry[f"{key}_g6"] = columns[key][index]["igd_g6"] + entry[f"{key}_n"] = columns[key][index]["n"] + per_seed.append(entry) + + blocks[name] = { + "protocol": { + "problem": proto.name, + "partitions": proto.partitions, + "pop": proto.pop, + "gens": proto.gens, + "seeds": seeds, + "igd": "Unsga3 PerformanceIndicators.InvertedGenerationalDistance on the feasible non-dominated front", + "reference": "ParetoFronts library sampler (ZDT 500 points; DTLZ2 Das-Dennis sphere at the same partitions)", + }, + "seeds": per_seed, + "summaries": { + key: summarize([float(row["igd"]) for row in columns[key]]) + for key, _, _ in MATING_CONFIGS + }, + "versus_default": comparisons, + "published_g6_checks": published_checks, + } + return blocks + + +def summary_cells(summary: dict | None) -> str: + if summary is None: + return "skip: no paired feasible seeds" + return ( + f"n={summary['n']} median={fmt(summary['median'])} " + f"mean={fmt(summary['mean'])} std={fmt(summary['std'])} " + f"min={fmt(summary['min'])} max={fmt(summary['max'])}" + ) + + +def markdown_report(payload: dict) -> str: + lines = [ + "# New surfaces: constrained problems and mating opt-ins", + "", + f"Measured {payload['measured_on']}. pymoo {payload['pymoo_version']}. " + f".NET SDK {payload['dotnet_version']}.", + "", + "Generated by `tools/oracle/run_new_surfaces.py`. " + "Every numeric cell is that run. A cell that was not measured is `skip:`. " + "This file does not replace [`WILCOXON-RESULTS.md`](WILCOXON-RESULTS.md).", + "", + "## Command", + "", + "```bash", + payload["command"], + "```", + "", + "C# side, one config, all seeds in one process:", + "", + "```bash", + "dotnet run --project tools/OracleCompare -c Release -- \\", + " --problem zdt1 --partitions 12 --pop 52 --gens 100 --seed 1 --seed-count 15 \\", + " --tournament rank-niche --mating independent", + "dotnet run --project tools/OracleCompare -c Release -- \\", + " --problem osy --partitions 12 --pop 52 --gens 250 --seed 1 --seed-count 15", + "dotnet run --project tools/OracleCompare -c Release -- \\", + " --problem c1dtlz1 --partitions 12 --pop 92 --gens 150 --seed 1 --seed-count 15", + "```", + "", + "Omitted `--gens` on `osy` and `tnk` is 250. Omitted `--gens` on `c1dtlz1` is 150. " + "Omitted `--tournament` is `RankNicheDistance`. Omitted `--mating` is `IndependentWithReplacement`.", + "", + "## Protocol", + "", + "Operators on both sides: SBX η=30, p_c=1, PM η=20, p_m=1/n, duplicate elimination on, " + "single-layer Das–Dennis. Seeds are `1 .. n` (published ZDT1 / DTLZ2 indexes). " + "C# and pymoo do not share an RNG stream. Pairing is by seed index.", + "", + "### Constrained", + "", + "| Problem | M | n_var | Partitions | Refs | Pop | Gens | C# tournament | C# mating | pymoo |", + "|---------|---|------:|------------|-----:|----:|-----:|---------------|-----------|-------|", + ] + for name, block in payload["constrained"].items(): + p = block["protocol"] + lines.append( + f"| {name} | {p['n_obj']} | {p['n_var']} | {p['partitions']} | {p['n_ref_dirs']} | " + f"{p['pop']} | {p['gens']} | {p['csharp_tournament']} | {p['csharp_mating']} | " + f"{p['pymoo_algorithm']} defaults |" + ) + lines += [ + "", + "Scored set for the comparison: feasible non-dominated subset of the **final population**, " + "both solvers, one `pymoo.indicators.igd.IGD` call, one shared reference front. " + "OSY and TNK use pymoo's shipped Pareto-front file (`pareto_front()`). " + "C1-DTLZ1 uses `pareto_front(das-dennis partitions=12)` at **n_var=7** (k=5). " + "That front is the DTLZ1 simplex (0.5 × the directions). The C1 inequality does not remove it. " + "pymoo's default `get_problem(\"c1dtlz1\")` is n_var=12 and is not this protocol.", + "", + "C# constructor defaults are not pymoo's mating. pymoo `UNSGA3` uses " + "`comp_by_rank_and_ref_line_dist` and two shuffled tournament passes. " + "The table is default Unsga3 against default pymoo UNSGA3 on a shared front definition. " + "It is not a mating-matched ablation. The mating section is that ablation, on unconstrained problems. " + "Each pymoo seed is run twice. The table is the first run. A second value is listed only when it differs. " + "The two runs are not averaged.", + "", + "`pymoo res.F` is reported in its own column. It is the survival result, not the feasible " + "non-dominated population, and it is not an input to the hypothesis tests.", + "", + "A seed with no feasible point on either side is `skip: no feasible points` and is left out of " + "the median and the tests. Tests run only when at least five seeds are paired.", + "", + "### Mating opt-ins", + "", + "| Problem | Partitions | Pop | Gens | Reference front | Seeds |", + "|---------|------------|----:|-----:|-----------------|-------|", + ] + for name, block in payload["mating"].items(): + p = block["protocol"] + lines.append( + f"| {name} | {p['partitions']} | {p['pop']} | {p['gens']} | library sampler | " + f"{p['seeds'][0]}..{p['seeds'][-1]} |" + ) + lines += [ + "", + "IGD is the C# feasible non-dominated front against `ParetoFronts` " + "(ZDT: 500 points; DTLZ2: Das–Dennis sphere at the same partitions). " + "That is the definition of the published ZDT1 C# column.", + "", + "| Column | Tournament | Mating pool |", + "|--------|------------|-------------|", + "| default | RankNicheDistance | IndependentWithReplacement |", + "| two_shuffled | RankNicheDistance | TwoShuffledPasses |", + "| algorithm2 | Algorithm2 | IndependentWithReplacement |", + "| both | Algorithm2 | TwoShuffledPasses |", + "", + "ZDT2 uses gens=250. Its default column is `RankNicheDistance`, not `PymooCompatible`. " + "There is no published ZDT2 Wilcoxon column. " + "DTLZ2's default column is also `RankNicheDistance`. " + "The published DTLZ2 column is `PymooCompatible` and is not rewritten here.", + "", + "Paired tests compare each opt-in with the default column (same seed index). " + "Lower IGD is better. p-values are two-sided and not multiplicity-adjusted.", + "", + ] + + lines += ["## Constrained results", ""] + if not payload["constrained"]: + lines += ["skip: constrained sweep not requested", ""] + for name, block in payload["constrained"].items(): + lines += [f"### {name.upper()}", ""] + lines.append(f"- C#: {summary_cells(block['csharp_summary'])}") + lines.append(f"- pymoo feasible ND: {summary_cells(block['pymoo_summary'])}") + lines.append(f"- Unpaired skips: {block['n_skipped_pairs']}") + p = block["protocol"] + lines.append( + "- Reference front: " + f"n={p['pf_points']} min={fmt_vec(p['pf_min'])} max={fmt_vec(p['pf_max'])}" + ) + tests = block["tests"] + if tests is None: + lines.append("- Hypothesis tests: skip: fewer than 5 paired feasible seeds") + else: + lines += [ + "", + "| Test | Statistic | p-value | Verdict (α=0.05) |", + "|------|-----------|---------|------------------|", + f"| Mann–Whitney U (U₁) | {fmt(tests['mannwhitney_U'])} | {fmt(tests['mannwhitney_p'])} | {tests['mannwhitney_verdict']} |", + f"| Wilcoxon signed-rank | {fmt(tests['wilcoxon_W'])} | {fmt(tests['wilcoxon_p'])} | {tests['wilcoxon_verdict']} |", + "", + f"- Median ratio (C# / pymoo), feasible ND: **{fmt(tests['median_ratio'])}**", + f"- Better median: **{tests['better_median']}**", + f"- n paired: **{tests['n']}**", + ] + lines += [ + "", + "| Seed | C# feasible ND IGD | pymoo feasible ND IGD | pymoo res.F IGD | C# n | pymoo n |", + "|------|--------------------|-----------------------|-----------------|-----:|--------:|", + ] + for row in block["seeds"]: + lines.append( + f"| {row['seed']} | {cell(row['csharp_igd'])} | {cell(row['pymoo_igd'])} | " + f"{cell(row['pymoo_resf_igd'])} | {row['csharp_n']} | {row['pymoo_n']} |" + ) + lines.append("") + diffs = [row for row in block["seeds"] if not row.get("pymoo_repeat_match", True)] + if diffs: + lines.append( + "A second pymoo run of the same seed is not always bit-identical on this host " + "(SciPy OpenBLAS, dynamic threads). The table cell is the first run. " + "The second run differed for:" + ) + for row in diffs: + lines.append( + f"- seed {row['seed']} feasible ND: table {cell(row['pymoo_igd'])}, " + f"repeat {cell(row['pymoo_igd_repeat'])}; " + f"res.F: table {cell(row['pymoo_resf_igd'])}, " + f"repeat {cell(row['pymoo_resf_igd_repeat'])}" + ) + lines.append("The two runs are not averaged. Hypothesis tests use the table column.") + else: + lines.append( + "A second pymoo run of each seed matched the table cell (feasible ND and res.F)." + ) + lines.append("") + + lines += ["## Mating opt-in results", ""] + if not payload["mating"]: + lines += ["skip: mating sweep not requested", ""] + for name, block in payload["mating"].items(): + lines += [f"### {name.upper()}", ""] + lines += [ + "| Config | n | median | mean | std | min | max |", + "|--------|--:|--------|------|-----|-----|-----|", + ] + for key, _, _ in MATING_CONFIGS: + summary = block["summaries"][key] + lines.append( + f"| {key} | {summary['n']} | {fmt(summary['median'])} | {fmt(summary['mean'])} | " + f"{fmt(summary['std'])} | {fmt(summary['min'])} | {fmt(summary['max'])} |" + ) + lines.append("") + if block["versus_default"]: + lines += [ + "| Opt-in vs default | U₁ | MWU p | W | signed-rank p | median ratio (opt-in / default) | better median |", + "|-------------------|----|-------|---|----------------|---------------------------------|---------------|", + ] + for key, tests in block["versus_default"].items(): + lines.append( + f"| {key} | {fmt(tests['mannwhitney_U'])} | {fmt(tests['mannwhitney_p'])} | " + f"{fmt(tests['wilcoxon_W'])} | {fmt(tests['wilcoxon_p'])} | " + f"{fmt(tests['median_ratio'])} | {tests['better_median']} |" + ) + lines.append("") + for check in block["published_g6_checks"]: + if check["mismatches"]: + lines.append( + f"Published G6 check ({check['compared_to']}): " + f"{check['n_match']}/{check['n_seeds']} match. " + "Mismatches (measured, published): " + + ", ".join( + f"seed {item['seed']} {item['measured_g6']} vs {item['published_g6']}" + for item in check["mismatches"] + ) + + ". The published file was not rewritten." + ) + else: + lines.append( + f"Published G6 check ({check['compared_to']}): " + f"{check['n_match']}/{check['n_seeds']} match. The published file was not rewritten." + ) + lines.append("") + header = "| Seed | " + " | ".join(key for key, _, _ in MATING_CONFIGS) + " |" + rule = "|------|" + "|".join("------" for _ in MATING_CONFIGS) + "|" + lines += [header, rule] + for row in block["seeds"]: + cells = " | ".join(fmt(row[key]) for key, _, _ in MATING_CONFIGS) + lines.append(f"| {row['seed']} | {cells} |") + lines.append("") + + lines += [ + "## Notes", + "", + "- Package defaults are unchanged: `RankNicheDistance`, `IndependentWithReplacement`, single-layer Das–Dennis.", + "- Front files under `tools/oracle/out/` are gitignored.", + "- Hypothesis tests reuse `mannwhitney_u` and `wilcoxon_signed_rank` from `run_multiseed_wilcoxon.py`.", + "", + ] + return "\n".join(lines) + "\n" + + +def fmt_vec(values: list[float]) -> str: + return "[" + ", ".join(fmt(v) for v in values) + "]" + + +def cell(value: object) -> str: + if value is None: + return "skip: no feasible points" + return fmt(float(value)) # type: ignore[arg-type] + + +def main() -> int: + parser = argparse.ArgumentParser() + parser.add_argument("--seeds", type=int, default=15) + parser.add_argument("--seed-start", type=int, default=1) + parser.add_argument( + "--constrained", + nargs="*", + default=["osy", "tnk", "c1dtlz1"], + choices=sorted(CONSTRAINED), + ) + parser.add_argument( + "--mating", + nargs="*", + default=["zdt1", "zdt2", "dtlz2"], + choices=sorted(MATING), + ) + parser.add_argument("--skip-constrained", action="store_true") + parser.add_argument("--skip-mating", action="store_true") + parser.add_argument("--skip-build", action="store_true") + parser.add_argument("--report", type=Path, default=ROOT / "docs" / "NEW-SURFACES-RESULTS.md") + args = parser.parse_args() + if args.seeds < 1: + print("--seeds must be at least 1", file=sys.stderr) + return 1 + + constrained = [] if args.skip_constrained else list(args.constrained) + mating = [] if args.skip_mating else list(args.mating) + seeds = list(range(args.seed_start, args.seed_start + args.seeds)) + OUT.mkdir(parents=True, exist_ok=True) + + print("Checking constrained formulations against the C# fixtures...", flush=True) + assert_formulations() + + if not args.skip_build: + print("Building OracleCompare (Release)...", flush=True) + build = subprocess.run( + ["dotnet", "build", str(ORACLE_COMPARE), "-c", "Release", "--nologo", "-v", "q"], + cwd=str(ROOT), + ) + if build.returncode != 0: + return build.returncode + + dotnet_version = subprocess.check_output(["dotnet", "--version"], text=True, cwd=str(ROOT)).strip() + + command = "python tools/oracle/run_new_surfaces.py" + if constrained != ["osy", "tnk", "c1dtlz1"] or mating != ["zdt1", "zdt2", "dtlz2"] or args.seeds != 15 or args.seed_start != 1: + command = "python tools/oracle/run_new_surfaces.py " + " ".join(sys.argv[1:]) + + payload: dict = { + "measured_on": date.today().isoformat(), + "pymoo_version": pymoo.__version__, + "dotnet_version": dotnet_version, + "command": command, + "n_seeds": len(seeds), + "seed_start": args.seed_start, + "constrained": measure_constrained(constrained, seeds) if constrained else {}, + "mating": measure_mating(mating, seeds) if mating else {}, + } + + json_path = OUT / "new_surfaces_results.json" + json_path.write_text(json.dumps(payload, indent=2), encoding="utf-8") + report = markdown_report(payload) + args.report.parent.mkdir(parents=True, exist_ok=True) + args.report.write_text(report, encoding="utf-8") + print(f"\nwrote {json_path}") + print(f"wrote {args.report}") + return 0 + + +if __name__ == "__main__": + raise SystemExit(main())