diff --git a/Cargo.lock b/Cargo.lock
index 0d0aba83d..ae259c8b2 100644
--- a/Cargo.lock
+++ b/Cargo.lock
@@ -449,6 +449,7 @@ dependencies = [
  "asap-types",
  "asap_sketchlib",
  "good_lp",
+ "highs",
  "serde",
  "serde_json",
  "thiserror 2.0.18",
diff --git a/crates/devtools/tests/stage_pipeline.rs b/crates/devtools/tests/stage_pipeline.rs
index 56dca67cd..1eef50103 100644
--- a/crates/devtools/tests/stage_pipeline.rs
+++ b/crates/devtools/tests/stage_pipeline.rs
@@ -1748,7 +1748,9 @@ fn search_flag_picks_the_search_algorithm() {
     assert!(report.get("routed_to_legacy").is_none(), "{report}");
     let milp_report = &report["milp"];
     assert!(milp_report["rounds"].as_u64().unwrap() >= 1, "{report}");
-    assert!(milp_report["gap"].as_f64().unwrap() <= 1e-9, "{report}");
+    // Certified: the gap is within the search's tolerance (1e-4), which
+    // also absorbs the solver's and the incumbent cutoff's tolerances.
+    assert!(milp_report["gap"].as_f64().unwrap() <= 1e-4, "{report}");
     assert!(
         milp_report["lower_bound"].as_f64().unwrap() <= price(&milp) * (1.0 + 1e-9),
         "{report}"
@@ -1789,8 +1791,9 @@ fn search_flag_is_checked() {
     assert!(stderr.contains("drop it with --exhaustive"), "{stderr}");
 }
 
-/// The MILP search stopped by its time limit before any valid plan says
-/// so (never a plan called optimal), and the limit needs the MILP search.
+/// The MILP search stopped by its time limit before its first round
+/// returns the warm-start plan (each query planned alone) with its gap,
+/// never called optimal, and the limit needs the MILP search.
 #[test]
 fn milp_time_limit_is_reported_and_needs_the_milp_search() {
     let document = document_of(&[
@@ -1802,11 +1805,15 @@ fn milp_time_limit_is_reported_and_needs_the_milp_search() {
         "0",
     ]);
     let selection = &document["stage3_plan_selection"];
-    assert!(selection["selected"].is_null(), "{selection}");
+    assert!(selection["selected"].is_string(), "{selection}");
+    let milp = &selection["search"]["milp"];
+    assert_eq!(milp["time_limit_reached"], true, "{selection}");
+    assert_eq!(milp["rounds"], 0, "{selection}");
+    let gap = milp["gap"].as_f64().unwrap();
+    assert!(gap > 1e-4 && gap <= 1.0, "{selection}");
     assert!(
-        selection["error"]
-            .as_str()
-            .is_some_and(|e| e.contains("reached its time limit")),
+        milp["lower_bound"].as_f64().unwrap()
+            <= selection["search"]["incumbent_cost"].as_f64().unwrap(),
         "{selection}"
     );
     let (_, stderr) = failure(&[
diff --git a/crates/integration-tests/tests/joint_search_equivalence.rs b/crates/integration-tests/tests/joint_search_equivalence.rs
index 2baf63a3f..196052bc4 100644
--- a/crates/integration-tests/tests/joint_search_equivalence.rs
+++ b/crates/integration-tests/tests/joint_search_equivalence.rs
@@ -23,6 +23,10 @@ use asap_types::workload::{
 };
 use planner_layering_common::*;
 
+/// The largest gap a certified search reports (the planner's
+/// `GAP_TOLERANCE`): the solver's own tolerances.
+const CERTIFIED_GAP: f64 = 1e-4;
+
 /// Every cost model, each objective of cost model 3 included, with
 /// `memory_gb` of cluster memory: small enough that some candidates do not
 /// fit, so the limits are exercised too.
@@ -1486,7 +1490,9 @@ fn three_way(
         .as_ref()
         .and_then(|milp| milp.gap)
         .expect("the MILP reports its gap");
-    if gap.abs() > 1e-9 {
+    // Certified: within the solver's tolerances (the planner's
+    // GAP_TOLERANCE), never reported as exactly 0 unless proven.
+    if gap.abs() > CERTIFIED_GAP {
         return Err(format!("{label}: MILP gap {gap}"));
     }
     Ok(true)
@@ -1573,6 +1579,63 @@ fn assert_three_way_on_seeds(
     assert!(mismatches.is_empty(), "{}", mismatches.join("\n"));
 }
 
+/// Queries over two gauges, two per gauge, are two independent clusters of
+/// the memo, each searched with its own MILP and lazy loop: with raw data
+/// kept and without, the MILP reports both clusters, each certified, and
+/// the three searches select plans of one price.
+#[test]
+fn three_way_on_two_independent_clusters() {
+    let mut workload = random_workload_over(3, &["m", "n"]);
+    let queries = [
+        "quantile_over_time(0.9, m[5m])",
+        "sum_over_time(m[1m])",
+        "max_over_time(n[5m])",
+        "quantile_over_time(0.99, n[10m])",
+    ];
+    workload.query_workload.repeating_queries = workload_of(
+        &queries.map(|q| (q, 60_000)),
+        AccuracyTarget::max_error(0.05),
+    )
+    .query_workload
+    .repeating_queries;
+    let roots = lower_promql_roots(&workload).unwrap();
+    let entries: Vec<_> = workload.query_workload.entries().collect();
+    let data = workload.data_workload.clone().unwrap();
+    for deployment in [
+        asap_logical_asap_optimization::DeploymentCapabilities::default(),
+        no_raw_data(),
+    ] {
+        let models = PlanningModels::builtin().with_deployment_capabilities(&deployment);
+        let label = format!("raw data retained {}", deployment.raw_data_retained);
+        assert!(
+            assert_three_way(&label, &roots, &entries, &data, models),
+            "{label}"
+        );
+        let run = run_stage1_to_stage3_with(
+            roots.iter().cloned().enumerate().collect(),
+            &entries,
+            &data,
+            models,
+            Stage1ToStage3Search::with_algorithm(SearchAlgorithm::MemoMilp),
+        )
+        .unwrap();
+        let milp = run.stage3_plan_selection.unwrap().search.milp.unwrap();
+        let clusters: Vec<_> = milp.clusters.iter().map(|c| c.queries.clone()).collect();
+        assert_eq!(clusters, [vec![0, 1], vec![2, 3]], "{label}");
+        assert!(
+            milp.clusters
+                .iter()
+                .all(|c| c.gap.is_some_and(|gap| gap <= CERTIFIED_GAP)),
+            "{label}: {milp:?}"
+        );
+        assert_eq!(
+            milp.rounds,
+            milp.clusters.iter().map(|c| c.rounds).sum::<u32>(),
+            "{label}"
+        );
+    }
+}
+
 /// The three searches agree on the random workloads under cost model 1.
 #[test]
 #[ignore = "the exhaustive search and the brute force over every seed take minutes in a debug build; run with --release --ignored"]
@@ -1843,6 +1906,76 @@ fn a_sketch_rollup_beside_other_quantiles_selects_the_exhaustive_price() {
     assert_memo_selects_the_exhaustive_price("Q2, Q3, Q4", &workload(&[q2, q3, q4]));
 }
 
+/// #839's workload (two 10-minute sums every 20 s, two 100-second
+/// quantiles of one metric at different label sets every second): the
+/// cheapest plan widens the quantiles' scans to one label set and shares
+/// one scan. With its term cuts (a source's scans apart or widened, held
+/// memory, ingestion copies) the memo MILP certifies the memo brute
+/// force's price.
+#[test]
+fn term_cuts_certify_the_brute_force_price_on_839() {
+    let entries: Vec<serde_json::Value> = [
+        (
+            "sum by (region, service, job) (sum_over_time(m[10m]))",
+            20000,
+        ),
+        ("sum by (region, service) (sum_over_time(m[10m]))", 20000),
+        (
+            "avg by (region, service, job) (quantile_over_time(0.99, m[100s]))",
+            1000,
+        ),
+        ("avg by (region) (quantile_over_time(0.95, m[100s]))", 1000),
+    ]
+    .iter()
+    .map(|(q, interval)| {
+        serde_json::json!({
+            "query": q,
+            "demand": {"fixed_interval": interval},
+            "requirements": {"accuracy": {"max_error": 0.05, "failure_probability": 0.3}},
+            "predictability": {"predictable": {"known_at": null}}
+        })
+    })
+    .collect();
+    let workload: PlanningWorkload = serde_json::from_value(serde_json::json!({
+        "query_workload": {"language": "promql", "query_batch": null, "repeating_queries": entries},
+        "data_workload": {
+            "arrival": {"ingesting": {"rows_per_second": 100000.0}},
+            "sources": {"m": {"series_count": 10000, "metric_type": "gauge", "values_non_negative": true, "sample_interval": 100}}
+        }
+    }))
+    .unwrap();
+    let roots = lower_promql_roots(&workload).unwrap();
+    let entries: Vec<_> = workload.query_workload.entries().collect();
+    let data = workload.data_workload.clone().unwrap();
+    let models = PlanningModels::builtin();
+    let run = |search| {
+        run_stage1_to_stage3_with(
+            roots.iter().cloned().enumerate().collect(),
+            &entries,
+            &data,
+            models,
+            search,
+        )
+        .unwrap()
+        .stage3_plan_selection
+        .unwrap()
+    };
+    let brute_force = run(Stage1ToStage3Search::with_algorithm(
+        SearchAlgorithm::MemoBruteForce,
+    ));
+    let milp = run(Stage1ToStage3Search::with_algorithm(
+        SearchAlgorithm::MemoMilp,
+    ));
+    let price = |s: &asap_plan_selection::Selection| s.costs[&s.selected].cost_per_second.unwrap();
+    let (b, m) = (price(&brute_force), price(&milp));
+    assert!(
+        (b - m).abs() <= 1e-9 * b,
+        "memo brute force {b}, memo MILP {m}"
+    );
+    let gap = milp.search.milp.and_then(|milp| milp.gap);
+    assert!(gap.is_some_and(|gap| gap <= CERTIFIED_GAP), "gap {gap:?}");
+}
+
 /// Debugging aid: the random workload of seed `JP_SEED` (over the gauge
 /// `m`) as stage_pipeline's `--workload` JSON, at `JP_SEED_OUT`.
 #[test]
@@ -1853,3 +1986,62 @@ fn dump_random_workload() {
     let workload = random_workload(seed);
     std::fs::write(&out, serde_json::to_string_pretty(&workload).unwrap()).unwrap();
 }
+
+/// With no time left for a single MILP round, memo-milp still returns a
+/// plan: the warm-start incumbent, each query planned alone by the legacy
+/// search with no sharing, offered before round 1. Under cost model 1 its
+/// cost is the sum of the queries' own plans' costs, and the run reports
+/// the time limit and a gap to the MILP lower bound.
+#[test]
+fn a_time_limited_memo_milp_run_returns_the_warm_start_plan() {
+    let queries = [
+        "quantile_over_time(0.9, m[5m])",
+        "sum_over_time(m[1m])",
+        "quantile_over_time(0.99, m[10m])",
+    ];
+    let workload = workload_of(
+        &queries.map(|q| (q, 60_000)),
+        AccuracyTarget::max_error(0.05),
+    );
+    let roots = lower_promql_roots(&workload).unwrap();
+    let entries: Vec<_> = workload.query_workload.entries().collect();
+    let data = workload.data_workload.clone().unwrap();
+    let models = PlanningModels::builtin();
+    let alone: f64 = roots
+        .iter()
+        .zip(&entries)
+        .map(|(root, entry)| {
+            let run = run_stage1_to_stage3_with(
+                vec![(0, root.clone())],
+                std::slice::from_ref(entry),
+                &data,
+                models,
+                Stage1ToStage3Search::default(),
+            )
+            .unwrap();
+            let selection = run.stage3_plan_selection.unwrap();
+            selection.costs[&selection.selected]
+                .cost_per_second
+                .unwrap()
+        })
+        .sum();
+    let run = run_stage1_to_stage3_with(
+        roots.iter().cloned().enumerate().collect(),
+        &entries,
+        &data,
+        models,
+        Stage1ToStage3Search {
+            milp_time_limit_ms: Some(0),
+            ..Stage1ToStage3Search::with_algorithm(SearchAlgorithm::MemoMilp)
+        },
+    )
+    .unwrap();
+    let selection = run.stage3_plan_selection.unwrap();
+    assert_eq!(selection.selected, "P0");
+    let cost = selection.search.incumbent_cost.unwrap();
+    assert!((cost - alone).abs() <= 1e-9 * alone, "{cost} vs {alone}");
+    let milp = selection.search.milp.unwrap();
+    assert!(milp.time_limit_reached, "{milp:?}");
+    assert!(milp.gap.is_some(), "{milp:?}");
+    assert!(run.stage3_selected_logical.is_some());
+}
diff --git a/crates/integration-tests/tests/tumbling_pane_merge.rs b/crates/integration-tests/tests/tumbling_pane_merge.rs
index 922222d7f..ffe767796 100644
--- a/crates/integration-tests/tests/tumbling_pane_merge.rs
+++ b/crates/integration-tests/tests/tumbling_pane_merge.rs
@@ -15,6 +15,7 @@ use asap_integration_tests::fixtures::lower_promql;
 use asap_integration_tests::planning::selected_dag;
 use asap_logical_asap_optimization::asap_aware_cse::{
     pane_candidates, window_composition_pass2_candidates,
+    window_composition_shared_chain_candidates,
 };
 use asap_types::ir::operator::operator_properties::TimeShift;
 use asap_types::ir::operator::pane_chain::expand_pane_chains;
@@ -854,6 +855,11 @@ fn queries_at_different_intervals_read_one_chain_of_the_finer_panes() {
         .into_iter()
         .find(|c| c.pane_candidates.len() == 2 && !c.shared_chains.is_empty())
         .expect("one shared chain");
+    // The whole set's shared chains, built alone, are the same candidates.
+    let alone = window_composition_shared_chain_candidates(&roots, &rule_candidates);
+    let first = alone.first().expect("one shared chain alone");
+    assert_eq!(first.shared_chains, shared.shared_chains);
+    assert_eq!(first.roots.as_ref().ok(), shared.roots.as_ref().ok());
     assert_eq!(
         (
             shared.shared_chains[0].pane_width_ms,
diff --git a/crates/logical-asap-optimization/src/asap_aware_cse.rs b/crates/logical-asap-optimization/src/asap_aware_cse.rs
index 417421154..24d7f9b95 100644
--- a/crates/logical-asap-optimization/src/asap_aware_cse.rs
+++ b/crates/logical-asap-optimization/src/asap_aware_cse.rs
@@ -502,6 +502,18 @@ pub fn identical_expression_pass2_candidates(
         .collect()
 }
 
+/// `roots` with the sharing candidates `subset` selects made one node each
+/// (and the ones that leaves one node): one identical-expression Pass 2
+/// candidate, for a search that decides which nodes to share rather than
+/// listing every closed set ([`identical_expression_pass2_candidates`]).
+pub fn with_shared_candidates(
+    roots: &[QueryRoot],
+    candidates: &[SharingCandidate],
+    subset: &[bool],
+) -> Vec<QueryRoot> {
+    apply(roots, candidates, subset).0
+}
+
 /// `roots` with the sharing candidates in `subset` shared, and the ones
 /// that leaves shared: those whose every node still read became one node (a
 /// node below a shared node other than the one kept is no longer read).
@@ -1261,29 +1273,54 @@ pub fn window_composition_pass2_candidates(
             pane_candidates: chosen.clone(),
             shared_chains: Vec::new(),
         });
-        // Shared chains on the members' own finest panes, then, when it
-        // groups them differently, on the gcd of their widths.
-        let on_member_panes = nested_window_groups(&chosen);
-        let on_gcd_panes = gcd_window_groups(&chosen);
-        push_shared_chains(roots, &chosen, &on_member_panes, &mut pass2_candidates);
-        if on_gcd_panes != on_member_panes {
-            push_shared_chains(roots, &chosen, &on_gcd_panes, &mut pass2_candidates);
-        }
-        // Coarser groupings read a finer grouping's chain, rolled up.
         let chosen_rollable: Vec<bool> = (0..candidates.len())
             .filter(|&i| mask >> i & 1 == 1)
             .map(|i| rollable[i])
             .collect();
-        push_shared_chains(
-            roots,
-            &chosen,
-            &rollup_window_groups(&chosen, &chosen_rollable),
-            &mut pass2_candidates,
-        );
+        push_every_shared_chain(roots, &chosen, &chosen_rollable, &mut pass2_candidates);
     }
     Ok(pass2_candidates)
 }
 
+/// The Pass 2 candidates of [`window_composition_pass2_candidates`] in
+/// which every one of `candidates` is in panes and reads a shared chain:
+/// that function's candidates of the whole set with shared chains, in its
+/// order, without building the other subsets'.
+pub fn window_composition_shared_chain_candidates(
+    roots: &[QueryRoot],
+    candidates: &[PaneCandidate],
+) -> Vec<WindowCompositionPass2Candidate> {
+    let rollable = rollable_pane_candidates(roots, candidates);
+    let mut pass2_candidates = Vec::new();
+    push_every_shared_chain(roots, candidates, &rollable, &mut pass2_candidates);
+    pass2_candidates
+}
+
+/// Push the shared-chain Pass 2 candidates of the subset `chosen` (each
+/// rollable or not, `chosen_rollable`) of [`window_composition_pass2_candidates`].
+fn push_every_shared_chain(
+    roots: &[QueryRoot],
+    chosen: &[PaneCandidate],
+    chosen_rollable: &[bool],
+    pass2_candidates: &mut Vec<WindowCompositionPass2Candidate>,
+) {
+    // Shared chains on the members' own finest panes, then, when it
+    // groups them differently, on the gcd of their widths.
+    let on_member_panes = nested_window_groups(chosen);
+    let on_gcd_panes = gcd_window_groups(chosen);
+    push_shared_chains(roots, chosen, &on_member_panes, pass2_candidates);
+    if on_gcd_panes != on_member_panes {
+        push_shared_chains(roots, chosen, &on_gcd_panes, pass2_candidates);
+    }
+    // Coarser groupings read a finer grouping's chain, rolled up.
+    push_shared_chains(
+        roots,
+        chosen,
+        &rollup_window_groups(chosen, chosen_rollable),
+        pass2_candidates,
+    );
+}
+
 /// Push the Pass 2 candidate in which each of `groups` reads its one shared
 /// chain: the groups whose members have a shared size, and the error of
 /// each group whose sizing fails. Nothing when no group has a shared size.
@@ -1976,6 +2013,37 @@ pub fn rollup_kind(family: &FieldDataType) -> Option<ExactKind> {
 /// read only by estimates, with the finest build of its algorithm over the
 /// same rows whose keys include its own ([`sketch_rollup_candidates`]).
 pub fn rollup_candidates(roots: &[QueryRoot]) -> Vec<RollupCandidate> {
+    rollup_candidates_of(roots, None)
+}
+
+/// The rollup candidate of coarse build `coarse` among
+/// [`rollup_candidates`]`(roots)`, the first when several, computed for that
+/// build alone (each coarse build's candidate depends only on `roots`, not
+/// on the other coarse builds).
+pub fn rollup_candidate_for(
+    roots: &[QueryRoot],
+    coarse: *const OperatorNode,
+) -> Option<RollupCandidate> {
+    rollup_candidates_of(roots, Some(coarse)).into_iter().next()
+}
+
+/// Whether adding `root` to a candidate's roots may change the rollup
+/// candidate of its coarse build `coarse`: `root` reaches `coarse` itself
+/// or a build that could answer it ([`finer_keys`]). Otherwise `root` adds
+/// no reader of `coarse` or of any build that could answer it (a reader
+/// reaches what it reads), no such build, and leaves the walk order of the
+/// other nodes as it was, so [`rollup_candidate_for`] is unchanged.
+pub fn may_change_rollup_candidate(coarse: &Rc<OperatorNode>, root: &Rc<OperatorNode>) -> bool {
+    OperatorNode::reachable(root)
+        .iter()
+        .any(|node| Rc::ptr_eq(node, coarse) || finer_keys(coarse, node).is_some())
+}
+
+/// [`rollup_candidates`], over the coarse builds that are `only` when set.
+fn rollup_candidates_of(
+    roots: &[QueryRoot],
+    only: Option<*const OperatorNode>,
+) -> Vec<RollupCandidate> {
     let mut nodes = Vec::new();
     let mut seen = std::collections::HashSet::new();
     for node in roots.iter().flat_map(root_operators) {
@@ -1999,8 +2067,9 @@ pub fn rollup_candidates(roots: &[QueryRoot]) -> Vec<RollupCandidate> {
             })
     };
     let builds: Vec<&Rc<OperatorNode>> = nodes.iter().filter(|n| grouped_exact(n)).collect();
+    let wanted = |coarse: &Rc<OperatorNode>| only.is_none_or(|only| Rc::as_ptr(coarse) == only);
     let mut candidates = Vec::new();
-    for &coarse in &builds {
+    for &coarse in builds.iter().filter(|coarse| wanted(coarse)) {
         let finals = finalized(coarse);
         // Every reader of the coarse state reads its finalized value.
         if finals.is_empty() || finals.len() != readers[&Rc::as_ptr(coarse)].len() {
@@ -2030,7 +2099,7 @@ pub fn rollup_candidates(roots: &[QueryRoot]) -> Vec<RollupCandidate> {
         }
         candidates.extend(best.map(|(_, candidate)| candidate));
     }
-    candidates.extend(sketch_rollup_candidates(&nodes, &readers));
+    candidates.extend(sketch_rollup_candidates(&nodes, &readers, only));
     candidates
 }
 
@@ -2046,6 +2115,7 @@ pub fn rollup_candidates(roots: &[QueryRoot]) -> Vec<RollupCandidate> {
 fn sketch_rollup_candidates(
     nodes: &[Rc<OperatorNode>],
     readers: &HashMap<*const OperatorNode, Vec<&Rc<OperatorNode>>>,
+    only: Option<*const OperatorNode>,
 ) -> Vec<RollupCandidate> {
     let builds: Vec<&Rc<OperatorNode>> = nodes.iter().filter(|n| regroupable(n)).collect();
     // The finer build with the most keys; among several (groupings that
@@ -2074,7 +2144,10 @@ fn sketch_rollup_candidates(
             .map(|(_, fine)| fine)
     };
     let mut candidates = Vec::new();
-    for &coarse in &builds {
+    for &coarse in builds
+        .iter()
+        .filter(|coarse| only.is_none_or(|only| Rc::as_ptr(coarse) == only))
+    {
         let Some(ASAPOp::SummaryAgg { family, .. }) = coarse.asap() else {
             continue;
         };
@@ -3925,6 +3998,32 @@ mod tests {
             .collect()
     }
 
+    /// One coarse build's rollup candidate, computed alone, is the one
+    /// [`rollup_candidates`] lists for it; a query reaching no build that
+    /// could answer it may not change it, one reaching a finer build may.
+    #[test]
+    fn one_coarse_builds_rollup_candidate_is_the_listed_one() {
+        let roots = exact(&[
+            "sum by (job) (m)",
+            "sum by (job, instance) (m)",
+            "sum by (job, instance, zone) (m)",
+        ]);
+        let listed = rollup_candidates(&roots);
+        assert!(!listed.is_empty());
+        for candidate in &listed {
+            let alone = rollup_candidate_for(&roots, Rc::as_ptr(&candidate.coarse)).unwrap();
+            assert!(Rc::ptr_eq(&alone.finer, &candidate.finer));
+        }
+        let coarse = &builds(&roots[0])[0];
+        assert!(rollup_candidate_for(&roots[..1], Rc::as_ptr(coarse)).is_none());
+        assert!(may_change_rollup_candidate(coarse, operator(&roots[1])));
+        let other_metric = exact(&["sum by (job, instance) (other)"]);
+        assert!(!may_change_rollup_candidate(
+            coarse,
+            operator(&other_metric[0])
+        ));
+    }
+
     /// `sum by (job)` is answered from `sum by (job, instance)`, which the
     /// workload computes anyway: the coarse query reads the finer query's
     /// finalized sums, summed by `job`, so it builds nothing over the raw
diff --git a/crates/physical-asap-optimization/src/lib.rs b/crates/physical-asap-optimization/src/lib.rs
index 1cf83984c..77ee7aac6 100644
--- a/crates/physical-asap-optimization/src/lib.rs
+++ b/crates/physical-asap-optimization/src/lib.rs
@@ -481,6 +481,57 @@ pub fn optimize_physical_asap_assigned(
     assigned_candidate(from_logical, &implemented, &mat_candidates, assignment).map(Some)
 }
 
+/// [`optimize_physical_asap_chosen`]'s choice per materialization
+/// candidate: from the logical node it implements, its readers and its
+/// options.
+pub type ChooseMaterialization<'a> =
+    dyn FnMut(*const OperatorNode, &[usize], &[Materialization]) -> Materialization + 'a;
+
+/// Stage 2 for one CandidateLogicalASAPDAG and one materialization
+/// assignment chosen node by node: `choose(logical, readers, options)`
+/// gives the option for each materialization candidate, `logical` being
+/// the address of the logical node of `roots` it implements (null when it
+/// implements none, e.g. a node the implementation adds) and `readers` the
+/// queries that read it. `None` when a chosen option is not one of the
+/// candidate's, or the assignment materializes a build at ingestion time
+/// above one at query time. The search that decodes one plan from a model
+/// solution uses it to build that plan alone, not every assignment.
+pub fn optimize_physical_asap_chosen(
+    from_logical: &str,
+    roots: &[QueryRoot],
+    entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
+    deployment: &DeploymentCapabilities,
+    choose: &mut ChooseMaterialization<'_>,
+) -> Result<Option<PhysicalCandidate>, PhysicalASAPOptimizationError> {
+    let mut memo = HashMap::new();
+    let implemented = map_roots(roots, &mut |node| implement(node, &mut memo))?;
+    let logical_of: HashMap<*const OperatorNode, *const OperatorNode> = memo
+        .iter()
+        .map(|(logical, implemented)| (Rc::as_ptr(implemented), *logical))
+        .collect();
+    let mat_candidates = materialization_candidates(&implemented, entries, data, deployment);
+    check_panes_on_range_steps(&mat_candidates, entries)?;
+    let mut assignment = Vec::with_capacity(mat_candidates.len());
+    for mat_candidate in &mat_candidates {
+        let logical = logical_of
+            .get(&Rc::as_ptr(&mat_candidate.node))
+            .copied()
+            .unwrap_or(std::ptr::null());
+        let options = mat_candidate.options();
+        let m = choose(logical, &mat_candidate.readers, &options);
+        if !options.contains(&m) {
+            return Ok(None);
+        }
+        assignment.push(m);
+    }
+    let nodes: Vec<_> = mat_candidates.iter().map(|c| Rc::clone(&c.node)).collect();
+    if !down_closed(&mat_candidates, &candidates_below(&nodes), &assignment) {
+        return Ok(None);
+    }
+    assigned_candidate(from_logical, &implemented, &mat_candidates, assignment).map(Some)
+}
+
 /// The physical candidate of the implemented roots `implemented` with
 /// `assignment[c]` chosen for `mat_candidates[c]`.
 fn assigned_candidate(
diff --git a/crates/plan-selection/Cargo.toml b/crates/plan-selection/Cargo.toml
index af7854e82..a4c1b6d30 100644
--- a/crates/plan-selection/Cargo.toml
+++ b/crates/plan-selection/Cargo.toml
@@ -16,6 +16,7 @@ serde = { version = "1", features = ["derive"] }
 serde_json = "1"
 
 good_lp = { version = "=1.15.3", default-features = false, features = ["highs"] }
+highs = "=2.4.0"
 
 [dev-dependencies]
 asap-logical-asap-optimization = { path = "../logical-asap-optimization", features = ["test-support"] }
diff --git a/crates/plan-selection/src/lib.rs b/crates/plan-selection/src/lib.rs
index f7274e9d5..05b83afc5 100644
--- a/crates/plan-selection/src/lib.rs
+++ b/crates/plan-selection/src/lib.rs
@@ -1074,6 +1074,42 @@ pub struct MilpReport {
     /// solver limit (iterations, memory, interrupt): the plan is then
     /// proven optimal only if `gap` is 0.
     pub time_limit_reached: bool,
+    /// Solutions solution decoding could not compose into a plan: each
+    /// keeps its objective as a lower bound, so the gap stays honest.
+    #[serde(default)]
+    pub decode_failures: u32,
+    /// Each independent cluster's search, ordered by its first query: one
+    /// cluster when the memo was searched as one. The fields above are
+    /// their totals: bounds, rounds, model sizes and solver time summed,
+    /// the gap from the summed bound and incumbent cost.
+    #[serde(default, skip_serializing_if = "Vec::is_empty")]
+    pub clusters: Vec<MilpClusterReport>,
+}
+
+/// One independent cluster of the MILP search: queries whose plans share
+/// no cost with other clusters' under cost model 1, searched with their
+/// own model and lazy loop.
+#[derive(Debug, Clone, Default, PartialEq, serde::Serialize, serde::Deserialize)]
+pub struct MilpClusterReport {
+    /// Its queries, by index in the workload.
+    pub queries: Vec<usize>,
+    /// The memo groups of its own workload.
+    pub groups: usize,
+    /// Its lazy loop's rounds.
+    pub rounds: u32,
+    /// (incumbent cost − lower bound) / incumbent cost, at least 0.
+    pub gap: Option<f64>,
+    /// Its MILP lower bound at the end.
+    pub lower_bound: Option<f64>,
+    /// Its selected plan's cost per second.
+    pub incumbent_cost: Option<f64>,
+    /// The wall time in milliseconds its solves took, together.
+    pub solver_time_ms: u64,
+    /// Whether it stopped at the search's time limit.
+    pub time_limit_reached: bool,
+    /// Solutions solution decoding could not compose into a plan.
+    #[serde(default)]
+    pub decode_failures: u32,
 }
 
 /// The size of a sub-expression memo (planner-layering design, "The
@@ -2923,6 +2959,28 @@ pub fn resize_for_failure_probability(
     with_resized_summaries(roots, &resize)
 }
 
+/// `roots` with each summary at the most any target can make
+/// [`resize_for_failure_probability`] enlarge it: classic HLL at precision
+/// 18 (its failure probability only falls with precision), every other
+/// summary sized for a vanishing per-value failure probability. A failure
+/// probability this plan misses, it misses at every size Stage 3 can give.
+pub(crate) fn resized_to_the_largest(
+    roots: &[QueryRoot],
+) -> Result<Vec<QueryRoot>, SchemaDerivationError> {
+    with_resized_summaries(roots, &|_, params| {
+        Some(match params {
+            SketchParams::Hll {
+                precision,
+                relative_error,
+            } => SketchParams::Hll {
+                precision: (*precision).max(18),
+                relative_error: Some(hll::classic_relative_error(*precision, *relative_error)),
+            },
+            params => with_failure_probability(params, 1e-300),
+        })
+    })
+}
+
 /// What Stage 3 accepted for each query ([`AcceptedAccuracy`]), or the
 /// first query whose output is not shown to meet its target, as a reason. An output meets its target when its theoretical root guarantee
 /// (composed in Stage 1 with the built-in model) satisfies it, or when the
diff --git a/crates/plan-selection/src/search/legacy/components.rs b/crates/plan-selection/src/search/legacy/components.rs
index 7d8caf43d..5c314ccd0 100644
--- a/crates/plan-selection/src/search/legacy/components.rs
+++ b/crates/plan-selection/src/search/legacy/components.rs
@@ -658,7 +658,7 @@ fn combined_pass2<Id>(runs: &[Stage1ToStage3Run<Id>]) -> Option<Option<Pass2Cand
 /// `rejection` of a component's search with its query numbers (`q<n>` at
 /// the start of its id or reason, counted within the component) renamed to
 /// the workload's.
-fn renumber(rejection: Rejection, component: &[usize]) -> Rejection {
+pub(crate) fn renumber(rejection: Rejection, component: &[usize]) -> Rejection {
     let rename = |text: &str| -> String {
         let Some(rest) = text.strip_prefix('q') else {
             return text.to_string();
diff --git a/crates/plan-selection/src/search/memo/brute_force.rs b/crates/plan-selection/src/search/memo/brute_force.rs
index f1f208795..c29652166 100644
--- a/crates/plan-selection/src/search/memo/brute_force.rs
+++ b/crates/plan-selection/src/search/memo/brute_force.rs
@@ -436,7 +436,7 @@ fn build(
 }
 
 /// Mark `g` and every group below it through the groups' children.
-fn mark_below(memo: &Memo, g: usize, below: &mut [bool]) {
+pub(crate) fn mark_below(memo: &Memo, g: usize, below: &mut [bool]) {
     if below[g] {
         return;
     }
@@ -448,7 +448,7 @@ fn mark_below(memo: &Memo, g: usize, below: &mut [bool]) {
 
 /// `node` with every node `built` maps (by address) replaced by its built
 /// node, rebuilt above them; `rebuilt` memoizes by address.
-fn substituted(
+pub(crate) fn substituted(
     node: &Rc<OperatorNode>,
     built: &HashMap<*const OperatorNode, Rc<OperatorNode>>,
     rebuilt: &mut HashMap<*const OperatorNode, Rc<OperatorNode>>,
diff --git a/crates/plan-selection/src/search/memo/mod.rs b/crates/plan-selection/src/search/memo/mod.rs
index 83fff6af2..5f3e71bfc 100644
--- a/crates/plan-selection/src/search/memo/mod.rs
+++ b/crates/plan-selection/src/search/memo/mod.rs
@@ -39,6 +39,9 @@
 //! materialization and priced by Stage 3. [`coefficients`](crate::search::milp::coefficients) turns Stage 3's per-node
 //! prices into the coefficients of [`milp`](crate::search::milp::model), the memo solved as a
 //! mixed-integer linear program ([`crate::SearchAlgorithm::MemoMilp`]).
+//! [`clusters`](crate::search::milp::clusters) splits the memo into independent clusters (queries no
+//! plan's cost links), each solved with its own model and lazy loop, their
+//! plans combined into the workload's.
 
 #![allow(dead_code)]
 
diff --git a/crates/plan-selection/src/search/memo/rules/scans.rs b/crates/plan-selection/src/search/memo/rules/scans.rs
index 38f99f75a..e99932f84 100644
--- a/crates/plan-selection/src/search/memo/rules/scans.rs
+++ b/crates/plan-selection/src/search/memo/rules/scans.rs
@@ -351,6 +351,58 @@ fn site_scan_readers(memo: &Memo, site: &CouplingSite) -> Vec<ScanReader> {
     out
 }
 
+/// The labels of each metric's scans `candidate` reads (its own node, its
+/// private nodes and its child groups' nodes), by the raw-horizon key of
+/// the metric's source ([`raw_horizon_key`]).
+pub(crate) fn scan_labels_read(
+    memo: &Memo,
+    candidate: &ImplementationCandidate,
+) -> Vec<(u64, BTreeSet<String>)> {
+    scan_reader(&candidate.node)
+        .into_iter()
+        .chain(reads(memo, candidate, scan_reader))
+        .filter_map(|reader| match reader.scan.non_asap() {
+            Some(NonASAPOp::Scan { source, .. }) => Some((raw_horizon_key(source), reader.labels)),
+            _ => None,
+        })
+        .collect()
+}
+
+/// The labels of each metric's scans among `candidate`'s own node and its
+/// private nodes (not its child groups' nodes), by the raw-horizon key of
+/// the metric's source.
+pub(crate) fn private_scan_labels(
+    memo: &Memo,
+    candidate: &ImplementationCandidate,
+) -> Vec<(u64, BTreeSet<String>)> {
+    let mut boundary = Boundary::of(memo, candidate);
+    let mut out = Vec::new();
+    let mut seen = HashSet::new();
+    let mut stack = vec![Rc::clone(&candidate.node)];
+    while let Some(node) = stack.pop() {
+        if !seen.insert(Rc::as_ptr(&node)) {
+            continue;
+        }
+        if let Some(reader) = scan_reader(&node) {
+            if let Some(NonASAPOp::Scan { source, .. }) = reader.scan.non_asap() {
+                out.push((raw_horizon_key(source), reader.labels));
+            }
+        }
+        for child in node.children() {
+            if boundary.group(child).is_none() {
+                stack.push(Rc::clone(child));
+            }
+        }
+    }
+    out
+}
+
+/// The raw-horizon key ([`raw_horizon_key`]) of a scan-labels site's
+/// metric.
+pub(crate) fn site_source_key(memo: &Memo, site: &CouplingSite) -> Option<u64> {
+    metric_of_site(memo, site).map(|source| raw_horizon_key(&source))
+}
+
 /// The metric of a scan-labels site.
 fn metric_of_site(memo: &Memo, site: &CouplingSite) -> Option<Source> {
     site_scan_readers(memo, site)
diff --git a/crates/plan-selection/src/search/memo/rules/windows.rs b/crates/plan-selection/src/search/memo/rules/windows.rs
index e1e38d818..2c30b65f7 100644
--- a/crates/plan-selection/src/search/memo/rules/windows.rs
+++ b/crates/plan-selection/src/search/memo/rules/windows.rs
@@ -23,6 +23,7 @@
 //! inputs where a rule reads exact values, and each rollup reading the
 //! build the joint rule would choose.
 
+use std::cell::RefCell;
 use std::collections::hash_map::DefaultHasher;
 use std::collections::{BTreeSet, HashMap, HashSet};
 use std::hash::{Hash, Hasher};
@@ -30,8 +31,10 @@ use std::rc::Rc;
 
 use asap_logical_asap_optimization::accuracy::estimators::size_for_all_readers;
 use asap_logical_asap_optimization::asap_aware_cse::{
-    pane_candidates, rollup_candidates, rollup_pass2_candidates, topk_reuse_candidates,
-    window_composition_pass2_candidates, with_pane_chain, with_resized_build, RollupCandidate,
+    may_change_rollup_candidate, pane_candidates, rollup_candidate_for, rollup_candidates,
+    rollup_pass2_candidates, topk_reuse_candidates, window_composition_pass2_candidates,
+    window_composition_shared_chain_candidates, with_pane_chain, with_resized_build,
+    RollupCandidate,
 };
 use asap_logical_asap_optimization::per_computation_candidates::PerComputationCandidateError;
 use asap_types::ir::flat::flatten;
@@ -464,10 +467,9 @@ pub(crate) fn add_rollup_windows(memo: &mut Memo, entries: &[QueryWorkloadEntry]
                     if panes.len() != 2 {
                         continue;
                     }
-                    let Ok(pass2_candidates) = window_composition_pass2_candidates(&roots, &panes)
-                    else {
-                        continue;
-                    };
+                    // Only the pair in panes, reading a shared chain.
+                    let pass2_candidates =
+                        window_composition_shared_chain_candidates(&roots, &panes);
                     for pass2_candidate in pass2_candidates {
                         if pass2_candidate.pane_candidates.len() != 2
                             || pass2_candidate.shared_chains.is_empty()
@@ -995,26 +997,55 @@ pub(crate) fn rollup_rules(memo: &Memo, present: &BTreeSet<(usize, usize)>) -> V
         Some(read) => Rc::clone(&read.coarse_candidate),
         None => Rc::clone(&candidate(key).node),
     };
-    // The joint rule's source of each coarse build over `nodes`, in group
-    // order.
-    let sources =
-        |nodes: &[(usize, Rc<OperatorNode>)]| -> Vec<(*const OperatorNode, Rc<OperatorNode>)> {
-            let mut nodes = nodes.to_vec();
-            nodes.sort_by_key(|(g, _)| *g);
-            let roots: Vec<QueryRoot> = nodes
-                .into_iter()
-                .map(|(_, node)| QueryRoot::Operator(node))
-                .collect();
-            rollup_candidates(&roots)
-                .into_iter()
-                .map(|r| (Rc::as_ptr(&r.coarse), r.finer))
-                .collect()
-        };
-    let source_of = |coarse: *const OperatorNode, nodes: &[(usize, Rc<OperatorNode>)]| {
-        sources(nodes)
+    // The joint rule's plan of `nodes`, in group order (one node per group),
+    // and its key: the same plan is asked about many times.
+    type PlanKey = Vec<(usize, *const OperatorNode)>;
+    let plan = |nodes: &[(usize, Rc<OperatorNode>)]| -> (PlanKey, Vec<QueryRoot>) {
+        let mut nodes = nodes.to_vec();
+        nodes.sort_by_key(|(g, _)| *g);
+        let key = nodes.iter().map(|(g, n)| (*g, Rc::as_ptr(n))).collect();
+        let roots = nodes
             .into_iter()
-            .find(|(c, _)| *c == coarse)
-            .map(|(_, finer)| finer)
+            .map(|(_, node)| QueryRoot::Operator(node))
+            .collect();
+        (key, roots)
+    };
+    // The joint rule's source of each coarse build over `nodes`.
+    type Sources = Rc<Vec<(Rc<OperatorNode>, Rc<OperatorNode>)>>;
+    let all_sources: RefCell<HashMap<PlanKey, Sources>> = RefCell::default();
+    let sources = |nodes: &[(usize, Rc<OperatorNode>)]| -> Sources {
+        let (key, roots) = plan(nodes);
+        Rc::clone(all_sources.borrow_mut().entry(key).or_insert_with(|| {
+            Rc::new(
+                rollup_candidates(&roots)
+                    .into_iter()
+                    .map(|r| (r.coarse, r.finer))
+                    .collect(),
+            )
+        }))
+    };
+    // The source of one coarse build over `nodes`, computed for that build
+    // alone.
+    type OneSourceKey = (*const OperatorNode, PlanKey);
+    let one_source: RefCell<HashMap<OneSourceKey, Option<Rc<OperatorNode>>>> = RefCell::default();
+    let source_of = |coarse: *const OperatorNode, nodes: &[(usize, Rc<OperatorNode>)]| {
+        let (key, roots) = plan(nodes);
+        one_source
+            .borrow_mut()
+            .entry((coarse, key))
+            .or_insert_with(|| rollup_candidate_for(&roots, coarse).map(|r| r.finer))
+            .clone()
+    };
+    // Whether adding candidate `key`'s plan node may change the source of
+    // coarse build `coarse` ([`may_change_rollup_candidate`]); when not, the
+    // source over a plan with it is the source over the plan without it.
+    type MayChangeKey = (*const OperatorNode, (usize, usize));
+    let may_change: RefCell<HashMap<MayChangeKey, bool>> = RefCell::default();
+    let may_change_source = |coarse: &Rc<OperatorNode>, key: (usize, usize)| {
+        *may_change
+            .borrow_mut()
+            .entry((Rc::as_ptr(coarse), key))
+            .or_insert_with(|| may_change_rollup_candidate(coarse, &plan_node(key)))
     };
     let mut rules = Vec::new();
     for &rollup in present {
@@ -1055,7 +1086,7 @@ pub(crate) fn rollup_rules(memo: &Memo, present: &BTreeSet<(usize, usize)>) -> V
                     .filter(move |(g, _)| *g != coarse_group && *g != finer_group)
             };
             // A third group's build that beats the finer build for `C`.
-            for other in others() {
+            for other in others().filter(|&o| may_change_source(&read.coarse_build, o)) {
                 let nodes = [r_node.clone(), f_node.clone(), (other.0, plan_node(other))];
                 if source_of(coarse, &nodes).is_none_or(|b| !Rc::ptr_eq(&b, &build)) {
                     rules.push(RollupRule {
@@ -1085,9 +1116,10 @@ pub(crate) fn rollup_rules(memo: &Memo, present: &BTreeSet<(usize, usize)>) -> V
                     Some(h) => plan_node(h),
                 };
                 let in_holder = reached(&holder_node);
-                for (other_coarse, finer_build) in sources(&nodes) {
+                for (other_coarse_build, finer_build) in sources(&nodes).iter() {
+                    let other_coarse = Rc::as_ptr(other_coarse_build);
                     if other_coarse == coarse
-                        || !Rc::ptr_eq(&finer_build, &build)
+                        || !Rc::ptr_eq(finer_build, &build)
                         || !in_holder.contains(&other_coarse)
                     {
                         continue;
@@ -1102,6 +1134,10 @@ pub(crate) fn rollup_rules(memo: &Memo, present: &BTreeSet<(usize, usize)>) -> V
                         .copied()
                         .filter(|(g, _)| !groups.contains(g))
                         .filter(|&z| {
+                            // `nodes`' source is `build`.
+                            if !may_change_source(other_coarse_build, z) {
+                                return false;
+                            }
                             let mut with = nodes.clone();
                             with.push((z.0, plan_node(z)));
                             source_of(other_coarse, &with).is_none_or(|b| !Rc::ptr_eq(&b, &build))
diff --git a/crates/plan-selection/src/search/milp/clusters.rs b/crates/plan-selection/src/search/milp/clusters.rs
new file mode 100644
index 000000000..a9b1584cb
--- /dev/null
+++ b/crates/plan-selection/src/search/milp/clusters.rs
@@ -0,0 +1,356 @@
+//! Independent clusters of the memo: queries whose plans share no cost
+//! under cost model 1, searched as separate workloads by [`crate::search::milp::model`]
+//! and combined into one plan.
+
+use asap_physical_asap_optimization::{optimize_physical_asap_assigned, PhysicalCandidate};
+use asap_types::ir::QueryRoot;
+use asap_types::physical::Materialization;
+use asap_types::workload::{DataWorkload, QueryWorkloadEntry};
+
+use crate::search::memo::{CouplingSite, GroupId, Memo};
+use crate::{CostModel, PlanSelector, PlanningModels, QuerySearchReport, Selection};
+
+/// The workload's queries (`roots`, in workload order) split into
+/// clusters no plan's cost links, each in workload order, ordered by their
+/// first query. Two queries are in one cluster when:
+///
+/// - they reach a common group;
+/// - a candidate of a group one reaches reads, requires or rolls up from a
+///   group the other reaches (children, exact children, required sites,
+///   a rollup's finer group);
+/// - candidates of groups each reaches carry one coupling site (shared
+///   build size, pane chain, scan labels, raw horizon);
+/// - they read one source while raw data is kept and priced on disk: the
+///   source's raw data is charged once, for its longest query-time scan.
+///
+/// Every query is in one cluster when the deployment calibrated operation
+/// costs over the whole workload: its prices are not a sum over clusters.
+pub(crate) fn split(
+    memo: &Memo,
+    roots: &[QueryRoot],
+    models: &PlanningModels<'_>,
+) -> Vec<Vec<usize>> {
+    let count = roots.len();
+    if !models
+        .plan_selection
+        .measured_operation_costs
+        .operators
+        .is_empty()
+    {
+        return vec![(0..count).collect()];
+    }
+    let mut parent: Vec<usize> = (0..count).collect();
+    fn find(parent: &mut [usize], i: usize) -> usize {
+        let mut root = i;
+        while parent[root] != root {
+            root = parent[root];
+        }
+        let mut at = i;
+        while parent[at] != root {
+            let next = parent[at];
+            parent[at] = root;
+            at = next;
+        }
+        root
+    }
+    fn union(parent: &mut [usize], a: usize, b: usize) {
+        let (ra, rb) = (find(parent, a), find(parent, b));
+        parent[ra.max(rb)] = ra.min(rb);
+    }
+    let first_query = |g: GroupId| memo.groups[g.0 as usize].queries.first().copied();
+    let mut site_query: std::collections::HashMap<CouplingSite, usize> =
+        std::collections::HashMap::new();
+    for (g, group) in memo.groups.iter().enumerate() {
+        for pair in group.queries.windows(2) {
+            union(&mut parent, pair[0], pair[1]);
+        }
+        let Some(own) = group.queries.first().copied() else {
+            continue;
+        };
+        for candidate in &memo.candidates[g] {
+            let linked = candidate
+                .implements
+                .iter()
+                .chain(&candidate.children)
+                .chain(&candidate.requires_exact_children)
+                .chain(candidate.requires_sites.iter().map(|(group, _)| group))
+                .chain(candidate.rollup_of.iter().map(|r| &r.finer_group));
+            for &other in linked {
+                if let Some(q) = first_query(other) {
+                    union(&mut parent, own, q);
+                }
+            }
+            for site in candidate
+                .coupling_sites
+                .iter()
+                .chain(candidate.requires_sites.iter().map(|(_, site)| site))
+            {
+                let q = *site_query.entry(*site).or_insert(own);
+                union(&mut parent, own, q);
+            }
+        }
+    }
+    // A query no group reaches (a scalar) shares nothing: it joins the
+    // first cluster rather than being a workload of its own.
+    let mut reached = vec![false; count];
+    for group in &memo.groups {
+        for &q in &group.queries {
+            reached[q] = true;
+        }
+    }
+    for (q, _) in reached.iter().enumerate().skip(1).filter(|(_, r)| !**r) {
+        union(&mut parent, 0, q);
+    }
+    if count > 0 && !reached[0] {
+        for q in 1..count {
+            union(&mut parent, 0, q);
+        }
+    }
+    let disk_priced = match &models.plan_selection.cost_model {
+        CostModel::UsageOverTime(prices) => prices.cost_per_disk_byte_second != 0.0,
+        _ => true,
+    };
+    if disk_priced && models.deployment_capabilities.raw_data_retained {
+        for queries in crate::search::legacy::components::source_components(roots) {
+            for pair in queries.windows(2) {
+                union(&mut parent, pair[0], pair[1]);
+            }
+        }
+    }
+    let mut clusters: Vec<Vec<usize>> = Vec::new();
+    let mut cluster_of_root: std::collections::HashMap<usize, usize> =
+        std::collections::HashMap::new();
+    for q in 0..count {
+        let root = find(&mut parent, q);
+        let c = *cluster_of_root.entry(root).or_insert_with(|| {
+            clusters.push(Vec::new());
+            clusters.len() - 1
+        });
+        clusters[c].push(q);
+    }
+    clusters
+}
+
+/// The workload's plan from its `clusters`' plans (`found[i]` the
+/// selection, physical candidate, logical roots and report of
+/// `clusters[i]`): the clusters' roots in workload order, each node
+/// materialized as its cluster's plan materializes it, checked and costed
+/// by Stage 3 over the whole workload. `None` when that cost is not the
+/// sum of the clusters' costs.
+pub(crate) fn combine(
+    found: Vec<(
+        Selection,
+        PhysicalCandidate,
+        Vec<QueryRoot>,
+        crate::MilpClusterReport,
+    )>,
+    clusters: &[Vec<usize>],
+    entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
+    models: PlanningModels<'_>,
+) -> Option<(Selection, PhysicalCandidate, Vec<QueryRoot>)> {
+    let count = entries.len();
+    let mut cluster_of = vec![0; count];
+    let mut combined: Vec<Option<QueryRoot>> = vec![None; count];
+    for (c, (queries, (_, _, roots, _))) in clusters.iter().zip(&found).enumerate() {
+        for (&q, root) in queries.iter().zip(roots) {
+            combined[q] = Some(root.clone());
+            cluster_of[q] = c;
+        }
+    }
+    let combined: Vec<QueryRoot> = combined.into_iter().collect::<Option<_>>()?;
+    let mut assignments: Vec<std::vec::IntoIter<Materialization>> = found
+        .iter()
+        .map(|(_, selected, _, _)| selected.assignment.clone().into_iter())
+        .collect();
+    let mut choose = |readers: &[usize]| assignments[cluster_of[*readers.first()?]].next();
+    let mut candidate = optimize_physical_asap_assigned(
+        "L1",
+        &combined,
+        entries,
+        data,
+        models.deployment_capabilities,
+        &mut choose,
+    )
+    .ok()??;
+    if assignments.iter_mut().any(|rest| rest.next().is_some()) {
+        return None;
+    }
+    candidate.id = found
+        .iter()
+        .map(|(s, _, _, _)| s.selected.as_str())
+        .collect::<Vec<_>>()
+        .join("+");
+    // Cost models 2 and 3 do not run the MILP yet: their peak cores and
+    // cluster capacity are not sums over clusters, and would be checked
+    // here on the combined plan.
+    let mut selector = PlanSelector::new(entries, data, models, false).ok()?;
+    selector.offer(&candidate).ok()?;
+    let (mut selection, selected) = selector.finish().ok()?;
+    let price = selection.costs[&selection.selected].cost_per_second?;
+    let costs: Vec<f64> = found
+        .iter()
+        .map(|(_, _, _, report)| report.incumbent_cost)
+        .collect::<Option<_>>()?;
+    let total: f64 = costs.iter().sum();
+    if (price - total).abs() > 1e-9 * total.abs() {
+        return None;
+    }
+    let mut queries = vec![QuerySearchReport::default(); count];
+    let mut rejected = Vec::new();
+    for (cluster, (s, _, _, _)) in clusters.iter().zip(&found) {
+        for (&q, report) in cluster.iter().zip(&s.search.queries) {
+            queries[q] = report.clone();
+        }
+        rejected.extend(
+            s.rejected
+                .iter()
+                .cloned()
+                .map(|r| crate::search::legacy::components::renumber(r, cluster)),
+        );
+    }
+    selection.rejected = rejected;
+    selection.search.queries = queries;
+    selection.search.candidates_checked += found
+        .iter()
+        .map(|(s, _, _, _)| s.search.candidates_checked)
+        .sum::<u64>();
+    selection.search.raw_scan_assignments_skipped = found
+        .iter()
+        .map(|(s, _, _, _)| s.search.raw_scan_assignments_skipped)
+        .sum();
+    selection.search.incumbent_cost = Some(price);
+    selection.search.milp = Some(milp_report(&found));
+    let roots = combined;
+    Some((selection, selected, roots))
+}
+
+/// The workload's MILP report from its clusters': bounds, rounds, model
+/// sizes and solver time summed, the gap from the summed bound and cost.
+fn milp_report(
+    found: &[(
+        Selection,
+        PhysicalCandidate,
+        Vec<QueryRoot>,
+        crate::MilpClusterReport,
+    )],
+) -> crate::MilpReport {
+    let reports: Vec<&crate::MilpReport> = found
+        .iter()
+        .filter_map(|(s, _, _, _)| s.search.milp.as_ref())
+        .collect();
+    let sum = |value: fn(&crate::MilpReport) -> Option<f64>| -> Option<f64> {
+        reports.iter().map(|r| value(r)).sum()
+    };
+    let lower_bound = sum(|r| r.lower_bound);
+    let incumbent_cost: Option<f64> = found.iter().map(|(_, _, _, c)| c.incumbent_cost).sum();
+    crate::MilpReport {
+        lower_bound,
+        first_lower_bound: sum(|r| r.first_lower_bound),
+        gap: crate::search::milp::model::gap(incumbent_cost, lower_bound),
+        rounds: reports.iter().map(|r| r.rounds).sum(),
+        variable_count: reports.iter().map(|r| r.variable_count).sum(),
+        constraint_count: reports.iter().map(|r| r.constraint_count).sum(),
+        solver_time_ms: reports.iter().map(|r| r.solver_time_ms).sum(),
+        time_limit_reached: reports.iter().any(|r| r.time_limit_reached),
+        decode_failures: reports.iter().map(|r| r.decode_failures).sum(),
+        clusters: found.iter().map(|(_, _, _, c)| c.clone()).collect(),
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::search::memo::{groups, CouplingKind, ImplementationCandidate};
+    use asap_logical_asap_optimization::per_computation_candidates::generate_per_computation_candidates;
+    use asap_logical_asap_optimization::test_support::lower_promql;
+    use asap_types::ir::properties::AccuracyTarget;
+
+    /// The memo of `queries` (Pass 1 groups only) and their roots.
+    fn memo(queries: &[&str]) -> (Memo, Vec<QueryRoot>) {
+        let roots: Vec<QueryRoot> = queries
+            .iter()
+            .map(|q| QueryRoot::Operator(lower_promql(q, AccuracyTarget::max_error(0.05))))
+            .collect();
+        let per_computation =
+            generate_per_computation_candidates(roots.iter().cloned().enumerate().collect())
+                .unwrap();
+        (groups::memo_of(&per_computation), roots)
+    }
+
+    /// The first candidate of the group whose root is query `q`.
+    fn root_candidate(memo: &mut Memo, q: usize) -> &mut ImplementationCandidate {
+        let g = memo
+            .groups
+            .iter()
+            .position(|g| g.roots_of.contains(&q))
+            .unwrap();
+        &mut memo.candidates[g][0]
+    }
+
+    /// Queries over two gauges are two clusters, those over one gauge
+    /// reaching its scan group; one coupling site (a shared build's size,
+    /// a source's raw horizon) carried by candidates of both joins them.
+    #[test]
+    fn a_shared_group_or_coupling_site_joins_clusters() {
+        let models = PlanningModels::builtin();
+        let (memo, roots) = memo(&[
+            "sum_over_time(m[5m])",
+            "max_over_time(n[5m])",
+            "quantile_over_time(0.9, m[1m])",
+        ]);
+        assert_eq!(split(&memo, &roots, &models), [vec![0, 2], vec![1]]);
+        for kind in [CouplingKind::SharedBuildSize, CouplingKind::RawHorizon] {
+            let mut linked = memo.clone();
+            let site = CouplingSite { key: 7, kind };
+            root_candidate(&mut linked, 1).coupling_sites.push(site);
+            assert_eq!(split(&linked, &roots, &models), [vec![0, 2], vec![1]]);
+            root_candidate(&mut linked, 2).coupling_sites.push(site);
+            assert_eq!(split(&linked, &roots, &models), [vec![0, 1, 2]], "{kind:?}");
+        }
+    }
+
+    /// Queries over one source through different scans share no group;
+    /// with raw data kept and disk priced, the source's raw data is charged
+    /// once for its longest scan, which joins them. Without a disk price,
+    /// or without raw data, they stay apart.
+    #[test]
+    fn a_source_s_charged_raw_data_joins_its_queries() {
+        let (memo, roots) = memo(&[
+            r#"sum_over_time(m{a="1"}[5m])"#,
+            r#"max_over_time(m{a="2"}[10m])"#,
+        ]);
+        let unpriced = PlanningModels::builtin();
+        assert_eq!(split(&memo, &roots, &unpriced), [vec![0], vec![1]]);
+        let mut settings = crate::PlanSelectionSettings::default();
+        let CostModel::UsageOverTime(prices) = &mut settings.cost_model else {
+            unreachable!("cost model 1 by default")
+        };
+        prices.cost_per_disk_byte_second = 1e-12;
+        let priced = PlanningModels::builtin().with_plan_selection(&settings);
+        assert_eq!(split(&memo, &roots, &priced), [vec![0, 1]]);
+        let no_raw = asap_logical_asap_optimization::DeploymentCapabilities {
+            raw_data_retained: false,
+            ..Default::default()
+        };
+        let priced_no_raw = priced.with_deployment_capabilities(&no_raw);
+        assert_eq!(split(&memo, &roots, &priced_no_raw), [vec![0], vec![1]]);
+    }
+
+    /// A candidate reading another query's group (a rollup or a top-k
+    /// reuse) joins the two queries.
+    #[test]
+    fn a_candidate_reading_another_query_s_group_joins_them() {
+        let models = PlanningModels::builtin();
+        let (mut memo, roots) = memo(&["sum_over_time(m[5m])", "max_over_time(n[5m])"]);
+        let other = memo
+            .groups
+            .iter()
+            .position(|g| g.roots_of.contains(&1))
+            .unwrap();
+        root_candidate(&mut memo, 0)
+            .children
+            .push(GroupId(other as u32));
+        assert_eq!(split(&memo, &roots, &models), [vec![0, 1]]);
+    }
+}
diff --git a/crates/plan-selection/src/search/milp/coefficients.rs b/crates/plan-selection/src/search/milp/coefficients.rs
index f2d18f208..849f65d37 100644
--- a/crates/plan-selection/src/search/milp/coefficients.rs
+++ b/crates/plan-selection/src/search/milp/coefficients.rs
@@ -133,6 +133,9 @@ pub(crate) struct ChoiceParts {
     pub build_keys: Vec<u32>,
     /// Per part (direct, then each slot), its nodes' keys.
     pub part_keys: Vec<HashSet<u32>>,
+    /// The candidate at the choice's site thresholds, which solution
+    /// decoding composes.
+    pub version: ImplementationCandidate,
 }
 
 pub(crate) fn choice_parts(
@@ -200,6 +203,32 @@ pub(crate) fn choice_parts(
             ) else {
                 continue;
             };
+            // Built for one query twice, a build may be built once only
+            // because of the second copy: the keys of the builds Stage 2
+            // also builds once for the one query alone.
+            let built_once_alone: Option<HashSet<u32>> = match built_for.as_slice() {
+                [a, b] if a == b => {
+                    let alone = optimize_physical_asap(
+                        "S",
+                        &roots[..1],
+                        &part_entries[..1],
+                        data,
+                        models.deployment_capabilities,
+                    )
+                    .unwrap_or_default();
+                    let mut keys = HashSet::new();
+                    for physical in &alone {
+                        let all = keys_of(&physical.dag, cache.keys_mut());
+                        for (i, node) in physical.dag.nodes.iter().enumerate() {
+                            if matches!(node.materialization, Materialization::MatAtQuery(_)) {
+                                keys.insert(all[i]);
+                            }
+                        }
+                    }
+                    Some(keys)
+                }
+                _ => None,
+            };
             let child_groups = children_of(memo, group, candidate);
             let child_key_sets: Vec<Rc<crate::cost::per_node::ChildKeys>> = child_groups
                 .iter()
@@ -323,6 +352,29 @@ pub(crate) fn choice_parts(
                             .iter()
                             .map(|&i| physical.dag.nodes[i].materialization)
                             .collect();
+                        let built_once_keys: Vec<u32> = builds
+                            .iter()
+                            .filter(|&&i| {
+                                matches!(
+                                    physical.dag.nodes[i].materialization,
+                                    Materialization::MatAtQuery(_)
+                                )
+                            })
+                            .map(|&i| widened[i])
+                            .collect();
+                        let built_once_with_other_group: Vec<u32> = match &built_once_alone {
+                            Some(alone) => builds
+                                .iter()
+                                .filter(|&&i| {
+                                    matches!(
+                                        physical.dag.nodes[i].materialization,
+                                        Materialization::MatAtQuery(_)
+                                    ) && !alone.contains(&keys[i])
+                                })
+                                .map(|&i| widened[i])
+                                .collect(),
+                            None => Vec::new(),
+                        };
                         if seen.contains(&assignment) {
                             continue;
                         }
@@ -375,6 +427,8 @@ pub(crate) fn choice_parts(
                                     _ => None,
                                 },
                                 scan_project_cost_per_unit_rate: parts.scan_project_cost,
+                                built_once_with_other_group,
+                                built_once_keys,
                             },
                             lazy: LazyTerms {
                                 thresholds: thresholds.clone(),
@@ -386,6 +440,7 @@ pub(crate) fn choice_parts(
                             slot_keys: parts.slot_keys,
                             build_keys: builds.iter().map(|&i| shallow[i]).collect(),
                             part_keys: parts.part_keys,
+                            version: candidate.clone(),
                         });
                     }
                 }
@@ -400,7 +455,7 @@ pub(crate) fn choice_parts(
 /// candidates, with the index of each threshold. A raw-horizon site has
 /// no threshold here: scan sharing is part of each choice
 /// ([`Choice::scan_children`]). A candidate without sites is itself.
-fn at_site_thresholds(
+pub(super) fn at_site_thresholds(
     memo: &Memo,
     candidate: &ImplementationCandidate,
 ) -> Vec<(ImplementationCandidate, Vec<(CouplingSite, usize)>)> {
@@ -903,6 +958,22 @@ fn coefficients_of(
         );
     }
     let mut part_keys = vec![direct.iter().map(|id| keys[position[id]]).collect()];
+    // Each query-time node a reader computes, keyed as one node in the
+    // plan: a scan by its source (scans of one source may be one widened
+    // scan), any other node by its key up to scan widening.
+    let node_work_key = |id: PhysicalASAPNodeId| -> u64 {
+        match &dag.nodes[position[&id]].payload {
+            Operator::NonASAP(NonASAPOp::Scan { source, .. }) => {
+                (1 << 63) | crate::search::memo::rules::scans::raw_horizon_key(source)
+            }
+            _ => u64::from(widened[position[&id]]),
+        }
+    };
+    let mut query_work: Vec<(u64, f64)> = direct
+        .iter()
+        .filter(|&&id| query_time(id))
+        .map(|&id| (node_work_key(id), priced[&id].work_ms))
+        .collect();
     let mut slots = Vec::new();
     let mut slot_keys = Vec::new();
     let mut slot_reads = Vec::new();
@@ -914,6 +985,13 @@ fn coefficients_of(
         let via = Via::Slot(slots.len());
         lazy_of(&part, via);
         let at_ingestion = matches!(node.materialization, Materialization::MatAtIngestion(_));
+        if !at_ingestion {
+            query_work.extend(
+                part.iter()
+                    .filter(|&&id| query_time(id))
+                    .map(|&id| (node_work_key(id), priced[&id].work_ms)),
+            );
+        }
         // A query-time slot's raw scans are charged as scans.
         let scans: Vec<(u64, f64, f64)> = match at_ingestion {
             true => Vec::new(),
@@ -975,28 +1053,36 @@ fn coefficients_of(
     for node in &mut at_ingestion.nodes {
         node.output_state.timing = ExecutionTiming::IngestionTime;
     }
-    let direct_at_ingestion = per_computation(
+    // `None` only when Stage 3 cannot place it at ingestion time (a
+    // grouping with no window below it); any other pricing error leaves
+    // its price unknown, charged 0 (a lower bound).
+    let direct_at_ingestion = match per_computation(
         &at_ingestion,
         entries,
         scan_lookback,
         horizons,
         data,
         models,
-    )
-    .ok()
-    .map(|ingestion_priced| PartCoefficients {
-        cost_per_unit_rate: direct
-            .iter()
-            .map(|id| {
-                let p = &ingestion_priced[id];
-                p.cpu_ops_per_second * prices.cost_per_cpu_op
-                    + p.scanned_bytes_per_second * prices.cost_per_scanned_byte
-                    + p.average_memory_bytes * prices.cost_per_memory_byte_second
-                    + p.disk_bytes as f64 * prices.cost_per_disk_byte_second
-            })
-            .sum(),
-        work_ms: 0.0,
-    });
+    ) {
+        Ok(ingestion_priced) => Some(PartCoefficients {
+            cost_per_unit_rate: direct
+                .iter()
+                .map(|id| {
+                    let p = &ingestion_priced[id];
+                    p.cpu_ops_per_second * prices.cost_per_cpu_op
+                        + p.scanned_bytes_per_second * prices.cost_per_scanned_byte
+                        + p.average_memory_bytes * prices.cost_per_memory_byte_second
+                        + p.disk_bytes as f64 * prices.cost_per_disk_byte_second
+                })
+                .sum(),
+            work_ms: 0.0,
+        }),
+        Err(SelectionError::Cost {
+            error: crate::cost::analytical::AnalyticalCostError::IngestionGroupingWithoutWindow,
+            ..
+        }) => None,
+        Err(_) => Some(PartCoefficients::default()),
+    };
     // The group's own node built once at query time, in memory: the state
     // it holds instead of its working memory.
     let mut built_once = dag.nodes[position[&root]].clone();
@@ -1032,6 +1118,7 @@ fn coefficients_of(
                     )
                 })
                 .collect(),
+            query_work,
         },
         lazy: LazyTerms {
             held,
@@ -1848,8 +1935,8 @@ mod tests {
         intervals_ms: &[u32],
         targets: &[AccuracyTarget],
     ) {
-        assert_model_is_a_lower_bound_in(queries, intervals_ms, targets, true);
-        assert_model_is_a_lower_bound_in(queries, intervals_ms, targets, false);
+        assert_model_is_a_lower_bound_in(queries, intervals_ms, targets, true, &data(), 0);
+        assert_model_is_a_lower_bound_in(queries, intervals_ms, targets, false, &data(), 0);
     }
 
     /// [`assert_model_is_a_lower_bound`] in a deployment that keeps raw
@@ -1859,7 +1946,9 @@ mod tests {
         intervals_ms: &[u32],
         targets: &[AccuracyTarget],
         raw_data_retained: bool,
-    ) {
+        data: &DataWorkload,
+        term_cuts_from_rounds: usize,
+    ) -> BTreeSet<&'static str> {
         let roots: Vec<(usize, QueryRoot)> = queries
             .iter()
             .zip(targets)
@@ -1871,7 +1960,6 @@ mod tests {
             .collect();
         let per_computation = generate_per_computation_candidates(roots).unwrap();
         let entries = entries(intervals_ms, targets);
-        let data = data();
         let deployment =
             asap_logical_asap_optimization::deployment_capabilities::DeploymentCapabilities {
                 raw_data_retained,
@@ -1880,14 +1968,20 @@ mod tests {
             };
         let models = PlanningModels::builtin().with_deployment_capabilities(&deployment);
         let memo = memo(&per_computation, &entries);
-        let input =
-            crate::search::milp::model::model_input_of(&memo, &entries, &data, models).unwrap();
+        let input = crate::search::milp::model::model_input_of(&memo, &entries, data, models)
+            .unwrap()
+            .0;
+        // With the term cuts the lazy loop may add from its first rounds.
+        let input = match term_cuts_from_rounds {
+            0 => input,
+            rounds => crate::search::milp::model::with_every_term_cut(&input, rounds),
+        };
         let checked = crate::search::milp::model::lower_bound_check(
             &memo,
             &per_computation,
             &input,
             &entries,
-            &data,
+            data,
             models,
         )
         .unwrap();
@@ -1900,6 +1994,115 @@ mod tests {
                 "{queries:?} {group_choice:?}: model {value} over exact {exact}"
             );
         }
+        // Each cut reads model variables that exist: a query's reach or
+        // held charge of a choice only of a group the query reads.
+        let reads = |q: usize, c: usize| {
+            input.group_queries[input.choices[c].group.0 as usize].contains(&q)
+        };
+        for cut in &input.cuts.cuts {
+            assert!(
+                cut.reached.iter().all(|&(q, c)| reads(q, c)),
+                "{queries:?}: {cut:?}"
+            );
+            if let crate::search::milp::cuts::Term::HeldPerReader(q, c) = cut.term {
+                assert!(reads(q, c), "{queries:?}: {cut:?}");
+            }
+        }
+        input.cuts.cuts.iter().map(|cut| cut.term.kind()).collect()
+    }
+
+    /// An exact top-k of 5-minute sums every 5 minutes over 1,000 series
+    /// sampled every 15 s, without raw data: only a summary maintained at
+    /// ingestion time is valid. Its scan and window are read only through
+    /// that ingestion-time slot, so no query-time work or held memory of
+    /// theirs may be charged (a random three-way workload's cluster).
+    #[test]
+    fn the_model_is_a_lower_bound_below_an_ingestion_time_slot() {
+        let data = DataWorkload {
+            arrival: DataArrival::Ingesting {
+                rows_per_second: 1_000.0 / 15.0,
+            },
+            sources: [(
+                "m".to_string(),
+                SourceFacts {
+                    series_count: Some(1_000),
+                    sample_interval: Some(asap_types::workload::DurationMs(15_000)),
+                    metric_type: Some(asap_types::workload::MetricType::Gauge),
+                    ..SourceFacts::default()
+                },
+            )]
+            .into_iter()
+            .collect(),
+        };
+        assert_model_is_a_lower_bound_in(
+            &["topk by (job) (3, sum_over_time(m[5m]))"],
+            &[300_000],
+            &[AccuracyTarget::exact()],
+            false,
+            &data,
+            0,
+        );
+    }
+
+    /// Two queries whose children are read by two ingestion-time slots,
+    /// without raw data: Stage 3 runs each shared child as one ingestion
+    /// copy, so the model charges its ingestion-time part once, not once
+    /// per slot reading it. Every 10 s, the two quantiles' query-time
+    /// node over their shared build is one node read by both queries of
+    /// one class: its held memory is charged at their union rate. Every
+    /// 10 s and every 20 s, its readers are in two classes no group's
+    /// class set joins: their union rate is at least the faster's.
+    #[test]
+    fn the_model_is_a_lower_bound_for_a_child_of_two_ingestion_time_slots() {
+        let within = AccuracyTarget::max_error(0.01);
+        for intervals in [
+            [300_000, 300_000],
+            [60_000, 300_000],
+            [10_000, 10_000],
+            [10_000, 20_000],
+        ] {
+            for window in ["5m", "1m"] {
+                assert_model_is_a_lower_bound_in(
+                    &[
+                        &format!("quantile_over_time(0.9, m[{window}])"),
+                        &format!("quantile_over_time(0.5, m[{window}])"),
+                    ],
+                    &intervals,
+                    &[within.clone(), within.clone()],
+                    false,
+                    &data(),
+                    0,
+                );
+            }
+            assert_model_is_a_lower_bound_in(
+                &["quantile_over_time(0.9, m[5m])", "sum_over_time(m[5m])"],
+                &intervals,
+                &[within.clone(), AccuracyTarget::exact()],
+                false,
+                &data(),
+                0,
+            );
+            assert_model_is_a_lower_bound_in(
+                &["quantile_over_time(0.9, m[5m])", "max_over_time(m[1m])"],
+                &intervals,
+                &[within.clone(), AccuracyTarget::exact()],
+                false,
+                &data(),
+                0,
+            );
+            let exact = AccuracyTarget::exact();
+            assert_model_is_a_lower_bound_in(
+                &[
+                    "max by (job) (sum_over_time(m[5m]))",
+                    "sum_over_time(m[5m])",
+                ],
+                &intervals,
+                &[exact.clone(), exact.clone()],
+                false,
+                &data(),
+                0,
+            );
+        }
     }
 
     /// Three quantiles of one metric over a 1-minute, a 5-minute and an
@@ -1947,6 +2150,83 @@ mod tests {
         );
     }
 
+    /// #839's shape with the term cuts of the lazy loop's first 30 rounds
+    /// (a source's scans apart or widened, held memory, ingestion copies):
+    /// each cut stays a lower bound of every plan it applies to.
+    #[test]
+    fn the_model_with_term_cuts_is_a_lower_bound() {
+        let target = AccuracyTarget::max_error(0.05).failure_probability(0.3);
+        let mut kinds = BTreeSet::new();
+        for raw_data_retained in [true, false] {
+            kinds.extend(assert_model_is_a_lower_bound_in(
+                &[
+                    "sum by (job) (sum_over_time(m[10m]))",
+                    "sum(sum_over_time(m[10m]))",
+                    "avg by (job) (quantile_over_time(0.99, m[100s]))",
+                    "avg(quantile_over_time(0.95, m[100s]))",
+                ],
+                &[20_000, 20_000, 1_000, 1_000],
+                &[
+                    target.clone(),
+                    target.clone(),
+                    target.clone(),
+                    target.clone(),
+                ],
+                raw_data_retained,
+                &data(),
+                30,
+            ));
+        }
+        // A fused choice: the top-k's heap implements the sum, the other
+        // query's root, so that query reaches a choice of a group it does
+        // not read.
+        let exact = AccuracyTarget::exact();
+        for raw_data_retained in [true, false] {
+            kinds.extend(assert_model_is_a_lower_bound_in(
+                &[
+                    "topk(5, sum by (job) (sum_over_time(m[5m])))",
+                    "sum by (job) (sum_over_time(m[5m]))",
+                ],
+                &[60_000, 60_000],
+                &[exact.clone(), exact.clone()],
+                raw_data_retained,
+                &data(),
+                30,
+            ));
+        }
+        // Ingestion copies (a pane build at ingestion beside a sum reading
+        // the scan built once) and held memory over the longest reader.
+        let within = AccuracyTarget::max_error(0.01);
+        kinds.extend(assert_model_is_a_lower_bound_in(
+            &["quantile_over_time(0.9, m[5m])", "sum_over_time(m[5m])"],
+            &[60_000, 300_000],
+            &[within, exact.clone()],
+            true,
+            &data(),
+            30,
+        ));
+        kinds.extend(assert_model_is_a_lower_bound_in(
+            &[
+                "max_over_time(m[5m]) - min_over_time(m[5m])",
+                "max_over_time(m[5m])",
+            ],
+            &[60_000, 60_000],
+            &[exact.clone(), exact],
+            true,
+            &data(),
+            30,
+        ));
+        for kind in [
+            "source scans",
+            "held per reader",
+            "held longest",
+            "shared build held",
+            "ingestion copy",
+        ] {
+            assert!(kinds.contains(kind), "no {kind} cut in {kinds:?}");
+        }
+    }
+
     /// One query reading two groups over one metric, whose scans may be one
     /// node: its latency and held memory stay lower bounds.
     #[test]
diff --git a/crates/plan-selection/src/search/milp/cuts.rs b/crates/plan-selection/src/search/milp/cuts.rs
new file mode 100644
index 000000000..6501e3081
--- /dev/null
+++ b/crates/plan-selection/src/search/milp/cuts.rs
@@ -0,0 +1,662 @@
+//! Term cuts (joint-planning-search design, "How we solve it", stronger
+//! cuts): after a round, each term of the solution's plan that the model
+//! underestimates gets a cut conditioned only on the choices that term
+//! depends on, so it corrects every plan sharing those choices, not just
+//! the one evaluated. Each cut's bound is a lower bound of the term in
+//! every plan whose solution makes those choices; the no-good cut of the
+//! solution stays as the fallback that ends the loop.
+//!
+//! The terms:
+//!
+//! - **A source's query-time scans.** The model charges the scans of one
+//!   source as one widened scan may cost (the dearest, plus a pass per
+//!   other scan). A plan reads them either apart, each at its own price,
+//!   or widened to one label set, each widened scan at least as dear as a
+//!   reader's scan widened to labels its plan's scans cover (Stage 3's
+//!   price of the [`crate::search::memo::CandidateOrigin::WidenedScans`]
+//!   twin at that label set). The cut charges the smaller of the two.
+//! - **A query's held memory** (latency sum × memory sum): the choices a
+//!   query reaches at query time fix its latency from below (each node key
+//!   once, at the least work any choice the plan may have gives it) and
+//!   the memory its parts hold; each held charge of those parts is at least
+//!   its bytes times that latency. Builds several groups share (one node,
+//!   left out of the model's held charges) get their held memory here.
+//! - **The ingestion-time copy of a group built once**: a group built once
+//!   per evaluation for its query-time readers that an ingestion-time slot
+//!   also reads is computed twice, once at each timing (Stage 2 builds one
+//!   copy per timing); the model charges only the query-time one.
+
+use std::collections::{BTreeMap, BTreeSet};
+
+use asap_types::physical::Materialization;
+
+use crate::search::milp::model::ModelInput;
+
+/// What the term cuts read about each choice of the model.
+#[derive(Debug, Clone, Default, PartialEq)]
+pub(crate) struct TermCuts {
+    /// Per choice, the labels of each source's scans it reads (its own
+    /// node, its private nodes and its child groups' nodes), by the
+    /// source's raw-horizon key.
+    pub scan_labels: Vec<Vec<(u64, BTreeSet<String>)>>,
+    /// Per choice, the labels of each source's scans it runs at query
+    /// time whenever it runs at query time: its own node's and its private
+    /// nodes' (child groups' scans are their own choices'), none for a
+    /// choice with an ingestion-time slot (whose scans below may run at
+    /// ingestion time only).
+    pub query_time_scan_labels: Vec<Vec<(u64, BTreeSet<String>)>>,
+    /// Per choice, its scans widened: (source key, label set, the price per
+    /// unit of rate of its query-time scans of that source widened to the
+    /// label set), from its candidate's widened-scans twins.
+    pub widened_scans: Vec<Vec<WidenedScan>>,
+    /// Per choice, whether every consumer that may compute its direct part
+    /// does so at query time (no ingestion-time slot above it): its
+    /// query-time scans then run whenever it is chosen.
+    pub always_at_query_time: Vec<bool>,
+    /// Per choice, the slowest rate any consumer computes its direct part
+    /// at: the least rate of its group's queries (a node built once runs at
+    /// the union rate of the queries reaching it, at least each one's).
+    pub rate_floor: Vec<f64>,
+    /// Per choice, the held memory of its query-time builds several
+    /// groups' choices share (one node, left out of
+    /// [`ModelInput::held`]): (the build's shared key, price per ms of
+    /// latency per unit of rate, whether each reader holds its own copy).
+    pub shared_slot_held: Vec<Vec<SharedBuildHeld>>,
+    /// Per choice, its direct part's price per second computed at
+    /// ingestion time (under an ingestion-time slot; 0 when it cannot run
+    /// there).
+    pub at_ingestion_cost: Vec<f64>,
+    /// The cuts so far.
+    pub cuts: Vec<TermCut>,
+}
+
+/// A choice's query-time scans of one source widened to one label set.
+#[derive(Debug, Clone, PartialEq)]
+pub(crate) struct WidenedScan {
+    /// The source's raw-horizon key.
+    pub source: u64,
+    pub labels: BTreeSet<String>,
+    pub cost_per_unit_rate: f64,
+}
+
+/// The memory a query-time build several groups' choices share holds
+/// while an evaluation runs.
+#[derive(Debug, Clone, Copy, PartialEq)]
+pub(crate) struct SharedBuildHeld {
+    /// The build's shared key ([`crate::search::milp::model::SlotCoefficients::shared_site`]).
+    pub key: (u64, usize),
+    /// Price per ms of latency per unit of rate.
+    pub cost_per_ms: f64,
+    /// Whether each reader holds its own copy.
+    pub per_reader: bool,
+}
+
+/// What makes a choice run at query time in a plan: other choices, and
+/// the model's query-time reach (query, choice).
+type Witness = (Vec<usize>, Vec<(usize, usize)>);
+
+/// One term cut: when every choice of `choices` is chosen and each query
+/// of `reached` reaches its choice at query time, the model's charge for
+/// `term` is at least `bound`.
+#[derive(Debug, Clone, PartialEq)]
+pub(crate) struct TermCut {
+    pub term: Term,
+    pub choices: Vec<usize>,
+    /// (query, choice): the model's query-time reach of the choice.
+    pub reached: Vec<(usize, usize)>,
+    /// Groups built once per evaluation.
+    pub built_once: Vec<usize>,
+    /// Nodes several groups' choices have, by key, whose charges the term
+    /// includes ([`crate::search::milp::model::ChoiceCoefficients::cross_shared`]).
+    pub nodes: Vec<u64>,
+    pub bound: f64,
+}
+
+/// A term of the plan's price the model charges through its own variables.
+#[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord)]
+pub(crate) enum Term {
+    /// The query-time scans of the source with this raw-horizon key and
+    /// the cut's nodes several groups' choices have: the model's
+    /// per-source scan charge and its passes (a pass stands for a reader's
+    /// node over a widened scan) and those nodes' charges.
+    SourceScans(u64),
+    /// The held memory of query `.0`'s own copy of choice `.1`'s parts.
+    HeldPerReader(usize, usize),
+    /// The held memory of choice `.0`'s parts over its longest reader.
+    HeldLongest(usize),
+    /// The held memory of a query-time build several groups' choices
+    /// share, by its shared key (a variable of its own, charged only
+    /// through these cuts).
+    SharedBuildHeld((u64, usize)),
+    /// The ingestion-time copy of a group built once per evaluation that an
+    /// ingestion-time slot reads, by group (a variable of its own, charged
+    /// only through these cuts). Its cut also holds only when some chosen
+    /// choice's ingestion-time slot reads the group.
+    IngestionCopy(usize),
+}
+
+impl Term {
+    /// The term's kind, for reports and tests.
+    pub(crate) fn kind(&self) -> &'static str {
+        match self {
+            Term::SourceScans(_) => "source scans",
+            Term::HeldPerReader(..) => "held per reader",
+            Term::HeldLongest(_) => "held longest",
+            Term::SharedBuildHeld(_) => "shared build held",
+            Term::IngestionCopy(_) => "ingestion copy",
+        }
+    }
+}
+
+/// The per-choice facts [`TermCuts`] reads that the model's input already
+/// has: rate floors and whether a choice always runs at query time.
+pub(crate) fn fill_from_input(input: &ModelInput, cuts: &mut TermCuts) {
+    let ingestion_slot_above = groups_below_ingestion_slots(input);
+    cuts.always_at_query_time = input
+        .choices
+        .iter()
+        .map(|choice| !ingestion_slot_above.contains(&(choice.group.0 as usize)))
+        .collect();
+    cuts.rate_floor = input
+        .choices
+        .iter()
+        .map(|choice| {
+            input.group_queries[choice.group.0 as usize]
+                .iter()
+                .map(|&q| input.queries[q].rate)
+                .fold(f64::INFINITY, f64::min)
+        })
+        .map(|rate| if rate.is_finite() { rate } else { 0.0 })
+        .collect();
+}
+
+/// The groups an ingestion-time slot of a choice may compute: its child
+/// groups through that slot, and theirs read directly, transitively.
+fn groups_below_ingestion_slots(input: &ModelInput) -> BTreeSet<usize> {
+    let mut below: BTreeSet<usize> = BTreeSet::new();
+    let mut stack: Vec<usize> = Vec::new();
+    for choice in &input.choices {
+        for (child, via) in choice.children.iter().zip(&choice.coefficients.child_via) {
+            if let crate::search::milp::model::Via::Slot(i) = via {
+                if matches!(
+                    choice.coefficients.slots[*i].materialization,
+                    Materialization::MatAtIngestion(_)
+                ) {
+                    stack.push(child.0 as usize);
+                }
+            }
+        }
+    }
+    while let Some(g) = stack.pop() {
+        if !below.insert(g) {
+            continue;
+        }
+        for choice in input.choices.iter().filter(|c| c.group.0 as usize == g) {
+            stack.extend(choice.children.iter().map(|h| h.0 as usize));
+        }
+    }
+    below
+}
+
+/// The term cuts of the solution with `chosen` choices that are not yet in
+/// `input`'s cuts and raise the model's charge.
+pub(crate) fn term_cuts_of(input: &ModelInput, chosen: &[usize]) -> Vec<TermCut> {
+    let facts = &input.cuts;
+    // Each chosen choice that runs at query time in every plan with the
+    // cut's choices: one that always does, or one a query reaches at query
+    // time through choices the cut then includes (`witness`).
+    let reach: Vec<Vec<usize>> = (0..input.queries.len())
+        .map(|q| query_time_reach(input, chosen, q))
+        .collect();
+    // A witness: nothing, the model's query-time reach of the choice by a
+    // query (when the model has that variable: the choice has query-time
+    // work), or a query's whole query-time reach set.
+    let has_reach_variable = |c: usize| {
+        input.choices[c]
+            .coefficients
+            .query_work
+            .iter()
+            .any(|&(_, work)| work > 0.0)
+    };
+    let witness = |c: usize| -> Option<Witness> {
+        if facts.always_at_query_time[c] {
+            return Some((Vec::new(), Vec::new()));
+        }
+        let q = reach.iter().position(|r| r.contains(&c))?;
+        // The model has the query's reach of the choice only for a query
+        // of the choice's group (a fused choice of another group may be
+        // reached through a group it implements).
+        let reads = input.group_queries[input.choices[c].group.0 as usize].contains(&q);
+        Some(match has_reach_variable(c) && reads {
+            true => (Vec::new(), vec![(q, c)]),
+            false => (reach[q].clone(), Vec::new()),
+        })
+    };
+    let at_query_time: Vec<usize> = chosen
+        .iter()
+        .copied()
+        .filter(|&c| witness(c).is_some())
+        .collect();
+    // Per source: its query-time scans in the plan, (choice, price per
+    // unit of rate, pass per unit of rate).
+    let mut scans: BTreeMap<u64, Vec<(usize, f64, f64)>> = BTreeMap::new();
+    for &c in &at_query_time {
+        let choice = &input.choices[c];
+        for &(source, cost, pass) in &choice.coefficients.absorbed_scans {
+            scans
+                .entry(source)
+                .or_default()
+                .push((c, cost, pass.min(cost)));
+        }
+        let total = choice.coefficients.direct.cost_per_unit_rate;
+        if let (Some(source), true) = (choice.scan_source, total > 0.0) {
+            scans.entry(source).or_default().push((
+                c,
+                total,
+                choice.scan_project_cost_per_unit_rate.min(total),
+            ));
+        }
+    }
+    let rate = |c: usize| facts.rate_floor[c];
+    let has_ingestion_slot = |c: usize| {
+        input.choices[c]
+            .coefficients
+            .slots
+            .iter()
+            .any(|slot| matches!(slot.materialization, Materialization::MatAtIngestion(_)))
+    };
+    let mut out = Vec::new();
+    for (source, members) in scans {
+        // The labels the plan's query-time scans of the source cover, at
+        // least: a widened plan's scans cover them.
+        let mut labels: BTreeSet<String> = BTreeSet::new();
+        for &c in &at_query_time {
+            for (s, read) in &facts.query_time_scan_labels[c] {
+                if *s == source {
+                    labels.extend(read.iter().cloned());
+                }
+            }
+        }
+        // Apart: each scan its own node at its own price (a plan makes
+        // scans one node only when it widens them to one label set, or as
+        // a scan group's one node, a single member).
+        let apart: f64 = members.iter().map(|&(c, cost, _)| rate(c) * cost).sum();
+        // Widened: at least the dearest scan, and at least any query-time
+        // choice's scans widened to a label set the plan's scans cover.
+        let mut widened = (0.0, None);
+        for &(c, cost, _) in &members {
+            if rate(c) * cost > widened.0 {
+                widened = (rate(c) * cost, Some(c));
+            }
+        }
+        for &c in &at_query_time {
+            if has_ingestion_slot(c) {
+                // Its widened scans may run at ingestion time only.
+                continue;
+            }
+            for scan in &facts.widened_scans[c] {
+                let charge = rate(c) * scan.cost_per_unit_rate;
+                if scan.source == source && scan.labels.is_subset(&labels) && charge > widened.0 {
+                    widened = (charge, Some(c));
+                }
+            }
+        }
+        // The nodes several groups' choices have among the choices reading
+        // the source: in the plan at least once each, at least the least
+        // charge of the choices having it (the model may charge a reader's
+        // node as a pass of the scan instead).
+        let mut nodes: BTreeMap<u64, f64> = BTreeMap::new();
+        let mut readers: Vec<usize> = Vec::new();
+        for &c in &at_query_time {
+            let reads = members.iter().any(|m| m.0 == c)
+                || facts.scan_labels[c].iter().any(|(s, _)| *s == source);
+            if !reads {
+                continue;
+            }
+            for &(key, cost) in &input.choices[c].coefficients.cross_shared {
+                let charge = nodes.entry(key).or_insert(f64::INFINITY);
+                *charge = charge.min(rate(c) * cost);
+                readers.push(c);
+            }
+        }
+        let bound = apart.min(widened.0) + nodes.values().sum::<f64>();
+        // What the model charges these scans at least.
+        let mut passes: BTreeMap<u64, f64> = BTreeMap::new();
+        let mut excess: f64 = 0.0;
+        for &(c, cost, pass) in &members {
+            let p = passes.entry(cost.to_bits()).or_insert(f64::INFINITY);
+            *p = p.min(rate(c) * pass);
+            excess = excess.max(rate(c) * (cost - pass));
+        }
+        let charged = excess + passes.values().sum::<f64>();
+        if bound <= charged * (1.0 + 1e-6) {
+            continue;
+        }
+        let mut choices: Vec<usize> = members.iter().map(|&(c, _, _)| c).collect();
+        choices.extend(widened.1);
+        choices.extend(readers);
+        let mut reached = Vec::new();
+        for c in choices.clone() {
+            let (more, by) = witness(c).expect("a query-time choice");
+            choices.extend(more);
+            reached.extend(by);
+        }
+        choices.sort_unstable();
+        choices.dedup();
+        reached.sort_unstable();
+        reached.dedup();
+        let cut = TermCut {
+            term: Term::SourceScans(source),
+            choices,
+            reached,
+            built_once: Vec::new(),
+            nodes: nodes.keys().copied().collect(),
+            bound,
+        };
+        if !facts.cuts.iter().any(|old| {
+            old.term == cut.term
+                && old.choices == cut.choices
+                && old.reached == cut.reached
+                && old.bound >= bound
+        }) {
+            out.push(cut);
+        }
+    }
+    out
+}
+
+/// The choices query `q` reaches at query time in the plan of `chosen`:
+/// from its roots through the chosen choices' children read directly or by
+/// a query-time slot.
+fn query_time_reach(input: &ModelInput, chosen: &[usize], q: usize) -> Vec<usize> {
+    let mut choice_of: BTreeMap<usize, usize> = BTreeMap::new();
+    for &c in chosen {
+        for g in &input.choices[c].implements {
+            choice_of.insert(g.0 as usize, c);
+        }
+    }
+    let mut reached: BTreeSet<usize> = BTreeSet::new();
+    let mut seen: BTreeSet<usize> = BTreeSet::new();
+    let mut stack: Vec<usize> = input.queries[q]
+        .roots
+        .iter()
+        .map(|g| g.0 as usize)
+        .collect();
+    while let Some(g) = stack.pop() {
+        if !seen.insert(g) {
+            continue;
+        }
+        let Some(&c) = choice_of.get(&g) else {
+            continue;
+        };
+        reached.insert(c);
+        let choice = &input.choices[c];
+        for (child, via) in choice.children.iter().zip(&choice.coefficients.child_via) {
+            let at_query = match via {
+                crate::search::milp::model::Via::Direct => true,
+                crate::search::milp::model::Via::Slot(i) => !matches!(
+                    choice.coefficients.slots[*i].materialization,
+                    Materialization::MatAtIngestion(_)
+                ),
+            };
+            if at_query {
+                stack.push(child.0 as usize);
+            }
+        }
+    }
+    reached.into_iter().collect()
+}
+
+/// The held-memory cuts of the solution with `chosen`
+/// choices.
+pub(crate) fn held_cuts_of(input: &ModelInput, chosen: &[usize]) -> Vec<TermCut> {
+    let facts = &input.cuts;
+    let mut out: Vec<TermCut> = Vec::new();
+    let group_of = |c: usize| input.choices[c].group.0 as usize;
+    let cut = |term: Term, choices: &[usize], bound: f64| TermCut {
+        term,
+        choices: choices.to_vec(),
+        reached: Vec::new(),
+        built_once: Vec::new(),
+        nodes: Vec::new(),
+        bound,
+    };
+    // Every query may read a part (nodes several groups share are one).
+    let all_queries: Vec<usize> = (0..input.queries.len()).collect();
+    // Per query: the choices it reaches at query time and the latency
+    // they give it at least.
+    let mut latency_of: Vec<(Vec<usize>, f64)> = Vec::new();
+    for q in 0..input.queries.len() {
+        let reached = query_time_reach(input, chosen, q);
+        let groups: BTreeSet<usize> = reached.iter().map(|&c| group_of(c)).collect();
+        // Each node key the query computes, at the least work a choice it
+        // reaches gives it, or any choice of a group outside the reached
+        // ones (whose node a plan may make one with it).
+        let mut keys: BTreeMap<u64, f64> = BTreeMap::new();
+        for &c in &reached {
+            for &(key, work) in &input.choices[c].coefficients.query_work {
+                if work > 0.0 {
+                    let least = keys.entry(key).or_insert(f64::INFINITY);
+                    *least = least.min(work);
+                }
+            }
+        }
+        let mut global: BTreeMap<u64, f64> = keys.clone();
+        for (c, choice) in input.choices.iter().enumerate() {
+            for &(key, work) in &choice.coefficients.query_work {
+                if let Some(least) = global.get_mut(&key) {
+                    *least = least.min(work);
+                    if !groups.contains(&group_of(c)) {
+                        let own = keys.get_mut(&key).expect("same keys");
+                        *own = own.min(work);
+                    }
+                }
+            }
+        }
+        let latency: f64 = keys.values().sum();
+        let modeled: f64 = global.values().sum();
+        let rate_q = input.queries[q].rate;
+        // The model's own held charges already use the least latency per
+        // key: a cut helps only when the reached choices give more.
+        if latency > modeled * (1.0 + 1e-9) {
+            for &c in &reached {
+                let held = input.held.get(c).copied().unwrap_or_default();
+                // The model's held charge of a query's copy exists only for
+                // a query of the choice's group.
+                let reads = input.group_queries[group_of(c)].contains(&q);
+                if reads && held.per_reader_per_ms > 0.0 {
+                    let share = reader_share(input, &all_queries);
+                    out.push(cut(
+                        Term::HeldPerReader(q, c),
+                        &reached,
+                        held.per_reader_per_ms * rate_q * share * latency,
+                    ));
+                }
+                if held.longest_per_ms > 0.0 {
+                    out.push(cut(
+                        Term::HeldLongest(c),
+                        &reached,
+                        held.longest_per_ms * facts.rate_floor[c] * latency,
+                    ));
+                }
+            }
+        }
+        latency_of.push((reached, latency));
+    }
+    // A build several groups' choices share: held while its readers run,
+    // at their union rate, for their rate-weighted mean latency (each
+    // reader its own copy) or the longest. Readers a plan with more choices
+    // may add lower the mean at most to the share of every query that may
+    // read it. Its copies may also stay apart, each read by some of the
+    // readers: then the copy a reader `q` reads holds at least at `q`'s
+    // rate for `q`'s share of the mean (or `q`'s latency).
+    let mut builds: BTreeMap<(u64, usize), (f64, bool, BTreeSet<usize>)> = BTreeMap::new();
+    for (q, (reached, _)) in latency_of.iter().enumerate() {
+        for &c in reached {
+            for held in &facts.shared_slot_held[c] {
+                let entry =
+                    builds
+                        .entry(held.key)
+                        .or_insert((f64::INFINITY, false, BTreeSet::new()));
+                entry.0 = entry.0.min(held.cost_per_ms);
+                entry.1 |= held.per_reader;
+                entry.2.insert(q);
+            }
+        }
+    }
+    for (key, (per_ms, per_reader, readers)) in builds {
+        let mut possible: BTreeSet<usize> = BTreeSet::new();
+        for (c, held) in facts.shared_slot_held.iter().enumerate() {
+            if held.iter().any(|h| h.key == key) {
+                possible.extend(input.group_queries[group_of(c)].iter().copied());
+            }
+        }
+        let rate_all: f64 = possible.iter().map(|&q| input.queries[q].rate).sum();
+        let mut classes: Vec<usize> = readers.iter().map(|&q| input.schedule[q]).collect();
+        classes.sort_unstable();
+        classes.dedup();
+        let union = input.union_rate.get(&classes).copied().unwrap_or_else(|| {
+            readers
+                .iter()
+                .map(|&q| input.queries[q].rate)
+                .fold(0.0, f64::max)
+        });
+        let longest = readers.iter().map(|&q| latency_of[q].1).fold(0.0, f64::max);
+        let weighted = match rate_all > 0.0 {
+            true => {
+                readers
+                    .iter()
+                    .map(|&q| input.queries[q].rate * latency_of[q].1)
+                    .sum::<f64>()
+                    / rate_all
+            }
+            false => 0.0,
+        };
+        let one_node = per_ms
+            * union
+            * match per_reader {
+                true => weighted.min(longest),
+                false => longest,
+            };
+        let apart = readers
+            .iter()
+            .map(|&q| {
+                let (rate, latency) = (input.queries[q].rate, latency_of[q].1);
+                per_ms
+                    * rate
+                    * match (per_reader, rate_all > 0.0) {
+                        (true, true) => (rate * latency / rate_all).min(latency),
+                        (true, false) => 0.0,
+                        (false, _) => latency,
+                    }
+            })
+            .fold(0.0, f64::max);
+        let bound = one_node.min(apart);
+        if bound > 0.0 {
+            let mut choices: Vec<usize> = readers
+                .iter()
+                .flat_map(|&q| latency_of[q].0.iter().copied())
+                .collect();
+            choices.sort_unstable();
+            choices.dedup();
+            out.push(cut(Term::SharedBuildHeld(key), &choices, bound));
+        }
+    }
+    out.retain(|cut| {
+        !facts.cuts.iter().any(|old| {
+            old.term == cut.term
+                && old.choices == cut.choices
+                && old.reached == cut.reached
+                && old.bound >= cut.bound
+        })
+    });
+    out
+}
+
+/// The least share of a part's held memory a reader among `queries` is
+/// charged per unit of its rate: the part runs at the union rate of its
+/// readers' classes, not their rates' sum, so a reader's share is `union /
+/// sum` of its rate (Stage 3's rate-weighted mean latency at the part's
+/// rate); the least ratio over the class sets of `queries` bounds it, and
+/// readers in classes no union rate joins run at least at their fastest's
+/// rate (`fastest / sum`).
+pub(crate) fn reader_share(input: &ModelInput, queries: &[usize]) -> f64 {
+    let fastest_over_sum = queries.iter().map(|&q| {
+        let fastest = input.queries[q].rate;
+        let sum: f64 = queries
+            .iter()
+            .map(|&r| input.queries[r].rate)
+            .filter(|&rate| rate <= fastest)
+            .sum();
+        match sum > 0.0 {
+            true => (fastest / sum).min(1.0),
+            false => 1.0,
+        }
+    });
+    input
+        .union_rate
+        .iter()
+        .filter_map(|(set, &union)| {
+            let rate_sum: f64 = queries
+                .iter()
+                .filter(|&&q| set.contains(&input.schedule[q]))
+                .map(|&q| input.queries[q].rate)
+                .sum();
+            let covered = set
+                .iter()
+                .all(|class| queries.iter().any(|&q| input.schedule[q] == *class));
+            (covered && rate_sum > 0.0).then_some((union / rate_sum).min(1.0))
+        })
+        .chain(fastest_over_sum)
+        .fold(1.0, f64::min)
+}
+
+/// The ingestion-time copy cuts of the solution with `chosen` choices and
+/// `built_once` groups: each group built once that a chosen choice's
+/// ingestion-time slot reads costs at least its chosen choice's direct
+/// part at ingestion time.
+pub(crate) fn ingestion_copy_cuts_of(
+    input: &ModelInput,
+    chosen: &[usize],
+    built_once: &BTreeSet<crate::search::memo::GroupId>,
+) -> Vec<TermCut> {
+    let facts = &input.cuts;
+    let mut choice_of: BTreeMap<usize, usize> = BTreeMap::new();
+    for &c in chosen {
+        choice_of.insert(input.choices[c].group.0 as usize, c);
+    }
+    let mut out: Vec<TermCut> = Vec::new();
+    for &p in chosen {
+        let parent = &input.choices[p];
+        for (child, via) in parent.children.iter().zip(&parent.coefficients.child_via) {
+            let crate::search::milp::model::Via::Slot(i) = via else {
+                continue;
+            };
+            if !matches!(
+                parent.coefficients.slots[*i].materialization,
+                Materialization::MatAtIngestion(_)
+            ) || !built_once.contains(child)
+            {
+                continue;
+            }
+            let h = child.0 as usize;
+            let Some(&c) = choice_of.get(&h) else {
+                continue;
+            };
+            let bound = facts.at_ingestion_cost.get(c).copied().unwrap_or(0.0);
+            if bound > 0.0 {
+                out.push(TermCut {
+                    term: Term::IngestionCopy(h),
+                    choices: vec![c],
+                    reached: Vec::new(),
+                    built_once: vec![h],
+                    nodes: Vec::new(),
+                    bound,
+                });
+            }
+        }
+    }
+    out.retain(|cut| !facts.cuts.contains(cut));
+    out.dedup();
+    out
+}
diff --git a/crates/plan-selection/src/search/milp/decode.rs b/crates/plan-selection/src/search/milp/decode.rs
new file mode 100644
index 000000000..71792af66
--- /dev/null
+++ b/crates/plan-selection/src/search/milp/decode.rs
@@ -0,0 +1,667 @@
+//! Solution decoding (joint-planning-search design, "One solution, one
+//! plan"): a MILP solution — the chosen choices `x` (a candidate at its
+//! site thresholds `u`, with its private slots' materializations), the
+//! groups built once `m`, the scans shared at a horizon — maps to the one
+//! physical plan it describes. Incumbent evaluation then runs Stage 3 once
+//! on that plan.
+//!
+//! What a solution decodes to:
+//! - each chosen choice's candidate at its site thresholds, its private
+//!   builds materialized as its assignment;
+//! - a group built once (`m`) or a scan shared at a horizon is one node
+//!   for every query; a scan a choice reads apart is that reader's own;
+//! - every other identical node is one node (the identical-expression
+//!   rule): not built once, Stage 3 costs it as its copies, recomputed per
+//!   consumer, so a plan keeping two identical non-scan nodes apart costs
+//!   the same or more.
+//!
+//! Known gaps, where a brute-force plan is not the decoding of any
+//! solution (the model still bounds it from below):
+//! - a group read apart whose shareable inner nodes are one node built
+//!   once (the model's per-choice `shared_nodes` choice);
+//! - scan-label widening (every metric's scans widened to one label set,
+//!   the identical-expression rule's variant) is no model variable yet:
+//!   the model charges a metric's scans as one widened scan may cost, so a
+//!   solution also decodes to its widened plan, and incumbent evaluation
+//!   runs on both.
+
+use std::collections::{BTreeSet, HashMap};
+use std::rc::Rc;
+
+use asap_logical_asap_optimization::asap_aware_cse::{
+    sharing_candidates, with_shared_candidates, with_union_scan_labels,
+};
+use asap_logical_asap_optimization::per_computation_candidates::{
+    recompose_guarantees, PerComputationCandidates,
+};
+use asap_physical_asap_optimization::{optimize_physical_asap_chosen, PhysicalCandidate};
+use asap_types::ir::{ASAPOp, NonASAPOp, OperatorNode, QueryRoot};
+use asap_types::physical::{Materialization, Medium};
+use asap_types::workload::{DataWorkload, QueryWorkloadEntry};
+
+use crate::search::memo::{GroupId, ImplementationCandidate, Memo};
+use crate::search::milp::model::Choice;
+use crate::PlanningModels;
+
+/// A decoded plan: its logical roots and its one physical candidate.
+pub(super) struct Decoded {
+    pub roots: Vec<QueryRoot>,
+    pub physical: PhysicalCandidate,
+}
+
+/// A solution's decoded plans.
+pub(super) struct Decoding {
+    pub plans: Vec<Decoded>,
+    /// Whether every combination of the free nodes' Stage 2 options was
+    /// evaluated: otherwise the choice set has plans incumbent evaluation
+    /// did not see, and the solution's objective stays a lower bound.
+    pub every_combination: bool,
+}
+
+/// The plan the solution with `chosen` choices (indices into `choices`,
+/// whose candidates at their thresholds are `versions`) and `built_once`
+/// groups describes, named `id`, with up to `max_free_combinations`
+/// combinations of its free nodes' options; `None` when it cannot be
+/// composed.
+#[allow(clippy::too_many_arguments)]
+pub(super) fn decode<Id>(
+    memo: &Memo,
+    per_computation: &PerComputationCandidates<Id>,
+    choices: &[Choice],
+    versions: &[ImplementationCandidate],
+    chosen: &[usize],
+    built_once: &BTreeSet<GroupId>,
+    entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
+    models: PlanningModels<'_>,
+    id: &str,
+    max_free_combinations: usize,
+) -> Option<Decoding> {
+    let n = memo.groups.len();
+    let mut choice_of: Vec<Option<usize>> = vec![None; n];
+    for &c in chosen {
+        choice_of[choices[c].group.0 as usize] = Some(c);
+    }
+    // Groups made one node for every query: built once, or a scan shared
+    // at a horizon.
+    let mut one_node: Vec<bool> = vec![false; n];
+    for g in built_once {
+        one_node[g.0 as usize] = true;
+    }
+    // Scans a choice reads apart: its own copy, not the scan group's node.
+    let mut apart_scans: BTreeSet<(usize, usize)> = BTreeSet::new();
+    for &c in chosen {
+        for (scan, horizon) in &choices[c].scan_children {
+            match horizon {
+                Some(_) => one_node[scan.0 as usize] = true,
+                None => {
+                    apart_scans.insert((choices[c].group.0 as usize, scan.0 as usize));
+                }
+            }
+        }
+    }
+    let plan: Vec<Option<ImplementationCandidate>> = choice_of
+        .iter()
+        .map(|c| c.map(|c| versions[c].clone()))
+        .collect();
+    let composed = compose(memo, per_computation, &plan, &one_node, &apart_scans)?;
+    // Each node's guarantee composed again over the chosen implementations.
+    let recomposed = recompose_guarantees(&composed.roots).ok()?;
+    let mut table = composed.table;
+    remap(&mut table, &composed.roots, &recomposed);
+    // Shared builds resized to one size are one build.
+    let group_choice: crate::search::memo::brute_force::GroupChoice = choice_of
+        .iter()
+        .map(|c| c.map(|c| choices[c].candidate))
+        .collect();
+    let resized = crate::search::memo::brute_force::resized_builds(memo, &group_choice, &plan);
+    let resized_shared =
+        crate::search::memo::rules::capability::share_equal_builds(&recomposed, &resized);
+    remap(&mut table, &recomposed, &resized_shared);
+    // Identical nodes are one node (the identical-expression rule), but a
+    // scan read apart stays each reader's: a scan shared at a horizon is
+    // already one node. Shared and not built once, a node costs what its
+    // copies cost (recomputed per consumer); the model's shared slots and
+    // nodes assume one node.
+    let (mut plans, mut every_combination) = implemented(
+        &resized_shared,
+        &table,
+        Some(false),
+        &choice_of,
+        choices,
+        entries,
+        data,
+        models,
+        id,
+        max_free_combinations,
+    );
+    // Scan widening is no model variable yet (a decoding gap): the model
+    // charges a metric's scans as one widened scan may cost, so the plan
+    // with every metric's scans widened to one label set, and identical
+    // scans one node, is evaluated too.
+    if let Ok(Some(widened)) = with_union_scan_labels(&resized_shared) {
+        let mut widened_table = table.clone();
+        remap(&mut widened_table, &resized_shared, &widened);
+        for share_scans in [true, false] {
+            let (widened_plans, every) = implemented(
+                &widened,
+                &widened_table,
+                Some(share_scans),
+                &choice_of,
+                choices,
+                entries,
+                data,
+                models,
+                &format!("{id}w{}", u8::from(share_scans)),
+                max_free_combinations,
+            );
+            plans.extend(widened_plans);
+            every_combination &= every;
+        }
+    }
+    Some(Decoding {
+        plans,
+        every_combination,
+    })
+}
+
+/// The physical plans of `roots` with identical nodes one node (scans too
+/// when `share_scans`), each node materialized as `table` says the
+/// solution chose, and whether every combination of the free nodes'
+/// options is among them.
+#[allow(clippy::too_many_arguments)]
+fn implemented(
+    roots: &[QueryRoot],
+    table: &HashMap<*const OperatorNode, Role>,
+    share_scans: Option<bool>,
+    choice_of: &[Option<usize>],
+    choices: &[Choice],
+    entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
+    models: PlanningModels<'_>,
+    id: &str,
+    max_free_combinations: usize,
+) -> (Vec<Decoded>, bool) {
+    let shared = match share_scans {
+        Some(share_scans) => {
+            let candidates = sharing_candidates(roots);
+            let subset: Vec<bool> = candidates
+                .iter()
+                .map(|c| {
+                    share_scans || !matches!(c.nodes[0].non_asap(), Some(NonASAPOp::Scan { .. }))
+                })
+                .collect();
+            with_shared_candidates(roots, &candidates, &subset)
+        }
+        None => roots.to_vec(),
+    };
+    let mut table = table.clone();
+    remap(&mut table, roots, &shared);
+    let wanted_of = |logical: *const OperatorNode| match table.get(&logical) {
+        Some(Role::Slot { group, index }) => Some(
+            choice_of[*group]
+                .and_then(|c| choices[c].assignment.get(*index).copied())
+                .unwrap_or(Materialization::NotMaterialized),
+        ),
+        Some(Role::BuiltOnce) => Some(Materialization::MatAtQuery(Medium::Memory)),
+        None => None,
+    };
+    // The nodes the solution sets no materialization for (identical nodes
+    // made one, which the model charges at the cheaper of once and per
+    // consumer): Stage 2's options for each, in Stage 2's order.
+    let mut free: Vec<Vec<Materialization>> = Vec::new();
+    let _ = optimize_physical_asap_chosen(
+        id,
+        &shared,
+        entries,
+        data,
+        models.deployment_capabilities,
+        &mut |logical, _readers, options| match wanted_of(logical) {
+            Some(wanted) if options.contains(&wanted) => wanted,
+            Some(_) => Materialization::NotMaterialized,
+            None => {
+                if options.len() > 1 {
+                    free.push(options.to_vec());
+                }
+                options[0]
+            }
+        },
+    );
+    // Every combination of the free nodes' options, up to
+    // `max_free_combinations`; beyond, each free node recomputed or each
+    // built once.
+    let count = free
+        .iter()
+        .try_fold(1usize, |n, options| n.checked_mul(options.len()))
+        .filter(|&n| n <= max_free_combinations);
+    let every_combination = count.is_some();
+    let combinations: Vec<Vec<Materialization>> = match count {
+        Some(n) => (0..n)
+            .map(|mut k| {
+                free.iter()
+                    .map(|options| {
+                        let m = options[k % options.len()];
+                        k /= options.len();
+                        m
+                    })
+                    .collect()
+            })
+            .collect(),
+        None => {
+            let pick = |m: Materialization| -> Vec<Materialization> {
+                free.iter()
+                    .map(|options| if options.contains(&m) { m } else { options[0] })
+                    .collect()
+            };
+            vec![
+                pick(Materialization::NotMaterialized),
+                pick(Materialization::MatAtQuery(Medium::Memory)),
+            ]
+        }
+    };
+    let mut out = Vec::new();
+    for (k, combination) in combinations.into_iter().enumerate() {
+        let mut next_free = 0usize;
+        let physical = optimize_physical_asap_chosen(
+            &format!("{id}-{k}"),
+            &shared,
+            entries,
+            data,
+            models.deployment_capabilities,
+            &mut |logical, _readers, options| match wanted_of(logical) {
+                Some(wanted) if options.contains(&wanted) => wanted,
+                Some(_) => Materialization::NotMaterialized,
+                None if options.len() > 1 => {
+                    next_free += 1;
+                    combination[next_free - 1]
+                }
+                None => options[0],
+            },
+        );
+        if let Ok(Some(physical)) = physical {
+            out.push(Decoded {
+                roots: shared.clone(),
+                physical,
+            });
+        }
+    }
+    (out, every_combination)
+}
+
+/// The most combinations of Stage 2 options of the nodes no solution
+/// variable sets that incumbent evaluation runs on.
+pub(super) const MAX_FREE_COMBINATIONS: usize = 16;
+
+/// What a node of the decoded plan is to the solution.
+#[derive(Debug, Clone, Copy)]
+enum Role {
+    /// The `index`-th private build of `group`'s chosen choice (in Stage
+    /// 2's order: children first), materialized as the choice's
+    /// `assignment[index]`.
+    Slot { group: usize, index: usize },
+    /// A group's own node, built once at query time.
+    BuiltOnce,
+}
+
+struct Composed {
+    roots: Vec<QueryRoot>,
+    table: HashMap<*const OperatorNode, Role>,
+}
+
+/// The query roots of `plan`, as [`crate::search::memo::brute_force::composed_roots`] composes
+/// them before recomposing guarantees, but with each `one_node` group built
+/// once for every query, and the roles of the solution's nodes.
+fn compose<Id>(
+    memo: &Memo,
+    per_computation: &PerComputationCandidates<Id>,
+    plan: &[Option<ImplementationCandidate>],
+    one_node: &[bool],
+    apart_scans: &BTreeSet<(usize, usize)>,
+) -> Option<Composed> {
+    memo.topological_order()?;
+    let queries = per_computation.roots.len();
+    let mut below: Vec<Vec<bool>> = vec![vec![false; memo.groups.len()]; queries];
+    for (g, group) in memo.groups.iter().enumerate() {
+        for &q in &group.roots_of {
+            crate::search::memo::brute_force::mark_below(memo, g, &mut below[q]);
+        }
+    }
+    let scope_of = |scope: usize, c: usize| {
+        if one_node[c] {
+            queries
+        } else if scope < queries && below[scope][c] {
+            scope
+        } else {
+            (0..queries).find(|&q| below[q][c]).unwrap_or(queries)
+        }
+    };
+    let building: Vec<Option<(usize, &ImplementationCandidate)>> = (0..memo.groups.len())
+        .map(|g| {
+            plan[g].as_ref().map(|c| (g, c)).or_else(|| {
+                plan.iter()
+                    .enumerate()
+                    .filter_map(|(h, c)| c.as_ref().map(|c| (h, c)))
+                    .find(|(_, c)| c.implements.contains(&GroupId(g as u32)))
+            })
+        })
+        .collect();
+    let mut state = State {
+        built: vec![HashMap::new(); queries + 1],
+        table: HashMap::new(),
+        apart_scans,
+    };
+    let roots: Vec<QueryRoot> = per_computation
+        .roots
+        .iter()
+        .enumerate()
+        .map(|(q, (_, root))| {
+            let group_of = |node: &Rc<OperatorNode>, state: &mut State| {
+                let g = memo
+                    .groups
+                    .iter()
+                    .position(|group| group.roots_of.contains(&q) && *group.node == **node)?;
+                build(
+                    memo,
+                    &building,
+                    one_node,
+                    scope_of(q, g),
+                    g,
+                    &scope_of,
+                    state,
+                )
+            };
+            Some(match root {
+                QueryRoot::Operator(node) => QueryRoot::Operator(group_of(node, &mut state)?),
+                QueryRoot::Scalar(expr) => {
+                    let mut missing = false;
+                    let mapped = expr.map_operator_refs(&mut |node| {
+                        group_of(node, &mut state).unwrap_or_else(|| {
+                            missing = true;
+                            Rc::clone(node)
+                        })
+                    });
+                    if missing {
+                        return None;
+                    }
+                    QueryRoot::Scalar(mapped)
+                }
+            })
+        })
+        .collect::<Option<_>>()?;
+    Some(Composed {
+        roots,
+        table: state.table,
+    })
+}
+
+struct State<'a> {
+    built: Vec<HashMap<usize, Rc<OperatorNode>>>,
+    table: HashMap<*const OperatorNode, Role>,
+    /// (reading group, scan group) pairs whose scan is the reader's own copy.
+    apart_scans: &'a BTreeSet<(usize, usize)>,
+}
+
+/// Group `g` built in `scope`, as [`crate::search::memo::brute_force`]'s `build`, recording
+/// its private builds and, built once, its own node.
+fn build(
+    memo: &Memo,
+    building: &[Option<(usize, &ImplementationCandidate)>],
+    one_node: &[bool],
+    scope: usize,
+    g: usize,
+    scope_of: &dyn Fn(usize, usize) -> usize,
+    state: &mut State,
+) -> Option<Rc<OperatorNode>> {
+    if let Some(node) = state.built[scope].get(&g) {
+        return Some(Rc::clone(node));
+    }
+    let (owner, candidate) = building[g]?;
+    let mut replaced: HashMap<*const OperatorNode, Rc<OperatorNode>> = HashMap::new();
+    for child in &candidate.children {
+        let c = child.0 as usize;
+        if state.apart_scans.contains(&(owner, c)) {
+            // A scan read apart: a node of the reader's own.
+            let copy = Rc::new(OperatorNode::clone(&memo.groups[c].node));
+            replaced.insert(Rc::as_ptr(&memo.groups[c].node), copy);
+            continue;
+        }
+        let node = build(
+            memo,
+            building,
+            one_node,
+            scope_of(scope, c),
+            c,
+            scope_of,
+            state,
+        )?;
+        replaced.insert(Rc::as_ptr(&memo.groups[c].node), node);
+    }
+    let mut rebuilt = HashMap::new();
+    let node = Rc::new(OperatorNode::clone(
+        &crate::search::memo::brute_force::substituted(&candidate.node, &replaced, &mut rebuilt),
+    ));
+    // The candidate's private builds, children first, as Stage 2 orders
+    // the choice's own sub-DAG; the root is a fresh node of its own.
+    let children: BTreeSet<*const OperatorNode> = candidate
+        .children
+        .iter()
+        .map(|c| Rc::as_ptr(&memo.groups[c.0 as usize].node))
+        .collect();
+    let mut seen = BTreeSet::new();
+    let mut order = Vec::new();
+    private_postorder(&candidate.node, &children, &mut seen, &mut order);
+    let mut index = 0;
+    for original in order {
+        if !is_build(&original) {
+            continue;
+        }
+        let new = if Rc::ptr_eq(&original, &candidate.node) {
+            Rc::as_ptr(&node)
+        } else {
+            rebuilt
+                .get(&Rc::as_ptr(&original))
+                .map_or(Rc::as_ptr(&original), Rc::as_ptr)
+        };
+        state.table.insert(
+            new,
+            Role::Slot {
+                group: owner,
+                index,
+            },
+        );
+        index += 1;
+    }
+    // Built once (the model's `once[g]`) overrides the slot assignment of
+    // a group's own node that is a build. The model builds a group once
+    // only through its own choice (`owner == g`).
+    if one_node[g] && owner == g {
+        state.table.insert(Rc::as_ptr(&node), Role::BuiltOnce);
+    }
+    for implemented in &candidate.implements {
+        state.built[scope].insert(implemented.0 as usize, Rc::clone(&node));
+    }
+    Some(node)
+}
+
+/// `node`'s nodes, children first, without entering `children` (the child
+/// groups' nodes).
+fn private_postorder(
+    node: &Rc<OperatorNode>,
+    children: &BTreeSet<*const OperatorNode>,
+    seen: &mut BTreeSet<*const OperatorNode>,
+    order: &mut Vec<Rc<OperatorNode>>,
+) {
+    let ptr = Rc::as_ptr(node);
+    if children.contains(&ptr) || !seen.insert(ptr) {
+        return;
+    }
+    for child in node.children() {
+        private_postorder(child, children, seen, order);
+    }
+    order.push(Rc::clone(node));
+}
+
+fn is_build(node: &OperatorNode) -> bool {
+    matches!(
+        node.asap(),
+        Some(ASAPOp::SummaryAgg { .. } | ASAPOp::PaneChain { .. })
+    )
+}
+
+/// Moves `table`'s entries from `before`'s nodes to the nodes in the same
+/// places of `after` (the same roots rebuilt node by node).
+fn remap(
+    table: &mut HashMap<*const OperatorNode, Role>,
+    before: &[QueryRoot],
+    after: &[QueryRoot],
+) {
+    let mut moved: HashMap<*const OperatorNode, Role> = HashMap::new();
+    let mut seen: BTreeSet<*const OperatorNode> = BTreeSet::new();
+    fn walk(
+        a: &Rc<OperatorNode>,
+        b: &Rc<OperatorNode>,
+        table: &HashMap<*const OperatorNode, Role>,
+        moved: &mut HashMap<*const OperatorNode, Role>,
+        seen: &mut BTreeSet<*const OperatorNode>,
+    ) {
+        if !seen.insert(Rc::as_ptr(a)) {
+            return;
+        }
+        // A `Project` scan widening put above a node: the node is below it.
+        let is_project = |n: &OperatorNode| matches!(n.non_asap(), Some(NonASAPOp::Project { .. }));
+        if is_project(b) && !is_project(a) {
+            if let [inner] = b.children().as_slice() {
+                seen.remove(&Rc::as_ptr(a));
+                return walk(a, inner, table, moved, seen);
+            }
+        }
+        if let Some(role) = table.get(&Rc::as_ptr(a)) {
+            moved.entry(Rc::as_ptr(b)).or_insert(*role);
+        }
+        for (x, y) in a.children().into_iter().zip(b.children()) {
+            walk(x, y, table, moved, seen);
+        }
+    }
+    for (a, b) in before.iter().zip(after) {
+        let pairs: Vec<(&Rc<OperatorNode>, &Rc<OperatorNode>)> = match (a, b) {
+            (QueryRoot::Operator(x), QueryRoot::Operator(y)) => vec![(x, y)],
+            (QueryRoot::Scalar(x), QueryRoot::Scalar(y)) => x
+                .operator_refs()
+                .into_iter()
+                .zip(y.operator_refs())
+                .collect(),
+            _ => Vec::new(),
+        };
+        for (x, y) in pairs {
+            walk(x, y, table, &mut moved, &mut seen);
+        }
+    }
+    *table = moved;
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::search::milp::model::{model_input_of, solve, tests::workload};
+    use asap_types::ir::properties::AccuracyTarget;
+
+    const QUERIES: [&str; 2] = ["sum_over_time(m[5m]) * 2", "sum_over_time(m[5m]) / 3"];
+
+    fn targets() -> [AccuracyTarget; 2] {
+        [AccuracyTarget::exact(), AccuracyTarget::exact()]
+    }
+
+    /// The model's optimum over a shared sum decodes to plans, every
+    /// combination of its free nodes' options among them.
+    #[test]
+    fn a_solution_decodes_to_its_plans() {
+        let (per_computation, entries, data) = workload(&QUERIES, &targets(), None);
+        let models = PlanningModels::builtin();
+        let memo = crate::search::memo::brute_force::memo(&per_computation, &entries);
+        let (input, versions) = model_input_of(&memo, &entries, &data, models).unwrap();
+        let solution = solve(&input).unwrap().unwrap();
+        let decoding = decode(
+            &memo,
+            &per_computation,
+            &input.choices,
+            &versions,
+            &solution.chosen,
+            &solution.built_once,
+            &entries,
+            &data,
+            models,
+            "t",
+            MAX_FREE_COMBINATIONS,
+        )
+        .expect("the optimum decodes");
+        assert!(!decoding.plans.is_empty());
+        assert!(decoding.every_combination);
+        for plan in &decoding.plans {
+            assert_eq!(plan.physical.roots.len(), QUERIES.len());
+        }
+    }
+
+    /// A group built once whose own node is a build: its node is built
+    /// once (the model's `once`), not materialized as the choice's slot.
+    #[test]
+    fn a_group_built_once_whose_node_is_a_build_is_built_once() {
+        let (per_computation, entries, _) = workload(&QUERIES, &targets(), None);
+        let memo = crate::search::memo::brute_force::memo(&per_computation, &entries);
+        // The shared sum's group, and a candidate of it cut down to its
+        // build: a group whose own node is a build.
+        let g = memo
+            .groups
+            .iter()
+            .position(|group| group.parents.len() > 1)
+            .expect("a shared group");
+        fn build_below(
+            node: &Rc<OperatorNode>,
+            children: &[*const OperatorNode],
+        ) -> Option<Rc<OperatorNode>> {
+            if children.contains(&Rc::as_ptr(node)) {
+                return None;
+            }
+            if is_build(node) {
+                return Some(Rc::clone(node));
+            }
+            node.children()
+                .into_iter()
+                .find_map(|child| build_below(child, children))
+        }
+        let mut plan: Vec<Option<ImplementationCandidate>> = memo
+            .candidates
+            .iter()
+            .map(|candidates| candidates.first().cloned())
+            .collect();
+        let (candidate, build) = memo.candidates[g]
+            .iter()
+            .find_map(|candidate| {
+                let children: Vec<*const OperatorNode> = candidate
+                    .children
+                    .iter()
+                    .map(|c| Rc::as_ptr(&memo.groups[c.0 as usize].node))
+                    .collect();
+                Some((candidate, build_below(&candidate.node, &children)?))
+            })
+            .expect("a candidate with a build");
+        plan[g] = Some(ImplementationCandidate {
+            node: build,
+            ..candidate.clone()
+        });
+        let mut one_node = vec![false; memo.groups.len()];
+        one_node[g] = true;
+        let composed =
+            compose(&memo, &per_computation, &plan, &one_node, &BTreeSet::new()).unwrap();
+        let roles: Vec<Role> = composed.table.values().copied().collect();
+        assert!(
+            roles.iter().any(|r| matches!(r, Role::BuiltOnce)),
+            "{roles:?}"
+        );
+        assert!(
+            !roles
+                .iter()
+                .any(|r| matches!(r, Role::Slot { group, .. } if *group == g)),
+            "{roles:?}"
+        );
+    }
+}
diff --git a/crates/plan-selection/src/search/milp/mod.rs b/crates/plan-selection/src/search/milp/mod.rs
index b4b448054..80baf74b6 100644
--- a/crates/plan-selection/src/search/milp/mod.rs
+++ b/crates/plan-selection/src/search/milp/mod.rs
@@ -3,5 +3,8 @@
 
 #![allow(dead_code)]
 
+pub(crate) mod clusters;
 pub(crate) mod coefficients;
+pub(crate) mod cuts;
+pub(crate) mod decode;
 pub(crate) mod model;
diff --git a/crates/plan-selection/src/search/milp/model.rs b/crates/plan-selection/src/search/milp/model.rs
index 352154f34..0bb11654a 100644
--- a/crates/plan-selection/src/search/milp/model.rs
+++ b/crates/plan-selection/src/search/milp/model.rs
@@ -15,11 +15,12 @@ use std::collections::{BTreeMap, BTreeSet, HashMap};
 
 use good_lp::solvers::highs::highs;
 use good_lp::{
-    constraint, variable, Expression, ProblemVariables, Solution, SolverModel, Variable, WithMipGap,
+    constraint, variable, Expression, ProblemVariables, Solution, SolverModel, Variable,
 };
 
-use crate::search::memo::{GroupId, Memo};
+use crate::search::memo::{GroupId, ImplementationCandidate, Memo};
 use crate::search::milp::coefficients::LazyTerms;
+use crate::search::milp::cuts::{SharedBuildHeld, Term, TermCuts, WidenedScan};
 use crate::{CostModel, PlanSelector, PlanningModels, ResourceUnitCost, Selection, SelectionError};
 use asap_logical_asap_optimization::per_computation_candidates::PerComputationCandidates;
 use asap_physical_asap_optimization::PhysicalCandidate;
@@ -110,6 +111,13 @@ pub(crate) struct ChoiceCoefficients {
     /// expressions, scans widened to one label set): the model's input
     /// moves them to `cross_shared`.
     pub direct_nodes: Vec<(u64, f64)>,
+    /// Each node a query reaching this choice at query time computes (its
+    /// direct part and its query-time slots), as (node key, work in ms per
+    /// computation). Nodes of several choices with one key may be one node
+    /// in the plan (identical expressions, scans of one source widened),
+    /// which a query's latency counts once: the model counts each key
+    /// once per query, at the least work any choice gives it.
+    pub query_work: Vec<(u64, f64)>,
 }
 
 /// A shareable node of a choice's direct part ([`ChoiceCoefficients::shared_nodes`]).
@@ -162,6 +170,13 @@ pub(crate) struct Choice {
     /// pass once (readers of one scan share its node), plus the largest
     /// excess over it.
     pub scan_project_cost_per_unit_rate: f64,
+    /// The keys of its private builds built once at query time that only
+    /// another group's identical build makes so: alone, the build has one
+    /// reader and Stage 2 recomputes it. A plan with this choice also has
+    /// another group's choice building the same node once.
+    pub built_once_with_other_group: Vec<u32>,
+    /// The keys of all its private builds built once at query time.
+    pub built_once_keys: Vec<u32>,
 }
 
 /// A shared group's own node, built once per evaluation when chosen.
@@ -203,6 +218,14 @@ pub(crate) struct ModelInput {
     /// With the groups that plan does not need: a plan needing any of them
     /// is another plan, not cut.
     pub excluded: Vec<(Vec<usize>, Vec<usize>, f64)>,
+    /// Per cut of `excluded` (by position; missing for the cuts that set
+    /// none), the groups built once and the groups not built once in the
+    /// plan it excludes: a solution that differs in either describes
+    /// another plan.
+    pub excluded_once: Vec<(Vec<GroupId>, Vec<GroupId>)>,
+    /// Per choice, the threshold index it reads each coupling site at (by
+    /// site key): the choices of one plan read each site at one threshold.
+    pub site_thresholds: Vec<Vec<(u64, usize)>>,
     /// Lower bounds on the omitted cost, added as `extra ≥ bound` when all
     /// of a choice set is chosen (local cuts).
     pub extra_bounds: Vec<(Vec<usize>, f64)>,
@@ -226,14 +249,27 @@ pub(crate) struct ModelInput {
     /// Per choice, the memory its parts hold while a query's evaluation
     /// runs ([`HeldAtQueryTime`]).
     pub held: Vec<HeldAtQueryTime>,
+    /// Per choice, its direct part's price per second computed at
+    /// ingestion time, under an ingestion-time slot (scans read one second
+    /// of rows; [`LazyTerms::direct_at_ingestion`]); `None` when it cannot
+    /// run at ingestion time.
+    pub direct_at_ingestion_cost: Vec<Option<f64>>,
+    /// The term cuts so far and the facts they read
+    /// ([`crate::search::milp::cuts`]).
+    pub cuts: TermCuts,
     /// The solver's time limit in seconds, if any.
     pub time_limit_s: Option<f64>,
+    /// Only solutions cheaper than this (the incumbent's price), by more
+    /// than [`CUTOFF_TOLERANCE`], are wanted: the solver prunes every
+    /// branch whose bound reaches it, and the model is infeasible when
+    /// none is left.
+    pub objective_cutoff: Option<f64>,
     /// Solve the linear relaxation (every binary in [0, 1]): a lower bound
     /// on the model's optimum.
     pub relaxed: bool,
 }
 
-/// The memory a choice's parts hold while an evaluation runs, priced per
+/// The memory a choice's parts hold while an evaluation runs, charged per
 /// ms of latency: Stage 3 charges `bytes · rate · E / 1000` at the memory
 /// price, `E` the latency of the queries reaching the part.
 #[derive(Debug, Clone, Copy, Default, PartialEq)]
@@ -291,6 +327,34 @@ fn held_at_query_time(
     held
 }
 
+/// The memory `choice`'s query-time builds that several groups' choices
+/// may share hold while an evaluation runs (left out of
+/// [`held_at_query_time`]).
+fn shared_build_held(
+    choice: &Choice,
+    lazy: &LazyTerms,
+    prices: &ResourceUnitCost,
+    shared_keys: &BTreeSet<(u64, usize)>,
+) -> Vec<SharedBuildHeld> {
+    lazy.held
+        .iter()
+        .filter_map(|term| {
+            let Via::Slot(i) = term.via else {
+                return None;
+            };
+            let slot = &choice.coefficients.slots[i];
+            let key = slot.shared_site.filter(|key| shared_keys.contains(key))?;
+            (!matches!(slot.materialization, Materialization::MatAtIngestion(_))).then_some(
+                SharedBuildHeld {
+                    key,
+                    cost_per_ms: prices.cost_per_memory_byte_second * term.bytes as f64 / 1000.0,
+                    per_reader: term.per_reader,
+                },
+            )
+        })
+        .collect()
+}
+
 /// The model's optimum.
 #[derive(Debug, Clone, PartialEq)]
 pub(crate) struct ModelSolution {
@@ -303,11 +367,52 @@ pub(crate) struct ModelSolution {
     /// Whether the solver proved it optimal (not stopped at its time
     /// limit): only then is `objective` a lower bound.
     pub optimal: bool,
+    /// When stopped at its time limit, the solver's proven lower bound on
+    /// every plan's price (its dual bound), when it has one.
+    pub proven_bound: Option<f64>,
     /// The model's size.
     pub variables: u64,
     pub constraints: u64,
 }
 
+/// Why [`solve`] returned no solution and no proof of infeasibility.
+#[derive(Debug, Clone, PartialEq)]
+pub(crate) enum SolveError {
+    /// The solver stopped at its time limit before finding a solution; its
+    /// proven lower bound on every plan's price, when it has one.
+    StoppedWithoutSolution { proven_bound: Option<f64> },
+    /// The solver failed.
+    Failed(String),
+}
+
+impl std::fmt::Display for SolveError {
+    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
+        match self {
+            Self::StoppedWithoutSolution { .. } => {
+                write!(f, "the solver stopped at its time limit without a solution")
+            }
+            Self::Failed(error) => write!(f, "{error}"),
+        }
+    }
+}
+
+/// A solve's column values, as a [`good_lp::Solution`].
+struct Columns {
+    status: good_lp::solvers::SolutionStatus,
+    /// Each variable's column, in the order the problem added them.
+    column_of: HashMap<Variable, usize>,
+    values: Vec<f64>,
+}
+
+impl good_lp::Solution for Columns {
+    fn status(&self) -> good_lp::solvers::SolutionStatus {
+        self.status
+    }
+    fn value(&self, variable: Variable) -> f64 {
+        self.values[self.column_of[&variable]]
+    }
+}
+
 /// A consumer of the model: a query, a materialized private slot of a
 /// choice, or a shared group built once.
 #[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, PartialOrd, Ord)]
@@ -317,20 +422,39 @@ enum Consumer {
     Group(GroupId),
 }
 
-/// Build and solve the model. `None` when it is infeasible.
-pub(crate) fn solve(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
+/// Build and solve the model. `None` when it is infeasible (no plan below
+/// [`ModelInput::objective_cutoff`], when set).
+pub(crate) fn solve(input: &ModelInput) -> Result<Option<ModelSolution>, SolveError> {
     // Prices per second are tiny (1e-6 and less) next to the solver's
     // feasibility tolerance (1e-7): solve in units where the largest
     // per-choice charge is about 1, and convert the objective back.
     let scale = price_scale(input);
     let scaled = scaled(input, scale);
-    let solution = solve_scaled(&scaled)?;
+    let solution = match solve_scaled(&scaled) {
+        Err(SolveError::StoppedWithoutSolution { proven_bound }) => {
+            return Err(SolveError::StoppedWithoutSolution {
+                proven_bound: proven_bound.map(|bound| bound / scale),
+            })
+        }
+        solution => solution?,
+    };
     Ok(solution.map(|mut solution| {
         solution.objective /= scale;
+        solution.proven_bound = solution.proven_bound.map(|bound| bound / scale);
         solution
     }))
 }
 
+/// How far below the incumbent's price (relatively) a solution must be
+/// to be wanted: the solver's optimum is exact only to its own tolerances
+/// (about 1e-7 of the scaled objective), so plans the model prices within
+/// them of the incumbent tie with it.
+const CUTOFF_TOLERANCE: f64 = 1e-6;
+
+/// How far above the cutoff, in the scaled objective, the solver may
+/// prune a node (HiGHS's default MIP feasibility tolerance).
+const SOLVER_PRUNING_TOLERANCE: f64 = 1e-6;
+
 /// How far below its objective a solve stopped at the solver's gap
 /// tolerance (a relative gap of about 1e-6 at most) may have its optimum.
 const GAP_TOLERANCE: f64 = 1e-4;
@@ -393,13 +517,22 @@ fn scaled(input: &ModelInput, scale: f64) -> ModelInput {
         held.per_reader_per_ms *= scale;
         held.longest_per_ms *= scale;
     }
+    for cost in scaled.direct_at_ingestion_cost.iter_mut().flatten() {
+        *cost *= scale;
+    }
     for (_, bound) in &mut scaled.extra_bounds {
         *bound *= scale;
     }
+    if let Some(cutoff) = &mut scaled.objective_cutoff {
+        *cutoff *= scale;
+    }
+    for cut in &mut scaled.cuts.cuts {
+        cut.bound *= scale;
+    }
     scaled
 }
 
-fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
+fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, SolveError> {
     let mut vars = ProblemVariables::new();
     let n = input.group_count;
     let relaxed = input.relaxed;
@@ -561,9 +694,44 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
     // Each choice's direct charge per schedule class of its query
     // consumers: one node serves the class's queries.
     let mut class_charges: BTreeMap<(usize, usize), Vec<Expression>> = BTreeMap::new();
+    // Reaches of a part that cannot run at ingestion time by an
+    // ingestion-time slot: fixed at 0.
+    let mut not_at_ingestion: Vec<Variable> = Vec::new();
+    // Per ingestion-time slot node, each slot consumer's charges for the
+    // parts it reaches.
+    let mut ingestion_charges: BTreeMap<String, BTreeMap<(usize, usize), Expression>> =
+        BTreeMap::new();
     // Each cross-group node's charges per consumer, and their largest sum.
     let mut cross_per_consumer: BTreeMap<u64, (Expression, f64)> = BTreeMap::new();
     for (&(k, c), &v) in &y {
+        // Under an ingestion-time slot the direct part runs at ingestion
+        // time, once per second, as Stage 3 prices it there: none of the
+        // query-time charges below apply. Slots that may be one shared
+        // node (charged once, `shared_slots` below) read one ingestion
+        // copy of their inputs: charged at the dearest of their
+        // consumers' charges. A part that cannot run at ingestion time is
+        // not reached so.
+        if let Consumer::Slot { choice, slot } = k {
+            let s = &input.choices[choice].coefficients.slots[slot];
+            if matches!(s.materialization, Materialization::MatAtIngestion(_)) {
+                let Some(cost) = input.direct_at_ingestion_cost.get(c).copied().flatten() else {
+                    not_at_ingestion.push(v);
+                    continue;
+                };
+                let node = match s.shared_site {
+                    Some((site, threshold)) => {
+                        format!("{site} {threshold} {:?}", s.materialization)
+                    }
+                    None => format!("choice {choice} slot {slot}"),
+                };
+                *ingestion_charges
+                    .entry(node)
+                    .or_default()
+                    .entry((choice, slot))
+                    .or_insert_with(|| Expression::from(0.0)) += cost * v;
+                continue;
+            }
+        }
         let total = input.choices[c].coefficients.direct.cost_per_unit_rate;
         let rate_or_slowest = match rate(k) {
             Some(r) => r,
@@ -580,35 +748,61 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             entry.0 += rate_or_slowest * cost * v;
             entry.1 += rate_or_slowest * cost;
         }
+        // A scan's charge at its consumer's rate: for a consumer whose rate
+        // is the reaching queries' union rate (a node built once or a
+        // slot), the slowest class's rate plus the excess of the union
+        // its rate selector picks.
+        let scans_here = !input.choices[c].coefficients.absorbed_scans.is_empty()
+            || (total != 0.0 && input.choices[c].scan_source.is_some());
+        let mut above_slowest: Vec<(f64, Variable)> = Vec::new();
+        if scans_here && rate(k).is_none() {
+            let g = match k {
+                Consumer::Group(g) => g.0 as usize,
+                Consumer::Slot { choice, .. } => input.choices[choice].group.0 as usize,
+                Consumer::Query(_) => unreachable!("a query's rate is constant"),
+            };
+            for (set, selector) in &rate_of[&g] {
+                let excess = input.union_rate[set] - rate_or_slowest;
+                if excess > 0.0 {
+                    above_slowest.push((excess, product(&mut vars, *selector, v)));
+                }
+            }
+        }
+        let at_rate = |cost: f64| -> Expression {
+            let mut charge = rate_or_slowest * cost * v;
+            for &(excess, z) in &above_slowest {
+                charge += excess * cost * z;
+            }
+            charge
+        };
         for &(source, cost, pass) in &input.choices[c].coefficients.absorbed_scans {
             let pass = pass.min(cost);
             *pass_charges
                 .entry((source, cost.to_bits()))
                 .or_default()
                 .entry((c, 0))
-                .or_insert_with(|| Expression::from(0.0)) += rate_or_slowest * pass * v;
+                .or_insert_with(|| Expression::from(0.0)) += at_rate(pass);
             *scan_charges
                 .entry(source)
                 .or_default()
                 .entry((c, 0))
-                .or_insert_with(|| Expression::from(0.0)) += rate_or_slowest * (cost - pass) * v;
+                .or_insert_with(|| Expression::from(0.0)) += at_rate(cost - pass);
         }
         if total == 0.0 {
             continue;
         }
         if let Some(source) = input.choices[c].scan_source {
-            let rate = rate_or_slowest;
             let pass = input.choices[c].scan_project_cost_per_unit_rate.min(total);
             *pass_charges
                 .entry((source, total.to_bits()))
                 .or_default()
                 .entry((c, 0))
-                .or_insert_with(|| Expression::from(0.0)) += rate * pass * v;
+                .or_insert_with(|| Expression::from(0.0)) += at_rate(pass);
             *scan_charges
                 .entry(source)
                 .or_default()
                 .entry((c, 0))
-                .or_insert_with(|| Expression::from(0.0)) += rate * (total - pass) * v;
+                .or_insert_with(|| Expression::from(0.0)) += at_rate(total - pass);
             continue;
         }
         let shared = shared_cost[c].min(total);
@@ -654,7 +848,31 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
     // Query-time slots' scans, at the slot's group's slowest rate, each
     // slot's apart from its choice's other parts.
     for (c, choice) in input.choices.iter().enumerate() {
-        let slowest = slowest_rate(choice.group.0 as usize);
+        let g = choice.group.0 as usize;
+        let slowest = slowest_rate(g);
+        let with_scans = choice.coefficients.slots.iter().any(|slot| {
+            !slot.scans.is_empty()
+                && !matches!(slot.materialization, Materialization::MatAtIngestion(_))
+        });
+        if !with_scans {
+            continue;
+        }
+        // At the slot's rate: the slowest class's, plus the excess of the
+        // union its group's rate selector picks.
+        let mut above_slowest: Vec<(f64, Variable)> = Vec::new();
+        for (set, selector) in &rate_of[&g] {
+            let excess = input.union_rate[set] - slowest;
+            if excess > 0.0 {
+                above_slowest.push((excess, product(&mut vars, *selector, x[c])));
+            }
+        }
+        let at_rate = |cost: f64| -> Expression {
+            let mut charge = slowest * cost * x[c];
+            for &(excess, z) in &above_slowest {
+                charge += excess * cost * z;
+            }
+            charge
+        };
         for (i, slot) in choice.coefficients.slots.iter().enumerate() {
             if matches!(slot.materialization, Materialization::MatAtIngestion(_)) {
                 continue;
@@ -665,12 +883,12 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
                     .entry((source, cost.to_bits()))
                     .or_default()
                     .entry((c, i + 1))
-                    .or_insert_with(|| Expression::from(0.0)) += slowest * pass * x[c];
+                    .or_insert_with(|| Expression::from(0.0)) += at_rate(pass);
                 *scan_charges
                     .entry(source)
                     .or_default()
                     .entry((c, i + 1))
-                    .or_insert_with(|| Expression::from(0.0)) += slowest * (cost - pass) * x[c];
+                    .or_insert_with(|| Expression::from(0.0)) += at_rate(cost - pass);
             }
         }
     }
@@ -680,9 +898,26 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         objective += m;
         scan_maxima.push((m, charges));
     }
-    for charges in scan_charges.into_values().chain(pass_charges.into_values()) {
+    // Each source's scan charges, for its term cuts.
+    let mut source_scan_charges: BTreeMap<u64, Vec<Variable>> = BTreeMap::new();
+    let sources = scan_charges
+        .keys()
+        .copied()
+        .chain(pass_charges.keys().map(|&(source, _)| source))
+        .map(Some)
+        .chain(ingestion_charges.keys().map(|_| None))
+        .collect::<Vec<Option<u64>>>();
+    for (source, charges) in sources.into_iter().zip(
+        scan_charges
+            .into_values()
+            .chain(pass_charges.into_values())
+            .chain(ingestion_charges.into_values()),
+    ) {
         let m = vars.add(variable().min(0));
         objective += m;
+        if let Some(source) = source {
+            source_scan_charges.entry(source).or_default().push(m);
+        }
         scan_maxima.push((m, charges.into_values().collect()));
     }
     for (c, choice) in input.choices.iter().enumerate() {
@@ -746,6 +981,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         }
     }
     let mut cross_vars: Vec<CrossNodeVars> = Vec::new();
+    let mut cross_charged: BTreeMap<u64, Variable> = BTreeMap::new();
     for (key, members) in cross_once {
         let (per_consumer, per_consumer_max) = cross_per_consumer
             .remove(&key)
@@ -755,6 +991,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         let charged = vars.add(variable().min(0));
         let u = vars.add(binary());
         objective += charged;
+        cross_charged.insert(key, charged);
         cross_vars.push((once, charged, u, members, per_consumer, big_m));
     }
     let mut shared_slot_vars: Vec<(Variable, Vec<(usize, f64)>)> = Vec::new();
@@ -809,10 +1046,11 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             at_query_time.insert((q, g), vars.add(variable().min(0).max(1)));
         }
     }
-    // A choice's query-time work, without its nodes another group's
-    // choice may share (shareable nodes and shared slots): Stage 3 counts a
-    // shared node once per query, so they are left out to keep latency a
-    // lower bound.
+    // A choice's query-time work as (node key, work): each node a query
+    // reaching the choice computes ([`ChoiceCoefficients::query_work`]);
+    // without that list, its parts' total under a key of its own, its
+    // shareable nodes and shared slots left out (Stage 3 counts a shared
+    // node once per query).
     let query_time_work = |c: usize| -> f64 {
         let coefficients = &input.choices[c].coefficients;
         let shared: f64 = coefficients.shared_nodes.iter().map(|n| n.work_ms).sum();
@@ -827,18 +1065,49 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
                 .map(|s| s.work_ms)
                 .sum::<f64>()
     };
+    let work_of = |c: usize| -> Vec<(u64, f64)> {
+        let listed = &input.choices[c].coefficients.query_work;
+        match listed.is_empty() {
+            false => listed.clone(),
+            // A key no node key takes: node keys are 32-bit, scan keys
+            // have the top bit.
+            true => vec![((1 << 62) | c as u64, query_time_work(c))],
+        }
+    };
+    // Per node key, the least and the most work any choice gives it: a
+    // query's latency counts each key once, at least at the least.
+    let mut key_work: HashMap<u64, (f64, f64)> = HashMap::new();
+    for c in 0..input.choices.len() {
+        for (key, work) in work_of(c) {
+            let entry = key_work.entry(key).or_insert((work, work));
+            entry.0 = entry.0.min(work);
+            entry.1 = entry.1.max(work);
+        }
+    }
     let held_of = |c: usize| input.held.get(c).copied().unwrap_or_default();
     let mut work_reached: HashMap<(usize, usize), Variable> = HashMap::new();
     for (c, choice) in input.choices.iter().enumerate() {
         let held = held_of(c);
-        if query_time_work(c) <= 0.0 && held.per_reader_per_ms <= 0.0 && held.longest_per_ms <= 0.0
-        {
+        let works = work_of(c).iter().any(|&(_, work)| work > 0.0);
+        if !works && held.per_reader_per_ms <= 0.0 && held.longest_per_ms <= 0.0 {
             continue;
         }
         for &q in &input.group_queries[choice.group.0 as usize] {
             work_reached.insert((q, c), vars.add(variable().min(0).max(1)));
         }
     }
+    // `key_reached[q, key]`: query `q` computes the node of `key` at
+    // query time, through any choice it reaches.
+    let mut key_reached: HashMap<(usize, u64), Variable> = HashMap::new();
+    for &(q, c) in work_reached.keys() {
+        for (key, work) in work_of(c) {
+            if work > 0.0 {
+                key_reached
+                    .entry((q, key))
+                    .or_insert_with(|| vars.add(variable().min(0).max(1)));
+            }
+        }
+    }
     // Held memory per (choice, query) for a reader's own copy, and per
     // choice for one copy over the longest reader.
     let mut held_per_reader: HashMap<(usize, usize), Variable> = HashMap::new();
@@ -856,7 +1125,42 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             held_longest.insert(c, h);
         }
     }
+    // The held memory of builds several groups share, charged only
+    // through their term cuts.
+    let mut shared_build_held: BTreeMap<(u64, usize), Variable> = BTreeMap::new();
+    let mut ingestion_copy: BTreeMap<usize, (Variable, Variable)> = BTreeMap::new();
+    for cut in &input.cuts.cuts {
+        match cut.term {
+            Term::SharedBuildHeld(key) => {
+                shared_build_held.entry(key).or_insert_with(|| {
+                    let h = vars.add(variable().min(0));
+                    objective += h;
+                    h
+                });
+            }
+            Term::IngestionCopy(g) => {
+                ingestion_copy.entry(g).or_insert_with(|| {
+                    let v = vars.add(variable().min(0));
+                    objective += v;
+                    // Whether a chosen choice's ingestion-time slot reads
+                    // the group: at least each such choice's `x`.
+                    (v, vars.add(variable().min(0).max(1)))
+                });
+            }
+            _ => {}
+        }
+    }
     let vars_exclusive = (vars.add(binary()), vars.add(binary()));
+    // One variable per (site, threshold): the site read at that threshold.
+    let site_threshold_vars: BTreeMap<(u64, usize), Variable> = input
+        .site_thresholds
+        .iter()
+        .flatten()
+        .copied()
+        .collect::<BTreeSet<_>>()
+        .into_iter()
+        .map(|key| (key, vars.add(binary())))
+        .collect();
     let source_count = input
         .raw_scans
         .iter()
@@ -869,9 +1173,27 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
     for &v in &kept_on_disk {
         objective += v;
     }
+    // A query's reach passing into each choice of a group it reaches
+    // (at all, and at query time): see the reach constraints.
+    let mut flow_all: HashMap<(usize, usize), Variable> = HashMap::new();
+    let mut flow_query_time: HashMap<(usize, usize), Variable> = HashMap::new();
+    for (c, choice) in input.choices.iter().enumerate() {
+        for &q in &input.group_queries[choice.group.0 as usize] {
+            flow_all.insert((q, c), vars.add(variable().min(0).max(1)));
+            flow_query_time.insert((q, c), vars.add(variable().min(0).max(1)));
+        }
+    }
     let variable_count = vars.len() as u64;
+    let column_of: HashMap<Variable, usize> = vars
+        .iter_variables_with_def()
+        .enumerate()
+        .map(|(column, (variable, _))| (variable, column))
+        .collect();
     let model = vars.minimise(objective.clone()).using(highs);
     let mut constraints = Vec::new();
+    for &v in &not_at_ingestion {
+        constraints.push(constraint!(v == 0));
+    }
     // Needed groups: every query's roots; children of chosen choices.
     for query in &input.queries {
         for root in &query.roots {
@@ -894,8 +1216,65 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             constraints.push(constraint!(need[child.0 as usize] >= x[c]));
         }
     }
+    // A group is needed only as a root or as a chosen choice's child: a
+    // choice nothing reads is no part of the plan. And a scan shared by
+    // one reader is that reader's own scan (read apart), the same plan.
+    let mut readers: Vec<Expression> = (0..n).map(|_| Expression::from(0.0)).collect();
+    let mut shared_scan_readers: Vec<Expression> = (0..n).map(|_| Expression::from(0.0)).collect();
+    for (c, choice) in input.choices.iter().enumerate() {
+        for child in &choice.children {
+            readers[child.0 as usize] += x[c];
+        }
+        // A fused choice needs the other groups it implements.
+        for implemented in choice.implements.iter().skip(1) {
+            readers[implemented.0 as usize] += x[c];
+        }
+        for (scan, horizon) in &choice.scan_children {
+            if horizon.is_some() {
+                shared_scan_readers[scan.0 as usize] += x[c];
+            }
+        }
+    }
+    let is_root: Vec<bool> = (0..n)
+        .map(|g| {
+            input
+                .queries
+                .iter()
+                .any(|q| q.roots.iter().any(|r| r.0 as usize == g))
+        })
+        .collect();
+    let is_scan: Vec<bool> = (0..n)
+        .map(|g| {
+            input
+                .choices
+                .iter()
+                .any(|c| c.group.0 as usize == g && c.scan_source.is_some())
+        })
+        .collect();
+    for (g, read) in readers.into_iter().enumerate() {
+        if !is_root[g] {
+            constraints.push(constraint!(need[g] <= read));
+        }
+    }
+    for (g, read) in shared_scan_readers.into_iter().enumerate() {
+        if is_scan[g] && !is_root[g] {
+            constraints.push(constraint!(read >= 2.0 * need[g]));
+        }
+    }
+    let mut choices_of: Vec<Vec<usize>> = vec![Vec::new(); n];
+    for (c, choice) in input.choices.iter().enumerate() {
+        choices_of[choice.group.0 as usize].push(c);
+    }
     for (&g, &v) in &once {
         constraints.push(constraint!(v <= need[g.0 as usize]));
+        // Built once through its own choice: a group a fused choice of
+        // another group implements has no node of its own to build once
+        // (solution decoding builds the fused node).
+        let mut own = Expression::from(0.0);
+        for &c in &choices_of[g.0 as usize] {
+            own += x[c];
+        }
+        constraints.push(constraint!(v <= own));
         let slot = input.group_slots.iter().find(|s| s.group == g).unwrap();
         if !slot.may_recompute {
             constraints.push(constraint!(v == need[g.0 as usize]));
@@ -910,29 +1289,47 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             }
         }
     }
-    // Aggregated over a group's choices (one is chosen when it is
-    // needed), which keeps the LP relaxation tight.
-    let mut choices_of: Vec<Vec<usize>> = vec![Vec::new(); n];
+    // The choices of other groups that implement each group (fused).
+    let mut fused_into: Vec<Vec<usize>> = vec![Vec::new(); n];
     for (c, choice) in input.choices.iter().enumerate() {
-        choices_of[choice.group.0 as usize].push(c);
+        for i in &choice.implements {
+            if *i != choice.group {
+                fused_into[i.0 as usize].push(c);
+            }
+        }
     }
+    // As a flow, which keeps the LP relaxation tight: a query reaching a
+    // group passes into its chosen choice (at most that choice's `x`, at
+    // least the reach in total, unless a fused choice implements the
+    // group), and reaches each child of the choices it passes into.
     for (g, of_g) in choices_of.iter().enumerate() {
-        let mut children: BTreeSet<usize> = BTreeSet::new();
-        for &c in of_g {
-            children.extend(input.choices[c].children.iter().map(|h| h.0 as usize));
-        }
-        for h in children {
-            let mut with_child = Expression::from(0.0);
+        for &q in &input.group_queries[g] {
+            let Some(&above) = reaches_at_all.get(&(q, g)) else {
+                continue;
+            };
+            let mut into = Expression::from(0.0);
+            let mut below: BTreeMap<usize, Expression> = BTreeMap::new();
             for &c in of_g {
-                if input.choices[c].children.iter().any(|k| k.0 as usize == h) {
-                    with_child += x[c];
+                let w = flow_all[&(q, c)];
+                constraints.push(constraint!(w <= x[c]));
+                into += w;
+                // A child listed twice (identical subtrees) is one flow.
+                let children: BTreeSet<usize> = input.choices[c]
+                    .children
+                    .iter()
+                    .map(|h| h.0 as usize)
+                    .collect();
+                for h in children {
+                    *below.entry(h).or_insert_with(|| Expression::from(0.0)) += w;
                 }
             }
-            for &q in &input.group_queries[g] {
-                if let (Some(&above), Some(&below)) =
-                    (reaches_at_all.get(&(q, g)), reaches_at_all.get(&(q, h)))
-                {
-                    constraints.push(constraint!(below >= above + with_child.clone() - 1));
+            for &c in &fused_into[g] {
+                into += x[c];
+            }
+            constraints.push(constraint!(into >= above));
+            for (h, flow) in below {
+                if let Some(&reached) = reaches_at_all.get(&(q, h)) {
+                    constraints.push(constraint!(reached >= flow));
                 }
             }
         }
@@ -992,6 +1389,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
     for (&(k, c), &v) in &y {
         let g = input.choices[c].group.0 as usize;
         constraints.push(constraint!(v >= passes[&(k, g)] + x[c] - 1));
+        constraints.push(constraint!(v <= x[c]));
     }
     for (&(k, g), &p) in &passes {
         let mut into = Expression::from(0.0);
@@ -1117,38 +1515,47 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             }
         }
     }
+    // As a flow, as the reach at all.
     for (g, of_g) in choices_of.iter().enumerate() {
-        let mut children: BTreeSet<usize> = BTreeSet::new();
-        for &c in of_g {
-            children.extend(input.choices[c].children.iter().map(|h| h.0 as usize));
-        }
-        for h in children {
-            let mut with_child = Expression::from(0.0);
+        for &q in &input.group_queries[g] {
+            let Some(&above) = at_query_time.get(&(q, g)) else {
+                continue;
+            };
+            let mut into = Expression::from(0.0);
+            let mut below: BTreeMap<usize, Expression> = BTreeMap::new();
             for &c in of_g {
+                let w = flow_query_time[&(q, c)];
+                constraints.push(constraint!(w <= x[c]));
+                into += w;
+                if let Some(&reached) = work_reached.get(&(q, c)) {
+                    constraints.push(constraint!(reached >= w));
+                }
                 let choice = &input.choices[c];
-                let at_query = choice
-                    .children
-                    .iter()
-                    .zip(&choice.coefficients.child_via)
-                    .any(|(child, via)| {
-                        child.0 as usize == h
-                            && match via {
-                                Via::Direct => true,
-                                Via::Slot(i) => !matches!(
-                                    choice.coefficients.slots[*i].materialization,
-                                    Materialization::MatAtIngestion(_)
-                                ),
-                            }
-                    });
-                if at_query {
-                    with_child += x[c];
+                // A child listed twice (identical subtrees) is one flow.
+                let mut at_query_children = BTreeSet::new();
+                for (h, via) in choice.children.iter().zip(&choice.coefficients.child_via) {
+                    let at_query = match via {
+                        Via::Direct => true,
+                        Via::Slot(i) => !matches!(
+                            choice.coefficients.slots[*i].materialization,
+                            Materialization::MatAtIngestion(_)
+                        ),
+                    };
+                    if at_query {
+                        at_query_children.insert(h.0 as usize);
+                    }
+                }
+                for h in at_query_children {
+                    *below.entry(h).or_insert_with(|| Expression::from(0.0)) += w;
                 }
             }
-            for &q in &input.group_queries[g] {
-                if let (Some(&above), Some(&below)) =
-                    (at_query_time.get(&(q, g)), at_query_time.get(&(q, h)))
-                {
-                    constraints.push(constraint!(below >= above + with_child.clone() - 1));
+            for &c in &fused_into[g] {
+                into += x[c];
+            }
+            constraints.push(constraint!(into >= above));
+            for (h, flow) in below {
+                if let Some(&reached) = at_query_time.get(&(q, h)) {
+                    constraints.push(constraint!(reached >= flow));
                 }
             }
         }
@@ -1156,14 +1563,20 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
     for (&(q, c), &w) in &work_reached {
         let g = input.choices[c].group.0 as usize;
         constraints.push(constraint!(w >= at_query_time[&(q, g)] + x[c] - 1));
+        for (key, work) in work_of(c) {
+            if work > 0.0 {
+                constraints.push(constraint!(key_reached[&(q, key)] >= w));
+            }
+        }
     }
-    // Each query's latency: the query-time work of the choices it reaches
-    // (a lower bound: compaction is left out).
+    // Each query's latency: the work of each node it computes at query
+    // time, counted once per node key at the least work any choice gives
+    // it (a lower bound: compaction is left out).
     let latency_of = |q: usize| {
         let mut latency = Expression::from(0.0);
-        for (&(r, c), &w) in &work_reached {
+        for (&(r, key), &v) in &key_reached {
             if r == q {
-                latency += query_time_work(c) * w;
+                latency += key_work[&key].0 * v;
             }
         }
         latency
@@ -1177,23 +1590,62 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
     // L the longest the query's latency can be; exact when w is 0 or 1.
     let longest_latency: Vec<f64> = (0..input.queries.len())
         .map(|q| {
-            choices_of
-                .iter()
-                .enumerate()
-                .filter(|(g, _)| input.group_queries[*g].contains(&q))
-                .map(|(_, of_g)| of_g.iter().map(|&c| query_time_work(c)).fold(0.0, f64::max))
+            key_reached
+                .keys()
+                .filter(|&&(r, _)| r == q)
+                .map(|(_, key)| key_work[key].1)
                 .sum()
         })
         .collect();
     // A part read by several queries runs at their classes' union rate,
     // not their rates' sum: per reader, its held memory is charged at
     // `union / sum` of that reader's rate (Stage 3's rate-weighted mean
-    // latency at the part's rate). The least ratio over the group's class
-    // sets bounds it from below.
-    let per_reader_share: Vec<f64> = input
-        .group_queries
-        .iter()
-        .map(|queries_of| {
+    // latency at the part's rate). Its readers may be any query of a group
+    // with a choice computing one of its nodes (a node several groups
+    // share is one node); the least ratio over those queries' class sets
+    // bounds it from below.
+    let mut groups_of_key: HashMap<u64, BTreeSet<usize>> = HashMap::new();
+    for (c, choice) in input.choices.iter().enumerate() {
+        for (key, _) in work_of(c) {
+            groups_of_key
+                .entry(key)
+                .or_default()
+                .insert(choice.group.0 as usize);
+        }
+    }
+    // Readers from several groups may be in classes no `union_rate` key
+    // joins; their union rate is still at least their fastest's rate. So
+    // `fastest / sum` bounds the share from below for each possible
+    // fastest reader, the sum over the readers no faster than it.
+    let fastest_over_sum = |queries_of: &BTreeSet<usize>| -> Option<f64> {
+        queries_of
+            .iter()
+            .map(|&q| {
+                let fastest = input.queries[q].rate;
+                let sum: f64 = queries_of
+                    .iter()
+                    .map(|&r| input.queries[r].rate)
+                    .filter(|&rate| rate <= fastest)
+                    .sum();
+                match sum > 0.0 {
+                    true => (fastest / sum).min(1.0),
+                    false => 1.0,
+                }
+            })
+            .reduce(f64::min)
+    };
+    let per_reader_share: Vec<f64> = (0..input.choices.len())
+        .map(|c| {
+            let mut queries_of: BTreeSet<usize> = input.group_queries
+                [input.choices[c].group.0 as usize]
+                .iter()
+                .copied()
+                .collect();
+            for (key, _) in work_of(c) {
+                for &g in &groups_of_key[&key] {
+                    queries_of.extend(input.group_queries[g].iter().copied());
+                }
+            }
             input
                 .union_rate
                 .iter()
@@ -1208,6 +1660,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
                         .all(|class| queries_of.iter().any(|&q| input.schedule[q] == *class));
                     (covered && rate_sum > 0.0).then_some((union / rate_sum).min(1.0))
                 })
+                .chain(fastest_over_sum(&queries_of))
                 .fold(1.0, f64::min)
         })
         .collect();
@@ -1215,7 +1668,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         let held = held_of(c);
         let longest = longest_latency[q];
         if let Some(&h) = held_per_reader.get(&(q, c)) {
-            let share = per_reader_share[input.choices[c].group.0 as usize];
+            let share = per_reader_share[c];
             let coefficient = held.per_reader_per_ms * input.queries[q].rate * share;
             constraints.push(constraint!(
                 h >= coefficient * (latency_of(q) - longest * (1 - w))
@@ -1252,7 +1705,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         }
     }
     // Cuts.
-    for (set, outside, bound) in &input.excluded {
+    for (i, (set, outside, bound)) in input.excluded.iter().enumerate() {
         let mut chosen = Expression::from(0.0);
         for &c in set {
             chosen += x[c];
@@ -1260,8 +1713,137 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         for &g in outside {
             chosen -= need[g];
         }
+        if let Some((built, not_built)) = input.excluded_once.get(i) {
+            for g in built {
+                if let Some(&v) = once.get(g) {
+                    chosen += v;
+                }
+            }
+            for g in not_built {
+                if let Some(&v) = once.get(g) {
+                    chosen -= v;
+                }
+            }
+        }
         constraints.push(constraint!(chosen <= *bound));
     }
+    // Each coupling site at one threshold: `x[c] ≤ u[site, threshold of c]`,
+    // at most one threshold per site.
+    let mut site_choice: BTreeMap<(u64, usize), Vec<usize>> = BTreeMap::new();
+    for (c, sites) in input.site_thresholds.iter().enumerate() {
+        for &(site, threshold) in sites {
+            site_choice.entry((site, threshold)).or_default().push(c);
+        }
+    }
+    let mut per_site: BTreeMap<u64, Expression> = BTreeMap::new();
+    for (key @ (site, _), readers) in &site_choice {
+        let u = site_threshold_vars[key];
+        *per_site.entry(*site).or_default() += u;
+        for &c in readers {
+            constraints.push(constraint!(x[c] <= u));
+        }
+    }
+    for (_, any) in per_site {
+        constraints.push(constraint!(any <= 1));
+    }
+    // A build built once only with another group's identical build: some
+    // other group's choice builds it too.
+    let mut builders: HashMap<u32, Vec<usize>> = HashMap::new();
+    for (c, choice) in input.choices.iter().enumerate() {
+        for &key in &choice.built_once_keys {
+            builders.entry(key).or_default().push(c);
+        }
+    }
+    for (c, choice) in input.choices.iter().enumerate() {
+        for key in &choice.built_once_with_other_group {
+            let mut others = Expression::from(0.0);
+            for &o in &builders[key] {
+                if input.choices[o].group != choice.group {
+                    others += x[o];
+                }
+            }
+            constraints.push(constraint!(x[c] <= others));
+        }
+    }
+    for (&g, &(_, read)) in &ingestion_copy {
+        for (p, choice) in input.choices.iter().enumerate() {
+            let reads = choice
+                .children
+                .iter()
+                .zip(&choice.coefficients.child_via)
+                .any(|(child, via)| {
+                    child.0 as usize == g
+                        && matches!(via, Via::Slot(i) if matches!(
+                            choice.coefficients.slots[*i].materialization,
+                            Materialization::MatAtIngestion(_)
+                        ))
+                });
+            if reads {
+                constraints.push(constraint!(read >= x[p]));
+            }
+        }
+    }
+    // Term cuts: when all of a cut's choices are chosen, its term's
+    // charge is at least its bound.
+    for cut in &input.cuts.cuts {
+        // A cut whose term or condition has no variable in this model is
+        // left out (with a charge of 0 it would exclude its plans).
+        let mut charge = Expression::from(0.0);
+        let charged = match cut.term {
+            Term::SourceScans(source) => {
+                let mut any = false;
+                for &m in source_scan_charges.get(&source).into_iter().flatten() {
+                    charge += m;
+                    any = true;
+                }
+                for key in &cut.nodes {
+                    if let Some(&m) = cross_charged.get(key) {
+                        charge += m;
+                        any = true;
+                    }
+                }
+                any
+            }
+            Term::HeldPerReader(q, c) => {
+                held_per_reader.get(&(q, c)).map(|&h| charge += h).is_some()
+            }
+            Term::HeldLongest(c) => held_longest.get(&c).map(|&h| charge += h).is_some(),
+            Term::SharedBuildHeld(key) => {
+                charge += shared_build_held[&key];
+                true
+            }
+            Term::IngestionCopy(g) => {
+                charge += ingestion_copy[&g].0;
+                true
+            }
+        };
+        let conditioned = cut.reached.iter().all(|qc| work_reached.contains_key(qc))
+            && cut
+                .built_once
+                .iter()
+                .all(|g| once.contains_key(&GroupId(*g as u32)));
+        if !charged || !conditioned {
+            continue;
+        }
+        let mut chosen = Expression::from(0.0);
+        for &c in &cut.choices {
+            chosen += x[c];
+        }
+        for qc in &cut.reached {
+            chosen += work_reached[qc];
+        }
+        for g in &cut.built_once {
+            chosen += once[&GroupId(*g as u32)];
+        }
+        let mut conditions = (cut.choices.len() + cut.reached.len() + cut.built_once.len()) as f64;
+        if let Term::IngestionCopy(g) = cut.term {
+            chosen += ingestion_copy[&g].1;
+            conditions += 1.0;
+        }
+        constraints.push(constraint!(
+            charge >= cut.bound * (chosen - (conditions - 1.0))
+        ));
+    }
     for ((set, bound), &e) in input.extra_bounds.iter().zip(&extra) {
         let mut chosen = Expression::from(0.0);
         for &c in set {
@@ -1278,6 +1860,62 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         }
         constraints.push(constraint!(x[*c] <= any));
     }
+    // A consumer reading a child group through a choice with such a
+    // requirement passes into one of the required choices (unless the
+    // child is built once, whose own consumer then does): the flow, not
+    // only the choice, meets the requirement.
+    for (c, allowed) in &input.requires {
+        let Some(h) = allowed.first().map(|&a| input.choices[a].group) else {
+            continue;
+        };
+        let choice = &input.choices[*c];
+        let Some(via) = choice
+            .children
+            .iter()
+            .zip(&choice.coefficients.child_via)
+            .find(|(child, _)| **child == h)
+            .map(|(_, via)| *via)
+        else {
+            continue;
+        };
+        let into_allowed = |k: Consumer| -> Option<Expression> {
+            let mut sum = Expression::from(0.0);
+            for &a in allowed {
+                sum += *y.get(&(k, a))?;
+            }
+            Some(sum)
+        };
+        let built_once = once.get(&h).copied();
+        let mut readers: Vec<(Consumer, Expression)> = Vec::new();
+        match via {
+            Via::Direct => {
+                for (&(k, d), &v) in &y {
+                    if d == *c {
+                        readers.push((k, Expression::from(v)));
+                    }
+                }
+            }
+            Via::Slot(i) => readers.push((
+                Consumer::Slot {
+                    choice: *c,
+                    slot: i,
+                },
+                Expression::from(x[*c]),
+            )),
+        }
+        for (k, passing) in readers {
+            let Some(sum) = into_allowed(k) else { continue };
+            match built_once {
+                Some(o) => constraints.push(constraint!(sum >= passing - o)),
+                None => constraints.push(constraint!(sum >= passing)),
+            }
+        }
+        if let Some(o) = built_once {
+            if let Some(sum) = into_allowed(Consumer::Group(h)) {
+                constraints.push(constraint!(sum >= o + x[*c] - 1));
+            }
+        }
+    }
     for (all_of, one_of) in &input.rollup_rules {
         let mut chosen = Expression::from(0.0);
         for set in all_of {
@@ -1338,20 +1976,83 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
         }
     }
     let constraint_count = constraints.len() as u64;
-    let mut model = model
+    let mut highs_model = model
         .with_all(constraints)
-        .with_mip_gap(0.0)
-        .map_err(|e| format!("MIP gap: {e:?}"))?;
-    // No absolute gap either: prices are scaled to about 1.
-    model = model.set_option("mip_abs_gap", 0.0);
+        .try_into_inner()
+        .map_err(|e| SolveError::Failed(e.to_string()))?;
+    // No relative or absolute gap: prices are scaled to about 1.
+    highs_model.set_option("mip_rel_gap", 0.0);
+    highs_model.set_option("mip_abs_gap", 0.0);
     if let Some(limit) = input.time_limit_s {
-        model = model.set_time_limit(limit.max(0.001));
+        highs_model.set_option("time_limit", limit.max(0.001));
+    }
+    let cutoff = input
+        .objective_cutoff
+        .filter(|_| !relaxed)
+        .map(|cutoff| cutoff * (1.0 - CUTOFF_TOLERANCE));
+    if let Some(cutoff) = cutoff {
+        // HiGHS leaves out the objective's constant.
+        highs_model.set_option(
+            "objective_bound",
+            cutoff - good_lp::IntoAffineExpression::constant(&objective),
+        );
+    }
+    let solved = highs_model
+        .try_solve()
+        .map_err(|e| SolveError::Failed(format!("HiGHS error while solving model: {e:?}")))?;
+    use highs::HighsModelStatus as Status;
+    let stopped_early = match solved.status() {
+        Status::Infeasible | Status::UnboundedOrInfeasible => return Ok(None),
+        Status::Optimal | Status::ObjectiveBound | Status::ObjectiveTarget => false,
+        Status::ReachedTimeLimit
+        | Status::ReachedSolutionLimit
+        | Status::ReachedInterrupt
+        | Status::ReachedIterationLimit
+        | Status::ReachedMemoryLimit => true,
+        status => {
+            return Err(SolveError::Failed(format!(
+                "HiGHS stopped with status {status:?}"
+            )))
+        }
+    };
+    // The solver's proven lower bound when it stopped early (the
+    // relaxation's optimum in a relaxed solve): every plan costs at least
+    // this much, solution found or not.
+    let proven_bound = || {
+        if relaxed {
+            None
+        } else {
+            solved
+                .double_info_value(c"mip_dual_bound")
+                .ok()
+                .filter(|bound| bound.is_finite())
+                // HiGHS leaves out the objective's constant; lowered by
+                // the solver's tolerance, as a gap-limit objective is.
+                .map(|bound| bound + good_lp::IntoAffineExpression::constant(&objective))
+                .map(|bound| bound - bound.abs() * GAP_TOLERANCE)
+        }
+    };
+    if solved.primal_solution_status() != highs::HighsSolutionStatus::Feasible {
+        return Err(SolveError::StoppedWithoutSolution {
+            proven_bound: proven_bound(),
+        });
     }
-    let solution = match model.solve() {
-        Ok(solution) => solution,
-        Err(good_lp::ResolutionError::Infeasible) => return Ok(None),
-        Err(error) => return Err(error.to_string()),
+    let gap = solved.mip_gap();
+    let status = match stopped_early {
+        true => good_lp::solvers::SolutionStatus::TimeLimit,
+        false if gap.is_finite() && gap > 0.0 => good_lp::solvers::SolutionStatus::GapLimit,
+        false => good_lp::solvers::SolutionStatus::Optimal,
+    };
+    let solution = Columns {
+        status,
+        column_of,
+        values: solved.get_solution().columns().to_vec(),
     };
+    // Nothing below the cutoff: the solver may still report a plan above
+    // it.
+    if cutoff.is_some_and(|cutoff| !stopped_early && solution.eval(&objective) >= cutoff) {
+        return Ok(None);
+    }
     let chosen: Vec<usize> = x
         .iter()
         .enumerate()
@@ -1379,6 +2080,7 @@ fn solve_scaled(input: &ModelInput) -> Result<Option<ModelSolution>, String> {
             solution.status(),
             good_lp::solvers::SolutionStatus::Optimal | good_lp::solvers::SolutionStatus::GapLimit
         ),
+        proven_bound: if stopped_early { proven_bound() } else { None },
         variables: variable_count,
         constraints: constraint_count,
     }))
@@ -1440,15 +2142,15 @@ pub(crate) type MemoSearchResult =
 /// The cheapest valid plan of `per_computation`'s memo under cost model 1
 /// (the caller routes cost models 2 and 3 to the legacy search).
 ///
-/// Lazy loop: solve the model (a lower bound: it omits memory held while a
-/// query runs, and checks composed guarantees only on its incumbents);
-/// price the solution's group choice exactly, as the reference search
-/// prices one ([`offer_group_choice`]: every threshold, identical-expression
-/// variant and Stage 2 materialization), which updates the incumbent; add
-/// a no-good cut excluding that group choice; stop when the model's bound
-/// reaches the incumbent. Each group choice is priced at most once, so the
-/// loop ends, and the incumbent is the reference search's optimum.
-pub(crate) fn search<Id: Clone>(
+/// The memo is split into independent clusters ([`crate::search::milp::clusters::split`]):
+/// queries linked through no shared group, coupling site or charged raw
+/// horizon. Each cluster is its own workload with its own model and lazy
+/// loop ([`ClusterLoop`]), the clusters' solves running in parallel; each
+/// loop stops once its incumbent is certified. Under cost model 1 the
+/// workload's cost is the sum of the clusters' costs, so the clusters'
+/// plans together are the workload's cheapest plan
+/// ([`crate::search::milp::clusters::combine`], which checks the sum).
+pub(crate) fn search<Id: Clone + std::fmt::Debug>(
     per_computation: &PerComputationCandidates<Id>,
     entries: &[QueryWorkloadEntry],
     data: &DataWorkload,
@@ -1459,174 +2161,656 @@ pub(crate) fn search<Id: Clone>(
     let memo = crate::search::memo::brute_force::memo(per_computation, entries);
     let memo_stats =
         crate::search::memo::brute_force::stats(&memo, started.elapsed().as_secs_f64() * 1e3);
-    let mut input = model_input_of(&memo, entries, data, models)?;
-    let mut selector = PlanSelector::new(entries, data, models, false)?;
-    let mut logical: HashMap<String, Vec<QueryRoot>> = HashMap::new();
-    let mut number = 0usize;
-    let mut iterations = 0u32;
-    let mut bound = None;
-    let mut first_bound = None;
-    let mut solver_ms = 0.0;
-    let mut size = (0, 0);
-    let mut time_limit_reached = false;
-    loop {
-        if time_limit_reached {
-            break;
-        }
-        let remaining_s =
-            time_limit_ms.map(|limit| limit as f64 / 1e3 - started.elapsed().as_secs_f64());
-        if remaining_s.is_some_and(|s| s <= 0.0) {
-            time_limit_reached = true;
-            break;
+    let roots: Vec<QueryRoot> = per_computation
+        .roots
+        .iter()
+        .map(|(_, root)| root.clone())
+        .collect();
+    let clusters = crate::search::milp::clusters::split(&memo, &roots, &models);
+    if clusters.len() > 1 {
+        if let Some(found) = search_clusters(
+            per_computation,
+            &clusters,
+            entries,
+            data,
+            models,
+            started,
+            time_limit_ms,
+        ) {
+            let mut found = found?;
+            found.0.search.memo = Some(memo_stats);
+            return Ok(found);
         }
-        input.time_limit_s = remaining_s;
-        iterations += 1;
-        let solve_started = std::time::Instant::now();
-        let solved = match solve(&input) {
-            Ok(solved) => solved,
-            // The solver stopped at the time limit before finding any
-            // solution: the incumbent so far stands.
-            Err(error) if remaining_s.is_some() && error.contains("NoSolutionFound") => {
-                time_limit_reached = true;
-                solver_ms += solve_started.elapsed().as_secs_f64() * 1e3;
-                break;
-            }
-            Err(error) => return Err(SelectionError::Solver(error)),
-        };
-        solver_ms += solve_started.elapsed().as_secs_f64() * 1e3;
-        // The solver stopped at the time limit: its plan is no bound, but
-        // it is still priced, as a last round, and may become the
-        // incumbent.
+    }
+    // One cluster, or clusters whose plans do not add up: the whole memo,
+    // unless the clusters spent the time limit.
+    if clusters.len() > 1
+        && time_limit_ms.is_some_and(|limit| started.elapsed().as_millis() >= u128::from(limit))
+    {
+        return Err(SelectionError::Solver(
+            "the MILP search reached its time limit: the clusters' plans did not combine into one plan whose cost is their sum, and no time is left to search the memo as one".to_string(),
+        ));
+    }
+    let all: Vec<usize> = (0..entries.len()).collect();
+    let mut cluster = ClusterLoop::new(
+        memo,
+        per_computation,
+        &all,
+        entries,
+        data,
+        models,
+        started,
+        time_limit_ms,
+    )?;
+    run_rounds(std::slice::from_mut(&mut cluster))?;
+    let (mut selection, selected, roots, report) = cluster.finish()?;
+    selection.search.memo = Some(memo_stats);
+    if let Some(milp) = &mut selection.search.milp {
+        milp.clusters = vec![report];
+    }
+    Ok((selection, selected, roots))
+}
+
+/// [`search`] over `clusters` (each a list of the workload's queries), each
+/// searched as its own workload, their plans combined into the workload's.
+/// `None` when a cluster's Pass 1 fails or the combined plan's cost is not
+/// the sum of the clusters' (the caller then searches the memo as one).
+fn search_clusters<Id: Clone + std::fmt::Debug>(
+    per_computation: &PerComputationCandidates<Id>,
+    clusters: &[Vec<usize>],
+    entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
+    models: PlanningModels<'_>,
+    started: std::time::Instant,
+    time_limit_ms: Option<u64>,
+) -> Option<MemoSearchResult> {
+    // Each cluster's workload: its queries' roots through Pass 1, pruned
+    // for the deployment as the caller pruned the workload's.
+    let mut workloads = Vec::with_capacity(clusters.len());
+    for queries in clusters {
+        let roots: Vec<(Id, QueryRoot)> = queries
+            .iter()
+            .map(|&q| per_computation.roots[q].clone())
+            .collect();
+        let inventory =
+            asap_logical_asap_optimization::per_computation_candidates::generate_per_computation_candidates(roots)
+                .ok()?
+                .without_unsupported_summaries(models.deployment_capabilities);
+        let cluster_entries: Vec<QueryWorkloadEntry> =
+            queries.iter().map(|&q| entries[q].clone()).collect();
+        workloads.push((inventory, cluster_entries));
+    }
+    let mut loops = Vec::with_capacity(clusters.len());
+    for (queries, (inventory, cluster_entries)) in clusters.iter().zip(&workloads) {
+        let memo = crate::search::memo::brute_force::memo(inventory, cluster_entries);
+        match ClusterLoop::new(
+            memo,
+            inventory,
+            queries,
+            cluster_entries,
+            data,
+            models,
+            started,
+            time_limit_ms,
+        ) {
+            Ok(cluster) => loops.push(cluster),
+            Err(error) => return Some(Err(error)),
+        }
+    }
+    if let Err(error) = run_rounds(&mut loops) {
+        return Some(Err(error));
+    }
+    let mut found = Vec::with_capacity(loops.len());
+    for (queries, cluster) in clusters.iter().zip(loops) {
+        match cluster.finish() {
+            Ok(result) => found.push(result),
+            Err(SelectionError::NoValidCandidate(rejected)) => {
+                return Some(Err(SelectionError::NoValidCandidate(
+                    rejected
+                        .into_iter()
+                        .map(|r| crate::search::legacy::components::renumber(r, queries))
+                        .collect(),
+                )))
+            }
+            Err(error) => return Some(Err(error)),
+        }
+    }
+    crate::search::milp::clusters::combine(found, clusters, entries, data, models).map(Ok)
+}
+
+/// Solve every cluster's model until each loop stops: up to one solver
+/// thread per core, each solving one cluster's model (HiGHS keeps one
+/// state per model), while this thread runs the incumbent evaluation of
+/// the solutions that come back. A loop's next solve starts as soon as
+/// its last solution is evaluated.
+fn run_rounds<Id: Clone>(loops: &mut [ClusterLoop<'_, Id>]) -> Result<(), SelectionError> {
+    /// A solver thread's stack: HiGHS recurses in presolve.
+    const SOLVER_STACK_BYTES: usize = 256 << 20;
+    let threads = std::thread::available_parallelism()
+        .map_or(1, |n| n.get())
+        .min(loops.len())
+        .max(1);
+    std::thread::scope(|scope| {
+        let (sender, receiver) = std::sync::mpsc::channel();
+        let mut ready: std::collections::VecDeque<usize> = (0..loops.len())
+            .filter(|&k| loops[k].next_round())
+            .collect();
+        let mut solving = 0usize;
+        loop {
+            while solving < threads {
+                let Some(k) = ready.pop_front() else {
+                    break;
+                };
+                // A queued solve's time limit counts from when it starts.
+                if !loops[k].start_solve() {
+                    continue;
+                }
+                let input = loops[k]
+                    .input
+                    .take()
+                    .expect("a loop's input between solves");
+                let sender = sender.clone();
+                std::thread::Builder::new()
+                    .stack_size(SOLVER_STACK_BYTES)
+                    .spawn_scoped(scope, move || {
+                        let solve_started = std::time::Instant::now();
+                        let solved = std::panic::catch_unwind(std::panic::AssertUnwindSafe(|| {
+                            solve(&input)
+                        }))
+                        .unwrap_or_else(|_| {
+                            Err(SolveError::Failed("the solver panicked".to_string()))
+                        });
+                        // The receiver is gone only once the search failed.
+                        let _ = sender.send((k, input, solved, solve_started.elapsed()));
+                    })
+                    .map_err(|e| SelectionError::Solver(format!("no solver thread: {e}")))?;
+                solving += 1;
+            }
+            if solving == 0 {
+                return Ok(());
+            }
+            let (k, input, solved, solve_time) = receiver
+                .recv()
+                .expect("every solver thread sends its result");
+            solving -= 1;
+            loops[k].input = Some(input);
+            loops[k].solver_ms += solve_time.as_secs_f64() * 1e3;
+            loops[k].take_solution(solved)?;
+            if loops[k].next_round() {
+                ready.push_back(k);
+            }
+        }
+    })
+}
+
+/// One cluster's lazy loop (joint-planning-search design, "One solution,
+/// one plan"): solve the model, whose objective is the MILP lower bound
+/// (coefficient derivation gives lower bounds); decode the solution to the
+/// one plan it describes ([`crate::search::milp::decode`]); run incumbent
+/// evaluation (Stage 3 on that plan, and on the few variants the model has
+/// no variable for yet: its widened scans, and the materializations of the
+/// nodes no variable sets), which may improve the incumbent; cut exactly that solution's choice set; stop when the MILP
+/// lower bound reaches the incumbent cost. Each choice set is evaluated at
+/// most once, so the loop ends. The decoding gaps listed in
+/// [`crate::search::milp::decode`] are plans no solution decodes to. A
+/// choice set cut before all its plans were evaluated (too many free
+/// combinations, or a decoding failure) keeps its objective as a lower
+/// bound, so the loop never certifies past it.
+///
+/// A round is [`Self::next_round`], [`solve`] of [`Self::input`] (on a
+/// solver thread, [`run_rounds`]), then [`Self::take_solution`].
+struct ClusterLoop<'a, Id> {
+    memo: Memo,
+    per_computation: &'a PerComputationCandidates<Id>,
+    /// Its queries, by index in the whole workload; query `i` of
+    /// `entries` is `queries[i]`.
+    queries: &'a [usize],
+    entries: &'a [QueryWorkloadEntry],
+    data: &'a DataWorkload,
+    models: PlanningModels<'a>,
+    started: std::time::Instant,
+    time_limit_ms: Option<u64>,
+    /// The model; `None` while a solver thread has it.
+    input: Option<ModelInput>,
+    /// Per choice of the model, its candidate at its site thresholds,
+    /// which solution decoding composes.
+    versions: Vec<ImplementationCandidate>,
+    selector: PlanSelector<'a>,
+    logical: HashMap<String, Vec<QueryRoot>>,
+    number: usize,
+    iterations: u32,
+    bound: Option<f64>,
+    first_bound: Option<f64>,
+    solver_ms: f64,
+    size: (u64, u64),
+    time_limit_reached: bool,
+    /// The solver's proven lower bound from the solve stopped at the time
+    /// limit, on every plan not cut before it.
+    proven_bound: Option<f64>,
+    /// The loop has stopped: certified, out of plans, or out of time.
+    done: bool,
+    /// The most combinations of free nodes' options decoding evaluates
+    /// ([`crate::search::milp::decode::MAX_FREE_COMBINATIONS`]).
+    max_free_combinations: usize,
+    /// The smallest objective of a solution whose choice set was cut with
+    /// some of its plans not evaluated (too many free combinations, or a
+    /// decoding failure): the lower bound never exceeds it.
+    unevaluated_bound: Option<f64>,
+    /// Solutions solution decoding could not compose into a plan.
+    decode_failures: u32,
+}
+
+impl<'a, Id: Clone> ClusterLoop<'a, Id> {
+    #[allow(clippy::too_many_arguments)]
+    fn new(
+        memo: Memo,
+        per_computation: &'a PerComputationCandidates<Id>,
+        queries: &'a [usize],
+        entries: &'a [QueryWorkloadEntry],
+        data: &'a DataWorkload,
+        models: PlanningModels<'a>,
+        started: std::time::Instant,
+        time_limit_ms: Option<u64>,
+    ) -> Result<Self, SelectionError> {
+        let (input, versions) = model_input_of(&memo, entries, data, models)?;
+        Ok(Self {
+            memo,
+            per_computation,
+            queries,
+            entries,
+            data,
+            models,
+            started,
+            time_limit_ms,
+            input: Some(input),
+            versions,
+            selector: PlanSelector::new(entries, data, models, false)?,
+            logical: HashMap::new(),
+            number: 0,
+            iterations: 0,
+            bound: None,
+            first_bound: None,
+            solver_ms: 0.0,
+            size: (0, 0),
+            time_limit_reached: false,
+            proven_bound: None,
+            done: false,
+            max_free_combinations: crate::search::milp::decode::MAX_FREE_COMBINATIONS,
+            unevaluated_bound: None,
+            decode_failures: 0,
+        })
+        .map(|mut cluster| {
+            cluster.offer_independent_plan();
+            cluster
+        })
+    }
+
+    /// Offer the selector the warm-start incumbent before round 1: each
+    /// query planned alone by the legacy search, with no sharing between
+    /// queries, the plans side by side. It only sets the starting
+    /// incumbent: certification still compares the MILP lower bound with
+    /// the incumbent, and a time-limited loop has a plan to report. A
+    /// query the legacy search finds no plan for leaves no warm start.
+    fn offer_independent_plan(&mut self) {
+        const ID: &str = "L0";
+        let mut roots = Vec::with_capacity(self.entries.len());
+        let mut assignments = Vec::with_capacity(self.entries.len());
+        for (q, (_, root)) in self.per_computation.roots.iter().enumerate() {
+            let Ok(inventory) =
+                asap_logical_asap_optimization::per_computation_candidates::generate_per_computation_candidates(
+                    vec![(q, root.clone())],
+                )
+            else {
+                return;
+            };
+            let run = crate::search::legacy::stage1_to_stage3::run_stage1_to_stage3_from_pass1(
+                inventory,
+                &self.entries[q..=q],
+                self.data,
+                self.models,
+                crate::Stage1ToStage3Search::default(),
+            );
+            let (Some(selected), Some(Ok(alone))) =
+                (run.stage3_selected, run.stage1_selected.map(|l| l.roots))
+            else {
+                return;
+            };
+            roots.extend(alone);
+            assignments.push(selected.assignment.into_iter());
+        }
+        // Each build is read by the one query whose plan it is in, and
+        // takes that plan's materialization.
+        let mut choose = |readers: &[usize]| assignments[*readers.first()?].next();
+        let Ok(Some(mut candidate)) =
+            asap_physical_asap_optimization::optimize_physical_asap_assigned(
+                ID,
+                &roots,
+                self.entries,
+                self.data,
+                self.models.deployment_capabilities,
+                &mut choose,
+            )
+        else {
+            return;
+        };
+        // Every query's assignment taken whole, or the composed plan is not
+        // the per-query plans.
+        if assignments.iter_mut().any(|rest| rest.next().is_some()) {
+            return;
+        }
+        candidate.id = "P0".to_string();
+        if self.selector.offer(&candidate).is_ok() && self.selector.keeps_from_logical(ID) {
+            self.logical.insert(ID.to_string(), roots);
+        }
+    }
+
+    /// Whether the loop has another round; `false` when it has stopped.
+    fn next_round(&mut self) -> bool {
+        if self.done || self.time_limit_reached {
+            self.done = true;
+            return false;
+        }
+        true
+    }
+
+    /// Start the round's solve: set its time limit to the time left now;
+    /// `false` (the loop stops at its time limit) when none is left.
+    fn start_solve(&mut self) -> bool {
+        let remaining_s = self
+            .time_limit_ms
+            .map(|limit| limit as f64 / 1e3 - self.started.elapsed().as_secs_f64());
+        if remaining_s.is_some_and(|s| s <= 0.0) {
+            self.time_limit_reached = true;
+            self.done = true;
+            return false;
+        }
+        let input = self.input.as_mut().expect("the loop's input");
+        input.time_limit_s = remaining_s;
+        input.objective_cutoff = self.selector.cheapest_price();
+        self.iterations += 1;
+        true
+    }
+
+    /// The rest of a round, given its solve's result: the solution's
+    /// incumbent evaluation and its cut, or the end of the loop.
+    fn take_solution(
+        &mut self,
+        solved: Result<Option<ModelSolution>, SolveError>,
+    ) -> Result<(), SelectionError> {
+        let memo = &self.memo;
+        let input = self.input.as_mut().expect("the loop's input");
+        let solved = match solved {
+            Ok(solved) => solved,
+            // The solver stopped at the time limit before finding any
+            // solution: the incumbent so far stands, and the solver's
+            // proven bound bounds every plan not cut.
+            Err(SolveError::StoppedWithoutSolution { proven_bound })
+                if input.time_limit_s.is_some() =>
+            {
+                self.proven_bound = proven_bound;
+                self.time_limit_reached = true;
+                self.done = true;
+                return Ok(());
+            }
+            Err(error) => return Err(SelectionError::Solver(error.to_string())),
+        };
+        // The solver stopped at the time limit: its plan is no bound, but
+        // it is still evaluated, as a last round, and may become the
+        // incumbent; the solver's proven bound bounds every plan not cut.
         let last_round = solved.as_ref().is_some_and(|s| !s.optimal);
-        time_limit_reached |= last_round;
+        self.time_limit_reached |= last_round;
+        if last_round {
+            self.proven_bound = solved.as_ref().and_then(|s| s.proven_bound);
+        }
         let Some(solution) = solved else {
-            // No plan left below the incumbent: it is optimal.
-            if let Some(incumbent) = selector.cheapest_price() {
-                bound = Some(incumbent);
+            // No plan left below the incumbent: it is optimal, but under a
+            // cutoff only to the cutoff's and the solver's tolerances, which
+            // the bound keeps, so the gap shows them.
+            if let Some(incumbent) = self.selector.cheapest_price() {
+                self.bound = Some(match input.objective_cutoff {
+                    Some(_) => (incumbent * (1.0 - CUTOFF_TOLERANCE)
+                        - SOLVER_PRUNING_TOLERANCE / price_scale(input))
+                    .max(0.0),
+                    None => incumbent,
+                });
             }
-            break;
+            self.done = true;
+            return Ok(());
         };
         if !last_round {
-            bound = Some(solution.objective);
-            first_bound.get_or_insert(solution.objective);
+            self.bound = Some(solution.objective);
+            self.first_bound.get_or_insert(solution.objective);
         }
-        size = (solution.variables, solution.constraints);
-        if let Some(incumbent) = selector.cheapest_price() {
+        self.size = (solution.variables, solution.constraints);
+        if let Some(incumbent) = self.selector.cheapest_price() {
             if !last_round && solution.objective >= incumbent - 1e-9 * incumbent.abs() {
-                break;
+                self.done = true;
+                return Ok(());
             }
         }
-        // The solution's group choice: the candidate each chosen choice
-        // takes for its group.
-        let mut group_choice: crate::search::memo::brute_force::GroupChoice =
-            vec![None; memo.groups.len()];
-        for &c in &solution.chosen {
-            group_choice[input.choices[c].group.0 as usize] = Some(input.choices[c].candidate);
+        // Solution decoding, then incumbent evaluation: Stage 3 on the plan
+        // the solution describes and its variants.
+        self.number += 1;
+        let id = format!("L{}", self.number);
+        let decoded = crate::search::milp::decode::decode(
+            memo,
+            self.per_computation,
+            &input.choices,
+            &self.versions,
+            &solution.chosen,
+            &solution.built_once,
+            self.entries,
+            self.data,
+            self.models,
+            &id,
+            self.max_free_combinations,
+        );
+        debug_assert!(decoded.is_some(), "a solution decodes to a plan");
+        self.decode_failures += u32::from(decoded.is_none());
+        let every_plan_evaluated = decoded.as_ref().is_some_and(|d| d.every_combination);
+        let plans = decoded.map_or_else(Vec::new, |d| d.plans);
+        for (k, decoded) in plans.into_iter().enumerate() {
+            let logical_id = decoded.physical.from_logical.clone();
+            let mut physical = decoded.physical;
+            physical.id = format!("P{}-{k}", self.number);
+            self.selector.offer(&physical)?;
+            if self.selector.keeps_from_logical(&logical_id) {
+                self.logical.insert(logical_id, decoded.roots);
+            }
+            let selector = &self.selector;
+            self.logical.retain(|id, _| selector.keeps_from_logical(id));
         }
-        // A group choice that fails a requirement the model does not
-        // encode is offered as any other: the reference evaluation finds
-        // no valid plan in it, and the cut below removes exactly it.
-        crate::search::memo::brute_force::offer_group_choice(
-            &memo,
-            per_computation,
-            &group_choice,
-            entries,
-            data,
-            models,
-            &mut selector,
-            &mut logical,
-            &mut number,
-        )?;
-        // Exclude it: every choice of each chosen (group, candidate).
-        let groups = group_choice.iter().flatten().count();
-        let excluded: Vec<usize> = input
-            .choices
-            .iter()
-            .enumerate()
-            .filter(|(_, choice)| group_choice[choice.group.0 as usize] == Some(choice.candidate))
-            .map(|(c, _)| c)
-            .collect();
+        // A choice set with plans not evaluated keeps its objective as a
+        // lower bound once cut; at the last round (no bound) it stays
+        // uncut, for the relaxation's bound.
+        if !every_plan_evaluated {
+            if last_round {
+                return Ok(());
+            }
+            self.unevaluated_bound = Some(
+                self.unevaluated_bound
+                    .map_or(solution.objective, |b| b.min(solution.objective)),
+            );
+        }
+        // Term cuts: each underestimated term, conditioned on the choices
+        // it depends on.
+        let mut term_cuts = crate::search::milp::cuts::term_cuts_of(input, &solution.chosen);
+        term_cuts.extend(crate::search::milp::cuts::held_cuts_of(
+            input,
+            &solution.chosen,
+        ));
+        term_cuts.extend(crate::search::milp::cuts::ingestion_copy_cuts_of(
+            input,
+            &solution.chosen,
+            &solution.built_once,
+        ));
+        input.cuts.cuts.extend(term_cuts);
+        // Cut exactly the decoded choice set: these choices, these groups
+        // built once or not, no other group needed.
         let mut covered = vec![false; memo.groups.len()];
-        for (g, c) in group_choice.iter().enumerate() {
-            if let Some(c) = c {
-                for i in &memo.candidates[g][*c].implements {
-                    covered[i.0 as usize] = true;
-                }
+        for &c in &solution.chosen {
+            for i in &input.choices[c].implements {
+                covered[i.0 as usize] = true;
             }
         }
         let outside = (0..memo.groups.len()).filter(|&g| !covered[g]).collect();
+        let built: Vec<GroupId> = solution.built_once.iter().copied().collect();
+        let not_built: Vec<GroupId> = input
+            .group_slots
+            .iter()
+            .map(|slot| slot.group)
+            .filter(|g| !solution.built_once.contains(g))
+            .collect();
+        let bound = (solution.chosen.len() + built.len()) as f64 - 1.0;
+        while input.excluded_once.len() < input.excluded.len() {
+            input.excluded_once.push((Vec::new(), Vec::new()));
+        }
         input
             .excluded
-            .push((excluded, outside, groups as f64 - 1.0));
-    }
-    if time_limit_reached && bound.is_none() {
-        // No round solved to optimality: the relaxation's optimum (all cuts
-        // kept) is still a lower bound for the gap.
-        input.relaxed = true;
-        input.time_limit_s = None;
-        // A failed relaxation leaves the gap unknown, never the incumbent lost.
-        if let Ok(Some(relaxation)) = solve(&input) {
-            bound = Some(relaxation.objective);
-        }
-    }
-    if !time_limit_reached && iterations == 1 && selector.cheapest_price().is_none() {
-        // Infeasible from the first round: name the queries with no plan
-        // of their own.
-        let mut without_plan = Vec::new();
-        for (q, entry) in entries.iter().enumerate().take(input.queries.len()) {
-            let mut alone = model_input_clone(&input);
-            for (r, query) in alone.queries.iter_mut().enumerate() {
-                if r != q {
-                    query.roots.clear();
+            .push((solution.chosen.clone(), outside, bound));
+        input.excluded_once.push((built, not_built));
+        Ok(())
+    }
+
+    /// The stopped loop's selection, physical candidate and logical roots,
+    /// with its MILP report, and its cluster's report.
+    fn finish(
+        self,
+    ) -> Result<
+        (
+            Selection,
+            PhysicalCandidate,
+            Vec<QueryRoot>,
+            crate::MilpClusterReport,
+        ),
+        SelectionError,
+    > {
+        let Self {
+            memo,
+            per_computation: _,
+            queries,
+            entries,
+            data: _,
+            models,
+            input,
+            selector,
+            mut logical,
+            iterations,
+            mut bound,
+            first_bound,
+            solver_ms,
+            size,
+            time_limit_reached,
+            unevaluated_bound,
+            decode_failures,
+            proven_bound,
+            ..
+        } = self;
+        let mut input = input.expect("the loop's input");
+        // The last solve's proven bound covers every plan not cut before it,
+        // an earlier optimal round's every plan not cut before that: both
+        // hold for the plans not cut at the end, the rest are evaluated or
+        // floored below.
+        if let Some(proven) = proven_bound {
+            bound = Some(bound.map_or(proven, |bound| bound.max(proven)));
+        }
+        if time_limit_reached && bound.is_none() {
+            // No round solved to optimality: the relaxation's optimum (all
+            // cuts kept) is still a lower bound for the gap.
+            input.relaxed = true;
+            input.time_limit_s = None;
+            // A failed relaxation leaves the gap unknown, never the
+            // incumbent lost.
+            if let Ok(Some(relaxation)) = solve(&input) {
+                bound = Some(relaxation.objective);
+            }
+        }
+        if !time_limit_reached && iterations == 1 && selector.cheapest_price().is_none() {
+            // Infeasible from the first round: name the queries with no plan
+            // of their own.
+            let mut without_plan = Vec::new();
+            for (q, entry) in entries.iter().enumerate().take(input.queries.len()) {
+                let mut alone = model_input_clone(&input);
+                for (r, query) in alone.queries.iter_mut().enumerate() {
+                    if r != q {
+                        query.roots.clear();
+                    }
+                }
+                if solve(&alone)
+                    .map_err(|e| SelectionError::Solver(e.to_string()))?
+                    .is_none()
+                {
+                    without_plan.push(format!("query {} ({})", queries[q], entry.query.0));
                 }
             }
-            if solve(&alone).map_err(SelectionError::Solver)?.is_none() {
-                without_plan.push(format!("query {q} ({})", entry.query.0));
+            if !without_plan.is_empty() {
+                return Err(SelectionError::Solver(format!(
+                    "no valid plan for {} (raw data {})",
+                    without_plan.join(", "),
+                    match models.deployment_capabilities.raw_data_retained {
+                        true => "retained",
+                        false => "not retained",
+                    }
+                )));
             }
         }
-        if !without_plan.is_empty() {
+        if time_limit_reached && selector.cheapest_price().is_none() {
             return Err(SelectionError::Solver(format!(
-                "no valid plan for {} (raw data {})",
-                without_plan.join(", "),
-                match models.deployment_capabilities.raw_data_retained {
-                    true => "retained",
-                    false => "not retained",
-                }
+                "the MILP search reached its time limit after {iterations} rounds before finding a valid plan"
             )));
         }
+        let (mut selection, selected) = selector.finish()?;
+        let price = selection.costs[&selection.selected].cost_per_second;
+        // Once no plan below the incumbent is left, the bound is the
+        // incumbent.
+        let bound = match (bound, unevaluated_bound) {
+            (Some(bound), Some(unevaluated)) => Some(bound.min(unevaluated)),
+            (bound, _) => bound,
+        };
+        let lower_bound = bound.map(|bound| match price {
+            Some(price) => bound.min(price),
+            None => bound,
+        });
+        selection.search.incumbent_cost = price;
+        let report = crate::MilpReport {
+            lower_bound,
+            first_lower_bound: first_bound,
+            gap: gap(price, lower_bound),
+            rounds: iterations,
+            variable_count: size.0,
+            constraint_count: size.1,
+            solver_time_ms: solver_ms as u64,
+            time_limit_reached,
+            decode_failures,
+            clusters: Vec::new(),
+        };
+        let cluster = crate::MilpClusterReport {
+            queries: queries.to_vec(),
+            groups: memo.groups.len(),
+            rounds: report.rounds,
+            gap: report.gap,
+            lower_bound: report.lower_bound,
+            incumbent_cost: price,
+            solver_time_ms: report.solver_time_ms,
+            time_limit_reached,
+            decode_failures,
+        };
+        selection.search.milp = Some(report);
+        let roots = logical
+            .remove(&selected.from_logical)
+            .expect("the selection's logical roots are kept");
+        Ok((selection, selected, roots, cluster))
     }
-    if time_limit_reached && selector.cheapest_price().is_none() {
-        return Err(SelectionError::Solver(format!(
-            "the MILP search reached its time limit after {iterations} rounds before finding a valid plan"
-        )));
+}
+
+/// (incumbent cost − lower bound) / incumbent cost, at least 0; `None`
+/// without both.
+pub(super) fn gap(incumbent_cost: Option<f64>, lower_bound: Option<f64>) -> Option<f64> {
+    match (incumbent_cost, lower_bound) {
+        (Some(price), Some(bound)) if price > 0.0 => Some(((price - bound) / price).max(0.0)),
+        (Some(_), Some(_)) => Some(0.0),
+        _ => None,
     }
-    let (mut selection, selected) = selector.finish()?;
-    let price = selection.costs[&selection.selected].cost_per_second;
-    selection.search.memo = Some(memo_stats);
-    // Once no plan below the incumbent is left, the bound is the incumbent.
-    let lower_bound = bound.map(|bound| match price {
-        Some(price) => bound.min(price),
-        None => bound,
-    });
-    selection.search.incumbent_cost = price;
-    selection.search.milp = Some(crate::MilpReport {
-        lower_bound,
-        first_lower_bound: first_bound,
-        gap: match (price, lower_bound) {
-            (Some(price), Some(bound)) if price > 0.0 => Some(((price - bound) / price).max(0.0)),
-            (Some(_), Some(_)) => Some(0.0),
-            _ => None,
-        },
-        rounds: iterations,
-        variable_count: size.0,
-        constraint_count: size.1,
-        solver_time_ms: solver_ms as u64,
-        time_limit_reached,
-    });
-    let roots = logical
-        .remove(&selected.from_logical)
-        .expect("the selection's logical roots are kept");
-    Ok((selection, selected, roots))
 }
 
 /// The model's input for `memo`: each group's shared choices (copies and
@@ -1637,7 +2821,7 @@ pub(crate) fn model_input_of(
     entries: &[QueryWorkloadEntry],
     data: &DataWorkload,
     models: PlanningModels<'_>,
-) -> Result<ModelInput, SelectionError> {
+) -> Result<(ModelInput, Vec<ImplementationCandidate>), SelectionError> {
     let scale = crate::cost::per_node::frozen_calibration(memo, entries, data, models)?;
     let mut cache = crate::cost::per_node::PriceCache::new(scale);
     let mut choices: Vec<Choice> = Vec::new();
@@ -1646,8 +2830,21 @@ pub(crate) fn model_input_of(
     let CostModel::UsageOverTime(prices) = &models.plan_selection.cost_model else {
         unreachable!("the caller routes cost models 2 and 3 to the legacy search");
     };
+    // Per choice, its candidate at its site thresholds (for solution
+    // decoding) and the threshold it reads each site at.
+    let mut versions: Vec<ImplementationCandidate> = Vec::new();
+    let mut site_thresholds: Vec<Vec<(u64, usize)>> = Vec::new();
+    let mut kept = Vec::new();
+    // Per candidate node (its group and node), its widened-scans twins'
+    // scans: (source key, label set, price per unit of rate).
+    let mut widened_of: HashMap<(usize, *const asap_types::ir::OperatorNode), Vec<WidenedScan>> =
+        HashMap::new();
+    let mut label_sets: HashMap<
+        crate::search::memo::CouplingSite,
+        (Option<u64>, Vec<BTreeSet<String>>),
+    > = HashMap::new();
     for g in 0..memo.groups.len() {
-        let (group_choices, lazy) = crate::search::milp::coefficients::choices(
+        let parts = crate::search::milp::coefficients::choice_parts(
             memo,
             GroupId(g as u32),
             entries,
@@ -1661,25 +2858,121 @@ pub(crate) fn model_input_of(
         // and its builds are sized for fewer readers but built per copy. So
         // the model keeps the shared choices (a lower bound), and the exact
         // pricing of a group choice prices both.
-        // Widened scans are reached through the identical-expression
-        // variants of the exact pricing, never chosen directly.
-        for (choice, lazy) in group_choices.into_iter().zip(lazy) {
-            if choice.copy_of.is_some()
-                || memo.candidates[g][choice.candidate].origin
-                    == crate::search::memo::CandidateOrigin::WidenedScans
-            {
+        for part in &parts {
+            let candidate = &memo.candidates[g][part.choice.candidate];
+            if candidate.origin != crate::search::memo::CandidateOrigin::WidenedScans {
                 continue;
             }
-            lazies.push(lazy.clone());
-            raw_scans.extend(
-                lazy.raw_scans
+            for (site, t) in &part.lazy.thresholds {
+                if site.kind != crate::search::memo::CouplingKind::ScanLabels {
+                    continue;
+                }
+                let (source, sets) = label_sets.entry(*site).or_insert_with(|| {
+                    (
+                        crate::search::memo::rules::scans::site_source_key(memo, site),
+                        crate::search::memo::rules::scans::label_thresholds(memo, site)
+                            .into_iter()
+                            .map(|t| t.labels)
+                            .collect(),
+                    )
+                });
+                let (Some(source), Some(labels)) = (*source, sets.get(*t)) else {
+                    continue;
+                };
+                let cost = part
+                    .choice
+                    .coefficients
+                    .absorbed_scans
                     .iter()
-                    .map(|(source, bytes, via)| (choices.len(), source.clone(), *bytes, *via)),
-            );
-            choices.push(choice);
+                    .filter(|(s, _, _)| *s == source)
+                    .map(|(_, cost, _)| *cost)
+                    .fold(0.0, f64::max);
+                if cost > 0.0 {
+                    widened_of
+                        .entry((g, std::rc::Rc::as_ptr(&candidate.node)))
+                        .or_default()
+                        .push(WidenedScan {
+                            source,
+                            labels: labels.clone(),
+                            cost_per_unit_rate: cost,
+                        });
+                }
+            }
         }
+        // Widened scans are not model variables: a solution also decodes
+        // to its widened plan ([`crate::search::milp::decode`]).
+        kept.extend(parts.into_iter().filter(|part| {
+            part.choice.copy_of.is_none()
+                && memo.candidates[g][part.choice.candidate].origin
+                    != crate::search::memo::CandidateOrigin::WidenedScans
+        }));
+    }
+    // A choice whose build is built once only with another group's
+    // identical build needs such a choice: without any, it describes no plan.
+    loop {
+        let mut groups_of_key: HashMap<u32, BTreeSet<GroupId>> = HashMap::new();
+        for part in &kept {
+            for &key in &part.choice.built_once_keys {
+                groups_of_key
+                    .entry(key)
+                    .or_default()
+                    .insert(part.choice.group);
+            }
+        }
+        let before = kept.len();
+        kept.retain(|part| {
+            part.choice
+                .built_once_with_other_group
+                .iter()
+                .all(|key| groups_of_key[key].iter().any(|&g| g != part.choice.group))
+        });
+        if kept.len() == before {
+            break;
+        }
+    }
+    let mut scan_labels = Vec::new();
+    let mut query_time_scan_labels = Vec::new();
+    let mut widened_scans = Vec::new();
+    for part in kept {
+        let g = part.choice.group.0 as usize;
+        let candidate = &memo.candidates[g][part.choice.candidate];
+        scan_labels.push(crate::search::memo::rules::scans::scan_labels_read(
+            memo, candidate,
+        ));
+        let at_ingestion = part
+            .choice
+            .coefficients
+            .slots
+            .iter()
+            .any(|slot| matches!(slot.materialization, Materialization::MatAtIngestion(_)));
+        query_time_scan_labels.push(match at_ingestion {
+            true => Vec::new(),
+            false => crate::search::memo::rules::scans::private_scan_labels(memo, candidate),
+        });
+        widened_scans.push(
+            widened_of
+                .get(&(g, std::rc::Rc::as_ptr(&candidate.node)))
+                .cloned()
+                .unwrap_or_default(),
+        );
+        let (choice, lazy) = (part.choice, part.lazy);
+        versions.push(part.version);
+        site_thresholds.push(
+            lazy.thresholds
+                .iter()
+                .filter(|(site, _)| site.kind != crate::search::memo::CouplingKind::RawHorizon)
+                .map(|(site, t)| (site.key, *t))
+                .collect(),
+        );
+        lazies.push(lazy.clone());
+        raw_scans.extend(
+            lazy.raw_scans
+                .iter()
+                .map(|(source, bytes, via)| (choices.len(), source.clone(), *bytes, *via)),
+        );
+        choices.push(choice);
     }
-    let mut input = model_input(memo, choices, entries, models)?;
+    let mut input = model_input(memo, choices, entries, data, models)?;
     let mut sources: Vec<String> = Vec::new();
     for (c, source, bytes, via) in raw_scans {
         let index = match sources.iter().position(|s| *s == source) {
@@ -1751,7 +3044,31 @@ pub(crate) fn model_input_of(
         .zip(&lazies)
         .map(|(choice, lazy)| held_at_query_time(choice, lazy, prices, &shared_keys))
         .collect();
-    Ok(input)
+    input.direct_at_ingestion_cost = lazies
+        .iter()
+        .map(|lazy| lazy.direct_at_ingestion.map(|p| p.cost_per_unit_rate))
+        .collect();
+    input.site_thresholds = site_thresholds;
+    input.cuts.shared_slot_held = input
+        .choices
+        .iter()
+        .zip(&lazies)
+        .map(|(choice, lazy)| shared_build_held(choice, lazy, prices, &shared_keys))
+        .collect();
+    input.cuts.at_ingestion_cost = lazies
+        .iter()
+        .map(|lazy| {
+            lazy.direct_at_ingestion
+                .map_or(0.0, |part| part.cost_per_unit_rate)
+        })
+        .collect();
+    input.cuts.scan_labels = scan_labels;
+    input.cuts.query_time_scan_labels = query_time_scan_labels;
+    input.cuts.widened_scans = widened_scans;
+    let mut cuts = std::mem::take(&mut input.cuts);
+    crate::search::milp::cuts::fill_from_input(&input, &mut cuts);
+    input.cuts = cuts;
+    Ok((input, versions))
 }
 
 /// Every group choice of `memo` with its exact price (the cheapest plan
@@ -1794,6 +3111,8 @@ pub(crate) fn lower_bound_check<Id: Clone>(
         };
         let mut restricted = model_input_clone(input);
         restricted.excluded.clear();
+        restricted.objective_cutoff = None;
+        restricted.excluded_once.clear();
         let outside: Vec<usize> = (0..restricted.choices.len())
             .filter(|&c| {
                 group_choice[restricted.choices[c].group.0 as usize]
@@ -1802,13 +3121,54 @@ pub(crate) fn lower_bound_check<Id: Clone>(
             .collect();
         restricted.excluded.push((outside, Vec::new(), 0.0));
         let value = solve(&restricted)
-            .map_err(SelectionError::Solver)?
+            .map_err(|e| SelectionError::Solver(e.to_string()))?
             .map(|s| s.objective);
         out.push((group_choice, exact, value));
     }
     Ok(out)
 }
 
+/// `input` with every term cut ([`crate::search::milp::cuts`]) the lazy
+/// loop may add from any of its first `max_solutions` solutions (each
+/// excluded in turn by its no-good cut, as the loop does), the no-good
+/// cuts dropped: a model that must stay a lower bound of every plan.
+#[cfg(test)]
+pub(crate) fn with_every_term_cut(input: &ModelInput, max_solutions: usize) -> ModelInput {
+    let mut cut = model_input_clone(input);
+    for _ in 0..max_solutions {
+        let Some(solution) = solve(&cut).unwrap() else {
+            break;
+        };
+        let mut found = crate::search::milp::cuts::term_cuts_of(&cut, &solution.chosen);
+        found.extend(crate::search::milp::cuts::held_cuts_of(
+            &cut,
+            &solution.chosen,
+        ));
+        found.extend(crate::search::milp::cuts::ingestion_copy_cuts_of(
+            &cut,
+            &solution.chosen,
+            &solution.built_once,
+        ));
+        cut.cuts.cuts.extend(found);
+        let built: Vec<GroupId> = solution.built_once.iter().copied().collect();
+        let not_built: Vec<GroupId> = cut
+            .group_slots
+            .iter()
+            .map(|slot| slot.group)
+            .filter(|g| !solution.built_once.contains(g))
+            .collect();
+        let bound = (solution.chosen.len() + built.len()) as f64 - 1.0;
+        while cut.excluded_once.len() < cut.excluded.len() {
+            cut.excluded_once.push((Vec::new(), Vec::new()));
+        }
+        cut.excluded.push((solution.chosen, Vec::new(), bound));
+        cut.excluded_once.push((built, not_built));
+    }
+    cut.excluded.clear();
+    cut.excluded_once.clear();
+    cut
+}
+
 fn model_input_clone(input: &ModelInput) -> ModelInput {
     ModelInput {
         group_count: input.group_count,
@@ -1819,23 +3179,330 @@ fn model_input_clone(input: &ModelInput) -> ModelInput {
         schedule: input.schedule.clone(),
         union_rate: input.union_rate.clone(),
         excluded: input.excluded.clone(),
+        excluded_once: input.excluded_once.clone(),
+        site_thresholds: input.site_thresholds.clone(),
         extra_bounds: input.extra_bounds.clone(),
         requires: input.requires.clone(),
         exclusive: input.exclusive.clone(),
         raw_scans: input.raw_scans.clone(),
         raw_data_retained: input.raw_data_retained,
         held: input.held.clone(),
+        direct_at_ingestion_cost: input.direct_at_ingestion_cost.clone(),
+        cuts: input.cuts.clone(),
         time_limit_s: input.time_limit_s,
+        objective_cutoff: input.objective_cutoff,
         relaxed: input.relaxed,
         rollup_rules: input.rollup_rules.clone(),
     }
 }
 
+/// Which choices keep their queries' accuracy targets reachable, as
+/// requirements on child groups' choices ([`ModelInput::requires`]). A
+/// plan's failure probability only grows up the plan (composition adds
+/// its inputs' failure probabilities, and a root's guarantee over a whole
+/// result is at least its per-value one) unless a node is exact whatever
+/// its input (a count of rows), so a group's node composed over its
+/// candidate and one child candidate, every other child exact, bounds it
+/// from below. With its summaries at the largest size Stage 3 can give
+/// them ([`crate::resized_to_the_largest`]), its data-workload statistics
+/// bound and every other unknown term at zero (the optimistic floor), its
+/// failure probability above the largest target of the queries a plan may
+/// lead to the group rules the pair out, under every size threshold of the
+/// candidate's and the child's coupling sites. No guarantee at all (an
+/// approximate input that no rule composes) counts as failing surely. A
+/// candidate that fails with every child exact is ruled out alone. A group
+/// a plan may read through a node exact whatever its input gets no rule.
+/// Nothing is ruled out unless the built-in model's semantics hold: the
+/// deployment's accuracy model may accept an output whatever its guarantee
+/// (as when guarantees are ignored), or its evidence may supply one. A
+/// failure only deeper compositions show is left to the exact pricing (the
+/// no-good cut of its group choice).
+fn guarantee_requirements(
+    memo: &Memo,
+    choices: &[Choice],
+    entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
+    models: PlanningModels<'_>,
+) -> Vec<(usize, Vec<usize>)> {
+    use asap_logical_asap_optimization::per_computation_candidates::recompose_guarantees;
+    use asap_types::ir::properties::{AccuracyGuarantee, AccuracyMetric};
+    use asap_types::ir::OperatorNode;
+    use std::rc::Rc;
+    // Only the built-in model's semantics: a model that may accept an
+    // output whatever its guarantee (also the ignore setting), or evidence,
+    // may accept what the guarantee misses.
+    if models
+        .effective_accuracy()
+        .may_accept_regardless_of_guarantee()
+        || models.evidence.may_provide_evidence()
+    {
+        return Vec::new();
+    }
+    // `node`'s guarantee, composed again over its inputs; `Err` when that
+    // fails.
+    let guarantee_of = |roots: &[QueryRoot]| -> Result<Option<AccuracyGuarantee>, ()> {
+        let roots = recompose_guarantees(roots).map_err(|_| ())?;
+        match &roots[0] {
+            QueryRoot::Operator(root) => Ok(root.guarantee.clone()),
+            _ => Err(()),
+        }
+    };
+    let exact = |node: &Rc<OperatorNode>| {
+        guarantee_of(&[QueryRoot::Operator(Rc::clone(node))])
+            .is_ok_and(|g| g.is_some_and(|g| g.is_exact()))
+    };
+    // The optimistic floor of the failure probability of `node`'s
+    // guarantee with its summaries at the largest size Stage 3's resizing
+    // can give them; `Some(1.0)` for no guarantee at its own size, which
+    // no resizing or node above repairs (an approximate input no rule
+    // composes); `None` when unknown. A guarantee lost only at the largest
+    // size (a family whose size cap lets its error grow there, e.g. a
+    // UnivMon heap at a tiny failure probability) is unknown, never 1.0.
+    let floor = |node: &Rc<OperatorNode>| -> Option<f64> {
+        let largest =
+            crate::resized_to_the_largest(&[QueryRoot::Operator(Rc::clone(node))]).ok()?;
+        let Some(guarantee) = guarantee_of(&largest).ok()? else {
+            return guarantee_of(&[QueryRoot::Operator(Rc::clone(node))])
+                .ok()?
+                .is_none()
+                .then_some(1.0);
+        };
+        crate::bind_statistics(&guarantee, &|statistic| {
+            crate::workload_statistic(statistic, None, data)
+        })
+        .optimistic_floor()
+        .spec
+        .failure_probability
+        .evaluate()
+    };
+    // The queries a plan may lead to each group: its roots' queries, and
+    // those of every group with a candidate reading it. A group is free
+    // when a plan may read it with no target: a candidate's output is
+    // exact over one of its approximate candidates (a count of rows), or a
+    // free group reads it.
+    let mut reaching: Vec<BTreeSet<usize>> = memo
+        .groups
+        .iter()
+        .map(|group| group.roots_of.iter().copied().collect())
+        .collect();
+    let mut free = vec![false; memo.groups.len()];
+    for candidates in &memo.candidates {
+        for candidate in candidates {
+            for &GroupId(h) in &candidate.children {
+                let h = h as usize;
+                let child_node = Rc::as_ptr(&memo.groups[h].node);
+                free[h] |= memo.candidates[h].iter().any(|d| {
+                    let replaced = HashMap::from([(child_node, Rc::clone(&d.node))]);
+                    let composed = crate::search::memo::brute_force::substituted(
+                        &candidate.node,
+                        &replaced,
+                        &mut HashMap::new(),
+                    );
+                    !exact(&d.node) && exact(&composed)
+                });
+            }
+        }
+    }
+    loop {
+        let mut changed = false;
+        for (g, candidates) in memo.candidates.iter().enumerate() {
+            let above = reaching[g].clone();
+            let children = candidates
+                .iter()
+                .flat_map(|c| c.children.iter())
+                .chain(&memo.groups[g].children);
+            for child in children {
+                let h = child.0 as usize;
+                let before = (reaching[h].len(), free[h]);
+                reaching[h].extend(&above);
+                free[h] |= free[g];
+                changed |= (reaching[h].len(), free[h]) != before;
+            }
+        }
+        if !changed {
+            break;
+        }
+    }
+    // The largest target of those queries: a plan needing the group meets
+    // at least one of them through it. None when one has no target, or
+    // the group is free.
+    let target = |g: usize| -> Option<f64> {
+        if free[g] {
+            return None;
+        }
+        let targets: Option<Vec<f64>> = reaching[g]
+            .iter()
+            .map(|&q| entries[q].requirements.accuracy.failure_probability)
+            .collect();
+        targets?.into_iter().max_by(f64::total_cmp)
+    };
+    let misses = |fp: Option<f64>, target: f64| fp.is_some_and(|fp| fp > target * (1.0 + 1e-9));
+    // Each candidate with a choice, at each threshold of its sites.
+    let mut choices_of: BTreeMap<(usize, usize), Vec<usize>> = BTreeMap::new();
+    for (c, choice) in choices.iter().enumerate() {
+        choices_of
+            .entry((choice.group.0 as usize, choice.candidate))
+            .or_default()
+            .push(c);
+    }
+    let versions: HashMap<(usize, usize), Vec<Rc<OperatorNode>>> = choices_of
+        .keys()
+        .map(|&(g, c)| {
+            let nodes =
+                crate::search::milp::coefficients::at_site_thresholds(memo, &memo.candidates[g][c])
+                    .into_iter()
+                    .map(|(version, _)| version.node)
+                    .collect();
+            ((g, c), nodes)
+        })
+        .collect();
+    let mut requires = Vec::new();
+    for (&(g, c), own) in &choices_of {
+        let Some(target) = target(g) else {
+            continue;
+        };
+        let nodes = &versions[&(g, c)];
+        if nodes.iter().all(|node| misses(floor(node), target)) {
+            requires.extend(own.iter().map(|&choice| (choice, Vec::new())));
+            continue;
+        }
+        for child in &memo.candidates[g][c].children {
+            let h = child.0 as usize;
+            let child_node = Rc::as_ptr(&memo.groups[h].node);
+            let of_child: Vec<(usize, usize)> = choices_of
+                .keys()
+                .copied()
+                .filter(|&(k, _)| k == h)
+                .collect();
+            let allowed: Vec<(usize, usize)> = of_child
+                .iter()
+                .copied()
+                .filter(|d| {
+                    nodes.iter().any(|node| {
+                        versions[d].iter().any(|child_version| {
+                            let replaced = HashMap::from([(child_node, Rc::clone(child_version))]);
+                            let composed = crate::search::memo::brute_force::substituted(
+                                node,
+                                &replaced,
+                                &mut HashMap::new(),
+                            );
+                            !misses(floor(&composed), target)
+                        })
+                    })
+                })
+                .collect();
+            if allowed.len() == of_child.len() {
+                continue;
+            }
+            let allowed: Vec<usize> = allowed
+                .iter()
+                .flat_map(|d| choices_of[d].iter().copied())
+                .collect();
+            requires.extend(own.iter().map(|&choice| (choice, allowed.clone())));
+        }
+    }
+    // The accuracy metric: a query's root guarantee must be in its target's
+    // metric (or one comparable directly), whatever its size. A root
+    // group's candidate decides that metric when everything below it that
+    // may approximate is in it or in one child candidate whose children
+    // only compute exactly (a group whose candidates and their children's
+    // are all exact): that pair's guarantee is the root's.
+    let mut exact_only = vec![false; memo.groups.len()];
+    if let Some(order) = memo.topological_order() {
+        for GroupId(h) in order {
+            let h = h as usize;
+            exact_only[h] = memo.candidates[h].iter().all(|candidate| {
+                crate::search::memo::rules::windows::is_exact(candidate)
+                    && candidate.children.iter().all(|c| exact_only[c.0 as usize])
+            });
+        }
+    }
+    let metric_of = |node: &Rc<OperatorNode>| {
+        let guarantee = guarantee_of(&[QueryRoot::Operator(Rc::clone(node))]).ok()??;
+        (!guarantee.is_exact()).then_some(guarantee.spec.accuracy_metric)
+    };
+    let metric_meets = |metric: AccuracyMetric, target: AccuracyMetric| {
+        metric == target
+            || (metric == AccuracyMetric::ABSOLUTE_STREAM_L2
+                && target == AccuracyMetric::ABSOLUTE_STREAM_L1)
+    };
+    for (&(g, c), own) in &choices_of {
+        let targets: Vec<AccuracyMetric> = memo.groups[g]
+            .roots_of
+            .iter()
+            .filter_map(|&q| {
+                let target = &entries[q].requirements.accuracy;
+                (!target.is_exact())
+                    .then_some(target.accuracy_metric)
+                    .flatten()
+            })
+            .collect();
+        if targets.is_empty() {
+            continue;
+        }
+        let misses = |node: &Rc<OperatorNode>| {
+            metric_of(node).is_some_and(|m| !targets.iter().all(|&t| metric_meets(m, t)))
+        };
+        let candidate = &memo.candidates[g][c];
+        let open: Vec<usize> = candidate
+            .children
+            .iter()
+            .map(|h| h.0 as usize)
+            .filter(|&h| !exact_only[h])
+            .collect();
+        match open[..] {
+            [] => {
+                if misses(&candidate.node) {
+                    requires.extend(own.iter().map(|&choice| (choice, Vec::new())));
+                }
+            }
+            [h] => {
+                let child_node = Rc::as_ptr(&memo.groups[h].node);
+                let of_child: Vec<(usize, usize)> = choices_of
+                    .keys()
+                    .copied()
+                    .filter(|&(k, _)| k == h)
+                    .collect();
+                let allowed: Vec<(usize, usize)> = of_child
+                    .iter()
+                    .copied()
+                    .filter(|&(_, d)| {
+                        let child = &memo.candidates[h][d];
+                        // A child reading approximate groups decides nothing.
+                        if !child.children.iter().all(|k| exact_only[k.0 as usize]) {
+                            return true;
+                        }
+                        let replaced = HashMap::from([(child_node, Rc::clone(&child.node))]);
+                        let composed = crate::search::memo::brute_force::substituted(
+                            &candidate.node,
+                            &replaced,
+                            &mut HashMap::new(),
+                        );
+                        !misses(&composed)
+                    })
+                    .collect();
+                if allowed.len() < of_child.len() {
+                    let allowed: Vec<usize> = allowed
+                        .iter()
+                        .flat_map(|d| choices_of[d].iter().copied())
+                        .collect();
+                    requires.extend(own.iter().map(|&choice| (choice, allowed.clone())));
+                }
+            }
+            // Several approximate children: no pair decides the metric.
+            _ => {}
+        }
+    }
+
+    requires
+}
+
 /// The model's input for `memo`'s `choices`.
 fn model_input(
     memo: &Memo,
     choices: Vec<Choice>,
     entries: &[QueryWorkloadEntry],
+    data: &DataWorkload,
     models: PlanningModels<'_>,
 ) -> Result<ModelInput, SelectionError> {
     let settings = models.plan_selection;
@@ -1927,6 +3594,9 @@ fn model_input(
             requires.push((c, allowed));
         }
     }
+    requires.extend(guarantee_requirements(
+        memo, &choices, entries, data, models,
+    ));
     // A reader of a scan shared at horizon `h` reads the scan group's
     // choice at that horizon, so every reader of one scan agrees on it.
     for (c, choice) in choices.iter().enumerate() {
@@ -1992,6 +3662,8 @@ fn model_input(
         schedule,
         union_rate,
         excluded: Vec::new(),
+        excluded_once: Vec::new(),
+        site_thresholds: Vec::new(),
         extra_bounds: Vec::new(),
         requires,
         exclusive,
@@ -1999,13 +3671,16 @@ fn model_input(
         raw_scans: Vec::new(),
         raw_data_retained: true,
         held: Vec::new(),
+        direct_at_ingestion_cost: Vec::new(),
+        cuts: TermCuts::default(),
         time_limit_s: None,
+        objective_cutoff: None,
         relaxed: false,
     })
 }
 
 #[cfg(test)]
-mod tests {
+pub(crate) mod tests {
     use super::*;
     use asap_types::ir::properties::AccuracyTarget;
 
@@ -2027,12 +3702,15 @@ mod tests {
                 cross_shared: Vec::new(),
                 absorbed_scans: Vec::new(),
                 direct_nodes: Vec::new(),
+                query_work: Vec::new(),
             },
             copy_of: None,
             scan_children: Vec::new(),
             scan_horizon_ms: None,
             scan_source: None,
             scan_project_cost_per_unit_rate: 0.0,
+            built_once_with_other_group: Vec::new(),
+            built_once_keys: Vec::new(),
         }
     }
 
@@ -2068,6 +3746,8 @@ mod tests {
             schedule: vec![0, 0],
             union_rate: [(vec![0], 1.0)].into(),
             excluded: Vec::new(),
+            excluded_once: Vec::new(),
+            site_thresholds: Vec::new(),
             extra_bounds: Vec::new(),
             requires: Vec::new(),
             exclusive: (Vec::new(), Vec::new()),
@@ -2075,7 +3755,10 @@ mod tests {
             raw_scans: Vec::new(),
             raw_data_retained: true,
             held: Vec::new(),
+            direct_at_ingestion_cost: Vec::new(),
+            cuts: TermCuts::default(),
             time_limit_s: None,
+            objective_cutoff: None,
             relaxed: false,
         };
         let solution = solve(&input).unwrap().unwrap();
@@ -2098,19 +3781,143 @@ mod tests {
         );
     }
 
+    /// One node (work key 7) computed by two groups' choices, each read by
+    /// one query: every 10 s and every 20 s, two classes no group's class
+    /// set joins (`union_rate` has only `[0]` and `[1]`). Stage 3 charges
+    /// its held memory per reader at union / sum of their rates, the union
+    /// at least the faster's: 0.1 / 0.15. The model charges no more.
+    #[test]
+    fn held_memory_of_a_node_read_from_two_classes_is_shared() {
+        let shared = |group: u32| {
+            let mut c = choice(group, &[], 0.0);
+            c.coefficients.query_work = vec![(7, 1.0)];
+            c
+        };
+        let input = ModelInput {
+            group_count: 2,
+            group_queries: vec![vec![0], vec![1]],
+            choices: vec![shared(0), shared(1)],
+            group_slots: Vec::new(),
+            queries: vec![
+                ModelQuery {
+                    rate: 0.1,
+                    roots: vec![GroupId(0)],
+                    max_latency_ms: None,
+                },
+                ModelQuery {
+                    rate: 0.05,
+                    roots: vec![GroupId(1)],
+                    max_latency_ms: None,
+                },
+            ],
+            schedule: vec![0, 1],
+            union_rate: [(vec![0], 0.1), (vec![1], 0.05)].into(),
+            excluded: Vec::new(),
+            excluded_once: Vec::new(),
+            site_thresholds: Vec::new(),
+            extra_bounds: Vec::new(),
+            requires: Vec::new(),
+            exclusive: (Vec::new(), Vec::new()),
+            rollup_rules: Vec::new(),
+            raw_scans: Vec::new(),
+            raw_data_retained: true,
+            held: vec![
+                HeldAtQueryTime {
+                    per_reader_per_ms: 1.0,
+                    longest_per_ms: 0.0,
+                };
+                2
+            ],
+            direct_at_ingestion_cost: Vec::new(),
+            cuts: TermCuts::default(),
+            time_limit_s: None,
+            relaxed: false,
+            objective_cutoff: None,
+        };
+        // Each reader: 1.0 per ms × its rate × 1 ms × 0.1 / 0.15.
+        let stage_3 = (0.1 + 0.05) * (0.1 / 0.15);
+        let solution = solve(&input).unwrap().unwrap();
+        assert!(
+            solution.objective <= stage_3 * (1.0 + 1e-9),
+            "model {} over {stage_3}",
+            solution.objective
+        );
+    }
+
     /// A group choice: per group, its candidate.
     type Plan = Vec<Option<usize>>;
 
-    /// The plans of `queries`' memo (each query every minute, at `targets`)
-    /// the model may return, as group choices: the model of every choice
-    /// with no latency bound, solved with a no-good cut on each group
-    /// choice it returns until it is infeasible; and the group choices of
-    /// the same candidates the reference search accepts
-    /// ([`crate::search::memo::rules::windows::requirements_met`]).
-    fn model_and_reference_plans(
+    /// `queries`' memo (each query every minute, at `targets`, over a
+    /// gauge `m` of 10,000 series) and its model's input: every choice
+    /// but copies and widened scans.
+    fn memo_and_input(queries: &[&str], targets: &[AccuracyTarget]) -> (Memo, ModelInput) {
+        memo_and_input_over(queries, targets, None)
+    }
+
+    /// A query over two identical subtrees reads one child group twice:
+    /// the child is one reach flow, so its plans stay feasible.
+    #[test]
+    fn a_child_read_twice_by_one_choice_is_one_flow() {
+        let (_, input) = memo_and_input(
+            &["sum_over_time(m[1m]) / sum_over_time(m[1m])"],
+            &[AccuracyTarget::exact()],
+        );
+        assert!(
+            input.choices.iter().any(|choice| {
+                let children: BTreeSet<_> = choice.children.iter().collect();
+                children.len() < choice.children.len()
+            }),
+            "a choice lists one child twice"
+        );
+        assert!(solve(&input).unwrap().is_some(), "a plan");
+    }
+
+    /// [`memo_and_input`] with `m`'s distinct items per evaluation bound
+    /// to `max_distinct`.
+    fn memo_and_input_over(
         queries: &[&str],
         targets: &[AccuracyTarget],
-    ) -> (BTreeSet<Plan>, BTreeSet<Plan>) {
+        max_distinct: Option<u64>,
+    ) -> (Memo, ModelInput) {
+        let (per_computation, entries, data) = workload(queries, targets, max_distinct);
+        let models = PlanningModels::builtin();
+        let memo = crate::search::memo::brute_force::memo(&per_computation, &entries);
+        let scale =
+            crate::cost::per_node::frozen_calibration(&memo, &entries, &data, models).unwrap();
+        let mut cache = crate::cost::per_node::PriceCache::new(scale);
+        let mut choices = Vec::new();
+        for g in 0..memo.groups.len() {
+            let (group_choices, _) = crate::search::milp::coefficients::choices(
+                &memo,
+                GroupId(g as u32),
+                &entries,
+                &data,
+                models,
+                &mut cache,
+            )
+            .unwrap();
+            choices.extend(group_choices.into_iter().filter(|choice| {
+                choice.copy_of.is_none()
+                    && memo.candidates[g][choice.candidate].origin
+                        != crate::search::memo::CandidateOrigin::WidenedScans
+            }));
+        }
+        let input = model_input(&memo, choices, &entries, &data, models).unwrap();
+        (memo, input)
+    }
+
+    /// `queries`' workload: each query every minute, at `targets`, over a
+    /// gauge `m` of 10,000 series whose distinct items per evaluation are
+    /// bounded by `max_distinct`.
+    pub(crate) fn workload(
+        queries: &[&str],
+        targets: &[AccuracyTarget],
+        max_distinct: Option<u64>,
+    ) -> (
+        PerComputationCandidates<usize>,
+        Vec<QueryWorkloadEntry>,
+        DataWorkload,
+    ) {
         use asap_logical_asap_optimization::per_computation_candidates::generate_per_computation_candidates;
         use asap_logical_asap_optimization::test_support::{entry_with, lower_promql};
         use asap_types::ir::schema::schema_support::with_promql_series_identity;
@@ -2149,35 +3956,27 @@ mod tests {
                 "m".to_string(),
                 SourceFacts {
                     series_count: Some(10_000),
+                    max_distinct_per_evaluation: max_distinct,
                     ..SourceFacts::default()
                 },
             )]
             .into_iter()
             .collect(),
         };
-        let models = PlanningModels::builtin();
-        let memo = crate::search::memo::brute_force::memo(&per_computation, &entries);
-        let scale =
-            crate::cost::per_node::frozen_calibration(&memo, &entries, &data, models).unwrap();
-        let mut cache = crate::cost::per_node::PriceCache::new(scale);
-        let mut choices = Vec::new();
-        for g in 0..memo.groups.len() {
-            let (group_choices, _) = crate::search::milp::coefficients::choices(
-                &memo,
-                GroupId(g as u32),
-                &entries,
-                &data,
-                models,
-                &mut cache,
-            )
-            .unwrap();
-            choices.extend(group_choices.into_iter().filter(|choice| {
-                choice.copy_of.is_none()
-                    && memo.candidates[g][choice.candidate].origin
-                        != crate::search::memo::CandidateOrigin::WidenedScans
-            }));
-        }
-        let mut input = model_input(&memo, choices, &entries, models).unwrap();
+        (per_computation, entries, data)
+    }
+
+    /// The plans of `queries`' memo ([`memo_and_input`]) the model may
+    /// return, as group choices: the model of every choice with no latency
+    /// bound, solved with a no-good cut on each group choice it returns
+    /// until it is infeasible; and the group choices of the same
+    /// candidates the reference search accepts
+    /// ([`crate::search::memo::rules::windows::requirements_met`]).
+    fn model_and_reference_plans(
+        queries: &[&str],
+        targets: &[AccuracyTarget],
+    ) -> (BTreeSet<Plan>, BTreeSet<Plan>) {
+        let (memo, mut input) = memo_and_input(queries, targets);
         assert!(
             input
                 .choices
@@ -2283,4 +4082,453 @@ mod tests {
             reference.len()
         );
     }
+
+    /// The average of per-series medians over 10 minutes, within 1%: the
+    /// average's exact sum and count over the medians have a guarantee only
+    /// when a rule composes the medians' error: an exact median, or a
+    /// relative error (DDSketch), but not a KLL's rank error. So the model
+    /// requires the average's choices that read the medians' group to read
+    /// one of those, never a KLL one.
+    #[test]
+    fn a_choice_requires_child_choices_that_keep_its_guarantee() {
+        let (memo, input) = memo_and_input(
+            &["avg by (instance) (quantile_over_time(0.5, m[10m]))"],
+            &[AccuracyTarget::max_error(0.01).failure_probability(0.01)],
+        );
+        let candidate = |c: usize| {
+            &memo.candidates[input.choices[c].group.0 as usize][input.choices[c].candidate]
+        };
+        let builds = |c: usize, algorithm: &str| {
+            asap_types::ir::OperatorNode::reachable(&candidate(c).node)
+                .iter()
+                .any(|n| {
+                    matches!(n.asap(), Some(asap_types::ir::ASAPOp::SummaryAgg { family, .. })
+                        if format!("{family:?}").contains(algorithm))
+                })
+        };
+        let exact = |c: usize| crate::search::memo::rules::windows::is_exact(candidate(c));
+        let mut checked = 0;
+        for (c, allowed) in &input.requires {
+            let child_group = |d: usize| candidate(*c).children.contains(&input.choices[d].group);
+            let kll_children: Vec<usize> = (0..input.choices.len())
+                .filter(|&d| child_group(d) && builds(d, "Kll"))
+                .collect();
+            if kll_children.is_empty() || !exact(*c) {
+                continue;
+            }
+            checked += 1;
+            assert!(
+                allowed.iter().all(|d| !kll_children.contains(d)),
+                "{allowed:?}"
+            );
+            assert!(allowed.iter().any(|&d| exact(d)), "{allowed:?}");
+        }
+        assert!(
+            checked > 0,
+            "the average's choices require a composable median"
+        );
+        // Solved, no exact average reads a KLL median.
+        let solution = solve(&input).unwrap().unwrap();
+        for &c in &solution.chosen {
+            for &d in &solution.chosen {
+                let reads = candidate(c).children.contains(&input.choices[d].group);
+                assert!(!(reads && exact(c) && builds(d, "Kll")), "{c} reads {d}");
+            }
+        }
+    }
+
+    /// Classic HLL's failure probability over more distinct items than its
+    /// contract covers is 1 at every precision: no size threshold (a larger
+    /// shared precision) repairs it, so ruling such a candidate out under
+    /// every threshold of its sites is ruling it out for good.
+    #[test]
+    fn an_uncertified_distinct_count_fails_at_every_precision() {
+        use asap_logical_asap_optimization::accuracy::estimators::hll::classic_failure_probability;
+        for precision in 4..=18 {
+            assert_eq!(classic_failure_probability(precision, 0.05, 1_000_000), 1.0);
+        }
+        // Within the contract, a precision large enough certifies it.
+        assert_eq!(classic_failure_probability(4, 0.05, 1_000), 1.0);
+        assert!(classic_failure_probability(18, 0.05, 1_000) < 1.0);
+    }
+
+    /// A 0.9 quantile over 5 minutes within a rank error of 1%: a DDSketch
+    /// estimate guarantees a relative error, which no rule compares with a
+    /// rank target, whatever its size. The model never chooses a plan
+    /// whose root reads a DDSketch, and still has one reading a KLL.
+    #[test]
+    fn a_choice_whose_guarantee_is_in_another_metric_is_ruled_out() {
+        use asap_types::ir::properties::AccuracyMetric;
+        let (memo, mut input) = memo_and_input(
+            &["quantile_over_time(0.9, m[5m])"],
+            &[AccuracyTarget::max_error(0.01)
+                .failure_probability(0.01)
+                .accuracy_metric(AccuracyMetric::Rank)],
+        );
+        let builds = |c: usize, algorithm: &str| {
+            let choice = &input.choices[c];
+            asap_types::ir::OperatorNode::reachable(
+                &memo.candidates[choice.group.0 as usize][choice.candidate].node,
+            )
+            .iter()
+            .any(|n| {
+                matches!(n.asap(), Some(asap_types::ir::ASAPOp::SummaryAgg { family, .. })
+                    if format!("{family:?}").contains(algorithm))
+            })
+        };
+        assert!((0..input.choices.len()).any(|c| builds(c, "DDSketch")));
+        for query in &mut input.queries {
+            query.max_latency_ms = None;
+        }
+        // Every plan the model may return reads no DDSketch.
+        let mut plans = 0;
+        while let Some(solution) = solve(&input).unwrap() {
+            plans += 1;
+            assert!(
+                solution.chosen.iter().all(|&c| !builds(c, "DDSketch")),
+                "{:?}",
+                solution.chosen
+            );
+            let groups = solution
+                .chosen
+                .iter()
+                .map(|&c| input.choices[c].group)
+                .collect::<BTreeSet<_>>()
+                .len();
+            input
+                .excluded
+                .push((solution.chosen.clone(), Vec::new(), groups as f64 - 1.0));
+        }
+        assert!(plans > 0);
+    }
+    /// The classic HLL `params` of every summary `c`'s candidate builds.
+    fn hll_precisions(memo: &Memo, input: &ModelInput, c: usize) -> Vec<(u8, Option<f64>)> {
+        let choice = &input.choices[c];
+        asap_types::ir::OperatorNode::reachable(
+            &memo.candidates[choice.group.0 as usize][choice.candidate].node,
+        )
+        .iter()
+        .filter_map(|n| match n.asap() {
+            Some(asap_types::ir::ASAPOp::SummaryAgg {
+                family: asap_types::ir::schema::FieldDataType::Sketch(kind, _),
+                ..
+            }) => match kind.params() {
+                asap_types::ir::schema::SketchParams::Hll {
+                    precision,
+                    relative_error,
+                } => Some((*precision, *relative_error)),
+                _ => None,
+            },
+            _ => None,
+        })
+        .collect()
+    }
+
+    /// Distinct items per minute within 5% with probability 0.99, over at
+    /// most 2,000 distinct items: classic HLL at its Pass 1 precision does
+    /// not certify that, and Stage 3 enlarges its precision until it does.
+    /// So no HLL choice is ruled out.
+    #[test]
+    fn a_summary_stage_3_enlarges_to_certify_is_kept() {
+        use asap_logical_asap_optimization::accuracy::estimators::hll;
+        let (memo, input) = memo_and_input_over(
+            &["distinct_over_time(m[1m])"],
+            &[AccuracyTarget::max_error(0.05).failure_probability(0.01)],
+            Some(2_000),
+        );
+        let hll_choices: Vec<usize> = (0..input.choices.len())
+            .filter(|&c| !hll_precisions(&memo, &input, c).is_empty())
+            .collect();
+        assert!(!hll_choices.is_empty());
+        // Its Pass 1 precision misses the target; Stage 3's meets it.
+        assert!(hll_choices.iter().any(|&c| {
+            hll_precisions(&memo, &input, c)
+                .iter()
+                .any(|&(precision, relative_error)| {
+                    let relative_error = hll::classic_relative_error(precision, relative_error);
+                    hll::classic_failure_probability(precision, relative_error, 2_000) > 0.01
+                        && hll::classic_failure_probability(18, relative_error, 2_000) <= 0.01
+                })
+        }));
+        for (c, allowed) in &input.requires {
+            assert!(
+                !(hll_choices.contains(c) && allowed.is_empty()),
+                "HLL choice {c} ruled out"
+            );
+        }
+    }
+
+    /// The number of series over 1,000,000 distinct items each: a count
+    /// of rows is exact whatever the per-series distinct estimates are, so
+    /// its target asks nothing of them, and an HLL that certifies no
+    /// failure probability over that many items is not ruled out.
+    #[test]
+    fn a_count_of_rows_asks_nothing_of_its_input() {
+        let (memo, input) = memo_and_input_over(
+            &["count(distinct_over_time(m[1m]))"],
+            &[AccuracyTarget::max_error(0.05).failure_probability(0.01)],
+            Some(1_000_000),
+        );
+        let hll_choices: Vec<usize> = (0..input.choices.len())
+            .filter(|&c| !hll_precisions(&memo, &input, c).is_empty())
+            .collect();
+        assert!(!hll_choices.is_empty());
+        for (c, allowed) in &input.requires {
+            assert!(
+                !(hll_choices.contains(c) && allowed.is_empty()),
+                "HLL choice {c} ruled out"
+            );
+        }
+    }
+    /// Samples per distinct value per series, within 5% with probability
+    /// 0.99: the distinct count is the divisor. A UnivMon divisor sized for
+    /// a vanishing failure probability hits its heap cap, and its relative
+    /// error grows past what a division composes, so the quotient has no
+    /// guarantee at the largest size only. That is unknown, not a sure
+    /// failure: no choice building a UnivMon is ruled out.
+    #[test]
+    fn a_guarantee_lost_only_at_the_largest_size_rules_nothing_out() {
+        let (memo, input) = memo_and_input(
+            &["count_over_time(m[1m]) / distinct_over_time(m[1m])"],
+            &[AccuracyTarget::max_error(0.05).failure_probability(0.01)],
+        );
+        let univmon = |c: usize| {
+            let choice = &input.choices[c];
+            asap_types::ir::OperatorNode::reachable(
+                &memo.candidates[choice.group.0 as usize][choice.candidate].node,
+            )
+            .iter()
+            .any(|n| {
+                matches!(n.asap(), Some(asap_types::ir::ASAPOp::SummaryAgg { family, .. })
+                    if format!("{family:?}").contains("UnivMon"))
+            })
+        };
+        let univmons: Vec<usize> = (0..input.choices.len()).filter(|&c| univmon(c)).collect();
+        assert!(!univmons.is_empty());
+        for (c, allowed) in &input.requires {
+            assert!(
+                !(univmons.contains(c) && allowed.is_empty()),
+                "UnivMon choice {c} ruled out"
+            );
+            // A requirement on a group with UnivMon choices allows them all.
+            for &d in &univmons {
+                let same_group = allowed
+                    .iter()
+                    .any(|&a| input.choices[a].group == input.choices[d].group);
+                assert!(
+                    !same_group || allowed.contains(&d),
+                    "choice {c} rules out UnivMon choice {d}"
+                );
+            }
+        }
+    }
+
+    /// A model of one query per root group, each every second in one
+    /// schedule class, with `choices`.
+    fn small_input(
+        group_queries: Vec<Vec<usize>>,
+        roots: &[u32],
+        choices: Vec<Choice>,
+    ) -> ModelInput {
+        ModelInput {
+            group_count: group_queries.len(),
+            group_queries,
+            choices,
+            group_slots: Vec::new(),
+            queries: roots
+                .iter()
+                .map(|&root| ModelQuery {
+                    rate: 1.0,
+                    roots: vec![GroupId(root)],
+                    max_latency_ms: None,
+                })
+                .collect(),
+            schedule: vec![0; roots.len()],
+            union_rate: [(vec![0], 1.0)].into(),
+            excluded: Vec::new(),
+            excluded_once: Vec::new(),
+            site_thresholds: Vec::new(),
+            extra_bounds: Vec::new(),
+            requires: Vec::new(),
+            exclusive: (Vec::new(), Vec::new()),
+            rollup_rules: Vec::new(),
+            raw_scans: Vec::new(),
+            raw_data_retained: true,
+            held: Vec::new(),
+            direct_at_ingestion_cost: Vec::new(),
+            cuts: TermCuts::default(),
+            time_limit_s: None,
+            objective_cutoff: None,
+            relaxed: false,
+        }
+    }
+
+    /// Two groups each cheapest at a different threshold of one coupling
+    /// site: one plan reads the site at one threshold, so the model pays
+    /// the dearer choice of one group (6), not both cheap ones (2).
+    #[test]
+    fn a_coupling_site_is_read_at_one_threshold() {
+        let mut input = small_input(
+            vec![vec![0], vec![1]],
+            &[0, 1],
+            vec![
+                choice(0, &[], 1.0),
+                choice(0, &[], 5.0),
+                choice(1, &[], 5.0),
+                choice(1, &[], 1.0),
+            ],
+        );
+        input.site_thresholds = vec![vec![(7, 0)], vec![(7, 1)], vec![(7, 0)], vec![(7, 1)]];
+        let solution = solve(&input).unwrap().unwrap();
+        assert!(
+            (solution.objective - 6.0).abs() < 1e-9,
+            "{}",
+            solution.objective
+        );
+        let thresholds: BTreeSet<usize> = solution
+            .chosen
+            .iter()
+            .map(|&c| input.site_thresholds[c][0].1)
+            .collect();
+        assert_eq!(thresholds.len(), 1, "{:?}", solution.chosen);
+    }
+
+    /// A group no chosen choice reads is not needed: a cut that frees the
+    /// groups its plan does not need cannot be escaped by choosing such a
+    /// group (free here) next to the cut plan.
+    #[test]
+    fn a_group_nothing_reads_is_not_needed() {
+        let mut input = small_input(
+            vec![vec![0], vec![0]],
+            &[0],
+            vec![
+                choice(0, &[], 1.0),
+                choice(0, &[], 3.0),
+                choice(1, &[], 0.0),
+            ],
+        );
+        let first = solve(&input).unwrap().unwrap();
+        assert_eq!(first.chosen, vec![0]);
+        input.excluded.push((vec![0], vec![1], 0.0));
+        let second = solve(&input).unwrap().unwrap();
+        assert_eq!(second.chosen, vec![1], "the cut plan with a dangling group");
+    }
+
+    /// A scan shared at a horizon by one reader is that reader's own scan
+    /// read apart: the model does not offer it shared, even when cheaper.
+    #[test]
+    fn a_scan_shared_by_one_reader_is_read_apart() {
+        let mut scan = choice(0, &[], 1.0);
+        scan.scan_source = Some(1);
+        scan.scan_horizon_ms = Some(60_000);
+        let mut shared = choice(1, &[0], 1.0);
+        shared.scan_children = vec![(GroupId(0), Some(60_000))];
+        let mut apart = choice(1, &[], 5.0);
+        apart.scan_children = vec![(GroupId(0), None)];
+        let input = small_input(vec![vec![0], vec![0]], &[1], vec![scan, shared, apart]);
+        let solution = solve(&input).unwrap().unwrap();
+        assert_eq!(solution.chosen, vec![2]);
+    }
+
+    /// A cut excludes its plan's choices with its groups built once or
+    /// not: the same choices with a group recomputed instead are another
+    /// plan, cut by a second cut.
+    #[test]
+    fn a_cut_keeps_the_same_choices_with_a_group_built_otherwise() {
+        let mut input = small_input(
+            vec![vec![0, 1], vec![0], vec![1]],
+            &[1, 2],
+            vec![
+                choice(0, &[], 10.0),
+                choice(1, &[0], 1.0),
+                choice(2, &[0], 1.0),
+            ],
+        );
+        input.group_slots = vec![GroupSlot {
+            group: GroupId(0),
+            may_recompute: true,
+        }];
+        let first = solve(&input).unwrap().unwrap();
+        assert_eq!(first.built_once, [GroupId(0)].into());
+        input.excluded.push((first.chosen.clone(), Vec::new(), 3.0));
+        input.excluded_once.push((vec![GroupId(0)], Vec::new()));
+        let second = solve(&input).unwrap().unwrap();
+        assert_eq!(second.chosen, first.chosen);
+        assert!(second.built_once.is_empty());
+        input.excluded.push((first.chosen.clone(), Vec::new(), 2.0));
+        input.excluded_once.push((Vec::new(), vec![GroupId(0)]));
+        assert_eq!(solve(&input).unwrap(), None);
+    }
+
+    /// A group a fused choice of another group implements has no node of
+    /// its own to build once: a group that must be built once when needed
+    /// is implemented by its own choice, though fusing is cheaper.
+    #[test]
+    fn a_group_built_once_is_implemented_by_its_own_choice() {
+        let mut fused = choice(1, &[], 1.0);
+        fused.implements = vec![GroupId(1), GroupId(0)];
+        let mut input = small_input(
+            vec![vec![0, 1], vec![0], vec![1]],
+            &[1, 2],
+            vec![
+                choice(0, &[], 10.0),
+                fused,
+                choice(1, &[0], 1.0),
+                choice(2, &[0], 1.0),
+            ],
+        );
+        input.group_slots = vec![GroupSlot {
+            group: GroupId(0),
+            may_recompute: false,
+        }];
+        let solution = solve(&input).unwrap().unwrap();
+        assert_eq!(solution.chosen, vec![0, 2, 3]);
+        assert_eq!(solution.built_once, [GroupId(0)].into());
+    }
+
+    /// A choice set cut with some of its plans not evaluated (here every
+    /// one: no combination of free nodes' options allowed) keeps its
+    /// objective as the lower bound, so the search is not certified; with
+    /// every combination evaluated it is.
+    #[test]
+    fn a_choice_set_with_unevaluated_plans_keeps_its_objective_as_the_bound() {
+        let queries = [
+            "sum by (job) (m)",
+            "sum by (job, instance) (m)",
+            "sum by (job, instance, zone) (m)",
+        ];
+        let targets = [
+            AccuracyTarget::exact(),
+            AccuracyTarget::exact(),
+            AccuracyTarget::exact(),
+        ];
+        let (per_computation, entries, data) = workload(&queries, &targets, None);
+        let all: Vec<usize> = (0..entries.len()).collect();
+        let run = |max_free_combinations: usize| {
+            let memo = crate::search::memo::brute_force::memo(&per_computation, &entries);
+            let mut cluster = ClusterLoop::new(
+                memo,
+                &per_computation,
+                &all,
+                &entries,
+                &data,
+                PlanningModels::builtin(),
+                std::time::Instant::now(),
+                None,
+            )
+            .unwrap();
+            cluster.max_free_combinations = max_free_combinations;
+            run_rounds(std::slice::from_mut(&mut cluster)).unwrap();
+            cluster.finish().unwrap().0.search.milp.unwrap()
+        };
+        let every = run(crate::search::milp::decode::MAX_FREE_COMBINATIONS);
+        assert!(
+            every.gap.is_some_and(|gap| gap <= GAP_TOLERANCE),
+            "{every:?}"
+        );
+        let none = run(0);
+        assert!(none.rounds > 1);
+        assert_eq!(none.lower_bound, none.first_lower_bound);
+        assert!(none.gap.unwrap() > 0.0, "{none:?}");
+    }
 }