SystemSimulationApp/skills/system-simulation/scripts/optimization_skill.py

#!/usr/bin/env python3
"""Deterministic single-objective optimizer for the system-simulation skill.

The optimizer treats the existing FastAPI backend as the authority for model
compilation and simulation.  It never rewrites the source project in place and
only changes explicitly selected continuous SI parameters in derived
candidates after the user confirms the complete plan.
"""

from __future__ import annotations

import argparse
import bisect
import copy
import csv
import hashlib
import io
import json
import math
import os
import platform
import random
import re
import sys
import tempfile
import time
import uuid
import xml.etree.ElementTree as ET
from dataclasses import dataclass, field
from pathlib import Path
from typing import Iterable, Mapping, Sequence


SCRIPT_DIRECTORY = Path(__file__).resolve().parent
if str(SCRIPT_DIRECTORY) not in sys.path:
    sys.path.insert(0, str(SCRIPT_DIRECTORY))

import simulation_skill as simulation  # noqa: E402


OPTIMIZATION_SCHEMA_VERSION = 1
OPTIMIZATION_RESULT_SCHEMA_VERSION = 2
STATISTIC_IMPLEMENTATION_VERSION = 1
CONTINUITY_POLICY = "user-assertion-explicit-false-veto-v1"
SEARCH_POLICY = "de-rand-1-bin-deferred-reflection-1d-endpoints-stagnation-v3"
NO_NEW_SUBMISSION_GENERATION_LIMIT = 3
TERMINAL_TREND_DIAGNOSTIC_VERSION = 1
TERMINAL_TREND_FRACTION = 0.05
TERMINAL_TREND_MAX_FRACTION = 0.20
TERMINAL_TREND_MIN_SAMPLES = 6
TERMINAL_TREND_RELATIVE_CHANGE_THRESHOLD = 0.01
TERMINAL_TREND_DIRECTIONAL_CONSISTENCY_THRESHOLD = 0.80
TERMINAL_TREND_RANGE_THRESHOLD = 0.02
MAX_DESIGN_VARIABLES = 16
MAX_RESPONSE_CONSTRAINTS = 16
MAX_SIMULATION_RUNS = 200
MAX_WALL_SECONDS = 7 * 24 * 60 * 60
RECOMMENDED_ALGORITHM_SEED = 0
RECOMMENDED_POPULATION_SIZE = 8
RECOMMENDED_MUTATION_FACTOR = 0.8
RECOMMENDED_CROSSOVER_PROBABILITY = 0.7
RECOMMENDED_VALIDATION_RELATIVE_TOLERANCE = 1e-8
RECOMMENDED_VALIDATION_ABSOLUTE_TOLERANCE = 0.0
OPTIMIZATION_ID_PATTERN = re.compile(r"^[A-Za-z0-9._-]{1,96}$")
IDENTIFIER_PATTERN = re.compile(r"^[A-Za-z][A-Za-z0-9._-]{0,63}$")
SUPPORTED_STATISTICS = {
    "final",
    "minimum",
    "maximum",
    "timeMean",
    "rms",
    "integral",
    "absoluteIntegral",
    "peakAbsolute",
}


def _metric_unit(series_unit: str, statistic_kind: str) -> str:
    if statistic_kind in {"integral", "absoluteIntegral"}:
        return f"({series_unit})*s" if series_unit else "s"
    return series_unit


class OptimizationError(simulation.SkillCliError):
    def __init__(
        self,
        code: str,
        message: str,
        details: object | None = None,
    ) -> None:
        super().__init__(code, message, exit_code=4, details=details)


class SearchStop(Exception):
    def __init__(self, reason: str) -> None:
        super().__init__(reason)
        self.reason = reason


@dataclass(frozen=True)
class TimeWindow:
    start: float
    end: float

    def as_dict(self) -> dict[str, float]:
        return {"start": self.start, "end": self.end}


@dataclass(frozen=True)
class StatisticSpec:
    kind: str
    window: TimeWindow | None

    def as_dict(self) -> dict[str, object]:
        return {
            "kind": self.kind,
            "window": self.window.as_dict() if self.window is not None else None,
        }


@dataclass(frozen=True)
class GoalSpec:
    kind: str
    value: float | None

    def as_dict(self) -> dict[str, object]:
        payload: dict[str, object] = {"kind": self.kind}
        if self.value is not None:
            payload["value"] = self.value
        return payload


@dataclass(frozen=True)
class ObjectiveSpec:
    result_key: str
    expected_unit: str
    statistic: StatisticSpec
    goal: GoalSpec

    def as_dict(self) -> dict[str, object]:
        return {
            "resultKey": self.result_key,
            "expectedUnit": self.expected_unit,
            "metricUnit": _metric_unit(self.expected_unit, self.statistic.kind),
            "statistic": self.statistic.as_dict(),
            "goal": self.goal.as_dict(),
        }


@dataclass(frozen=True)
class DesignVariableSpec:
    id: str
    component_id: str
    parameter: str
    unit: str
    lower: float
    upper: float

    def as_dict(self) -> dict[str, object]:
        return {
            "id": self.id,
            "componentId": self.component_id,
            "parameter": self.parameter,
            "unit": self.unit,
            "lower": self.lower,
            "upper": self.upper,
        }


@dataclass(frozen=True)
class ResponseConstraintSpec:
    id: str
    result_key: str
    expected_unit: str
    statistic: StatisticSpec
    lower: float | None
    upper: float | None
    tolerance: float
    scale: float

    def as_dict(self) -> dict[str, object]:
        return {
            "id": self.id,
            "resultKey": self.result_key,
            "expectedUnit": self.expected_unit,
            "metricUnit": _metric_unit(self.expected_unit, self.statistic.kind),
            "statistic": self.statistic.as_dict(),
            "lower": self.lower,
            "upper": self.upper,
            "tolerance": self.tolerance,
            "scale": self.scale,
        }


@dataclass(frozen=True)
class AlgorithmSpec:
    name: str
    seed: int
    population_size: int
    mutation_factor: float
    crossover_probability: float

    def as_dict(self) -> dict[str, object]:
        return {
            "name": self.name,
            "seed": self.seed,
            "populationSize": self.population_size,
            "mutationFactor": self.mutation_factor,
            "crossoverProbability": self.crossover_probability,
            "workers": 1,
            "strategy": "DE/rand/1/bin",
            "searchPolicy": SEARCH_POLICY,
            "updating": "deferred",
            "boundaryHandling": "reflection",
            "oneDimensionalEndpointSeedingPolicy": (
                "exactBoundsWhenOneDesignVariable"
            ),
            "noNewSubmissionGenerationLimit": (
                NO_NEW_SUBMISSION_GENERATION_LIMIT
            ),
        }


@dataclass(frozen=True)
class BudgetSpec:
    max_simulation_runs: int
    max_wall_seconds: float

    def as_dict(self) -> dict[str, object]:
        return {
            "maxSimulationRuns": self.max_simulation_runs,
            "maxWallSeconds": self.max_wall_seconds,
            "reservedFreshVerificationRuns": 1,
        }


@dataclass(frozen=True)
class ValidationSpec:
    relative_tolerance: float
    absolute_tolerance: float

    def as_dict(self) -> dict[str, object]:
        return {
            "freshRuns": 1,
            "relativeTolerance": self.relative_tolerance,
            "absoluteTolerance": self.absolute_tolerance,
        }


@dataclass(frozen=True)
class OptimizationSpec:
    objective: ObjectiveSpec
    design_variables: tuple[DesignVariableSpec, ...]
    constraints: tuple[ResponseConstraintSpec, ...]
    algorithm: AlgorithmSpec
    budget: BudgetSpec
    validation: ValidationSpec

    def as_dict(self) -> dict[str, object]:
        return {
            "optimizationSchemaVersion": OPTIMIZATION_SCHEMA_VERSION,
            "objective": self.objective.as_dict(),
            "designVariables": [item.as_dict() for item in self.design_variables],
            "constraints": [item.as_dict() for item in self.constraints],
            "algorithm": self.algorithm.as_dict(),
            "budget": self.budget.as_dict(),
            "validation": self.validation.as_dict(),
        }


@dataclass(frozen=True)
class ResolvedDesignVariable:
    spec: DesignVariableSpec
    label: str
    quantity: str
    initial: float
    was_explicit: bool
    original_value: object
    catalog_optimization_eligible: bool | None

    def as_dict(self) -> dict[str, object]:
        return {
            **self.spec.as_dict(),
            "label": self.label,
            "quantity": self.quantity,
            "initial": self.initial,
            "wasExplicit": self.was_explicit,
            "originalValue": self.original_value,
            "continuity": {
                "machineVerified": False,
                "userAssertionRequired": True,
                "catalogOptimizationEligible": (
                    self.catalog_optimization_eligible
                ),
                "editorAbsent": True,
                "optionsAbsent": True,
            },
        }


@dataclass(frozen=True)
class RuntimePlan:
    source: simulation.SourceFile
    spec_source: simulation.SourceFile
    spec: OptimizationSpec
    inspection: dict[str, object]
    baseline_xml: bytes
    variables: dict[str, dict[str, object]]
    resolved_design_variables: tuple[ResolvedDesignVariable, ...]
    output_directory: Path
    base_url: str
    timeout: float
    confirmation_token: str

    def presentation_dict(self) -> dict[str, object]:
        """Return the strict allowlist used for an ordinary user-facing plan."""
        objective = self.spec.objective
        objective_metadata = self.variables[objective.result_key]
        payload: dict[str, object] = {
            "objective": {
                "resultKey": objective.result_key,
                "componentId": objective_metadata.get("componentId"),
                "label": objective_metadata.get("label"),
                "quantity": objective_metadata.get("quantity"),
                "statistic": objective.statistic.as_dict(),
                "goal": objective.goal.as_dict(),
                "seriesUnit": objective.expected_unit,
                "metricUnit": _metric_unit(
                    objective.expected_unit,
                    objective.statistic.kind,
                ),
            },
            "designVariables": [
                {
                    "id": item.spec.id,
                    "componentId": item.spec.component_id,
                    "parameter": item.spec.parameter,
                    "label": item.label,
                    "quantity": item.quantity,
                    "current": item.initial,
                    "lower": item.spec.lower,
                    "upper": item.spec.upper,
                    "unit": item.spec.unit,
                }
                for item in self.resolved_design_variables
            ],
            "constraints": [
                {
                    "id": item.id,
                    "resultKey": item.result_key,
                    "componentId": self.variables[item.result_key].get(
                        "componentId"
                    ),
                    "label": self.variables[item.result_key].get("label"),
                    "quantity": self.variables[item.result_key].get("quantity"),
                    "statistic": item.statistic.as_dict(),
                    "lower": item.lower,
                    "upper": item.upper,
                    "tolerance": item.tolerance,
                    "seriesUnit": item.expected_unit,
                    "metricUnit": _metric_unit(
                        item.expected_unit,
                        item.statistic.kind,
                    ),
                }
                for item in self.spec.constraints
            ],
            "runLimits": {
                "maxSimulationBudgetSlots": (
                    self.spec.budget.max_simulation_runs
                ),
                "searchLaunchTimeLimitSeconds": (
                    self.spec.budget.max_wall_seconds
                ),
                "inFlightSimulationMayFinishAfterLimit": True,
                "freshVerificationSlotsReserved": 1,
                "verificationRunsOnlyIfFeasibleCandidateFound": True,
            },
            "assumptions": {
                "continuousLinearSiDesignVariableIds": [
                    item.spec.id for item in self.resolved_design_variables
                ],
                "portsTopologyModesUnchanged": True,
            },
            "outputDirectory": str(self.output_directory),
            "sourceWillBeOverwritten": False,
        }
        non_default_settings: dict[str, object] = {}
        search_settings: dict[str, object] = {}
        algorithm = self.spec.algorithm
        if algorithm.seed != RECOMMENDED_ALGORITHM_SEED:
            search_settings["randomSeed"] = algorithm.seed
        if algorithm.population_size != RECOMMENDED_POPULATION_SIZE:
            search_settings["populationSize"] = algorithm.population_size
        if algorithm.mutation_factor != RECOMMENDED_MUTATION_FACTOR:
            search_settings["mutationFactor"] = algorithm.mutation_factor
        if (
            algorithm.crossover_probability
            != RECOMMENDED_CROSSOVER_PROBABILITY
        ):
            search_settings["crossoverProbability"] = (
                algorithm.crossover_probability
            )
        if search_settings:
            non_default_settings["search"] = search_settings

        verification_settings: dict[str, object] = {}
        validation = self.spec.validation
        if (
            validation.relative_tolerance
            != RECOMMENDED_VALIDATION_RELATIVE_TOLERANCE
        ):
            verification_settings["relativeTolerance"] = (
                validation.relative_tolerance
            )
        if (
            validation.absolute_tolerance
            != RECOMMENDED_VALIDATION_ABSOLUTE_TOLERANCE
        ):
            verification_settings["absoluteTolerance"] = (
                validation.absolute_tolerance
            )
        if verification_settings:
            non_default_settings["verification"] = verification_settings

        constraint_scales = []
        for constraint in self.spec.constraints:
            bounds = [
                abs(bound)
                for bound in (constraint.lower, constraint.upper)
                if bound is not None
            ]
            recommended_scale = max(bounds, default=0.0) or 1.0
            if constraint.scale != recommended_scale:
                constraint_scales.append(
                    {"id": constraint.id, "scale": constraint.scale}
                )
        if constraint_scales:
            non_default_settings["constraintRankingScales"] = (
                constraint_scales
            )
        if non_default_settings:
            payload["nonDefaultSettings"] = non_default_settings

        search_run_limit = self.spec.budget.max_simulation_runs - 1
        population_size = self.spec.algorithm.population_size
        if (search_run_limit - population_size) // population_size < 1:
            payload["attention"] = [
                (
                    "The budget can cover the initial candidate set and reserve "
                    "fresh verification, but it cannot cover one complete search "
                    "update cycle if every candidate is unique."
                )
            ]
        return payload

    def public_dict(self) -> dict[str, object]:
        search_run_limit = self.spec.budget.max_simulation_runs - 1
        population_size = self.spec.algorithm.population_size
        full_generations_with_unique_candidates = max(
            0, (search_run_limit - population_size) // population_size
        )
        warnings: list[dict[str, object]] = []
        design_variable_ids = [
            item.spec.id for item in self.resolved_design_variables
        ]
        if full_generations_with_unique_candidates == 0:
            warnings.append(
                {
                    "code": "OPTIMIZATION_BUDGET_INITIAL_POPULATION_ONLY",
                    "message": (
                        "The simulation budget covers the initial population and "
                        "fresh verification, but no complete DE generation if every "
                        "candidate is unique."
                    ),
                }
            )
        return {
            "ok": True,
            "command": "optimization-plan",
            "confirmationRequired": True,
            "source": {
                "path": str(self.source.path),
                "sha256": self.source.sha256,
                "format": self.source.format,
            },
            "spec": {
                "path": str(self.spec_source.path),
                "sha256": self.spec_source.sha256,
                "resolved": self.spec.as_dict(),
            },
            "baselineSystemXmlSha256": _sha256(self.baseline_xml),
            "objectiveMetadata": self.variables[self.spec.objective.result_key],
            "objectiveMetricUnit": _metric_unit(
                self.spec.objective.expected_unit,
                self.spec.objective.statistic.kind,
            ),
            "constraintMetadata": [
                self.variables[item.result_key] for item in self.spec.constraints
            ],
            "designVariables": [
                item.as_dict() for item in self.resolved_design_variables
            ],
            "parameterContinuity": {
                "policy": CONTINUITY_POLICY,
                "machineVerified": False,
                "backendDeclarationRequired": False,
                "userAssertionRequiredFor": design_variable_ids,
            },
            "requiredAssertions": [
                {
                    "code": "OPTIMIZATION_CONTINUITY_USER_ASSERTION",
                    "designVariableIds": design_variable_ids,
                    "message": (
                        "Confirm that every listed design variable is continuous, "
                        "uses a linear SI scale, and cannot change ports, topology, "
                        "modes, or other discrete model structure. The backend does "
                        "not verify this property."
                    ),
                }
            ],
            "execution": {
                "algorithm": self.spec.algorithm.as_dict(),
                "searchRunLimit": search_run_limit,
                "totalRunLimit": self.spec.budget.max_simulation_runs,
                "fullGenerationsWithUniqueCandidates": (
                    full_generations_with_unique_candidates
                ),
                "sequential": True,
                "backendBaseUrl": self.base_url,
                "readTimeoutSeconds": self.timeout,
                "outputDirectory": str(self.output_directory),
                "sourceWillBeOverwritten": False,
            },
            "warnings": warnings,
            "planHash": self.confirmation_token,
            "confirmationToken": self.confirmation_token,
        }


@dataclass
class Trial:
    evaluation_id: int
    stage: str
    simulation_id: str
    parameters: dict[str, float]
    status: str
    duration_seconds: float
    objective_value: float | None = None
    objective_loss: float | None = None
    constraints: list[dict[str, object]] = field(default_factory=list)
    feasible: bool = False
    total_constraint_violation: float | None = None
    failure_code: str | None = None
    failure_message: str | None = None
    cache_reuse_count: int = 0

    def as_dict(self) -> dict[str, object]:
        return {
            "evaluationId": self.evaluation_id,
            "stage": self.stage,
            "simulationId": self.simulation_id,
            "parameters": self.parameters,
            "status": self.status,
            "durationSeconds": self.duration_seconds,
            "objectiveValue": self.objective_value,
            "objectiveLoss": self.objective_loss,
            "constraints": self.constraints,
            "feasible": self.feasible,
            "totalConstraintViolation": self.total_constraint_violation,
            "failureCode": self.failure_code,
            "failureMessage": self.failure_message,
            "cacheReuseCount": self.cache_reuse_count,
        }


def _sha256(data: bytes) -> str:
    return hashlib.sha256(data).hexdigest()


def _mapping(value: object, path: str) -> dict[str, object]:
    if not isinstance(value, dict):
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must be an object.",
        )
    return value


def _list(value: object, path: str) -> list[object]:
    if not isinstance(value, list):
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must be an array.",
        )
    return value


def _keys(
    value: Mapping[str, object],
    *,
    path: str,
    required: set[str],
    optional: set[str] = frozenset(),
) -> None:
    missing = sorted(required - set(value))
    unknown = sorted(set(value) - required - optional)
    if missing or unknown:
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} contains missing or unsupported fields.",
            {"path": path, "missing": missing, "unknown": unknown},
        )


def _text(value: object, path: str, *, identifier: bool = False) -> str:
    if not isinstance(value, str) or not value:
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must be a non-empty string.",
        )
    if identifier and IDENTIFIER_PATTERN.fullmatch(value) is None:
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must match {IDENTIFIER_PATTERN.pattern}.",
        )
    return value


def _number(value: object, path: str) -> float:
    if not isinstance(value, (int, float)) or isinstance(value, bool):
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must be a finite number.",
        )
    parsed = float(value)
    if not math.isfinite(parsed):
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must be a finite number.",
        )
    return parsed


def _integer(value: object, path: str) -> int:
    if not isinstance(value, int) or isinstance(value, bool):
        raise simulation.InputError(
            "OPTIMIZATION_SPEC_INVALID",
            f"{path} must be an integer.",
        )
    return value


def _optional_number(value: object, path: str) -> float | None:
    return None if value is None else _number(value, path)


def _parse_statistic(value: object, path: str) -> StatisticSpec:
    payload = _mapping(value, path)
    _keys(payload, path=path, required={"kind", "window"})
    kind = _text(payload["kind"], f"{path}.kind")
    if kind not in SUPPORTED_STATISTICS:
        raise simulation.InputError(
            "OPTIMIZATION_STATISTIC_UNSUPPORTED",
            f"{path}.kind is not supported.",
            {"supported": sorted(SUPPORTED_STATISTICS)},
        )
    raw_window = payload["window"]
    if raw_window is None:
        window = None
    else:
        window_payload = _mapping(raw_window, f"{path}.window")
        _keys(
            window_payload,
            path=f"{path}.window",
            required={"start", "end"},
        )
        start = _number(window_payload["start"], f"{path}.window.start")
        end = _number(window_payload["end"], f"{path}.window.end")
        if not start < end or not math.isfinite(end - start):
            raise simulation.InputError(
                "OPTIMIZATION_WINDOW_INVALID",
                f"{path}.window requires start < end with a finite numeric span.",
            )
        window = TimeWindow(start, end)
    return StatisticSpec(kind, window)


def _parse_goal(value: object, path: str) -> GoalSpec:
    payload = _mapping(value, path)
    _keys(payload, path=path, required={"kind"}, optional={"value"})
    kind = _text(payload["kind"], f"{path}.kind")
    if kind not in {"minimize", "maximize", "target"}:
        raise simulation.InputError(
            "OPTIMIZATION_GOAL_UNSUPPORTED",
            f"{path}.kind must be minimize, maximize, or target.",
        )
    if kind == "target":
        if "value" not in payload:
            raise simulation.InputError(
                "OPTIMIZATION_SPEC_INVALID",
                f"{path}.value is required for a target goal.",
            )
        target = _number(payload["value"], f"{path}.value")
    else:
        if "value" in payload:
            raise simulation.InputError(
                "OPTIMIZATION_SPEC_INVALID",
                f"{path}.value is only valid for a target goal.",
            )
        target = None
    return GoalSpec(kind, target)


def parse_optimization_spec(value: object) -> OptimizationSpec:
    payload = _mapping(value, "optimization spec")
    _keys(
        payload,
        path="optimization spec",
        required={
            "optimizationSchemaVersion",
            "objective",
            "designVariables",
            "constraints",
            "algorithm",
            "budget",
            "validation",
        },
    )
    schema_version = _integer(
        payload["optimizationSchemaVersion"],
        "optimizationSchemaVersion",
    )
    if schema_version != OPTIMIZATION_SCHEMA_VERSION:
        raise simulation.InputError(
            "OPTIMIZATION_SCHEMA_VERSION_UNSUPPORTED",
            f"optimizationSchemaVersion must be {OPTIMIZATION_SCHEMA_VERSION}.",
        )

    objective_payload = _mapping(payload["objective"], "objective")
    _keys(
        objective_payload,
        path="objective",
        required={"resultKey", "expectedUnit", "statistic", "goal"},
    )
    objective = ObjectiveSpec(
        result_key=_text(objective_payload["resultKey"], "objective.resultKey"),
        expected_unit=(
            objective_payload["expectedUnit"]
            if isinstance(objective_payload["expectedUnit"], str)
            else _text(objective_payload["expectedUnit"], "objective.expectedUnit")
        ),
        statistic=_parse_statistic(objective_payload["statistic"], "objective.statistic"),
        goal=_parse_goal(objective_payload["goal"], "objective.goal"),
    )

    raw_design_variables = _list(payload["designVariables"], "designVariables")
    if not 1 <= len(raw_design_variables) <= MAX_DESIGN_VARIABLES:
        raise simulation.InputError(
            "OPTIMIZATION_DESIGN_VARIABLE_COUNT_INVALID",
            f"designVariables must contain 1 to {MAX_DESIGN_VARIABLES} entries.",
        )
    design_variables: list[DesignVariableSpec] = []
    for index, raw_item in enumerate(raw_design_variables):
        path = f"designVariables[{index}]"
        item = _mapping(raw_item, path)
        _keys(
            item,
            path=path,
            required={"id", "componentId", "parameter", "unit", "lower", "upper"},
        )
        unit = item["unit"]
        if not isinstance(unit, str):
            raise simulation.InputError(
                "OPTIMIZATION_SPEC_INVALID",
                f"{path}.unit must be a string.",
            )
        lower = _number(item["lower"], f"{path}.lower")
        upper = _number(item["upper"], f"{path}.upper")
        if not lower < upper or not math.isfinite(upper - lower):
            raise simulation.InputError(
                "OPTIMIZATION_BOUNDS_INVALID",
                f"{path} requires lower < upper with a finite numeric span.",
            )
        design_variables.append(
            DesignVariableSpec(
                id=_text(item["id"], f"{path}.id", identifier=True),
                component_id=_text(item["componentId"], f"{path}.componentId"),
                parameter=_text(item["parameter"], f"{path}.parameter"),
                unit=unit,
                lower=lower,
                upper=upper,
            )
        )
    ids = [item.id for item in design_variables]
    targets = [(item.component_id, item.parameter) for item in design_variables]
    if len(ids) != len(set(ids)) or len(targets) != len(set(targets)):
        raise simulation.InputError(
            "OPTIMIZATION_DESIGN_VARIABLE_DUPLICATE",
            "Design-variable IDs and component/parameter targets must be unique.",
        )

    raw_constraints = _list(payload["constraints"], "constraints")
    if len(raw_constraints) > MAX_RESPONSE_CONSTRAINTS:
        raise simulation.InputError(
            "OPTIMIZATION_CONSTRAINT_COUNT_INVALID",
            f"constraints may contain at most {MAX_RESPONSE_CONSTRAINTS} entries.",
        )
    constraints: list[ResponseConstraintSpec] = []
    for index, raw_item in enumerate(raw_constraints):
        path = f"constraints[{index}]"
        item = _mapping(raw_item, path)
        _keys(
            item,
            path=path,
            required={
                "id",
                "resultKey",
                "expectedUnit",
                "statistic",
                "lower",
                "upper",
                "tolerance",
                "scale",
            },
        )
        unit = item["expectedUnit"]
        if not isinstance(unit, str):
            raise simulation.InputError(
                "OPTIMIZATION_SPEC_INVALID",
                f"{path}.expectedUnit must be a string.",
            )
        lower = _optional_number(item["lower"], f"{path}.lower")
        upper = _optional_number(item["upper"], f"{path}.upper")
        tolerance = _number(item["tolerance"], f"{path}.tolerance")
        scale = _number(item["scale"], f"{path}.scale")
        if lower is None and upper is None:
            raise simulation.InputError(
                "OPTIMIZATION_CONSTRAINT_INVALID",
                f"{path} requires at least one finite bound.",
            )
        if lower is not None and upper is not None and lower > upper:
            raise simulation.InputError(
                "OPTIMIZATION_CONSTRAINT_INVALID",
                f"{path} requires lower <= upper.",
            )
        if tolerance < 0 or scale <= 0:
            raise simulation.InputError(
                "OPTIMIZATION_CONSTRAINT_INVALID",
                f"{path}.tolerance must be non-negative and scale must be positive.",
            )
        constraints.append(
            ResponseConstraintSpec(
                id=_text(item["id"], f"{path}.id", identifier=True),
                result_key=_text(item["resultKey"], f"{path}.resultKey"),
                expected_unit=unit,
                statistic=_parse_statistic(item["statistic"], f"{path}.statistic"),
                lower=lower,
                upper=upper,
                tolerance=tolerance,
                scale=scale,
            )
        )
    constraint_ids = [item.id for item in constraints]
    if len(constraint_ids) != len(set(constraint_ids)):
        raise simulation.InputError(
            "OPTIMIZATION_CONSTRAINT_DUPLICATE",
            "Constraint IDs must be unique.",
        )

    algorithm_payload = _mapping(payload["algorithm"], "algorithm")
    _keys(
        algorithm_payload,
        path="algorithm",
        required={
            "name",
            "seed",
            "populationSize",
            "mutationFactor",
            "crossoverProbability",
        },
    )
    algorithm_name = _text(algorithm_payload["name"], "algorithm.name")
    if algorithm_name != "differentialEvolution":
        raise simulation.InputError(
            "OPTIMIZATION_ALGORITHM_UNSUPPORTED",
            "Only differentialEvolution is supported in optimization schema 1.",
        )
    seed = _integer(algorithm_payload["seed"], "algorithm.seed")
    population_size = _integer(
        algorithm_payload["populationSize"], "algorithm.populationSize"
    )
    mutation_factor = _number(
        algorithm_payload["mutationFactor"], "algorithm.mutationFactor"
    )
    crossover_probability = _number(
        algorithm_payload["crossoverProbability"],
        "algorithm.crossoverProbability",
    )
    if not 0 <= seed <= 2**32 - 1:
        raise simulation.InputError(
            "OPTIMIZATION_ALGORITHM_INVALID", "algorithm.seed is out of range."
        )
    if not 4 <= population_size <= 50:
        raise simulation.InputError(
            "OPTIMIZATION_ALGORITHM_INVALID",
            "algorithm.populationSize must be between 4 and 50.",
        )
    if not 0 < mutation_factor <= 2 or not 0 <= crossover_probability <= 1:
        raise simulation.InputError(
            "OPTIMIZATION_ALGORITHM_INVALID",
            "Mutation must be in (0, 2] and crossover probability in [0, 1].",
        )
    algorithm = AlgorithmSpec(
        algorithm_name,
        seed,
        population_size,
        mutation_factor,
        crossover_probability,
    )

    budget_payload = _mapping(payload["budget"], "budget")
    _keys(
        budget_payload,
        path="budget",
        required={"maxSimulationRuns", "maxWallSeconds"},
    )
    max_runs = _integer(
        budget_payload["maxSimulationRuns"], "budget.maxSimulationRuns"
    )
    max_wall_seconds = _number(
        budget_payload["maxWallSeconds"], "budget.maxWallSeconds"
    )
    if not population_size + 1 <= max_runs <= MAX_SIMULATION_RUNS:
        raise simulation.InputError(
            "OPTIMIZATION_BUDGET_INVALID",
            "maxSimulationRuns must fit the initial population plus one fresh "
            f"verification run and may not exceed {MAX_SIMULATION_RUNS}.",
        )
    if not 10 <= max_wall_seconds <= MAX_WALL_SECONDS:
        raise simulation.InputError(
            "OPTIMIZATION_BUDGET_INVALID",
            f"maxWallSeconds must be between 10 and {MAX_WALL_SECONDS}.",
        )
    budget = BudgetSpec(max_runs, max_wall_seconds)

    validation_payload = _mapping(payload["validation"], "validation")
    _keys(
        validation_payload,
        path="validation",
        required={"relativeTolerance", "absoluteTolerance"},
    )
    relative_tolerance = _number(
        validation_payload["relativeTolerance"],
        "validation.relativeTolerance",
    )
    absolute_tolerance = _number(
        validation_payload["absoluteTolerance"],
        "validation.absoluteTolerance",
    )
    if relative_tolerance < 0 or absolute_tolerance < 0:
        raise simulation.InputError(
            "OPTIMIZATION_VALIDATION_INVALID",
            "Validation tolerances must be non-negative.",
        )
    validation = ValidationSpec(relative_tolerance, absolute_tolerance)

    return OptimizationSpec(
        objective=objective,
        design_variables=tuple(design_variables),
        constraints=tuple(constraints),
        algorithm=algorithm,
        budget=budget,
        validation=validation,
    )


def _validate_output_target(path_text: str) -> Path:
    path = Path(path_text).expanduser().resolve(strict=False)
    if path.exists():
        if not path.is_dir():
            raise simulation.InputError(
                "OUTPUT_DIR_NOT_DIRECTORY",
                f"Output path is not a directory: {path}",
            )
        try:
            if any(path.iterdir()):
                raise simulation.InputError(
                    "OUTPUT_DIR_NOT_EMPTY",
                    "Optimization output directory must be new or empty.",
                    {"path": str(path)},
                )
        except OSError as exc:
            raise simulation.ArtifactError(
                "OUTPUT_DIR_UNREADABLE", str(exc), {"path": str(path)}
            ) from exc
    return path


def _validate_plan_file_target(
    path_text: str,
    output_directory_text: str,
) -> Path:
    requested_path = Path(path_text).expanduser()
    if requested_path.exists() or requested_path.is_symlink():
        raise simulation.InputError(
            "OPTIMIZATION_PLAN_FILE_ALREADY_EXISTS",
            "The saved optimization plan path must be new.",
            {"path": str(requested_path.absolute())},
        )
    try:
        path = requested_path.resolve(strict=False)
        output_directory = Path(output_directory_text).expanduser().resolve(
            strict=False
        )
    except (OSError, RuntimeError) as exc:
        raise simulation.InputError(
            "OPTIMIZATION_PLAN_FILE_PATH_INVALID",
            "The saved plan file path could not be resolved safely.",
            {"path": str(requested_path.absolute())},
        ) from exc
    if path.exists():
        raise simulation.InputError(
            "OPTIMIZATION_PLAN_FILE_ALREADY_EXISTS",
            "The saved optimization plan path must be new.",
            {"path": str(path)},
        )
    if (
        path == output_directory
        or path.is_relative_to(output_directory)
        or output_directory.is_relative_to(path)
    ):
        raise simulation.InputError(
            "OPTIMIZATION_PLAN_FILE_OUTPUT_CONFLICT",
            "The saved plan file and optimization output directory must be separate sibling paths.",
            {
                "planFile": str(path),
                "outputDirectory": str(output_directory),
            },
        )
    return path


def _write_new_private_file(path: Path, data: bytes) -> None:
    """Exclusively create a plan receipt readable only by its owner."""
    created = False
    descriptor: int | None = None
    try:
        path.parent.mkdir(parents=True, exist_ok=True)
        flags = os.O_WRONLY | os.O_CREAT | os.O_EXCL
        if hasattr(os, "O_BINARY"):
            flags |= os.O_BINARY
        if hasattr(os, "O_NOFOLLOW"):
            flags |= os.O_NOFOLLOW
        descriptor = os.open(path, flags, 0o600)
        created = True
        if hasattr(os, "fchmod"):
            os.fchmod(descriptor, 0o600)
        with os.fdopen(descriptor, "wb") as handle:
            descriptor = None
            handle.write(data)
            handle.flush()
            os.fsync(handle.fileno())
    except OSError as exc:
        if descriptor is not None:
            try:
                os.close(descriptor)
            except OSError:
                pass
        try:
            if created and path.exists():
                path.unlink()
        except OSError:
            pass
        raise simulation.ArtifactError(
            "OUTPUT_WRITE_FAILED", str(exc), {"path": str(path)}
        ) from exc


def _parameter_contracts(
    inspection: Mapping[str, object],
) -> dict[tuple[str, str], tuple[dict[str, object], dict[str, object] | None]]:
    system = inspection.get("system")
    if not isinstance(system, Mapping):
        return {}
    details = system.get("componentDetails")
    if not isinstance(details, list):
        return {}
    contracts: dict[
        tuple[str, str], tuple[dict[str, object], dict[str, object] | None]
    ] = {}
    for detail in details:
        if not isinstance(detail, dict):
            continue
        component_id = detail.get("id")
        compiled = detail.get("compiled")
        if not isinstance(component_id, str) or not isinstance(compiled, dict):
            continue
        source_component = detail.get("source")
        source_mapping = source_component if isinstance(source_component, dict) else None
        parameters = compiled.get("parameters")
        if not isinstance(parameters, list):
            continue
        for parameter in parameters:
            if isinstance(parameter, dict) and isinstance(parameter.get("name"), str):
                contracts[(component_id, str(parameter["name"]))] = (
                    parameter,
                    source_mapping,
                )
    return contracts


def _source_parameter(
    source_component: Mapping[str, object] | None,
    parameter_name: str,
) -> tuple[bool, object]:
    if not isinstance(source_component, Mapping):
        return False, None
    data = source_component.get("data")
    if not isinstance(data, Mapping):
        return False, None
    parameters = data.get("parameters")
    if not isinstance(parameters, Mapping) or parameter_name not in parameters:
        return False, None
    return True, parameters[parameter_name]


def _resolve_design_variables(
    spec: OptimizationSpec,
    inspection: Mapping[str, object],
) -> tuple[ResolvedDesignVariable, ...]:
    contracts = _parameter_contracts(inspection)
    resolved: list[ResolvedDesignVariable] = []
    for variable in spec.design_variables:
        target = (variable.component_id, variable.parameter)
        contract_pair = contracts.get(target)
        if contract_pair is None:
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_UNKNOWN",
                "A design variable does not identify a compiled component parameter.",
                {"componentId": variable.component_id, "parameter": variable.parameter},
            )
        contract, source_component = contract_pair
        catalog_optimization_eligible: bool | None = None
        if "optimizationEligible" in contract:
            optimization_eligible = contract["optimizationEligible"]
            if type(optimization_eligible) is not bool:
                raise simulation.BackendError(
                    "BACKEND_PARAMETER_CONTRACT_INVALID",
                    "The optional optimizationEligible parameter field must be a boolean.",
                    {
                        "componentId": variable.component_id,
                        "parameter": variable.parameter,
                        "optimizationEligible": optimization_eligible,
                    },
                )
            if not optimization_eligible:
                raise simulation.InputError(
                    "OPTIMIZATION_PARAMETER_INELIGIBLE",
                    "The backend contract explicitly marks this parameter as "
                    "ineligible for continuous optimization.",
                    {
                        "componentId": variable.component_id,
                        "parameter": variable.parameter,
                        "optimizationEligible": False,
                    },
                )
            catalog_optimization_eligible = True
        if "editor" in contract:
            editor = contract["editor"]
            if not isinstance(editor, str) or not editor:
                raise simulation.BackendError(
                    "BACKEND_PARAMETER_CONTRACT_INVALID",
                    "The optional editor parameter field must be a non-empty string.",
                    {
                        "componentId": variable.component_id,
                        "parameter": variable.parameter,
                        "editor": editor,
                    },
                )
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_DISCRETE_METADATA",
                "A parameter with an editor cannot be used for continuous optimization.",
                {
                    "componentId": variable.component_id,
                    "parameter": variable.parameter,
                    "optimizationEligible": contract.get("optimizationEligible"),
                    "editor": editor,
                },
            )
        if "options" in contract:
            options = contract["options"]
            if not isinstance(options, list) or not options:
                raise simulation.BackendError(
                    "BACKEND_PARAMETER_CONTRACT_INVALID",
                    "The optional options parameter field must be a non-empty array.",
                    {
                        "componentId": variable.component_id,
                        "parameter": variable.parameter,
                        "options": options,
                    },
                )
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_DISCRETE_METADATA",
                "A parameter with discrete options cannot be used for continuous optimization.",
                {
                    "componentId": variable.component_id,
                    "parameter": variable.parameter,
                    "optimizationEligible": contract.get("optimizationEligible"),
                    "options": options,
                },
            )
        contract_unit = contract.get("unit")
        if not isinstance(contract_unit, str) or contract_unit != variable.unit:
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_UNIT_MISMATCH",
                "The design-variable unit must exactly match the backend SI contract.",
                {
                    "id": variable.id,
                    "expected": contract_unit,
                    "received": variable.unit,
                },
            )
        current = contract.get("value")
        initial = _number(
            current,
            f"compiled parameter {variable.component_id}.{variable.parameter}.value",
        )
        if not variable.lower <= initial <= variable.upper:
            raise simulation.InputError(
                "OPTIMIZATION_INITIAL_OUTSIDE_BOUNDS",
                "The current SI parameter value must lie inside the declared search bounds.",
                {"id": variable.id, "initial": initial},
            )
        minimum = contract.get("minimum")
        maximum = contract.get("maximum")
        if isinstance(minimum, (int, float)) and not isinstance(minimum, bool):
            minimum_value = float(minimum)
            minimum_exclusive = contract.get("minimumExclusive") is True
            if variable.lower < minimum_value or (
                minimum_exclusive and variable.lower <= minimum_value
            ):
                raise simulation.InputError(
                    "OPTIMIZATION_BOUNDS_OUTSIDE_CONTRACT",
                    "A lower search bound violates the registered parameter contract.",
                    {
                        "id": variable.id,
                        "minimum": minimum_value,
                        "minimumExclusive": minimum_exclusive,
                    },
                )
        if isinstance(maximum, (int, float)) and not isinstance(maximum, bool):
            maximum_value = float(maximum)
            if variable.upper > maximum_value:
                raise simulation.InputError(
                    "OPTIMIZATION_BOUNDS_OUTSIDE_CONTRACT",
                    "An upper search bound violates the registered parameter contract.",
                    {"id": variable.id, "maximum": maximum_value},
                )
        was_explicit, original_value = _source_parameter(
            source_component, variable.parameter
        )
        resolved.append(
            ResolvedDesignVariable(
                spec=variable,
                label=str(contract.get("label") or variable.parameter),
                quantity=str(contract.get("quantity") or "dimensionless"),
                initial=initial,
                was_explicit=was_explicit,
                original_value=original_value,
                catalog_optimization_eligible=catalog_optimization_eligible,
            )
        )
    return tuple(resolved)


def _validate_result_contracts(
    spec: OptimizationSpec,
    variables: Mapping[str, Mapping[str, object]],
) -> None:
    requested: list[tuple[str, str, str]] = [
        (
            "objective",
            spec.objective.result_key,
            spec.objective.expected_unit,
        )
    ]
    requested.extend(
        (f"constraint {constraint.id}", constraint.result_key, constraint.expected_unit)
        for constraint in spec.constraints
    )
    for label, result_key, expected_unit in requested:
        metadata = variables.get(result_key)
        if metadata is None:
            raise simulation.InputError(
                "OPTIMIZATION_RESULT_VARIABLE_UNKNOWN",
                f"The {label} result key is not declared by the compiled model.",
                {"resultKey": result_key},
            )
        actual_unit = metadata.get("unit")
        if not isinstance(actual_unit, str):
            actual_unit = ""
        if actual_unit != expected_unit:
            raise simulation.InputError(
                "OPTIMIZATION_RESULT_UNIT_MISMATCH",
                f"The {label} expected unit does not match result metadata.",
                {
                    "resultKey": result_key,
                    "expected": expected_unit,
                    "received": actual_unit,
                },
            )


def _confirmation_token(
    *,
    source_sha256: str,
    spec_sha256: str,
    baseline_xml_sha256: str,
    output_directory: Path,
    base_url: str,
    timeout: float,
) -> str:
    payload = {
        "contract": "system-simulation-optimization-confirmation-v3",
        "continuityPolicy": CONTINUITY_POLICY,
        "searchPolicy": SEARCH_POLICY,
        "sourceSha256": source_sha256,
        "specSha256": spec_sha256,
        "baselineSystemXmlSha256": baseline_xml_sha256,
        "outputDirectory": str(output_directory),
        "baseUrl": base_url,
        "readTimeoutSeconds": timeout,
    }
    canonical = json.dumps(
        payload,
        ensure_ascii=False,
        sort_keys=True,
        separators=(",", ":"),
        allow_nan=False,
    ).encode("utf-8")
    return _sha256(canonical)


def build_runtime_plan(
    input_path: str,
    spec_path: str,
    output_directory: str,
    *,
    base_url: str,
    timeout: float,
) -> RuntimePlan:
    source = simulation.load_source(input_path, "json")
    if source.format != "json":
        raise simulation.InputError(
            "OPTIMIZATION_SOURCE_FORMAT_UNSUPPORTED",
            "Optimization schema 1 requires an editable ReactFlow project JSON v1 source.",
        )
    spec_source = simulation.load_source(spec_path, "json")
    spec = parse_optimization_spec(spec_source.parsed)
    target_output = _validate_output_target(output_directory)
    component_ids = list(
        dict.fromkeys(item.component_id for item in spec.design_variables)
    )
    inspection = simulation.inspect_source(
        source,
        base_url=base_url,
        timeout=timeout,
        component_ids=component_ids,
        variable_limit=1_000_000,
    )
    variables = simulation._available_variables(inspection)
    _validate_result_contracts(spec, variables)
    resolved = _resolve_design_variables(spec, inspection)
    baseline_xml = simulation._json_to_xml(
        source,
        base_url=base_url,
        timeout=timeout,
    )
    token = _confirmation_token(
        source_sha256=source.sha256,
        spec_sha256=spec_source.sha256,
        baseline_xml_sha256=_sha256(baseline_xml),
        output_directory=target_output,
        base_url=base_url,
        timeout=timeout,
    )
    return RuntimePlan(
        source=source,
        spec_source=spec_source,
        spec=spec,
        inspection=inspection,
        baseline_xml=baseline_xml,
        variables=variables,
        resolved_design_variables=resolved,
        output_directory=target_output,
        base_url=base_url,
        timeout=timeout,
        confirmation_token=token,
    )


def command_plan(args: argparse.Namespace) -> int:
    present = bool(getattr(args, "present", False))
    plan_file_text = getattr(args, "plan_file", None)
    if present != bool(plan_file_text):
        raise simulation.InputError(
            "OPTIMIZATION_PLAN_PRESENTATION_OPTIONS_REQUIRED",
            "--present and --plan-file must be supplied together.",
        )
    plan_file = (
        _validate_plan_file_target(plan_file_text, args.output_dir)
        if plan_file_text
        else None
    )
    plan = build_runtime_plan(
        args.input,
        args.spec,
        args.output_dir,
        base_url=args.base_url,
        timeout=args.timeout,
    )
    audit_payload = plan.public_dict()
    if plan_file is not None:
        serialized = (
            json.dumps(
                audit_payload,
                ensure_ascii=False,
                allow_nan=False,
                indent=2,
            )
            + "\n"
        ).encode("utf-8")
        _write_new_private_file(plan_file, serialized)
    if present:
        simulation.emit_json(
            {
                "ok": True,
                "command": "optimization-plan",
                "view": "presentation",
                "confirmationRequired": True,
                "presentation": plan.presentation_dict(),
            }
        )
    else:
        simulation.emit_json(audit_payload)
    return 0


def _loss(value: float, goal: GoalSpec) -> float:
    if goal.kind == "minimize":
        return value
    if goal.kind == "maximize":
        return -value
    assert goal.value is not None
    return abs(value - goal.value)


def _interpolate(
    times: Sequence[float],
    values: Sequence[float],
    point: float,
) -> float:
    index = bisect.bisect_left(times, point)
    if index < len(times) and times[index] == point:
        return values[index]
    if index == 0 or index == len(times):
        raise OptimizationError(
            "OPTIMIZATION_WINDOW_NOT_COVERED",
            "A completed simulation does not cover a requested statistic window.",
            {"point": point, "start": times[0], "end": times[-1]},
        )
    left_time = times[index - 1]
    right_time = times[index]
    time_span = right_time - left_time
    if not math.isfinite(time_span):
        raise OptimizationError(
            "OPTIMIZATION_METRIC_OVERFLOW",
            "Interpolation time span overflowed the finite numeric contract.",
        )
    fraction = (point - left_time) / time_span
    interpolated = (
        (1.0 - fraction) * values[index - 1] + fraction * values[index]
    )
    if not math.isfinite(interpolated):
        raise OptimizationError(
            "OPTIMIZATION_METRIC_OVERFLOW",
            "Interpolated result value overflowed the finite numeric contract.",
        )
    return interpolated


def _windowed_series(
    times: Sequence[float],
    values: Sequence[float],
    window: TimeWindow | None,
) -> tuple[list[float], list[float]]:
    if len(times) != len(values) or not times:
        raise OptimizationError(
            "OPTIMIZATION_RESULT_SERIES_INVALID",
            "Objective and constraint series must be non-empty and match time length.",
        )
    converted_times: list[float] = []
    converted_values: list[float] = []
    previous: float | None = None
    for raw_time, raw_value in zip(times, values):
        if (
            not isinstance(raw_time, (int, float))
            or isinstance(raw_time, bool)
            or not math.isfinite(float(raw_time))
            or not isinstance(raw_value, (int, float))
            or isinstance(raw_value, bool)
            or not math.isfinite(float(raw_value))
        ):
            raise OptimizationError(
                "OPTIMIZATION_RESULT_SERIES_INVALID",
                "Objective and constraint series must contain finite numbers.",
            )
        current_time = float(raw_time)
        if previous is not None and current_time <= previous:
            raise OptimizationError(
                "OPTIMIZATION_RESULT_TIME_INVALID",
                "Simulation time must be strictly increasing for optimization metrics.",
            )
        converted_times.append(current_time)
        converted_values.append(float(raw_value))
        previous = current_time
    if window is None:
        return converted_times, converted_values
    if window.start < converted_times[0] or window.end > converted_times[-1]:
        raise OptimizationError(
            "OPTIMIZATION_WINDOW_NOT_COVERED",
            "A completed simulation does not cover a requested statistic window.",
            {
                "window": window.as_dict(),
                "resultStart": converted_times[0],
                "resultEnd": converted_times[-1],
            },
        )
    start_index = bisect.bisect_left(converted_times, window.start)
    end_index = bisect.bisect_right(converted_times, window.end)
    selected_times = converted_times[start_index:end_index]
    selected_values = converted_values[start_index:end_index]
    if not selected_times or selected_times[0] != window.start:
        selected_times.insert(0, window.start)
        selected_values.insert(
            0,
            _interpolate(converted_times, converted_values, window.start),
        )
    if selected_times[-1] != window.end:
        selected_times.append(window.end)
        selected_values.append(
            _interpolate(converted_times, converted_values, window.end)
        )
    return selected_times, selected_values


def statistic_value(
    times: Sequence[float],
    values: Sequence[float],
    statistic: StatisticSpec,
) -> float:
    selected_times, selected_values = _windowed_series(
        times, values, statistic.window
    )
    if statistic.kind == "final":
        return selected_values[-1]
    if statistic.kind == "minimum":
        return min(selected_values)
    if statistic.kind == "maximum":
        return max(selected_values)
    if statistic.kind == "peakAbsolute":
        return max(abs(value) for value in selected_values)
    if len(selected_times) < 2 or selected_times[-1] <= selected_times[0]:
        raise OptimizationError(
            "OPTIMIZATION_STATISTIC_WINDOW_TOO_SHORT",
            f"{statistic.kind} requires a positive-duration time window.",
        )
    duration = selected_times[-1] - selected_times[0]
    if not math.isfinite(duration):
        raise OptimizationError(
            "OPTIMIZATION_METRIC_OVERFLOW",
            "Statistic time span overflowed the finite numeric contract.",
        )

    def finite_sum(terms: Iterable[float]) -> float:
        try:
            value = math.fsum(terms)
        except (OverflowError, ValueError) as exc:
            raise OptimizationError(
                "OPTIMIZATION_METRIC_OVERFLOW",
                "Statistic integration overflowed the finite numeric contract.",
            ) from exc
        if not math.isfinite(value):
            raise OptimizationError(
                "OPTIMIZATION_METRIC_OVERFLOW",
                "Statistic integration overflowed the finite numeric contract.",
            )
        return value

    if statistic.kind == "integral":
        return finite_sum(
            0.5 * (left_value + right_value) * (right_time - left_time)
            for left_time, right_time, left_value, right_value in zip(
                selected_times,
                selected_times[1:],
                selected_values,
                selected_values[1:],
            )
        )
    if statistic.kind == "absoluteIntegral":
        return finite_sum(
            0.5 * (abs(left_value) + abs(right_value)) * (right_time - left_time)
            for left_time, right_time, left_value, right_value in zip(
                selected_times,
                selected_times[1:],
                selected_values,
                selected_values[1:],
            )
        )
    if statistic.kind == "timeMean":
        integral = statistic_value(
            selected_times,
            selected_values,
            StatisticSpec("integral", None),
        )
        mean = integral / duration
        if not math.isfinite(mean):
            raise OptimizationError(
                "OPTIMIZATION_METRIC_OVERFLOW",
                "Time-mean calculation overflowed the finite numeric contract.",
            )
        return mean
    if statistic.kind == "rms":
        square_integral = finite_sum(
            0.5 * (left_value**2 + right_value**2) * (right_time - left_time)
            for left_time, right_time, left_value, right_value in zip(
                selected_times,
                selected_times[1:],
                selected_values,
                selected_values[1:],
            )
        )
        rms = math.sqrt(max(square_integral / duration, 0.0))
        if not math.isfinite(rms):
            raise OptimizationError(
                "OPTIMIZATION_METRIC_OVERFLOW",
                "RMS calculation overflowed the finite numeric contract.",
            )
        return rms
    raise AssertionError(f"Unhandled statistic kind: {statistic.kind}")


def _planned_objective_end(plan: RuntimePlan) -> float | None:
    window = plan.spec.objective.statistic.window
    if window is not None:
        return window.end
    system = plan.inspection.get("system")
    simulation_settings = (
        system.get("simulation") if isinstance(system, Mapping) else None
    )
    if isinstance(simulation_settings, Mapping):
        for key in ("t_stop", "tStop", "endTime", "stopTime"):
            value = simulation_settings.get(key)
            if (
                isinstance(value, (int, float))
                and not isinstance(value, bool)
                and math.isfinite(float(value))
            ):
                return float(value)
    try:
        root = ET.fromstring(plan.baseline_xml)
    except ET.ParseError:
        return None
    simulation_element = root.find("./Simulation")
    if simulation_element is None:
        return None
    for attribute in ("tStop", "endTime", "stopTime"):
        raw_value = simulation_element.get(attribute)
        if raw_value is None:
            continue
        try:
            value = float(raw_value)
        except ValueError:
            continue
        if math.isfinite(value):
            return value
    return None


def _objective_endpoint_trend(
    result: Mapping[str, object],
    objective: ObjectiveSpec,
    *,
    expected_end: float | None = None,
) -> dict[str, object]:
    """Describe recent movement for a final-value objective without claiming steady state."""
    base: dict[str, object] = {
        "diagnosticVersion": TERMINAL_TREND_DIAGNOSTIC_VERSION,
        "applicable": objective.statistic.kind == "final",
        "steadyStateProven": False,
        "heuristic": {
            "preferredTailFraction": TERMINAL_TREND_FRACTION,
            "maximumTailFraction": TERMINAL_TREND_MAX_FRACTION,
            "minimumSamples": TERMINAL_TREND_MIN_SAMPLES,
            "relativeChangeThreshold": (
                TERMINAL_TREND_RELATIVE_CHANGE_THRESHOLD
            ),
            "directionalConsistencyThreshold": (
                TERMINAL_TREND_DIRECTIONAL_CONSISTENCY_THRESHOLD
            ),
            "relativeRangeThreshold": TERMINAL_TREND_RANGE_THRESHOLD,
        },
    }
    if objective.statistic.kind != "final":
        return {**base, "status": "notApplicable"}
    if result.get("status") != "completed" or result.get("success") is not True:
        return {
            **base,
            "status": "unavailable",
            "reason": "freshVerificationNotCompleted",
        }
    raw_series = result.get("series")
    if not isinstance(raw_series, Mapping):
        return {
            **base,
            "status": "unavailable",
            "reason": "freshVerificationSeriesMissing",
        }
    raw_times = raw_series.get("time")
    raw_values = raw_series.get(objective.result_key)
    if not isinstance(raw_times, list) or not isinstance(raw_values, list):
        return {
            **base,
            "status": "unavailable",
            "reason": "freshVerificationSeriesMissing",
        }
    try:
        times, values = _windowed_series(
            raw_times,
            raw_values,
            objective.statistic.window,
        )
    except OptimizationError as exc:
        return {
            **base,
            "status": "unavailable",
            "reason": "freshVerificationSeriesInvalid",
            "diagnosticErrorCode": exc.code,
        }
    required_end = (
        objective.statistic.window.end
        if objective.statistic.window is not None
        else expected_end
    )
    if required_end is None:
        return {
            **base,
            "status": "unavailable",
            "reason": "plannedEndpointUnavailable",
        }
    end_tolerance = max(1e-12, abs(required_end) * 1e-12)
    if not math.isclose(
        times[-1], required_end, rel_tol=1e-12, abs_tol=end_tolerance
    ):
        return {
            **base,
            "status": "unavailable",
            "reason": "plannedEndpointNotCovered",
            "expectedEnd": required_end,
            "actualSeriesEnd": times[-1],
        }
    if len(times) < TERMINAL_TREND_MIN_SAMPLES:
        return {
            **base,
            "status": "insufficientData",
            "availableSamples": len(times),
        }
    duration = times[-1] - times[0]
    if duration <= 0.0 or not math.isfinite(duration):
        return {
            **base,
            "status": "insufficientData",
            "availableSamples": len(times),
        }
    preferred_start = times[-1] - TERMINAL_TREND_FRACTION * duration
    start_index = bisect.bisect_left(times, preferred_start)
    start_index = min(start_index, len(times) - TERMINAL_TREND_MIN_SAMPLES)
    tail_duration = times[-1] - times[start_index]
    if tail_duration > TERMINAL_TREND_MAX_FRACTION * duration:
        return {
            **base,
            "status": "insufficientData",
            "availableSamples": len(times),
        }
    tail_values = values[start_index:]
    increments = [
        right - left for left, right in zip(tail_values, tail_values[1:])
    ]
    positive = sum(item > 0.0 for item in increments)
    negative = sum(item < 0.0 for item in increments)
    nonzero_increments = positive + negative
    directional_consistency = (
        max(positive, negative) / nonzero_increments
        if nonzero_increments
        else 0.0
    )
    delta = tail_values[-1] - tail_values[0]
    maximum_absolute_full = max(abs(item) for item in values)
    scale = max(
        max(abs(item) for item in tail_values),
        1e-6 * maximum_absolute_full,
        sys.float_info.min,
    )
    relative_change = delta / scale
    relative_range = (max(tail_values) - min(tail_values)) / scale
    average_slope = delta / tail_duration
    if any(
        not math.isfinite(item)
        for item in (delta, scale, relative_change, relative_range, average_slope)
    ):
        return {
            **base,
            "status": "unavailable",
            "reason": "terminalTrendNumericOverflow",
        }
    directional_change_detected = (
        abs(relative_change) >= TERMINAL_TREND_RELATIVE_CHANGE_THRESHOLD
        and directional_consistency
        >= TERMINAL_TREND_DIRECTIONAL_CONSISTENCY_THRESHOLD
    )
    tail_variability_detected = relative_range >= TERMINAL_TREND_RANGE_THRESHOLD
    material_change = directional_change_detected or tail_variability_detected
    direction = "flat"
    if material_change and not directional_change_detected:
        direction = "fluctuating"
    elif delta > 0.0:
        direction = "increasing"
    elif delta < 0.0:
        direction = "decreasing"
    detection_reasons: list[str] = []
    if directional_change_detected:
        detection_reasons.append("directionalChange")
    if tail_variability_detected:
        detection_reasons.append("tailVariability")
    return {
        **base,
        "status": (
            "materialChangeDetected"
            if material_change
            else "noMaterialChangeDetected"
        ),
        "tailStart": times[start_index],
        "tailEnd": times[-1],
        "samplesUsed": len(tail_values),
        "startValue": tail_values[0],
        "endValue": tail_values[-1],
        "change": delta,
        "averageSlopePerSecond": average_slope,
        "relativeChange": relative_change,
        "absoluteRelativeChange": abs(relative_change),
        "relativeRange": relative_range,
        "direction": direction,
        "directionalConsistency": directional_consistency,
        "nonzeroIncrementCount": nonzero_increments,
        "detectionReasons": detection_reasons,
        "directionalChangeDetected": directional_change_detected,
        "tailVariabilityDetected": tail_variability_detected,
        "materialChangeDetected": material_change,
    }


def _completed_metrics(
    result: Mapping[str, object],
    spec: OptimizationSpec,
) -> tuple[float, float, list[dict[str, object]], bool, float]:
    status = str(result.get("status") or "")
    if status != "completed" or result.get("success") is not True:
        raise OptimizationError(
            "OPTIMIZATION_SIMULATION_INCOMPLETE",
            "Only a complete successful simulation can be scored.",
            {"status": status, "success": result.get("success")},
        )
    raw_variables = result.get("variables")
    if not isinstance(raw_variables, list):
        raise OptimizationError(
            "OPTIMIZATION_RESULT_METADATA_INVALID",
            "A completed simulation does not contain result-variable metadata.",
        )
    variables: dict[str, dict[str, object]] = {}
    for item in raw_variables:
        if not isinstance(item, dict) or not isinstance(item.get("key"), str):
            raise OptimizationError(
                "OPTIMIZATION_RESULT_METADATA_INVALID",
                "Completed result-variable metadata is outside the backend contract.",
            )
        key = str(item["key"])
        if key in variables:
            raise OptimizationError(
                "OPTIMIZATION_RESULT_METADATA_INVALID",
                "Completed result-variable metadata contains a duplicate key.",
                {"resultKey": key},
            )
        variables[key] = item
    requested_metadata: list[tuple[str, str, str]] = [
        ("objective", spec.objective.result_key, spec.objective.expected_unit)
    ]
    requested_metadata.extend(
        (
            f"constraint {constraint.id}",
            constraint.result_key,
            constraint.expected_unit,
        )
        for constraint in spec.constraints
    )
    for label, result_key, expected_unit in requested_metadata:
        metadata = variables.get(result_key)
        if metadata is None:
            raise OptimizationError(
                "OPTIMIZATION_RESULT_METADATA_MISMATCH",
                f"The completed result does not declare the {label} result key.",
                {"resultKey": result_key},
            )
        actual_unit = metadata.get("unit")
        if not isinstance(actual_unit, str) or actual_unit != expected_unit:
            raise OptimizationError(
                "OPTIMIZATION_RESULT_METADATA_MISMATCH",
                f"The completed result changed the {label} unit after planning.",
                {
                    "resultKey": result_key,
                    "expected": expected_unit,
                    "received": actual_unit,
                },
            )
    series = result.get("series")
    if not isinstance(series, Mapping):
        raise OptimizationError(
            "OPTIMIZATION_RESULT_SERIES_INVALID",
            "Simulation result does not contain a series object.",
        )
    times = series.get("time")
    objective_values = series.get(spec.objective.result_key)
    if not isinstance(times, list) or not isinstance(objective_values, list):
        raise OptimizationError(
            "OPTIMIZATION_RESULT_VARIABLE_MISSING",
            "The completed result is missing time or the objective series.",
            {"resultKey": spec.objective.result_key},
        )
    objective_value = statistic_value(
        times, objective_values, spec.objective.statistic
    )
    objective_loss = _loss(objective_value, spec.objective.goal)
    if not math.isfinite(objective_loss):
        raise OptimizationError(
            "OPTIMIZATION_METRIC_OVERFLOW",
            "The objective loss overflowed the finite numeric contract.",
        )

    constraint_results: list[dict[str, object]] = []
    total_violation = 0.0
    feasible = True
    for constraint in spec.constraints:
        raw_values = series.get(constraint.result_key)
        if not isinstance(raw_values, list):
            raise OptimizationError(
                "OPTIMIZATION_RESULT_VARIABLE_MISSING",
                "The completed result is missing a constraint series.",
                {"constraint": constraint.id, "resultKey": constraint.result_key},
            )
        value = statistic_value(times, raw_values, constraint.statistic)
        lower_violation = (
            max(0.0, constraint.lower - constraint.tolerance - value)
            if constraint.lower is not None
            else 0.0
        )
        upper_violation = (
            max(0.0, value - constraint.upper - constraint.tolerance)
            if constraint.upper is not None
            else 0.0
        )
        raw_violation = max(lower_violation, upper_violation)
        normalized_violation = raw_violation / constraint.scale
        lower_margin = (
            value - constraint.lower if constraint.lower is not None else None
        )
        upper_margin = (
            constraint.upper - value if constraint.upper is not None else None
        )
        finite_margins = [
            margin for margin in (lower_margin, upper_margin) if margin is not None
        ]
        margin = min(finite_margins) if finite_margins else None
        derived_values = [
            lower_violation,
            upper_violation,
            raw_violation,
            normalized_violation,
            *finite_margins,
            total_violation + normalized_violation,
        ]
        if any(not math.isfinite(item) for item in derived_values):
            raise OptimizationError(
                "OPTIMIZATION_METRIC_OVERFLOW",
                "A response-constraint calculation overflowed the finite numeric contract.",
                {"constraint": constraint.id},
            )
        item_feasible = raw_violation == 0.0
        feasible = feasible and item_feasible
        total_violation += normalized_violation
        constraint_results.append(
            {
                "id": constraint.id,
                "resultKey": constraint.result_key,
                "value": value,
                "unit": _metric_unit(
                    constraint.expected_unit, constraint.statistic.kind
                ),
                "seriesUnit": constraint.expected_unit,
                "lower": constraint.lower,
                "upper": constraint.upper,
                "tolerance": constraint.tolerance,
                "margin": margin,
                "violation": raw_violation,
                "normalizedViolation": normalized_violation,
                "feasible": item_feasible,
            }
        )
    return (
        objective_value,
        objective_loss,
        constraint_results,
        feasible,
        total_violation,
    )


def _format_float(value: float) -> str:
    normalized = 0.0 if value == 0.0 else value
    return format(normalized, ".17g")


def _format_report_float(value: float) -> str:
    normalized = 0.0 if value == 0.0 else value
    return format(normalized, ".12g")


def _reflect_unit_interval(value: float) -> float:
    if not math.isfinite(value):
        raise AssertionError("Optimizer generated a non-finite coordinate.")
    wrapped = value % 2.0
    reflected = wrapped if wrapped <= 1.0 else 2.0 - wrapped
    return 0.0 if reflected == 0.0 else reflected


def _markdown_text(value: object) -> str:
    return (
        str(value)
        .replace("&", "&amp;")
        .replace("<", "&lt;")
        .replace(">", "&gt;")
        .replace("|", "&#124;")
        .replace("`", "&#96;")
        .replace("[", "&#91;")
        .replace("]", "&#93;")
        .replace("(", "&#40;")
        .replace(")", "&#41;")
        .replace("!", "&#33;")
        .replace("\r", " ")
        .replace("\n", " ")
    )


def _candidate_parameters(
    plan: RuntimePlan,
    normalized_vector: Sequence[float],
) -> dict[str, float]:
    if len(normalized_vector) != len(plan.resolved_design_variables):
        raise AssertionError("Candidate vector dimension does not match the plan.")
    values: dict[str, float] = {}
    for coordinate, resolved in zip(
        normalized_vector, plan.resolved_design_variables
    ):
        if not math.isfinite(float(coordinate)):
            raise AssertionError("Optimizer generated a non-finite coordinate.")
        clipped = min(1.0, max(0.0, float(coordinate)))
        spec = resolved.spec
        value = spec.lower + clipped * (spec.upper - spec.lower)
        if not math.isfinite(value):
            raise OptimizationError(
                "OPTIMIZATION_CANDIDATE_OVERFLOW",
                "A normalized candidate did not map to a finite SI parameter value.",
                {"id": spec.id},
            )
        values[spec.id] = 0.0 if value == 0.0 else value
    return values


def _normalized_initial(plan: RuntimePlan) -> list[float]:
    normalized: list[float] = []
    for resolved in plan.resolved_design_variables:
        coordinate = (resolved.initial - resolved.spec.lower) / (
            resolved.spec.upper - resolved.spec.lower
        )
        if not math.isfinite(coordinate):
            raise OptimizationError(
                "OPTIMIZATION_CANDIDATE_OVERFLOW",
                "The baseline parameter could not be normalized to a finite value.",
                {"id": resolved.spec.id},
            )
        normalized.append(coordinate)
    return normalized


def _candidate_xml(plan: RuntimePlan, parameters: Mapping[str, float]) -> bytes:
    try:
        root = ET.fromstring(plan.baseline_xml)
    except ET.ParseError as exc:  # pragma: no cover - backend contract guard
        raise simulation.BackendError(
            "BACKEND_BASELINE_XML_INVALID",
            "Backend-generated baseline System XML is not parseable.",
        ) from exc
    components_parent = root.find("./Components")
    if components_parent is None:
        raise simulation.BackendError(
            "BACKEND_BASELINE_XML_INVALID",
            "Backend-generated System XML has no Components element.",
        )
    components = {
        component.get("id"): component
        for component in components_parent.findall("./Component")
        if component.get("id") is not None
    }
    for resolved in plan.resolved_design_variables:
        component = components.get(resolved.spec.component_id)
        if component is None:
            raise simulation.BackendError(
                "BACKEND_BASELINE_XML_INVALID",
                "An optimized component is missing from baseline System XML.",
                {"componentId": resolved.spec.component_id},
            )
        parameter = next(
            (
                item
                for item in component.findall("./Parameter")
                if item.get("name") == resolved.spec.parameter
            ),
            None,
        )
        if parameter is None:
            raise simulation.BackendError(
                "BACKEND_BASELINE_XML_INVALID",
                "An optimized parameter is missing from baseline System XML.",
                {
                    "componentId": resolved.spec.component_id,
                    "parameter": resolved.spec.parameter,
                },
            )
        parameter.set("value", _format_float(parameters[resolved.spec.id]))
    return ET.tostring(root, encoding="utf-8", xml_declaration=True)


def _optimized_project(
    plan: RuntimePlan,
    parameters: Mapping[str, float],
) -> dict[str, object]:
    project = copy.deepcopy(plan.source.parsed)
    if not isinstance(project, dict):  # pragma: no cover - load_source guarantees this
        raise AssertionError("ReactFlow project root is not an object.")
    nodes = project.get("nodes")
    if not isinstance(nodes, list):
        raise simulation.InputError(
            "PROJECT_NODES_INVALID",
            "The source project no longer contains a nodes array.",
        )
    nodes_by_id = {
        node.get("id"): node
        for node in nodes
        if isinstance(node, dict) and isinstance(node.get("id"), str)
    }
    for resolved in plan.resolved_design_variables:
        node = nodes_by_id.get(resolved.spec.component_id)
        if not isinstance(node, dict):
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_UNKNOWN",
                "An optimized component is missing from the project copy.",
            )
        data = node.get("data")
        if not isinstance(data, dict):
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_UNKNOWN",
                "An optimized component has no data object.",
            )
        raw_parameters = data.get("parameters")
        if raw_parameters is None:
            raw_parameters = {}
            data["parameters"] = raw_parameters
        if not isinstance(raw_parameters, dict):
            raise simulation.InputError(
                "OPTIMIZATION_PARAMETER_UNKNOWN",
                "An optimized component has no parameter object.",
            )
        raw_parameters[resolved.spec.parameter] = parameters[resolved.spec.id]
    return project


def _candidate_cache_key(parameters: Mapping[str, float], order: Sequence[str]) -> bytes:
    return json.dumps(
        [_format_float(parameters[item]) for item in order],
        separators=(",", ":"),
    ).encode("ascii")


def _trial_rank(trial: Trial) -> tuple[float, ...]:
    if trial.status != "completed" or trial.objective_loss is None:
        return (2.0, math.inf, math.inf)
    if trial.feasible:
        return (0.0, trial.objective_loss, float(trial.evaluation_id))
    violation = (
        trial.total_constraint_violation
        if trial.total_constraint_violation is not None
        else math.inf
    )
    return (1.0, violation, trial.objective_loss, float(trial.evaluation_id))


def _is_close(first: float, second: float, validation: ValidationSpec) -> bool:
    return math.isclose(
        first,
        second,
        rel_tol=validation.relative_tolerance,
        abs_tol=validation.absolute_tolerance,
    )


def _atomic_write(path: Path, data: bytes) -> None:
    temporary = path.with_name(f".{path.name}.{uuid.uuid4().hex}.tmp")
    try:
        with temporary.open("xb") as stream:
            stream.write(data)
            stream.flush()
            os.fsync(stream.fileno())
        os.replace(temporary, path)
    except OSError as exc:
        try:
            temporary.unlink(missing_ok=True)
        except OSError:
            pass
        raise simulation.ArtifactError(
            "OPTIMIZATION_ARTIFACT_WRITE_FAILED",
            str(exc),
            {"path": str(path)},
        ) from exc


def _atomic_json(path: Path, payload: object) -> None:
    data = (
        json.dumps(payload, ensure_ascii=False, allow_nan=False, indent=2) + "\n"
    ).encode("utf-8")
    _atomic_write(path, data)


def _artifact(path: Path) -> dict[str, object]:
    data = path.read_bytes()
    return {
        "path": str(path),
        "sha256": _sha256(data),
        "sizeBytes": len(data),
    }


class OptimizationRunner:
    def __init__(self, plan: RuntimePlan, optimization_id: str) -> None:
        self.plan = plan
        self.optimization_id = optimization_id
        self.output_directory = plan.output_directory
        self.progress_directory = self.output_directory / "simulation-progress"
        self.events_path = self.output_directory / "optimization-events.jsonl"
        self.trials: list[Trial] = []
        self.cache: dict[bytes, Trial] = {}
        self.optimizer_calls = 0
        self.search_proposals = 0
        self.cache_hits = 0
        self.backend_submissions = 0
        self.search_backend_submissions = 0
        self.verification_backend_submissions = 0
        self.completed_generations = 0
        self.generations_with_new_backend_submissions = 0
        self.generations_without_new_backend_submissions = 0
        self.ending_no_new_submission_generation_streak = 0
        self.max_no_new_submission_generation_streak = 0
        self.final_population_unique_candidates: int | None = None
        self.optimizer_call_limit: int | None = None
        self.current_simulation_id: str | None = None
        self.started_monotonic = time.monotonic()
        self.started_unix = time.time()
        self.search_run_limit = plan.spec.budget.max_simulation_runs - 1
        self._events: io.TextIOWrapper | None = None

    @property
    def parameter_order(self) -> list[str]:
        return [item.spec.id for item in self.plan.resolved_design_variables]

    def _counts_payload(self) -> dict[str, object]:
        search_trials = [trial for trial in self.trials if trial.stage == "search"]
        verification_trials = [
            trial for trial in self.trials if trial.stage == "verification"
        ]
        return {
            "optimizerCalls": self.optimizer_calls,
            "candidateRequestsAllStages": self.optimizer_calls,
            "searchProposals": self.search_proposals,
            "backendSubmissions": self.backend_submissions,
            "backendSubmissionSlotsConsumed": self.backend_submissions,
            "searchBackendSubmissions": self.search_backend_submissions,
            "searchSubmissionSlotsConsumed": (
                self.search_backend_submissions
            ),
            "verificationBackendSubmissions": (
                self.verification_backend_submissions
            ),
            "verificationSubmissionSlotsConsumed": (
                self.verification_backend_submissions
            ),
            "cacheHits": self.cache_hits,
            "generations": self.completed_generations,
            "generationsWithNewBackendSubmissions": (
                self.generations_with_new_backend_submissions
            ),
            "generationsWithoutNewBackendSubmissions": (
                self.generations_without_new_backend_submissions
            ),
            "endingNoNewSubmissionGenerationStreak": (
                self.ending_no_new_submission_generation_streak
            ),
            "maxNoNewSubmissionGenerationStreak": (
                self.max_no_new_submission_generation_streak
            ),
            "finalPopulationUniqueCandidates": (
                self.final_population_unique_candidates
            ),
            "remainingSearchRunBudget": max(
                0, self.search_run_limit - self.search_backend_submissions
            ),
            "completedTrials": sum(
                trial.status == "completed" for trial in self.trials
            ),
            "failedTrials": sum(
                trial.status != "completed" for trial in self.trials
            ),
            "searchCompletedTrials": sum(
                trial.status == "completed" for trial in search_trials
            ),
            "searchFailedTrials": sum(
                trial.status != "completed" for trial in search_trials
            ),
            "verificationCompletedTrials": sum(
                trial.status == "completed" for trial in verification_trials
            ),
            "verificationFailedTrials": sum(
                trial.status != "completed" for trial in verification_trials
            ),
            "recordedTrials": len(self.trials),
            "unrecordedSubmissionSlotsConsumed": max(
                0, self.backend_submissions - len(self.trials)
            ),
            "searchUnrecordedSubmissionSlotsConsumed": max(
                0, self.search_backend_submissions - len(search_trials)
            ),
            "verificationUnrecordedSubmissionSlotsConsumed": max(
                0,
                self.verification_backend_submissions
                - len(verification_trials),
            ),
        }

    def _search_summary(self, termination_reason: str) -> dict[str, object]:
        categories = {
            "simulationBudgetExhausted": "budget",
            "searchWallTimeReached": "wallTime",
            "populationCollapsedAfterDuplicateStagnation": "stagnation",
            "duplicateProposalStagnation": "stagnation",
            "optimizerCallLimitReached": "safeguard",
            "userCancelled": "cancelled",
            "error": "error",
        }
        return {
            "terminationReason": termination_reason,
            "terminationCategory": categories.get(termination_reason, "other"),
            "searchConvergenceEstablished": False,
            "populationCollapsedToSingleCandidate": (
                termination_reason
                == "populationCollapsedAfterDuplicateStagnation"
            ),
            "candidateRequests": self.search_proposals,
            "backendSubmissions": self.search_backend_submissions,
            "submissionSlotsConsumed": self.search_backend_submissions,
            "submissionCounting": "reservedBeforeProgressLogAndBackendStream",
            "cacheHits": self.cache_hits,
            "cacheHitRate": (
                self.cache_hits / self.search_proposals
                if self.search_proposals
                else 0.0
            ),
            "budget": {
                "limit": self.search_run_limit,
                "used": self.search_backend_submissions,
                "unused": max(
                    0, self.search_run_limit - self.search_backend_submissions
                ),
                "exhausted": (
                    self.search_backend_submissions >= self.search_run_limit
                ),
            },
            "generations": {
                "completed": self.completed_generations,
                "withNewBackendSubmissions": (
                    self.generations_with_new_backend_submissions
                ),
                "withoutNewBackendSubmissions": (
                    self.generations_without_new_backend_submissions
                ),
                "endingNoNewSubmissionStreak": (
                    self.ending_no_new_submission_generation_streak
                ),
                "maximumNoNewSubmissionStreak": (
                    self.max_no_new_submission_generation_streak
                ),
            },
            "finalPopulationUniqueCandidates": (
                self.final_population_unique_candidates
            ),
            "noNewSubmissionGenerationLimit": (
                NO_NEW_SUBMISSION_GENERATION_LIMIT
            ),
            "optimizerCallLimit": self.optimizer_call_limit,
        }

    def _timing_payload(self) -> dict[str, object]:
        durations = [trial.duration_seconds for trial in self.trials]
        return {
            "wallClockSeconds": max(0.0, time.monotonic() - self.started_monotonic),
            "simulationRuns": len(durations),
            "simulationDurationTotalSeconds": (
                math.fsum(durations) if durations else 0.0
            ),
            "simulationDurationMinimumSeconds": (
                min(durations) if durations else None
            ),
            "simulationDurationMaximumSeconds": (
                max(durations) if durations else None
            ),
            "simulationDurationMeanSeconds": (
                math.fsum(durations) / len(durations) if durations else None
            ),
        }

    def _emit(self, payload: Mapping[str, object]) -> None:
        event = dict(payload)
        if self._events is not None:
            serialized = json.dumps(
                event, ensure_ascii=False, allow_nan=False, separators=(",", ":")
            )
            self._events.write(serialized + "\n")
            self._events.flush()
        simulation.emit_json(event)

    def _write_evaluations(self) -> None:
        columns = [
            "evaluationId",
            "stage",
            "status",
            "feasible",
            "objectiveValue",
            "objectiveLoss",
            "totalConstraintViolation",
            "simulationId",
            "durationSeconds",
            "failureCode",
            "cacheReuseCount",
        ]
        columns.extend(f"parameter.{item}" for item in self.parameter_order)
        for constraint in self.plan.spec.constraints:
            columns.extend(
                (
                    f"constraint.{constraint.id}.value",
                    f"constraint.{constraint.id}.margin",
                    f"constraint.{constraint.id}.violation",
                )
            )
        stream = io.StringIO(newline="")
        writer = csv.DictWriter(stream, fieldnames=columns, lineterminator="\n")
        writer.writeheader()
        for trial in self.trials:
            row: dict[str, object] = {
                "evaluationId": trial.evaluation_id,
                "stage": trial.stage,
                "status": trial.status,
                "feasible": str(trial.feasible).lower(),
                "objectiveValue": (
                    _format_float(trial.objective_value)
                    if trial.objective_value is not None
                    else ""
                ),
                "objectiveLoss": (
                    _format_float(trial.objective_loss)
                    if trial.objective_loss is not None
                    else ""
                ),
                "totalConstraintViolation": (
                    _format_float(trial.total_constraint_violation)
                    if trial.total_constraint_violation is not None
                    else ""
                ),
                "simulationId": trial.simulation_id,
                "durationSeconds": _format_float(trial.duration_seconds),
                "failureCode": trial.failure_code or "",
                "cacheReuseCount": trial.cache_reuse_count,
            }
            for key, value in trial.parameters.items():
                row[f"parameter.{key}"] = _format_float(value)
            constraints = {
                str(item.get("id")): item for item in trial.constraints
            }
            for constraint in self.plan.spec.constraints:
                item = constraints.get(constraint.id, {})
                for name in ("value", "margin", "violation"):
                    raw_value = item.get(name)
                    row[f"constraint.{constraint.id}.{name}"] = (
                        _format_float(float(raw_value))
                        if isinstance(raw_value, (int, float))
                        and not isinstance(raw_value, bool)
                        else ""
                    )
            writer.writerow(row)
        _atomic_write(
            self.output_directory / "evaluations.csv",
            stream.getvalue().encode("utf-8"),
        )

    def _best_trial(self, *, feasible_only: bool) -> Trial | None:
        candidates = [
            trial
            for trial in self.trials
            if trial.stage == "search"
            and trial.status == "completed"
            and trial.objective_loss is not None
            and (trial.feasible or not feasible_only)
        ]
        return min(candidates, key=_trial_rank) if candidates else None

    def _checkpoint(self, termination_reason: str | None = None) -> None:
        best = self._best_trial(feasible_only=True)
        payload: dict[str, object] = {
            "optimizationResultSchemaVersion": OPTIMIZATION_RESULT_SCHEMA_VERSION,
            "optimizationId": self.optimization_id,
            "sourceSha256": self.plan.source.sha256,
            "specSha256": self.plan.spec_source.sha256,
            "counts": self._counts_payload(),
            "bestFeasibleSearch": best.as_dict() if best is not None else None,
            "trials": [trial.as_dict() for trial in self.trials],
        }
        if termination_reason is not None:
            payload["terminationReason"] = termination_reason
        _atomic_json(self.output_directory / "checkpoint.json", payload)
        self._write_evaluations()

    def _inner_event_sink(self, evaluation_id: int, simulation_id: str):
        def sink(payload: object) -> None:
            if not isinstance(payload, Mapping):
                return
            if payload.get("event") != "progress":
                return
            self._emit(
                {
                    "event": "optimization-simulation-progress",
                    "optimizationId": self.optimization_id,
                    "evaluationId": evaluation_id,
                    "simulationId": simulation_id,
                    "searchRunLimit": self.search_run_limit,
                    "simulation": dict(payload),
                }
            )

        return sink

    def _evaluate(
        self,
        normalized_vector: Sequence[float],
        *,
        stage: str,
        bypass_cache: bool = False,
        retain_result: bool = False,
    ) -> tuple[Trial, dict[str, object] | None]:
        elapsed = time.monotonic() - self.started_monotonic
        if stage == "search":
            if self.search_backend_submissions >= self.search_run_limit:
                raise SearchStop("simulationBudgetExhausted")
            if elapsed >= self.plan.spec.budget.max_wall_seconds:
                raise SearchStop("searchWallTimeReached")
        elif self.backend_submissions >= self.plan.spec.budget.max_simulation_runs:
            raise SearchStop("simulationBudgetExhausted")
        self.optimizer_calls += 1
        if stage == "search":
            self.search_proposals += 1
        parameters = _candidate_parameters(self.plan, normalized_vector)
        cache_key = _candidate_cache_key(parameters, self.parameter_order)
        if not bypass_cache and cache_key in self.cache:
            cached = self.cache[cache_key]
            cached.cache_reuse_count += 1
            self.cache_hits += 1
            self._checkpoint()
            self._emit(
                {
                    "event": "optimization-cache-hit",
                    "optimizationId": self.optimization_id,
                    "evaluationId": cached.evaluation_id,
                    "parameters": parameters,
                }
            )
            return cached, None

        evaluation_id = self.backend_submissions + 1
        simulation_id = f"{self.optimization_id}.{evaluation_id:04d}"
        xml = _candidate_xml(self.plan, parameters)
        progress_path = self.progress_directory / f"evaluation-{evaluation_id:04d}.jsonl"
        started = time.monotonic()
        self.backend_submissions += 1
        if stage == "search":
            self.search_backend_submissions += 1
        else:
            self.verification_backend_submissions += 1
        self._emit(
            {
                "event": "optimization-evaluation-started",
                "optimizationId": self.optimization_id,
                "evaluationId": evaluation_id,
                "stage": stage,
                "simulationId": simulation_id,
                "parameters": parameters,
            }
        )
        result: dict[str, object] | None = None
        backend_error: dict[str, object] | None = None
        self.current_simulation_id = simulation_id
        try:
            result, backend_error = simulation._read_simulation_stream(
                self.plan.base_url,
                xml,
                simulation_id,
                self.plan.timeout,
                progress_path,
                event_sink=self._inner_event_sink(evaluation_id, simulation_id),
                full_result=None,
            )
        except BaseException:
            raise
        else:
            self.current_simulation_id = None
        duration = time.monotonic() - started

        contract_error: OptimizationError | None = None
        if (
            result is not None
            and str(result.get("status") or "") == "completed"
            and result.get("success") is True
        ):
            try:
                (
                    objective_value,
                    objective_loss,
                    constraints,
                    feasible,
                    total_violation,
                ) = _completed_metrics(result, self.plan.spec)
            except OptimizationError as exc:
                contract_error = exc
                trial = Trial(
                    evaluation_id=evaluation_id,
                    stage=stage,
                    simulation_id=simulation_id,
                    parameters=parameters,
                    status="failed",
                    duration_seconds=duration,
                    failure_code=exc.code,
                    failure_message=exc.message,
                )
            else:
                trial = Trial(
                    evaluation_id=evaluation_id,
                    stage=stage,
                    simulation_id=simulation_id,
                    parameters=parameters,
                    status="completed",
                    duration_seconds=duration,
                    objective_value=objective_value,
                    objective_loss=objective_loss,
                    constraints=constraints,
                    feasible=feasible,
                    total_constraint_violation=total_violation,
                )
        else:
            result_status = (
                str(result.get("status") or "failed")
                if isinstance(result, Mapping)
                else "failed"
            )
            if result_status == "completed":
                result_status = "failed"
            failure_message = (
                str(result.get("message") or "Simulation did not complete successfully.")
                if isinstance(result, Mapping)
                else "Simulation ended without a usable result."
            )
            trial = Trial(
                evaluation_id=evaluation_id,
                stage=stage,
                simulation_id=simulation_id,
                parameters=parameters,
                status=result_status,
                duration_seconds=duration,
                failure_code="OPTIMIZATION_SIMULATION_INCOMPLETE",
                failure_message=failure_message,
            )
            if backend_error is not None:
                trial.failure_message = json.dumps(
                    backend_error, ensure_ascii=False, allow_nan=False
                )
        self.trials.append(trial)
        if not bypass_cache:
            self.cache[cache_key] = trial
        self._checkpoint()
        self._emit(
            {
                "event": "optimization-evaluation-completed",
                "optimizationId": self.optimization_id,
                "evaluation": trial.as_dict(),
                "bestFeasibleSearch": (
                    self._best_trial(feasible_only=True).as_dict()
                    if self._best_trial(feasible_only=True) is not None
                    else None
                ),
            }
        )
        if contract_error is not None:
            raise contract_error
        return trial, result if retain_result else None

    def _initial_population(self, random_source: random.Random) -> list[list[float]]:
        population_size = self.plan.spec.algorithm.population_size
        dimension = len(self.plan.resolved_design_variables)
        population = [[0.0] * dimension for _ in range(population_size)]
        for coordinate in range(dimension):
            strata = list(range(population_size))
            random_source.shuffle(strata)
            for row, stratum in enumerate(strata):
                population[row][coordinate] = (
                    stratum + random_source.random()
                ) / population_size
        population[0] = _normalized_initial(self.plan)
        if dimension == 1:
            available_rows = iter(range(population_size - 1, 0, -1))
            for endpoint in (0.0, 1.0):
                if any(row[0] == endpoint for row in population):
                    continue
                population[next(available_rows)] = [endpoint]
        return population

    def _population_unique_candidate_count(
        self, population: Sequence[Sequence[float]]
    ) -> int:
        return len(
            {
                _candidate_cache_key(
                    _candidate_parameters(self.plan, vector),
                    self.parameter_order,
                )
                for vector in population
            }
        )

    def _search(self) -> tuple[str, int]:
        algorithm = self.plan.spec.algorithm
        random_source = random.Random(algorithm.seed)
        population = self._initial_population(random_source)
        trials: list[Trial] = []
        termination_reason = "simulationBudgetExhausted"
        maximum_optimizer_calls = max(100, self.search_run_limit * 20)
        self.optimizer_call_limit = maximum_optimizer_calls
        try:
            for index, vector in enumerate(population):
                trial, _ = self._evaluate(vector, stage="search")
                if index == 0 and trial.status != "completed":
                    raise OptimizationError(
                        "OPTIMIZATION_BASELINE_FAILED",
                        "The baseline model must complete successfully before search starts.",
                        trial.as_dict(),
                    )
                trials.append(trial)

            self.final_population_unique_candidates = (
                self._population_unique_candidate_count(population)
            )
            while (
                self.search_backend_submissions < self.search_run_limit
                and self.search_proposals < maximum_optimizer_calls
            ):
                submissions_before = self.search_backend_submissions
                cache_hits_before = self.cache_hits
                source_population = [list(vector) for vector in population]
                source_trials = list(trials)
                next_population = [list(vector) for vector in population]
                next_trials = list(trials)
                for target_index in range(len(population)):
                    if self.search_proposals >= maximum_optimizer_calls:
                        raise SearchStop("optimizerCallLimitReached")
                    available = [
                        index
                        for index in range(len(source_population))
                        if index != target_index
                    ]
                    first_index, second_index, third_index = random_source.sample(
                        available, 3
                    )
                    mutant = [
                        _reflect_unit_interval(
                            source_population[first_index][coordinate]
                            + algorithm.mutation_factor
                            * (
                                source_population[second_index][coordinate]
                                - source_population[third_index][coordinate]
                            )
                        )
                        for coordinate in range(len(source_population[target_index]))
                    ]
                    mandatory_coordinate = random_source.randrange(len(mutant))
                    candidate = [
                        (
                            mutant[coordinate]
                            if coordinate == mandatory_coordinate
                            or random_source.random()
                            <= algorithm.crossover_probability
                            else source_population[target_index][coordinate]
                        )
                        for coordinate in range(len(mutant))
                    ]
                    candidate_trial, _ = self._evaluate(
                        candidate, stage="search"
                    )
                    if _trial_rank(candidate_trial) < _trial_rank(
                        source_trials[target_index]
                    ):
                        next_population[target_index] = candidate
                        next_trials[target_index] = candidate_trial
                population = next_population
                trials = next_trials
                self.completed_generations += 1
                new_submissions = (
                    self.search_backend_submissions - submissions_before
                )
                new_cache_hits = self.cache_hits - cache_hits_before
                if new_submissions:
                    self.generations_with_new_backend_submissions += 1
                    self.ending_no_new_submission_generation_streak = 0
                else:
                    self.generations_without_new_backend_submissions += 1
                    self.ending_no_new_submission_generation_streak += 1
                    self.max_no_new_submission_generation_streak = max(
                        self.max_no_new_submission_generation_streak,
                        self.ending_no_new_submission_generation_streak,
                    )
                self.final_population_unique_candidates = (
                    self._population_unique_candidate_count(population)
                )
                self._emit(
                    {
                        "event": "optimization-generation-completed",
                        "optimizationId": self.optimization_id,
                        "generation": self.completed_generations,
                        "backendSubmissions": self.backend_submissions,
                        "searchRunLimit": self.search_run_limit,
                        "backendSubmissionsThisGeneration": new_submissions,
                        "cacheHitsThisGeneration": new_cache_hits,
                        "consecutiveGenerationsWithoutNewBackendSubmissions": (
                            self.ending_no_new_submission_generation_streak
                        ),
                        "populationUniqueCandidates": (
                            self.final_population_unique_candidates
                        ),
                    }
                )
                self._checkpoint()
                if (
                    self.ending_no_new_submission_generation_streak
                    >= NO_NEW_SUBMISSION_GENERATION_LIMIT
                ):
                    termination_reason = (
                        "populationCollapsedAfterDuplicateStagnation"
                        if self.final_population_unique_candidates == 1
                        else "duplicateProposalStagnation"
                    )
                    break
            if (
                termination_reason == "simulationBudgetExhausted"
                and self.search_proposals >= maximum_optimizer_calls
                and self.search_backend_submissions < self.search_run_limit
            ):
                termination_reason = "optimizerCallLimitReached"
        except SearchStop as stop:
            termination_reason = stop.reason
        self.final_population_unique_candidates = (
            self._population_unique_candidate_count(population)
        )
        return termination_reason, self.completed_generations

    def cancel_active(self) -> None:
        if self.current_simulation_id is None:
            return
        simulation_id = self.current_simulation_id
        try:
            simulation.http_json(
                self.plan.base_url,
                "/api/system-xml/simulations/"
                + simulation.urllib.parse.quote(simulation_id, safe="")
                + "/cancel",
                method="POST",
                payload={"reason": "user"},
                timeout=self.plan.timeout,
            )
        except simulation.SkillCliError:
            pass
        finally:
            self.current_simulation_id = None

    def _source_is_unchanged(self) -> bool:
        try:
            return _sha256(self.plan.source.path.read_bytes()) == self.plan.source.sha256
        except OSError:
            return False

    def _verification_comparison(
        self,
        search: Trial,
        verification: Trial,
    ) -> dict[str, object]:
        validation = self.plan.spec.validation
        comparisons: list[dict[str, object]] = []
        objective_matches = (
            search.objective_value is not None
            and verification.objective_value is not None
            and _is_close(
                search.objective_value,
                verification.objective_value,
                validation,
            )
        )
        comparisons.append(
            {
                "id": "objective",
                "search": search.objective_value,
                "verification": verification.objective_value,
                "matches": objective_matches,
            }
        )
        search_constraints = {
            str(item.get("id")): item for item in search.constraints
        }
        verification_constraints = {
            str(item.get("id")): item for item in verification.constraints
        }
        for constraint in self.plan.spec.constraints:
            first = search_constraints.get(constraint.id, {}).get("value")
            second = verification_constraints.get(constraint.id, {}).get("value")
            matches = (
                isinstance(first, (int, float))
                and not isinstance(first, bool)
                and isinstance(second, (int, float))
                and not isinstance(second, bool)
                and _is_close(float(first), float(second), validation)
            )
            comparisons.append(
                {
                    "id": constraint.id,
                    "search": first,
                    "verification": second,
                    "matches": matches,
                }
            )
        source_unchanged = self._source_is_unchanged()
        return {
            "passed": (
                verification.status == "completed"
                and verification.feasible
                and source_unchanged
                and all(item["matches"] is True for item in comparisons)
            ),
            "sourceUnchanged": source_unchanged,
            "relativeTolerance": validation.relative_tolerance,
            "absoluteTolerance": validation.absolute_tolerance,
            "metrics": comparisons,
        }

    def _improvement(
        self, baseline: Trial | None, best: Trial | None
    ) -> dict[str, object] | None:
        if (
            baseline is None
            or best is None
            or baseline.objective_value is None
            or best.objective_value is None
        ):
            return None
        goal = self.plan.spec.objective.goal
        baseline_loss = _loss(baseline.objective_value, goal)
        best_loss = _loss(best.objective_value, goal)
        absolute = baseline_loss - best_loss
        if not math.isfinite(absolute):
            absolute = None
        relative = (
            absolute / abs(baseline_loss)
            if absolute is not None and baseline_loss != 0.0
            else None
        )
        if relative is not None and not math.isfinite(relative):
            relative = None
        return {
            "baselineValue": baseline.objective_value,
            "verifiedValue": best.objective_value,
            "baselineLoss": baseline_loss,
            "verifiedLoss": best_loss,
            "lossReduction": absolute,
            "relativeLossReduction": relative,
        }

    def _write_best_artifacts(
        self,
        best: Trial,
        result: Mapping[str, object],
    ) -> dict[str, object]:
        created_paths: list[Path] = []
        try:
            with tempfile.TemporaryDirectory(
                prefix=".best-artifacts-",
                dir=self.output_directory,
            ) as staging_text:
                staging = Path(staging_text)
                best_parameters_path = staging / "best-parameters.json"
                simulation._write_new_file(
                    best_parameters_path,
                    (
                        json.dumps(
                            {
                                "sourceSha256": self.plan.source.sha256,
                                "specSha256": self.plan.spec_source.sha256,
                                "parameters": [
                                    {
                                        **resolved.spec.as_dict(),
                                        "value": best.parameters[resolved.spec.id],
                                    }
                                    for resolved in self.plan.resolved_design_variables
                                ],
                            },
                            ensure_ascii=False,
                            allow_nan=False,
                            indent=2,
                        )
                        + "\n"
                    ).encode("utf-8"),
                )

                best_xml_path = staging / "best-system.xml"
                simulation._write_new_file(
                    best_xml_path,
                    _candidate_xml(self.plan, best.parameters),
                )

                best_project_path = staging / "best-project.json"
                best_project = _optimized_project(self.plan, best.parameters)
                simulation._write_new_file(
                    best_project_path,
                    (
                        json.dumps(
                            best_project,
                            ensure_ascii=False,
                            allow_nan=False,
                            indent=2,
                        )
                        + "\n"
                    ).encode("utf-8"),
                )

                result_path = staging / "result.json"
                simulation._write_result_json(result_path, result)
                variables, series = simulation._validate_result_shape(result)
                csv_path = staging / "results.csv"
                csv_bytes = simulation._download_csv(
                    self.plan.base_url,
                    self.plan.timeout,
                    simulation._project_name(self.plan.inspection),
                    variables,
                    series,
                )
                simulation._write_new_file(csv_path, csv_bytes)

                metadata = simulation._variable_metadata_by_key(variables)
                chart_keys = list(
                    dict.fromkeys(
                        [self.plan.spec.objective.result_key]
                        + [item.result_key for item in self.plan.spec.constraints]
                    )
                )
                missing = [
                    key
                    for key in chart_keys
                    if key not in metadata or key not in series
                ]
                if missing:
                    raise OptimizationError(
                        "OPTIMIZATION_RESULT_VARIABLE_MISSING",
                        "Verified result is missing a requested report curve.",
                        {"missing": missing},
                    )
                staged_chart_paths = [
                    Path(path)
                    for path in simulation._write_charts(
                        staging,
                        "separate",
                        chart_keys,
                        metadata,
                        series,
                    )
                ]

                staged_artifacts = {
                    "bestParameters": best_parameters_path,
                    "bestSystemXml": best_xml_path,
                    "bestProject": best_project_path,
                    "result": result_path,
                    "csv": csv_path,
                }
                final_artifacts: dict[str, object] = {}
                for name, staged_path in staged_artifacts.items():
                    final_path = self.output_directory / staged_path.name
                    simulation._write_new_file(final_path, staged_path.read_bytes())
                    created_paths.append(final_path)
                    final_artifacts[name] = _artifact(final_path)

                final_chart_paths: list[Path] = []
                for staged_path in staged_chart_paths:
                    final_path = self.output_directory / staged_path.name
                    simulation._write_new_file(final_path, staged_path.read_bytes())
                    created_paths.append(final_path)
                    final_chart_paths.append(final_path)
                final_artifacts["charts"] = [
                    _artifact(path) for path in final_chart_paths
                ]
                return final_artifacts
        except BaseException as exc:
            cleanup_failures: list[dict[str, str]] = []
            for path in reversed(created_paths):
                try:
                    path.unlink(missing_ok=True)
                except OSError as cleanup_exc:
                    cleanup_failures.append(
                        {"path": str(path), "error": str(cleanup_exc)}
                    )
            if cleanup_failures:
                raise simulation.ArtifactError(
                    "OPTIMIZATION_ARTIFACT_ROLLBACK_FAILED",
                    "Failed to remove partially committed best artifacts.",
                    {"failures": cleanup_failures},
                ) from exc
            raise

    def _discard_generated_artifacts(
        self, artifacts: Mapping[str, object]
    ) -> None:
        output_root = self.output_directory.resolve()
        paths: list[Path] = []
        for value in artifacts.values():
            items = value if isinstance(value, list) else [value]
            for item in items:
                if not isinstance(item, Mapping):
                    continue
                raw_path = item.get("path")
                if not isinstance(raw_path, str):
                    continue
                path = Path(raw_path).resolve()
                if not path.is_relative_to(output_root):
                    raise simulation.ArtifactError(
                        "OPTIMIZATION_ARTIFACT_PATH_INVALID",
                        "A generated artifact resolved outside the optimization output directory.",
                        {"path": str(path)},
                    )
                paths.append(path)
        for path in paths:
            try:
                path.unlink(missing_ok=True)
            except OSError as exc:
                raise simulation.ArtifactError(
                    "OPTIMIZATION_ARTIFACT_REMOVE_FAILED",
                    str(exc),
                    {"path": str(path)},
                ) from exc

    def _report_markdown(self, result: Mapping[str, object]) -> str:
        objective = self.plan.spec.objective
        objective_metric_unit = _metric_unit(
            objective.expected_unit, objective.statistic.kind
        )
        baseline = result.get("baseline")
        best_search = result.get("bestSearch")
        best_diagnostic_search = result.get("bestDiagnosticSearch")
        search_comparison = (
            best_search
            if isinstance(best_search, Mapping)
            else best_diagnostic_search
        )
        best = result.get("verifiedBest")
        counts = result.get("counts")
        count_values = counts if isinstance(counts, Mapping) else {}
        search = result.get("search")
        search_values = search if isinstance(search, Mapping) else {}
        search_generations = search_values.get("generations")
        generation_values = (
            search_generations
            if isinstance(search_generations, Mapping)
            else {}
        )
        timing = result.get("timing")
        timing_values = timing if isinstance(timing, Mapping) else {}
        goal_text = objective.goal.kind
        if objective.goal.value is not None:
            goal_text += f" {_format_report_float(objective.goal.value)}"
        window_text = (
            "完整返回时段"
            if objective.statistic.window is None
            else (
                f"[{_format_report_float(objective.statistic.window.start)}, "
                f"{_format_report_float(objective.statistic.window.end)}] s"
            )
        )

        def trial_number(trial: object, field_name: str) -> str:
            if not isinstance(trial, Mapping):
                return "—"
            value = trial.get(field_name)
            if not isinstance(value, (int, float)) or isinstance(value, bool):
                return "—"
            return _format_report_float(float(value))

        def constraint_map(trial: object) -> dict[str, Mapping[str, object]]:
            if not isinstance(trial, Mapping):
                return {}
            raw_constraints = trial.get("constraints")
            if not isinstance(raw_constraints, list):
                return {}
            return {
                str(item.get("id")): item
                for item in raw_constraints
                if isinstance(item, Mapping)
            }

        status = str(result.get("solutionStatus") or "")
        search_points = count_values.get("searchCompletedTrials")
        search_scope = (
            f"本次完整完成的 {search_points} 个搜索仿真"
            if isinstance(search_points, int) and not isinstance(search_points, bool)
            else "本次完整完成的搜索仿真"
        )
        conclusion = {
            "verified": (
                f"这是{search_scope}中表现最好的可行候选，"
                "且已通过一次独立复验；这不证明搜索收敛、系统达到稳态或全局最优。"
            ),
            "noFeasibleCandidate": (
                "预算内没有找到完整可行候选；所列诊断点不能称为优化方案。"
            ),
            "verificationFailed": (
                "搜索阶段找到了可行候选，但新鲜复验未通过；不交付已验证方案。"
            ),
            "aborted": "优化在完成可靠复验前中止；不交付已验证方案。",
        }.get(status, "本次运行没有形成可交付的已验证方案。")
        full_generations_with_unique_candidates = max(
            0,
            (
                self.search_run_limit
                - self.plan.spec.algorithm.population_size
            )
            // self.plan.spec.algorithm.population_size,
        )
        endpoint_seeding_text = (
            "- 单变量边界覆盖：初始种群明确包含上下界"
            if len(self.plan.resolved_design_variables) == 1
            else (
                "- 多变量边界覆盖：不枚举全部边界组合；初始种群使用"
                "逐坐标拉丁超立方分层"
            )
        )
        termination_reason = str(result.get("terminationReason") or "")
        termination_descriptions = {
            "simulationBudgetExhausted": "搜索仿真预算已用完",
            "searchWallTimeReached": "搜索达到墙钟时间上限",
            "populationCollapsedAfterDuplicateStagnation": (
                "连续多代只生成已缓存候选，且种群已塌缩到同一个精确候选；"
                "这是重复停滞状态，不是收敛判定"
            ),
            "duplicateProposalStagnation": (
                "连续多代只生成已缓存候选，但种群尚未合一；这是重复候选停滞，"
                "不是收敛判定"
            ),
            "optimizerCallLimitReached": (
                "候选请求触发内部防死循环上限；这是安全停止，不是收敛判定"
            ),
            "userCancelled": "用户取消了搜索",
            "error": "搜索因错误中止",
        }
        termination_description = termination_descriptions.get(
            termination_reason, "搜索已停止"
        )
        project_name = simulation._project_name(self.plan.inspection)
        status_label = {
            "verified": "候选已复验",
            "noFeasibleCandidate": "未找到可行候选",
            "verificationFailed": "候选复验失败",
            "aborted": "运行已中止",
        }.get(status, "没有已验证方案")
        parameter_section_title = (
            "最佳参数" if status == "verified" else "参数范围"
        )
        search_comparison_label = (
            "搜索阶段最佳可行点"
            if isinstance(best_search, Mapping)
            else "约束违反最小的完整诊断点（不可行）"
        )
        search_constraint_column = (
            "搜索最佳"
            if isinstance(best_search, Mapping)
            else "搜索诊断（不可行）"
        )
        lines = [
            f"# {_markdown_text(project_name)} 优化结果 — {status_label}",
            "",
            f"- 结论：{conclusion}",
            f"- 目标：`{goal_text}` `{_markdown_text(objective.result_key)}` / `{objective.statistic.kind}`（{_markdown_text(objective_metric_unit or 'dimensionless')}）",
            f"- 目标时间窗：{window_text}",
            f"- 候选复验状态：`{result.get('solutionStatus')}`",
            "",
            f"## {parameter_section_title}",
            "",
            (
                "- 连续性声明：后端不验证参数是否连续；以下变量由用户在确认"
                "完整计划时声明为线性 SI 连续参数，且不会改变端口、拓扑或模式。"
            ),
            "",
            "| ID | 组件参数 | 起点 | 下界 | 上界 | 验证值 | 单位 |",
            "| --- | --- | ---: | ---: | ---: | ---: | --- |",
        ]
        best_parameters = (
            best.get("parameters") if isinstance(best, Mapping) else {}
        )
        if not isinstance(best_parameters, Mapping):
            best_parameters = {}
        boundary_notes: list[str] = []
        for variable in self.plan.resolved_design_variables:
            verified = best_parameters.get(variable.spec.id)
            lines.append(
                "| {id} | `{component}.{parameter}` | {initial} | {lower} | {upper} | {best} | {unit} |".format(
                    id=variable.spec.id,
                    component=_markdown_text(variable.spec.component_id),
                    parameter=_markdown_text(variable.spec.parameter),
                    initial=_format_report_float(variable.initial),
                    lower=_format_report_float(variable.spec.lower),
                    upper=_format_report_float(variable.spec.upper),
                    best=(
                        _format_report_float(float(verified))
                        if isinstance(verified, (int, float))
                        and not isinstance(verified, bool)
                        else "—"
                    ),
                    unit=_markdown_text(variable.spec.unit or "dimensionless"),
                )
            )
            if isinstance(verified, (int, float)) and not isinstance(verified, bool):
                boundary_name: str | None = None
                boundary_tolerance = max(
                    math.ulp(variable.spec.lower),
                    math.ulp(variable.spec.upper),
                    sys.float_info.epsilon
                    * (variable.spec.upper - variable.spec.lower)
                    * 4.0,
                )
                if abs(float(verified) - variable.spec.lower) <= boundary_tolerance:
                    boundary_name = "下界"
                elif abs(float(verified) - variable.spec.upper) <= boundary_tolerance:
                    boundary_name = "上界"
                if boundary_name is not None:
                    boundary_notes.append(
                        "`{}` 位于当前确认范围的{}".format(
                            _markdown_text(variable.spec.id), boundary_name
                        )
                    )
        if boundary_notes:
            lines.extend(
                [
                    "",
                    "- 边界说明：{}。这只说明当前范围内得到的是边界候选；"
                    "在确认更宽范围符合物理、安全和组件合同前，不据此建议放宽边界。".format(
                        "；".join(boundary_notes)
                    ),
                ]
            )
        lines.extend(
            [
                "",
                "## 目标对比",
                "",
                "| 项目 | 数值 |",
                "| --- | ---: |",
                f"| 基线 | {trial_number(baseline, 'objectiveValue')} |",
                f"| {search_comparison_label} | {trial_number(search_comparison, 'objectiveValue')} |",
                f"| 新鲜复验 | {trial_number(best, 'objectiveValue')} |",
                "",
            ]
        )
        if isinstance(best_diagnostic_search, Mapping):
            lines.extend(
                [
                    "- 该诊断点的归一化约束违反总量：{}；它不是可行方案。".format(
                        trial_number(
                            best_diagnostic_search,
                            "totalConstraintViolation",
                        )
                    ),
                    "",
                ]
            )
        endpoint_trend = result.get("objectiveEndpointTrend")
        endpoint_lines: list[str] = []
        if isinstance(endpoint_trend, Mapping) and endpoint_trend.get("applicable"):
            endpoint_lines.extend(["## 终点趋势检查", ""])
            trend_status = endpoint_trend.get("status")
            if trend_status == "materialChangeDetected":
                relative_change = endpoint_trend.get("absoluteRelativeChange")
                relative_percent = (
                    _format_report_float(float(relative_change) * 100.0)
                    if isinstance(relative_change, (int, float))
                    and not isinstance(relative_change, bool)
                    else "—"
                )
                relative_range = endpoint_trend.get("relativeRange")
                relative_range_percent = (
                    _format_report_float(float(relative_range) * 100.0)
                    if isinstance(relative_range, (int, float))
                    and not isinstance(relative_range, bool)
                    else "—"
                )
                if endpoint_trend.get("direction") == "fluctuating":
                    endpoint_lines.append(
                        "- 复验曲线在末段 [{start}, {end}] s 仍有明显波动："
                        "首尾净变化约 {change}%，峰峰范围约 {range_}%；因此 `final` "
                        "只代表终点快照，不能据此声称达到稳态。".format(
                            start=trial_number(endpoint_trend, "tailStart"),
                            end=trial_number(endpoint_trend, "tailEnd"),
                            change=relative_percent,
                            range_=relative_range_percent,
                        )
                    )
                else:
                    endpoint_lines.append(
                        "- 复验曲线在末段 [{start}, {end}] s 仍{direction}，"
                        "首尾净变化约 {percent}%；因此 `final` 只代表终点快照，"
                        "不能据此声称达到稳态。".format(
                            start=trial_number(endpoint_trend, "tailStart"),
                            end=trial_number(endpoint_trend, "tailEnd"),
                            direction=(
                                "上升"
                                if endpoint_trend.get("direction") == "increasing"
                                else "下降"
                            ),
                            percent=relative_percent,
                        )
                    )
                endpoint_lines.append(
                    "- 末段平均变化率：{} {}/s；这是启发式风险提示，不改变复验状态。".format(
                        trial_number(endpoint_trend, "averageSlopePerSecond"),
                        _markdown_text(objective_metric_unit or "dimensionless"),
                    )
                )
            elif trend_status == "noMaterialChangeDetected":
                endpoint_lines.append(
                    "- 启发式检查未发现明显的末端变化；这仍不是稳态证明。"
                )
            elif trend_status == "insufficientData":
                endpoint_lines.append(
                    "- 完整复验曲线的末段样本不足，无法可靠检查终点附近的变化；"
                    "不能据此判断稳态。"
                )
            else:
                endpoint_lines.append(
                    "- 趋势诊断不可用（{}）；没有完整覆盖计划终点的新鲜序列时，"
                    "不推断终点趋势。".format(
                        _markdown_text(endpoint_trend.get("reason") or "原因未提供")
                    )
                )
            endpoint_lines.append("")
        lines.extend(["## 约束与验证", ""])
        baseline_constraints = constraint_map(baseline)
        search_constraints = constraint_map(search_comparison)
        verified_constraints = constraint_map(best)
        if self.plan.spec.constraints:
            lines.extend(
                [
                    f"| 约束 | 下界 | 上界 | 容差 | 基线 | {search_constraint_column} | 新鲜复验 | 复验余量 | 统计量单位 |",
                    "| --- | ---: | ---: | ---: | ---: | ---: | ---: | ---: | --- |",
                ]
            )
            for constraint in self.plan.spec.constraints:
                baseline_item = baseline_constraints.get(constraint.id, {})
                search_item = search_constraints.get(constraint.id, {})
                verified_item = verified_constraints.get(constraint.id, {})
                lines.append(
                    "| {id} | {lower} | {upper} | {tolerance} | {baseline} | {search} | {verified} | {margin} | {unit} |".format(
                        id=_markdown_text(constraint.id),
                        lower=(
                            _format_report_float(constraint.lower)
                            if constraint.lower is not None
                            else "—"
                        ),
                        upper=(
                            _format_report_float(constraint.upper)
                            if constraint.upper is not None
                            else "—"
                        ),
                        tolerance=_format_report_float(constraint.tolerance),
                        baseline=trial_number(baseline_item, "value"),
                        search=trial_number(search_item, "value"),
                        verified=trial_number(verified_item, "value"),
                        margin=trial_number(verified_item, "margin"),
                        unit=_markdown_text(
                            _metric_unit(
                                constraint.expected_unit,
                                constraint.statistic.kind,
                            )
                            or "dimensionless"
                        ),
                    )
                )
            lines.append("")
        else:
            lines.extend(["- 未声明响应约束。", ""])
        verification = result.get("verification")
        if isinstance(verification, Mapping):
            lines.append(
                "- 新鲜复验通过：{}".format(
                    "是" if verification.get("passed") is True else "否"
                )
            )
            lines.append(
                "- 复验结束时源工程未变化：{}".format(
                    "是" if verification.get("sourceUnchanged") is True else "否"
                )
            )
            lines.append(
                f"- 相对/绝对复验容差：`{verification.get('relativeTolerance')}` / `{verification.get('absoluteTolerance')}`"
            )
            raw_metrics = verification.get("metrics")
            if isinstance(raw_metrics, list):
                for metric in raw_metrics:
                    if isinstance(metric, Mapping):
                        lines.append(
                            "- `{id}`：搜索 `{search}`，复验 `{verification}`，匹配 {matches}".format(
                                id=_markdown_text(metric.get("id")),
                                search=trial_number(metric, "search"),
                                verification=trial_number(metric, "verification"),
                                matches=(
                                    "是" if metric.get("matches") is True else "否"
                                ),
                            )
                        )
        else:
            lines.append("- 未生成可验证的可行候选。")
        lines.append("")
        lines.extend(endpoint_lines)
        lines.extend(
            [
                "## 搜索过程",
                "",
                f"- 搜索停止：`{termination_reason}` — {termination_description}",
                f"- 搜索候选评估：{count_values.get('searchProposals', '—')}",
                f"- 搜索请求（计入仿真预算）：{count_values.get('searchSubmissionSlotsConsumed', '—')} / {self.search_run_limit}",
                f"- 搜索完成 / 失败或不完整：{count_values.get('searchCompletedTrials', '—')} / {count_values.get('searchFailedTrials', '—')}",
                f"- 独立复验请求（计入仿真预算）：{count_values.get('verificationSubmissionSlotsConsumed', '—')}；完成 / 失败或不完整：{count_values.get('verificationCompletedTrials', '—')} / {count_values.get('verificationFailedTrials', '—')}",
                f"- 计入仿真预算的请求总数：{count_values.get('backendSubmissionSlotsConsumed', '—')}",
                f"- 未形成试验记录的预算请求：{count_values.get('unrecordedSubmissionSlotsConsumed', '—')}",
                f"- 重复候选（缓存命中）：{count_values.get('cacheHits', '—')}",
                f"- 完整代：{generation_values.get('completed', '—')}；产生新仿真的代：{generation_values.get('withNewBackendSubmissions', '—')}；纯重复代：{generation_values.get('withoutNewBackendSubmissions', '—')}",
                f"- 未使用的搜索仿真额度：{count_values.get('remainingSearchRunBudget', '—')}",
                f"- 仿真耗时：总计 {trial_number(timing_values, 'simulationDurationTotalSeconds')} s；单次最短 / 最长 / 平均 {trial_number(timing_values, 'simulationDurationMinimumSeconds')} / {trial_number(timing_values, 'simulationDurationMaximumSeconds')} / {trial_number(timing_values, 'simulationDurationMeanSeconds')} s",
                "",
                "## 算法与预算",
                "",
                (
                    "- 算法：标准库实现的 `DE/rand/1/bin`，每代冻结供体种群、"
                    "越界坐标反射回可行域，顺序执行，不做局部抛光"
                ),
                endpoint_seeding_text,
                f"- 随机种子：`{self.plan.spec.algorithm.seed}`",
                f"- 种群大小：{self.plan.spec.algorithm.population_size}",
                f"- 变异因子 / 交叉概率：{_format_report_float(self.plan.spec.algorithm.mutation_factor)} / {_format_report_float(self.plan.spec.algorithm.crossover_probability)}",
                f"- 后端仿真上限：{self.plan.spec.budget.max_simulation_runs}（搜索最多 {self.search_run_limit}，预留 1 次新鲜复验）",
                f"- 候选都唯一时预算可覆盖的完整 DE 代数：{full_generations_with_unique_candidates}",
                f"- 搜索软墙钟上限：{_format_report_float(self.plan.spec.budget.max_wall_seconds)} s",
                "",
                "## 审计信息",
                "",
                f"- 源工程：`{_markdown_text(self.plan.source.path)}`",
                f"- 源 SHA-256：`{self.plan.source.sha256}`",
                f"- 计划哈希：`{self.plan.confirmation_token}`",
                f"- 搜索策略：`{SEARCH_POLICY}`",
                "",
                "## 产物",
                "",
                "- `optimization-result.json`：权威结构化汇总",
                "- `evaluations.csv`：逐次已记录后端提交尝试摘要",
                "- `optimization-events.jsonl`：优化与内层仿真进度",
                "- `checkpoint.json`：中断时的已完成试验快照",
            ]
        )
        if result.get("solutionStatus") == "verified":
            lines.extend(
                [
                    "- `best-project.json` / `best-system.xml`：验证后的派生模型；源文件未覆盖",
                    "- `result.json` / `results.csv` / `curve-*.svg`：最终新鲜复验结果",
                    "",
                    "派生工程会把被优化参数的原表达式替换为普通 SI 数值；其他工程字段保持不变。",
                ]
            )
        return "\n".join(lines) + "\n"

    def _result_payload(
        self,
        *,
        termination_reason: str,
        solution_status: str,
        baseline: Trial | None,
        best_search: Trial | None,
        best_diagnostic_search: Trial | None,
        verification_trial: Trial | None,
        verification: Mapping[str, object] | None,
        artifacts: Mapping[str, object],
        objective_endpoint_trend: Mapping[str, object] | None = None,
        source_unchanged: bool | None = None,
        error: Mapping[str, object] | None = None,
    ) -> dict[str, object]:
        source_status = (
            self._source_is_unchanged()
            if source_unchanged is None
            else source_unchanged
        )
        claim = {
            "verified": (
                "best verified feasible point among completed search simulations "
                "within the declared budget; "
                "not a proof of global optimality, search convergence, or steady state"
            ),
            "noFeasibleCandidate": (
                "no complete feasible candidate found within the declared budget"
            ),
            "verificationFailed": (
                "a feasible search candidate was found but fresh verification failed"
            ),
            "aborted": "optimization aborted before a verified solution was established",
        }.get(solution_status, "no verified solution was established")
        verification_payload = (
            dict(verification) if verification is not None else None
        )
        if verification_payload is not None:
            verification_payload["backendSubmissions"] = (
                self.verification_backend_submissions
            )
            verification_payload["submissionSlotsConsumed"] = (
                self.verification_backend_submissions
            )
            verification_payload["meaning"] = (
                "fresh candidate reproducibility and feasibility; not search convergence"
            )
        result: dict[str, object] = {
            "optimizationResultSchemaVersion": OPTIMIZATION_RESULT_SCHEMA_VERSION,
            "optimizationId": self.optimization_id,
            "solutionStatus": solution_status,
            "terminationReason": termination_reason,
            "search": self._search_summary(termination_reason),
            "claim": claim,
            "source": {
                "path": str(self.plan.source.path),
                "sha256": self.plan.source.sha256,
                "unchanged": source_status,
            },
            "spec": {
                "path": str(self.plan.spec_source.path),
                "sha256": self.plan.spec_source.sha256,
                "resolved": self.plan.spec.as_dict(),
            },
            "planHash": self.plan.confirmation_token,
            "continuityPolicy": {
                "id": CONTINUITY_POLICY,
                "machineVerified": False,
                "confirmedUserAssertion": True,
                "designVariableIds": [
                    item.spec.id
                    for item in self.plan.resolved_design_variables
                ],
            },
            "objectiveMetric": {
                "resultKey": self.plan.spec.objective.result_key,
                "statistic": self.plan.spec.objective.statistic.as_dict(),
                "seriesUnit": self.plan.spec.objective.expected_unit,
                "unit": _metric_unit(
                    self.plan.spec.objective.expected_unit,
                    self.plan.spec.objective.statistic.kind,
                ),
            },
            "environment": {
                "python": platform.python_version(),
                "platform": platform.platform(),
                "statisticImplementationVersion": STATISTIC_IMPLEMENTATION_VERSION,
                "backendBaseUrl": self.plan.base_url,
                "backendBuildFingerprint": None,
            },
            "counts": self._counts_payload(),
            "timing": self._timing_payload(),
            "baseline": baseline.as_dict() if baseline is not None else None,
            "bestSearch": best_search.as_dict() if best_search is not None else None,
            "bestDiagnosticSearch": (
                best_diagnostic_search.as_dict()
                if best_diagnostic_search is not None
                else None
            ),
            "verifiedBest": (
                verification_trial.as_dict()
                if solution_status == "verified" and verification_trial is not None
                else None
            ),
            "verification": verification_payload,
            "objectiveEndpointTrend": (
                dict(objective_endpoint_trend)
                if objective_endpoint_trend is not None
                else None
            ),
            "improvement": self._improvement(
                baseline,
                verification_trial if solution_status == "verified" else None,
            ),
            "trials": [trial.as_dict() for trial in self.trials],
            "artifacts": dict(artifacts),
        }
        if error is not None:
            result["error"] = dict(error)
        return result

    def _finalize_result(self, payload: dict[str, object]) -> dict[str, object]:
        report_path = self.output_directory / "report.md"
        _atomic_write(report_path, self._report_markdown(payload).encode("utf-8"))
        artifacts = payload.get("artifacts")
        if isinstance(artifacts, dict):
            supporting_artifacts = {
                "plan": self.output_directory / "optimization-plan.json",
                "report": report_path,
                "evaluations": self.output_directory / "evaluations.csv",
                "events": self.events_path,
                "checkpoint": self.output_directory / "checkpoint.json",
            }
            for name, path in supporting_artifacts.items():
                if path.is_file():
                    artifacts[name] = _artifact(path)
        result_path = self.output_directory / "optimization-result.json"
        _atomic_json(result_path, payload)
        return payload

    def run(self) -> dict[str, object]:
        simulation._prepare_output_directory(str(self.output_directory))
        try:
            self.progress_directory.mkdir()
            self._events = self.events_path.open(
                "x", encoding="utf-8", newline="\n"
            )
        except OSError as exc:
            raise simulation.ArtifactError(
                "OPTIMIZATION_OUTPUT_INITIALIZATION_FAILED",
                str(exc),
                {"path": str(self.output_directory)},
            ) from exc
        try:
            _atomic_json(
                self.output_directory / "optimization-plan.json",
                self.plan.public_dict(),
            )
            self._emit(
                {
                    "event": "optimization-started",
                    "optimizationId": self.optimization_id,
                    "source": str(self.plan.source.path),
                    "sourceSha256": self.plan.source.sha256,
                    "maxSimulationRuns": self.plan.spec.budget.max_simulation_runs,
                    "searchRunLimit": self.search_run_limit,
                    "outputDirectory": str(self.output_directory),
                }
            )
            termination_reason, _ = self._search()
            baseline = next(
                (
                    trial
                    for trial in self.trials
                    if trial.stage == "search" and trial.evaluation_id == 1
                ),
                None,
            )
            best_search = self._best_trial(feasible_only=True)
            diagnostic_best = self._best_trial(feasible_only=False)
            if best_search is None:
                self._checkpoint(termination_reason)
                payload = self._result_payload(
                    termination_reason=termination_reason,
                    solution_status="noFeasibleCandidate",
                    baseline=baseline,
                    best_search=None,
                    best_diagnostic_search=diagnostic_best,
                    verification_trial=None,
                    verification=None,
                    artifacts={},
                    objective_endpoint_trend=_objective_endpoint_trend(
                        {},
                        self.plan.spec.objective,
                        expected_end=_planned_objective_end(self.plan),
                    ),
                )
                self._emit(
                    {
                        "event": "optimization-finished",
                        "optimizationId": self.optimization_id,
                        "solutionStatus": payload["solutionStatus"],
                        "terminationReason": payload["terminationReason"],
                    }
                )
                return self._finalize_result(payload)
            if not self._source_is_unchanged():
                raise OptimizationError(
                    "OPTIMIZATION_SOURCE_CHANGED",
                    "The source project changed after planning; fresh verification was not started.",
                )

            best_vector = [
                (best_search.parameters[resolved.spec.id] - resolved.spec.lower)
                / (resolved.spec.upper - resolved.spec.lower)
                for resolved in self.plan.resolved_design_variables
            ]
            verification_trial, verification_result = self._evaluate(
                best_vector,
                stage="verification",
                bypass_cache=True,
                retain_result=True,
            )
            comparison = self._verification_comparison(
                best_search, verification_trial
            )
            objective_endpoint_trend = _objective_endpoint_trend(
                verification_result or {},
                self.plan.spec.objective,
                expected_end=_planned_objective_end(self.plan),
            )
            verified = comparison["passed"] is True and verification_result is not None
            solution_status = "verified" if verified else "verificationFailed"
            artifacts: dict[str, object] = {}
            if verified:
                assert verification_result is not None
                artifacts.update(
                    self._write_best_artifacts(
                        verification_trial, verification_result
                    )
                )
                final_source_unchanged = self._source_is_unchanged()
                if not final_source_unchanged:
                    self._discard_generated_artifacts(artifacts)
                    artifacts.clear()
                    comparison["sourceUnchanged"] = False
                    comparison["passed"] = False
                    verified = False
                    solution_status = "verificationFailed"
            else:
                final_source_unchanged = self._source_is_unchanged()
            self._checkpoint(termination_reason)
            payload = self._result_payload(
                termination_reason=termination_reason,
                solution_status=solution_status,
                baseline=baseline,
                best_search=best_search,
                best_diagnostic_search=None,
                verification_trial=verification_trial,
                verification=comparison,
                artifacts=artifacts,
                objective_endpoint_trend=objective_endpoint_trend,
                source_unchanged=final_source_unchanged,
            )
            self._emit(
                {
                    "event": "optimization-finished",
                    "optimizationId": self.optimization_id,
                    "solutionStatus": payload["solutionStatus"],
                    "terminationReason": payload["terminationReason"],
                }
            )
            return self._finalize_result(payload)
        finally:
            if self._events is not None:
                self._events.close()
                self._events = None

    def finalize_error(
        self,
        code: str,
        message: str,
        *,
        termination_reason: str,
    ) -> None:
        if not self.output_directory.exists():
            return
        try:
            self._checkpoint(termination_reason)
            baseline = self.trials[0] if self.trials else None
            best_feasible = self._best_trial(feasible_only=True)
            best_diagnostic = (
                None
                if best_feasible is not None
                else self._best_trial(feasible_only=False)
            )
            payload = self._result_payload(
                termination_reason=termination_reason,
                solution_status="aborted",
                baseline=baseline,
                best_search=best_feasible,
                best_diagnostic_search=best_diagnostic,
                verification_trial=None,
                verification=None,
                artifacts={},
                objective_endpoint_trend=_objective_endpoint_trend(
                    {},
                    self.plan.spec.objective,
                    expected_end=_planned_objective_end(self.plan),
                ),
                error={"code": code, "message": message},
            )
            self._finalize_result(payload)
        except (OSError, ValueError, simulation.SkillCliError):
            pass


def command_optimize(args: argparse.Namespace) -> int:
    if not args.confirmed:
        raise simulation.InputError(
            "OPTIMIZATION_CONFIRMATION_REQUIRED",
            "Run optimization-plan first and explicitly confirm that exact plan.",
        )
    plan = build_runtime_plan(
        args.input,
        args.spec,
        args.output_dir,
        base_url=args.base_url,
        timeout=args.timeout,
    )
    mismatches: dict[str, object] = {}
    if args.expected_source_sha256 != plan.source.sha256:
        mismatches["sourceSha256"] = {
            "expected": args.expected_source_sha256,
            "actual": plan.source.sha256,
        }
    if args.expected_spec_sha256 != plan.spec_source.sha256:
        mismatches["specSha256"] = {
            "expected": args.expected_spec_sha256,
            "actual": plan.spec_source.sha256,
        }
    if args.confirmation_token != plan.confirmation_token:
        mismatches["confirmationToken"] = "does not match the current plan"
    if mismatches:
        raise simulation.InputError(
            "OPTIMIZATION_CONFIRMATION_STALE",
            "Source, spec, backend conversion, or output path changed after planning.",
            mismatches,
        )
    optimization_id = args.optimization_id or f"opt-{uuid.uuid4().hex[:24]}"
    if OPTIMIZATION_ID_PATTERN.fullmatch(optimization_id) is None:
        raise simulation.InputError(
            "OPTIMIZATION_ID_INVALID",
            "Optimization ID must match [A-Za-z0-9._-]{1,96}.",
        )
    runner = OptimizationRunner(plan, optimization_id)
    try:
        result = runner.run()
    except KeyboardInterrupt:
        runner.cancel_active()
        runner.finalize_error(
            "OPTIMIZATION_CANCELLED_BY_USER",
            "Optimization was interrupted by the user.",
            termination_reason="userCancelled",
        )
        raise
    except simulation.SkillCliError as exc:
        runner.cancel_active()
        runner.finalize_error(
            exc.code,
            exc.message,
            termination_reason="error",
        )
        raise
    simulation.emit_json(
        {
            "event": "optimization-completed",
            "optimizationId": optimization_id,
            "solutionStatus": result.get("solutionStatus"),
            "terminationReason": result.get("terminationReason"),
            "result": str(plan.output_directory / "optimization-result.json"),
            "report": str(plan.output_directory / "report.md"),
            "outputDirectory": str(plan.output_directory),
        }
    )
    return 0 if result.get("solutionStatus") == "verified" else 4


def build_parser() -> argparse.ArgumentParser:
    parser = argparse.ArgumentParser(
        description=(
            "Plan and run bounded, single-objective optimization for a "
            "SystemSimulationApp ReactFlow project."
        )
    )
    parser.add_argument(
        "--base-url",
        default=simulation.DEFAULT_BASE_URL,
        help="Loopback FastAPI base URL (default: %(default)s).",
    )
    parser.add_argument(
        "--timeout",
        default=simulation.DEFAULT_TIMEOUT_SECONDS,
        type=simulation.positive_timeout,
        help="HTTP/heartbeat read timeout in seconds (default: %(default)s).",
    )
    subparsers = parser.add_subparsers(dest="command", required=True)

    plan_parser = subparsers.add_parser(
        "plan",
        help="Validate and preview an immutable optimization plan without running it.",
    )
    plan_parser.add_argument("input")
    plan_parser.add_argument("--spec", required=True)
    plan_parser.add_argument("--output-dir", required=True)
    plan_parser.add_argument(
        "--plan-file",
        help=(
            "Write the complete auditable plan to this new path outside the "
            "optimization output directory; requires --present."
        ),
    )
    plan_parser.add_argument(
        "--present",
        action="store_true",
        help=(
            "Emit only the user-facing plan allowlist; requires --plan-file."
        ),
    )
    plan_parser.set_defaults(handler=command_plan)

    optimize_parser = subparsers.add_parser(
        "optimize",
        help="Execute an explicitly confirmed optimization plan.",
    )
    optimize_parser.add_argument("input")
    optimize_parser.add_argument("--spec", required=True)
    optimize_parser.add_argument("--output-dir", required=True)
    optimize_parser.add_argument("--expected-source-sha256", required=True)
    optimize_parser.add_argument("--expected-spec-sha256", required=True)
    optimize_parser.add_argument("--confirmation-token", required=True)
    optimize_parser.add_argument("--confirmed", action="store_true")
    optimize_parser.add_argument("--optimization-id")
    optimize_parser.set_defaults(handler=command_optimize)
    return parser


def main(argv: Sequence[str] | None = None) -> int:
    simulation.configure_standard_streams()
    parser = build_parser()
    args = parser.parse_args(argv)
    try:
        args.base_url = simulation.validate_base_url(args.base_url)
        return int(args.handler(args))
    except simulation.SkillCliError as exc:
        simulation.emit_json(
            simulation.stable_error_payload(exc),
            stream=sys.stderr,
        )
        return exc.exit_code
    except KeyboardInterrupt:
        simulation.emit_json(
            simulation.stable_error_payload(
                OptimizationError(
                    "OPTIMIZATION_CANCELLED_BY_USER",
                    "Optimization was interrupted and no additional candidate will be started.",
                )
            ),
            stream=sys.stderr,
        )
        return 4


if __name__ == "__main__":
    raise SystemExit(main())