diff --git a/aqpbm-cli/src/cli.rs b/aqpbm-cli/src/cli.rs
index 1caa964a..613a075e 100644
--- a/aqpbm-cli/src/cli.rs
+++ b/aqpbm-cli/src/cli.rs
@@ -236,6 +236,17 @@ pub struct SketchbenchArgs {
     /// merge exactly, so a gap is a defect; for KLL it is the result.
     #[arg(long, default_value_t = 2, help_heading = "Measurement content")]
     pub merge_shards: usize,
+    /// Which label columns the `hydra-*` rows group by, comma-separated
+    /// indices into the columns before the value column: any non-empty subset
+    /// (`0`, `1`, `0,1`). The grid stores every subset, so this picks the one
+    /// asked and scored. Every other target ignores it.
+    #[arg(
+        long,
+        value_delimiter = ',',
+        default_value = "0",
+        help_heading = "Measurement content"
+    )]
+    pub group_columns: Vec<usize>,
     /// How the merge operation splits the stream into its shards: `contiguous`
     /// cuts it into stretches, `interleaved` deals it round-robin by record, so
     /// each shard is a sample of the whole stream. The shard sizes are the same
@@ -248,17 +259,6 @@ pub struct SketchbenchArgs {
         help_heading = "Measurement content"
     )]
     pub merge_split: MergeSplit,
-    /// Which label columns the `hydra-*` rows group by, comma-separated
-    /// indices into the columns before the value column: any non-empty subset
-    /// (`0`, `1`, `0,1`). The grid stores every subset, so this picks the one
-    /// asked and scored. Every other target ignores it.
-    #[arg(
-        long,
-        value_delimiter = ',',
-        default_value = "0",
-        help_heading = "Measurement content"
-    )]
-    pub group_columns: Vec<usize>,
 
     /// Path to append JSONL records to. `-` or omitted sends them to stdout.
     #[arg(long, help_heading = "Output")]
diff --git a/aqpbm-cli/src/external.rs b/aqpbm-cli/src/external.rs
index 3f65ca8d..3d897fb4 100644
--- a/aqpbm-cli/src/external.rs
+++ b/aqpbm-cli/src/external.rs
@@ -271,6 +271,8 @@ fn shape_column(data_type: &str) -> ColumnSpec {
         string: None,
         child_of: None,
         fan_out: None,
+        scale_by: None,
+        scale_range: None,
     }
 }
 
diff --git a/aqpbm-cli/src/flatten_record.rs b/aqpbm-cli/src/flatten_record.rs
index 921b24ec..5937e173 100644
--- a/aqpbm-cli/src/flatten_record.rs
+++ b/aqpbm-cli/src/flatten_record.rs
@@ -259,6 +259,8 @@ mod tests {
                 string: None,
                 child_of: None,
                 fan_out: None,
+                scale_by: None,
+                scale_range: None,
             },
             1000,
         );
diff --git a/aqpbm-cli/src/main.rs b/aqpbm-cli/src/main.rs
index 2cd85e4b..a70fd611 100644
--- a/aqpbm-cli/src/main.rs
+++ b/aqpbm-cli/src/main.rs
@@ -260,6 +260,8 @@ fn dataset_spec(args: &SketchbenchArgs) -> Result<InputDataSetSpec> {
             },
             child_of: None,
             fan_out: None,
+            scale_by: None,
+            scale_range: None,
         },
         args.size.expect("clap requires --size") as u64,
     )))
@@ -374,9 +376,9 @@ fn run_sketchbench(args: SketchbenchArgs) -> Result<()> {
         width,
         workers: args.workers.max(1),
         merge_shards: args.merge_shards,
-        merge_split: args.merge_split,
         group_columns,
         per_group_out: args.per_group_out.as_ref().map(std::path::PathBuf::from),
+        merge_split: args.merge_split,
         comparator: args.comparator.clone(),
         runs: args.runs,
         warmup_runs: args.warmup_runs,
diff --git a/aqpbm-cli/src/repeat.rs b/aqpbm-cli/src/repeat.rs
index d90aaf63..4b613028 100644
--- a/aqpbm-cli/src/repeat.rs
+++ b/aqpbm-cli/src/repeat.rs
@@ -193,6 +193,8 @@ mod tests {
                 string: None,
                 child_of: None,
                 fan_out: None,
+                scale_by: None,
+                scale_range: None,
             },
             1000,
         );
diff --git a/aqpbm-core/schema/merged_record.schema.json b/aqpbm-core/schema/merged_record.schema.json
index eeb150a3..a8b77e9d 100644
--- a/aqpbm-core/schema/merged_record.schema.json
+++ b/aqpbm-core/schema/merged_record.schema.json
@@ -431,6 +431,28 @@
           "format": "uint64",
           "minimum": 0
         },
+        "scale_by": {
+          "description": "Index of an earlier `string` column whose values scale this one. Each\ndistinct label value gets a fixed factor in `scale_range`, and every\nrow's value is multiplied by its label's factor. `f64` columns only.",
+          "type": [
+            "integer",
+            "null"
+          ],
+          "format": "uint",
+          "minimum": 0
+        },
+        "scale_range": {
+          "description": "`[lo, hi]` with `0 < lo <= hi`: the factors are log-uniform in it.\nRequired with `scale_by`.",
+          "type": [
+            "array",
+            "null"
+          ],
+          "items": {
+            "type": "number",
+            "format": "double"
+          },
+          "maxItems": 2,
+          "minItems": 2
+        },
         "shift": {
           "description": "Added to every value, for a column that should not start where its\ndistribution naturally does. Under the monotonic rule it is the value the\nseries starts from.",
           "type": [
diff --git a/aqpbm-core/src/atomic_costs.rs b/aqpbm-core/src/atomic_costs.rs
index 449b0d2f..11ce9cdf 100644
--- a/aqpbm-core/src/atomic_costs.rs
+++ b/aqpbm-core/src/atomic_costs.rs
@@ -74,6 +74,15 @@ pub struct MeasuredAt {
     pub merge_operand_items: Option<u64>,
     /// The value column's distribution as generated, e.g. Zipf with its skew.
     pub distribution: Option<DataDistribution>,
+    /// The datagen spec (yaml stem) a `hydra-*` row was measured on, named by
+    /// the driver: the record carries the description, not its name.
+    #[serde(default, skip_serializing_if = "Option::is_none")]
+    pub dataset: Option<String>,
+    /// A `hydra-*` row's label columns `d`: each record is inserted into
+    /// `2^d - 1` subsets, so its insert cost holds at this width only. Absent
+    /// for every other row.
+    #[serde(default, skip_serializing_if = "Option::is_none")]
+    pub schema_width: Option<u32>,
 }
 
 impl MeasuredAt {
@@ -275,7 +284,10 @@ pub(crate) fn reduce_scores(record: &MergedRecord) -> Result<AtomicCostEntry, Sk
         insert_cpu_secs,
         merge_cpu_secs,
         query_cpu_secs,
-        measured_at: measured_at(record, groups_measured),
+        measured_at: MeasuredAt {
+            schema_width: query_accuracy.get("schema_width").map(|&d| d as u32),
+            ..measured_at(record, groups_measured)
+        },
         query_accuracy,
         accuracy_metric: String::new(),
     })
@@ -320,6 +332,8 @@ fn measured_at(record: &MergedRecord, groups_measured: Option<f64>) -> MeasuredA
             .filter(|&m| m > 0)
             .map(|m| items / m as u64),
         distribution: value.map(|c| c.distribution.clone()),
+        dataset: None,
+        schema_width: None,
     }
 }
 
@@ -420,6 +434,8 @@ mod tests {
                 string: None,
                 child_of: None,
                 fan_out: None,
+                scale_by: None,
+                scale_range: None,
             },
             1_000_000,
         )
@@ -483,6 +499,25 @@ mod tests {
 
     /// A one-column table: its column gives the keys and the values, and the
     /// merge operand is the stream split `merge_shards` ways.
+    /// A Hydra row's insert is per record at its schema width, which the
+    /// row names; a row without the score names none.
+    #[test]
+    fn measured_at_carries_a_hydra_rows_schema_width() {
+        let mut record = full_record();
+        assert_eq!(
+            reduce_one(&record, METRIC)
+                .unwrap()
+                .measured_at
+                .schema_width,
+            None
+        );
+        record.query.accuracy =
+            Some(serde_json::json!({"relative_error_mean": 0.1, "schema_width": 4.0}));
+        let entry = reduce_one(&record, METRIC).expect("fully populated record");
+        assert_eq!(entry.measured_at.schema_width, Some(4));
+        assert!((entry.insert_cpu_secs - 0.5e-6).abs() < 1e-12);
+    }
+
     #[test]
     fn measured_at_reads_a_single_column_workload() {
         let mut record = full_record();
@@ -721,6 +756,8 @@ mod tests {
                 value_range: None,
                 merge_operand_items: None,
                 distribution: None,
+                dataset: None,
+                schema_width: None,
             },
         };
         let json = serde_json::to_string(&entry).unwrap();
diff --git a/aqpbm-core/src/benchmark_result/schema.rs b/aqpbm-core/src/benchmark_result/schema.rs
index 1ce5a7df..50a14022 100644
--- a/aqpbm-core/src/benchmark_result/schema.rs
+++ b/aqpbm-core/src/benchmark_result/schema.rs
@@ -525,6 +525,8 @@ mod tests {
                 string: None,
                 child_of: None,
                 fan_out: None,
+                scale_by: None,
+                scale_range: None,
             },
             1_000_000,
         )
diff --git a/aqpbm-core/src/erp.rs b/aqpbm-core/src/erp.rs
index 9bddc906..8c5d43ca 100644
--- a/aqpbm-core/src/erp.rs
+++ b/aqpbm-core/src/erp.rs
@@ -183,6 +183,8 @@ mod tests {
                     string: None,
                     child_of: None,
                     fan_out: None,
+                    scale_by: None,
+                    scale_range: None,
                 },
                 100_000,
             ),
diff --git a/aqpbm-datagen/src/column.rs b/aqpbm-datagen/src/column.rs
index 092bea87..072e0f91 100644
--- a/aqpbm-datagen/src/column.rs
+++ b/aqpbm-datagen/src/column.rs
@@ -45,6 +45,15 @@ pub struct ColumnSpec {
     /// checked against the distribution's domain like `cardinality` is.
     #[serde(default, skip_serializing_if = "Option::is_none")]
     pub fan_out: Option<u64>,
+    /// Index of an earlier `string` column whose values scale this one. Each
+    /// distinct label value gets a fixed factor in `scale_range`, and every
+    /// row's value is multiplied by its label's factor. `f64` columns only.
+    #[serde(default, skip_serializing_if = "Option::is_none")]
+    pub scale_by: Option<usize>,
+    /// `[lo, hi]` with `0 < lo <= hi`: the factors are log-uniform in it.
+    /// Required with `scale_by`.
+    #[serde(default, skip_serializing_if = "Option::is_none")]
+    pub scale_range: Option<[f64; 2]>,
 }
 
 impl ColumnSpec {
@@ -77,6 +86,35 @@ impl ColumnSpec {
             (Some(_), Some(fan_out)) => self.validate_child(fan_out)?,
         }
 
+        match (self.scale_by, self.scale_range) {
+            (None, None) => {}
+            (Some(_), None) => {
+                return Err(DataGenError::BadParam(
+                    "scale_by: needs scale_range, the [lo, hi] the factors are drawn from".into(),
+                ))
+            }
+            (None, Some(_)) => {
+                return Err(DataGenError::BadParam(
+                    "scale_range: applies only to a scale_by column".into(),
+                ))
+            }
+            (Some(_), Some([lo, hi])) => {
+                // An integer column would round the scaled value back down.
+                if self.data_type != "f64" {
+                    return Err(DataGenError::BadParam(format!(
+                        "scale_by: a scaled value is fractional, so data_type has to be \
+                         `f64`, not '{}'",
+                        self.data_type
+                    )));
+                }
+                if !(lo > 0.0 && hi >= lo && hi.is_finite()) {
+                    return Err(DataGenError::BadParam(format!(
+                        "scale_range [{lo}, {hi}] needs 0 < lo <= hi"
+                    )));
+                }
+            }
+        }
+
         // `cardinality` restates the domain. Silently preferring one over the
         // other would let a report cite a key space the column never had.
         if let Some(claimed) = self.cardinality {
@@ -246,6 +284,19 @@ impl ColumnSpec {
         }
     }
 
+    /// The fixed factor of one `scale_by` label value: a SplitMix64 hash of
+    /// the value's FNV-1a and this column's seed, as `u` in `[0, 1)`, mapped to
+    /// `lo·(hi/lo)^u`.
+    pub(crate) fn scale_factor(&self, label: &str) -> f64 {
+        let [lo, hi] = self.scale_range.expect("scale_range was validated above");
+        let fnv = label.bytes().fold(0xCBF2_9CE4_8422_2325u64, |h, b| {
+            (h ^ b as u64).wrapping_mul(0x0000_0100_0000_01B3)
+        });
+        let x = crate::table::mix_connected_bits(fnv ^ self.distribution.seed());
+        let u = (x >> 11) as f64 / (1u64 << 53) as f64;
+        lo * (hi / lo).powf(u)
+    }
+
     /// Render a `child_of` column: each row's parent value, a `.`, and the
     /// draw's 0-based rank in the domain as the child index.
     pub(crate) fn render_child(
diff --git a/aqpbm-datagen/src/lib.rs b/aqpbm-datagen/src/lib.rs
index e90ac559..91e63d6c 100644
--- a/aqpbm-datagen/src/lib.rs
+++ b/aqpbm-datagen/src/lib.rs
@@ -63,6 +63,8 @@ mod tests {
             string: None,
             child_of: None,
             fan_out: None,
+            scale_by: None,
+            scale_range: None,
         }
     }
 
@@ -662,6 +664,122 @@ mod tests {
         }
     }
 
+    // ---------- scale_by ----------
+
+    /// service (25 labels) and a latency scaled by it, next to the same
+    /// latency unscaled.
+    fn scaled_table() -> TableDescription {
+        TableDescription {
+            column_num: 3,
+            column_label: vec!["service".into(), "plain".into(), "scaled".into()],
+            column_spec: vec![
+                column(zipf(25, 1.1, 1), "string"),
+                column(pareto(2.0, 1000.0, 5), "f64"),
+                ColumnSpec {
+                    scale_by: Some(0),
+                    scale_range: Some([1.0, 10.0]),
+                    ..column(pareto(2.0, 1000.0, 5), "f64")
+                },
+            ],
+            column_connected: Vec::new(),
+            row_num: 20_000,
+        }
+    }
+
+    /// Every row is its unscaled draw times one factor per label, the factors
+    /// stay in the range and spread across it.
+    #[test]
+    fn each_label_scales_its_rows_by_one_factor_in_the_range() {
+        let d = scaled_table();
+        let t = d.generate().unwrap();
+        assert_eq!(t, d.generate().unwrap());
+        let service = t.data[0].as_string().unwrap();
+        let plain = t.data[1].as_f64().unwrap();
+        let scaled = t.data[2].as_f64().unwrap();
+        let mut factors: std::collections::HashMap<&str, f64> = Default::default();
+        for ((s, p), v) in service.iter().zip(plain).zip(scaled) {
+            let f = *factors
+                .entry(s)
+                .or_insert_with(|| d.column_spec[2].scale_factor(s));
+            assert_eq!(*v, p * f, "'{s}'");
+        }
+        assert_eq!(factors.len(), 25);
+        let (min, max) = factors
+            .values()
+            .fold((f64::MAX, f64::MIN), |(a, b), &f| (a.min(f), b.max(f)));
+        assert!(min >= 1.0 && max <= 10.0, "{min}..{max}");
+        // 25 log-uniform factors over a decade span most of it.
+        assert!(min < 2.0 && max > 5.0, "{min}..{max}");
+    }
+
+    /// The factor depends on the label and the seed only.
+    #[test]
+    fn the_scale_factor_is_fixed_by_label_and_seed() {
+        let spec = scaled_table().column_spec[2].clone();
+        assert_eq!(spec.scale_factor("abcd"), spec.scale_factor("abcd"));
+        assert_ne!(spec.scale_factor("abcd"), spec.scale_factor("abce"));
+        let reseeded = ColumnSpec {
+            distribution: pareto(2.0, 1000.0, 6),
+            ..spec.clone()
+        };
+        assert_ne!(spec.scale_factor("abcd"), reseeded.scale_factor("abcd"));
+        let flat = ColumnSpec {
+            scale_range: Some([3.0, 3.0]),
+            ..spec
+        };
+        assert_eq!(flat.scale_factor("abcd"), 3.0);
+    }
+
+    #[test]
+    fn bad_scaled_columns_are_refused_by_name() {
+        type Edit = fn(&mut TableDescription);
+        let cases: Vec<(&str, Edit)> = vec![
+            ("earlier column", |d| d.column_spec[2].scale_by = Some(2)),
+            ("earlier column", |d| d.column_spec[2].scale_by = Some(9)),
+            ("'plain', which is not a string", |d| {
+                d.column_spec[2].scale_by = Some(1)
+            }),
+            ("needs scale_range", |d| d.column_spec[2].scale_range = None),
+            ("only to a scale_by", |d| d.column_spec[2].scale_by = None),
+            ("has to be `f64`, not 'i64'", |d| {
+                d.column_spec[2].data_type = "i64".into()
+            }),
+            ("has to be `f64`, not 'u64'", |d| {
+                d.column_spec[2].data_type = "u64".into()
+            }),
+            ("0 < lo <= hi", |d| {
+                d.column_spec[2].scale_range = Some([0.0, 10.0])
+            }),
+            ("0 < lo <= hi", |d| {
+                d.column_spec[2].scale_range = Some([5.0, 2.0])
+            }),
+        ];
+        for (want, edit) in cases {
+            let mut d = scaled_table();
+            edit(&mut d);
+            let err = d.generate().unwrap_err().to_string();
+            assert!(err.contains(want), "wanted '{want}' in: {err}");
+        }
+    }
+
+    /// The eval's latency spec scales by service and leaves the labels alone.
+    #[test]
+    fn hydra_http_latency_scales_by_service() {
+        let path = std::path::Path::new("../configs/datagen/hydra_http_latency.yaml");
+        let d = TableDescription::from_path(path).unwrap();
+        assert_eq!(d.column_label[1], "service");
+        assert_eq!(d.column_spec[4].scale_by, Some(1));
+        let mut plain = d.clone();
+        plain.column_spec[4].scale_by = None;
+        plain.column_spec[4].scale_range = None;
+        plain.row_num = 1000;
+        let mut scaled = d;
+        scaled.row_num = 1000;
+        let (a, b) = (plain.generate().unwrap(), scaled.generate().unwrap());
+        assert_eq!(a.data[..4], b.data[..4]);
+        assert_ne!(a.data[4], b.data[4]);
+    }
+
     // ---------- serde ----------
 
     #[test]
diff --git a/aqpbm-datagen/src/table.rs b/aqpbm-datagen/src/table.rs
index cb64c73d..611c71dc 100644
--- a/aqpbm-datagen/src/table.rs
+++ b/aqpbm-datagen/src/table.rs
@@ -16,9 +16,9 @@ use crate::{inject_interval_bursts, BurstSpec};
 
 /// Scramble a latent rank before dividing it among connected columns. This
 /// SplitMix64 finalizer spreads nearby ranks across the word while keeping the
-/// mapping deterministic.
+/// mapping deterministic. Also hashes a `scale_by` label to its factor.
 #[inline]
-fn mix_connected_bits(mut x: u64) -> u64 {
+pub(crate) fn mix_connected_bits(mut x: u64) -> u64 {
     x = x.wrapping_add(0x9E37_79B9_7F4A_7C15);
     x = (x ^ (x >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9);
     x = (x ^ (x >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB);
@@ -140,6 +140,25 @@ impl TableDescription {
             }
         }
 
+        // Like a parent, the scaling label is generated first.
+        for (i, spec) in self.column_spec.iter().enumerate() {
+            let Some(by) = spec.scale_by else {
+                continue;
+            };
+            let label = &self.column_label[i];
+            if by >= i {
+                return Err(DataGenError::BadParam(format!(
+                    "column '{label}': scale_by {by} has to name an earlier column, one of 0..{i}"
+                )));
+            }
+            if self.column_spec[by].data_type != "string" {
+                return Err(DataGenError::BadParam(format!(
+                    "column '{label}': scale_by names '{}', which is not a string column",
+                    self.column_label[by]
+                )));
+            }
+        }
+
         let mut claimed: HashMap<&str, usize> = HashMap::new();
         for (g, group) in self.column_connected.iter().enumerate() {
             if group.len() < 2 {
@@ -265,10 +284,18 @@ impl TableDescription {
                     &own_draw
                 }
             };
-            let rendered = match spec.child_of {
+            let mut rendered = match spec.child_of {
                 Some(parent) => spec.render_child(raw, data[parent].as_string()?),
                 None => spec.render(raw),
             };
+            if let (Some(by), Ok(ColumnData::Float64(values))) = (spec.scale_by, &mut rendered) {
+                let mut factors: HashMap<&str, f64> = HashMap::new();
+                for (v, label) in values.iter_mut().zip(data[by].as_string()?) {
+                    *v *= *factors
+                        .entry(label)
+                        .or_insert_with(|| spec.scale_factor(label));
+                }
+            }
             data.push(rendered.map_err(|e| {
                 DataGenError::BadParam(format!("column '{}': {e}", self.column_label[i]))
             })?);
diff --git a/configs/datagen/hydra_flows.yaml b/configs/datagen/hydra_flows.yaml
new file mode 100644
index 00000000..5a365913
--- /dev/null
+++ b/configs/datagen/hydra_flows.yaml
@@ -0,0 +1,63 @@
+# The §6.3 eval's `flows` metric (export_autosketch_eval_table.py SCHEMAS,
+# #196), for `study_saturation.py --phase hydra`. Labels in schema order, then
+# the value:
+#
+#   dst_subnet 1000 values, Zipf 1.1
+#   dst_port 1000 values, Zipf 1.2
+#   proto 3 values, equally likely (#196's skew 0.0)
+#   src_ip uniform over 1e6 (distinct sources)
+#
+# Not generated: #196's DDoS burst, where the least popular subnet (rank 999)
+# gets 5% of the traffic from 1e4 distinct sources. datagen draws each column
+# independently, so it cannot tie a value's sources to one label value. The
+# burst subnet is one group of share >= 5% with few distinct sources; the
+# study's 5%-coverage column here reads the top Zipf subnets instead.
+#
+#   approxbench sketchbench --variant hydra-hll --library lib \
+#     --spec configs/datagen/hydra_flows.yaml --dtype i64 \
+#     --config "rows=3 cols=1024" --group-columns 0 \
+#     --operations query --metrics accuracy
+#
+# Three labels cost `2^3 - 1 = 7` cell insertions per record.
+
+column_num: 4
+column_label: [dst_subnet, dst_port, proto, src_ip]
+row_num: 1000000
+
+column_spec:
+  - data_type: string
+    string:
+      min_len: 4
+      max_len: 4
+    distribution:
+      kind: zipf
+      skewness: 1.1
+      population_size: 1000
+      seed: 1
+
+  - data_type: string
+    string:
+      min_len: 4
+      max_len: 4
+    distribution:
+      kind: zipf
+      skewness: 1.2
+      population_size: 1000
+      seed: 2
+
+  - data_type: string
+    string:
+      min_len: 2
+      max_len: 2
+    distribution:
+      kind: uniform
+      lower_bound: 0.0
+      upper_bound: 3.0
+      seed: 3
+
+  - data_type: i64
+    distribution:
+      kind: uniform
+      lower_bound: 0.0
+      upper_bound: 1000000.0
+      seed: 4
diff --git a/configs/datagen/hydra_http.yaml b/configs/datagen/hydra_http.yaml
new file mode 100644
index 00000000..d0ed960a
--- /dev/null
+++ b/configs/datagen/hydra_http.yaml
@@ -0,0 +1,69 @@
+# The §6.3 eval's `http` metric (export_autosketch_eval_table.py SCHEMAS, #196),
+# for `study_saturation.py --phase hydra`: the eval reads Hydra's accuracy on
+# data with its own label schema and skews, since a group's error depends on
+# its share of the fanned-out mass. Labels in schema order, then the value:
+#
+#   region 4 values, Zipf 0.5
+#   service 25 values, Zipf 1.1 (independent of region, as #196 has it)
+#   endpoint 25 per service (child_of service), Zipf 1.1
+#   status 4 values, Zipf 2.0
+#   user_id Zipf 0.8 over 1e6 (distinct users)
+#
+#   approxbench sketchbench --variant hydra-hll --library lib \
+#     --spec configs/datagen/hydra_http.yaml --dtype i64 \
+#     --config "rows=3 cols=1024" --group-columns 0,1 \
+#     --operations query --metrics accuracy
+#
+# Four labels cost `2^4 - 1 = 15` cell insertions per record. The study sets
+# `row_num` per N and offsets every seed per study seed.
+
+column_num: 5
+column_label: [region, service, endpoint, status, user_id]
+row_num: 1000000
+
+column_spec:
+  - data_type: string
+    string:
+      min_len: 4
+      max_len: 4
+    distribution:
+      kind: zipf
+      skewness: 0.5
+      population_size: 4
+      seed: 1
+
+  - data_type: string
+    string:
+      min_len: 4
+      max_len: 4
+    distribution:
+      kind: zipf
+      skewness: 1.1
+      population_size: 25
+      seed: 2
+
+  - data_type: string
+    child_of: 1
+    fan_out: 25
+    distribution:
+      kind: zipf
+      skewness: 1.1
+      population_size: 25
+      seed: 3
+
+  - data_type: string
+    string:
+      min_len: 3
+      max_len: 3
+    distribution:
+      kind: zipf
+      skewness: 2.0
+      population_size: 4
+      seed: 4
+
+  - data_type: i64
+    distribution:
+      kind: zipf
+      skewness: 0.8
+      population_size: 1000000
+      seed: 5
diff --git a/configs/datagen/hydra_http_latency.yaml b/configs/datagen/hydra_http_latency.yaml
new file mode 100644
index 00000000..147fbb7c
--- /dev/null
+++ b/configs/datagen/hydra_http_latency.yaml
@@ -0,0 +1,67 @@
+# `hydra_http.yaml` with the latency value instead of user_id, for the eval's
+# quantile RQEs (template 16, p99 latency): the same labels, skews and seeds,
+# and an f64 value drawn Pareto(a = 2, scale 1000), the classic set's tail
+# (COST_PARETO_ALPHA), times a per-service factor: each service gets a fixed
+# factor log-uniform in [1, 10] (`scale_by: 1`), so services differ in latency
+# scale and keep the tail index, as in #196's runner.
+#
+# Not modeled: #196 also makes one endpoint 10x slower; there is no injected
+# anomaly here.
+#
+#   approxbench sketchbench --variant hydra-kll --library lib \
+#     --spec configs/datagen/hydra_http_latency.yaml --dtype f64 \
+#     --config "rows=3 cols=1024 cell_k=200" --group-columns 1 \
+#     --operations query --metrics accuracy
+
+column_num: 5
+column_label: [region, service, endpoint, status, latency]
+row_num: 1000000
+
+column_spec:
+  - data_type: string
+    string:
+      min_len: 4
+      max_len: 4
+    distribution:
+      kind: zipf
+      skewness: 0.5
+      population_size: 4
+      seed: 1
+
+  - data_type: string
+    string:
+      min_len: 4
+      max_len: 4
+    distribution:
+      kind: zipf
+      skewness: 1.1
+      population_size: 25
+      seed: 2
+
+  - data_type: string
+    child_of: 1
+    fan_out: 25
+    distribution:
+      kind: zipf
+      skewness: 1.1
+      population_size: 25
+      seed: 3
+
+  - data_type: string
+    string:
+      min_len: 3
+      max_len: 3
+    distribution:
+      kind: zipf
+      skewness: 2.0
+      population_size: 4
+      seed: 4
+
+  - data_type: f64
+    scale_by: 1          # service
+    scale_range: [1.0, 10.0]
+    distribution:
+      kind: pareto
+      alpha: 2.0
+      scale: 1000.0
+      seed: 5
diff --git a/docs/aqpbm-datagen.md b/docs/aqpbm-datagen.md
index 3ff1f839..4c925d9f 100644
--- a/docs/aqpbm-datagen.md
+++ b/docs/aqpbm-datagen.md
@@ -37,6 +37,8 @@ pub struct ColumnSpec {
     pub string: Option<StringOpts>, // if this column is about string data, description about how string looks like
     pub child_of: Option<usize>, // index of an earlier string column this one nests under; see child_of below
     pub fan_out: Option<u64>, // children per parent value; required with child_of
+    pub scale_by: Option<usize>, // index of an earlier string column whose values scale this f64 column; see scale_by below
+    pub scale_range: Option<[f64; 2]>, // [lo, hi] the per-label factors are log-uniform in; required with scale_by
 }
 
 enum DataDistribution {
@@ -110,6 +112,31 @@ Each of the following is a description error:
 
 `configs/datagen/hydra_hier.yaml` uses it for a region → service → endpoint hierarchy.
 
+### scale_by
+
+A value column drawn independently of the labels gives every group the same distribution.
+A `data_type: f64` column with `scale_by: <index>` and `scale_range: [lo, hi]` multiplies each row's value
+by a factor fixed per distinct value of the earlier string column at that index.
+The factor is `lo·(hi/lo)^u`, log-uniform in `[lo, hi]`, with `u ∈ [0, 1)` a hash of the label value and the column's seed:
+the same label always gets the same factor, whatever the row count, and a new seed gives new factors.
+The draws themselves are unchanged, so a Pareto column keeps its tail index and each group's scale moves.
+
+```yaml
+- data_type: f64
+  scale_by: 1        # each service's latencies × its own factor in [1, 10]
+  scale_range: [1.0, 10.0]
+  distribution: {kind: pareto, alpha: 2.0, scale: 1000.0, seed: 5}
+```
+
+Each of the following is a description error:
+
+- `scale_by` without `scale_range`, or `scale_range` without `scale_by`
+- `scale_by` naming the column itself, a later column, or no column, or a column that is not `data_type: string`
+- a scaled column that is not `f64` (an integer column would round the scaled value)
+- a `scale_range` without `0 < lo <= hi`
+
+`configs/datagen/hydra_http_latency.yaml` scales latency by service.
+
 ## Output
 
 The output is a struct containing data:
diff --git a/docs/saturation_study.md b/docs/saturation_study.md
index a2ad6f7d..f0319949 100644
--- a/docs/saturation_study.md
+++ b/docs/saturation_study.md
@@ -264,6 +264,68 @@ point in `crossover.csv`. The plotting scripts draw one config per sketch
 (rows=3 cols=1024, lg_k=12, k=200, alpha=0.01, the configs their bounds are
 written for).
 
+## Hydra (`--phase hydra`)
+
+A Hydra grid's error for a group depends on the group's share of the
+fanned-out mass `N·(2^d − 1)`, not on a per-key shape, so it is measured on
+data with the §6.3 eval's own schemas rather than on the Zipf grid above
+(docs/rqe_optimizer_hydra.md). `HYDRA_SWEEPS` holds two sweeps:
+
+- `eval`: `hydra-hll` and `hydra-univmon-cardinality` on
+  `configs/datagen/hydra_http.yaml` and `hydra_flows.yaml`, `hydra-kll` on
+  `hydra_http_latency.yaml` (the specs follow `export_autosketch_eval_table.py`'s
+  `SCHEMAS`; their headers say what is not modelled). Latency is Pareto(a = 2,
+  scale 1000) times a fixed per-service factor, log-uniform in `[1, 10]`
+  (`scale_by`, docs/aqpbm-datagen.md; the factors follow the study seed); no
+  endpoint anomaly is injected. Grid `R = 3`,
+  `W ∈ {1024, 4096, 16384}`, every non-empty label subset (15 for http, 7 for
+  flows), `N ∈ {1e5, 1e6, 1e7}`, merged from `{1, 4, 16}` interleaved shards,
+  seeds 1..3.
+- `generic`: every statistic (`HYDRA_CELLS`) on `hydra_hier_d2.yaml`,
+  `hydra_hier_d3.yaml`, `hydra_hier.yaml` (`d` = 2, 3, 4), every subset,
+  `W ∈ {1024, 4096}`, `N = 1e6`, shards `{1, 4}`, seeds 1..2.
+
+Inner cells are the canonical ones, except UnivMon's (`HYDRA_UNIVMON_CELL`,
+3 × 256 × 8 layers, heap 64): the canonical 5 × 2048 cell makes a
+3 × 16384 grid 32 GB. A study seed `s` rewrites the spec's `row_num` and adds
+`1000·(s − 1)` to every seed, so seed 1 is the file as written.
+
+Each run asks one grouping with `--group-columns` and `--per-group-out`, and
+adds one row to `hydra_saturation.csv`: `variant, statistic, config, R, W,
+dataset, group_columns` (label names joined by `,`, in schema order; quoted),
+`schema_width, records, fanned_mass, merge_shards, merge_split` (empty for
+the plain query), `seed`, the record's `err_mean, err_p50, err_p90, err_max,
+groups_scored` (in the family's own metric), and `err_max_cov_0.01`,
+`err_max_cov_0.05`: the worst error over the groups holding at least that
+share `n_q / records` of the records, plus always the largest group (an
+unscored group has no error but can be the largest). Rows are written as runs
+finish; `--resume` keeps them and runs the rest, and a failed run is reported
+and rerun by the next `--resume`.
+
+`--jobs` runs in parallel, admitted while the running runs' estimated memory
+(`m` shard grids plus the data, `hydra_mem_bytes`) fits `--hydra-mem-gb`.
+`--hydra-shard i/n` takes every `n`-th run of the fixed run list from the
+`i`-th, so `n` machines split a sweep; concatenate their tables afterwards.
+`--hydra-ns`, `--hydra-ws`, `--hydra-seeds`, `--hydra-variants`,
+`--hydra-datasets` and `--hydra-groupings` narrow the sweeps.
+
+```sh
+# Node i of 5 (0-based), 48 parallel runs:
+python3 scripts/study_saturation.py --phase hydra --jobs 48 \
+    --hydra-shard i/5 --out out_hydra_i         # add --resume after an interruption
+# Then, in one place:
+awk 'FNR > 1 || NR == 1' out_hydra_*/hydra_saturation.csv > hydra_saturation.csv
+```
+
+`--phase optimizer-cost` also measures the eval variants at each eval `W` on
+their datasets (full schema, `COST_N` records, 3 runs after 1 warm-up, `--families`
+includes `hydra`):
+`insert_cpu_secs` per record at that schema width, `merge_cpu_secs` per grid
+merge, `query_cpu_secs` per probe (a group, or a group and quantile), memory
+per grid; `measured_at` names the `dataset` and `schema_width`. Their accuracy
+column is one pass's `err_max` at the default grouping; the optimizer reads
+Hydra accuracy from `hydra_saturation.csv`.
+
 ## Results (config grid, N up to 1e7)
 
 `--phase accuracy --n-max 1e7 --seeds 3 --jobs 16 --out out_grid_1e7`: 595
diff --git a/rqe-optimizer/src/lib.rs b/rqe-optimizer/src/lib.rs
index e7868d5d..a3b59131 100644
--- a/rqe-optimizer/src/lib.rs
+++ b/rqe-optimizer/src/lib.rs
@@ -628,6 +628,8 @@ pub(crate) mod test_support {
             value_range: None,
             merge_operand_items: None,
             distribution: None,
+            dataset: None,
+            schema_width: None,
         }
     }
 
diff --git a/rqe-optimizer/src/saturation.rs b/rqe-optimizer/src/saturation.rs
index bbf623ee..cea60c71 100644
--- a/rqe-optimizer/src/saturation.rs
+++ b/rqe-optimizer/src/saturation.rs
@@ -1320,6 +1320,8 @@ mod tests {
                     population_size: 1_000,
                     seed: 1,
                 })),
+                dataset: None,
+                schema_width: None,
             },
             ..deployment(TOPK, 1).config
         }
diff --git a/rqe-optimizer/tests/mixed_grouping.rs b/rqe-optimizer/tests/mixed_grouping.rs
index b33c564d..59c5f3d5 100644
--- a/rqe-optimizer/tests/mixed_grouping.rs
+++ b/rqe-optimizer/tests/mixed_grouping.rs
@@ -92,6 +92,8 @@ fn cost_row(sketch: &str, params: serde_json::Value) -> AtomicCostEntry {
             value_range: None,
             merge_operand_items: None,
             distribution: None,
+            dataset: None,
+            schema_width: None,
         },
     }
 }
diff --git a/scripts/study_saturation.py b/scripts/study_saturation.py
index 82ba8fad..fe8a1a4b 100644
--- a/scripts/study_saturation.py
+++ b/scripts/study_saturation.py
@@ -30,6 +30,8 @@
 The grid also holds every config at the optimizer-cost shape unless
 --no-cost-shape: θ = 1.1 and K = 1e4 join --thetas and --cardinalities (a full
 row and column, so the grid stays a full cross) and a = 2 joins --alphas.
+--phase hydra measures Hydra's per-group error on the eval's datasets and a
+generic hierarchy (HYDRA_SWEEPS) into hydra_saturation.csv; see `hydra`.
 --resume keeps the complete curves of an interrupted accuracy run;
 --cost-rows 3 times only the rows=3 Vector2D configs.
 
@@ -46,16 +48,22 @@
       seed-mean error per N of the sketch merged from m contiguous shards
       (m=1 is the plain single-sketch query)
   rqe_atomic_costs.json (+ _raw, _grid .jsonl)       --phase optimizer-cost only
+  hydra_saturation.csv (+ hydra_accuracy.jsonl)      --phase hydra only: one row
+      per run, HYDRA_COLUMNS
 """
 
 import argparse
 import csv
+import io
+import itertools
 import json
 import math
 import os
+import re
 import subprocess
 import sys
-from concurrent.futures import ThreadPoolExecutor
+import threading
+from concurrent.futures import ThreadPoolExecutor, as_completed
 
 FREQ_CONFIGS = [f"rows={r} cols={c}" for r in (3, 5) for c in (256, 1024, 4096, 16384)]
 
@@ -122,6 +130,81 @@
 ]
 EXACT_METRIC = "relative_error"
 
+# --phase hydra (docs/rqe_optimizer_hydra.md; ASAPQuery#777 §6.3): Hydra's
+# per-group error depends on the group's share of the fanned-out mass, so the
+# eval reads accuracy measured on datasets with its own schemas
+# (configs/datagen/hydra_http*.yaml, hydra_flows.yaml, matching
+# export_autosketch_eval_table.py's SCHEMAS). Each run asks one grouping (by
+# label name, every non-empty subset of the schema) at grid R x W, N records,
+# merged from m interleaved shards (1 = the plain query), one seed, and writes
+# its per-group errors (--per-group-out), from which the coverage columns are
+# read. The other statistics run on the generic hydra_hier*.yaml sweep
+# (schema width 2, 3, 4).
+HYDRA_R = 3
+# Inner cells: the canonical ones (docs/rqe_optimizer_hydra.md §1), except
+# UnivMon's, whose 5 x 2048 x 8-layer cell (~650 KB) makes a 3 x 16384 grid
+# 32 GB; this one is ~100 KB.
+HYDRA_UNIVMON_CELL = "cell_heap_size=64 cell_sketch_row=3 cell_sketch_col=256 cell_layer_size=8"
+HYDRA_CELLS = {
+    "hydra-cms": "cell_rows=3 cell_cols=512",
+    "hydra-cs": "cell_rows=3 cell_cols=512",
+    "hydra-hll": "",
+    "hydra-kll": "cell_k=200",
+    "hydra-kll-cdf": "cell_k=200",
+    **{f"hydra-univmon-{s}": HYDRA_UNIVMON_CELL
+       for s in ("cardinality", "l1-norm", "l2-norm", "entropy", "sum")},
+}
+# The statistic each variant answers, a column of hydra_saturation.csv.
+HYDRA_STATISTICS = {
+    "hydra-cms": "frequency", "hydra-cs": "frequency", "hydra-hll": "cardinality",
+    "hydra-kll": "quantile", "hydra-kll-cdf": "cdf",
+    "hydra-univmon-cardinality": "cardinality", "hydra-univmon-l1-norm": "l1_norm",
+    "hydra-univmon-l2-norm": "l2_norm", "hydra-univmon-entropy": "entropy",
+    "hydra-univmon-sum": "sum",
+}
+# Bytes per grid cell, measured at N = 1e6 (data and truth excluded), and per
+# record of generated data plus its exact per-group truth: --hydra-mem-gb
+# admits a run only while the running ones' estimates fit. A merge holds all
+# m shard grids at once.
+HYDRA_CELL_BYTES = {"hydra-cms": 6_500, "hydra-cs": 6_500, "hydra-hll": 17_000,
+                    "hydra-kll": 10_000, "hydra-kll-cdf": 10_000}
+HYDRA_UNIVMON_CELL_BYTES = 100_000
+HYDRA_RECORD_BYTES = 250
+# (sweep, datasets as (yaml stem, value dtype, variants), W, N, merge shards,
+# seeds). hydra_http_latency scales each service's latency by its own factor,
+# log-uniform in [1, 10] and redrawn per seed (hydra_spec offsets the seed);
+# no endpoint anomaly is injected.
+HYDRA_SWEEPS = {
+    "eval": ([("hydra_http", "i64", ["hydra-hll", "hydra-univmon-cardinality"]),
+              ("hydra_flows", "i64", ["hydra-hll", "hydra-univmon-cardinality"]),
+              ("hydra_http_latency", "f64", ["hydra-kll"])],
+             (1024, 4096, 16384), (100_000, 1_000_000, 10_000_000), (1, 4, 16), 3),
+    "generic": ([(stem, "i64", list(HYDRA_CELLS))
+                 for stem in ("hydra_hier_d2", "hydra_hier_d3", "hydra_hier")],
+                (1024, 4096), (1_000_000,), (1, 4), 2),
+}
+HYDRA_DATAGEN = os.path.join(os.path.dirname(os.path.abspath(__file__)), "..", "configs",
+                             "datagen")
+HYDRA_TABLE = "hydra_saturation.csv"
+# Coverage thresholds τ: err_max_cov_τ is the worst error over the groups of
+# share >= τ plus the largest group (the eval's covers_share rule).
+HYDRA_COVERS = (0.01, 0.05)
+HYDRA_COLUMNS = [
+    "variant", "statistic", "config", "R", "W", "dataset", "group_columns", "schema_width",
+    "records", "fanned_mass", "merge_shards", "merge_split", "seed", "err_mean", "err_p50",
+    "err_p90", "err_max", "groups_scored",
+] + [f"err_max_cov_{tau}" for tau in HYDRA_COVERS]
+# What makes a row the same measurement, for --resume.
+HYDRA_KEY = ["variant", "config", "dataset", "group_columns", "records", "merge_shards",
+             "seed"]
+# The optimizer-cost table's Hydra rows: the eval's variants on its datasets
+# at each eval W. Their accuracy column is the default grouping's err_max; the
+# optimizer reads Hydra accuracy from hydra_saturation.csv.
+HYDRA_COST_METRIC = "err_max"
+# Fewer passes than COST_RUNS: one pass inserts COST_N records into up to
+# 2^4 - 1 subsets each, and a UnivMon grid's merge pass rebuilds 16 of them.
+HYDRA_COST_RUNS = (3, 1)
+
 SUMMARY_COLUMNS = [
     "family", "sketch", "config", "dist", "param", "cardinality", "n_sat",
     "final_error", "error_metric", "insert_cpu_secs", "merge_cpu_secs",
@@ -333,6 +416,40 @@ def optimizer_cost(args, families):
     for variant, comparator, spec in EXACT_COST_ROWS:
         metrics[variant] = EXACT_METRIC
         rows.append((variant, "exact", [], comparator, ["--spec", spec, "--dtype", "i64"]))
+    # Hydra: each eval variant on its datasets at the full schema, at each
+    # eval W, COST_N records (the spec as written, seed 1). Insert is per
+    # record at that schema width, merge per grid, query per probe, memory
+    # per grid; measured_at names the dataset.
+    hydra_datasets = {}
+    if "hydra" in families:
+        datasets, ws = HYDRA_SWEEPS["eval"][0], HYDRA_SWEEPS["eval"][1]
+        for stem, dtype, variants in datasets:
+            with open(os.path.join(args.datagen, f"{stem}.yaml")) as f:
+                text = f.read()
+            spec = os.path.join(args.out, f"{stem}_n{COST_N}.yaml")
+            with open(spec, "w") as f:
+                f.write(hydra_spec(text, COST_N, 1))
+            for variant in variants:
+                # --hydra-variants / --hydra-ws narrow these rows as they do
+                # the hydra phase (e.g. to leave out a variant whose 16-shard
+                # merge at the largest W outgrows the machine).
+                if args.hydra_variants and variant not in args.hydra_variants.split(","):
+                    continue
+                metrics[variant] = HYDRA_COST_METRIC
+                hydra_ws = {int(x) for x in args.hydra_ws.split(",")} if args.hydra_ws else None
+                for w in ws[:1] if args.one_config else ws:
+                    if hydra_ws and w not in hydra_ws:
+                        continue
+                    config = hydra_config(variant, w)
+                    rows.append((variant, "lib", ["--config", config], None,
+                                 ["--spec", spec, "--dtype", dtype]))
+                    # The reducer reads the schema width off the scores; the
+                    # dataset is the one thing the record cannot name.
+                    key = (variant, w, len(hydra_labels(text)))
+                    if key in hydra_datasets:
+                        sys.exit(f"Hydra cost rows {key} on {hydra_datasets[key]} and {stem} "
+                                 "cannot be told apart")
+                    hydra_datasets[key] = stem
 
     def sketchbench(variant, library, config, data, passes):
         subprocess.run(
@@ -343,13 +460,15 @@ def sketchbench(variant, library, config, data, passes):
     for variant, library, config, comparator, data in rows:
         print(f"  {variant} {' '.join(config[1:])}", file=sys.stderr)
         # Cost first: flatten keeps the first pass's query timings.
+        runs, warmup = HYDRA_COST_RUNS if variant.startswith("hydra-") else \
+            (COST_RUNS, COST_WARMUP)
         sketchbench(variant, library, config, data, [
             "--operations", "insert,query,merge", "--metrics", "throughput,cpu,memory",
-            "--merge-shards", str(MERGE_SHARDS), "--runs", str(COST_RUNS),
-            "--warmup-runs", str(COST_WARMUP)])
+            "--merge-shards", str(MERGE_SHARDS), "--runs", str(runs),
+            "--warmup-runs", str(warmup)])
         sketchbench(variant, library, config, data, [
-            "--operations", "query", "--metrics", "accuracy", "--comparator", comparator,
-            "--runs", "1", "--warmup-runs", "0"])
+            "--operations", "query", "--metrics", "accuracy", "--runs", "1",
+            "--warmup-runs", "0"] + (["--comparator", comparator] if comparator else []))
 
     subprocess.run([args.binary, "flatten", raw, "--output", grid], check=True)
     reduce = subprocess.run(
@@ -367,6 +486,17 @@ def sketchbench(variant, library, config, data, passes):
     if summary != expected:
         sys.exit(f"atomic-costs did not keep all {len(rows)} rows: {summary!r}")
     seed_mean_accuracy(args, table, rows, metrics)
+    if hydra_datasets:
+        with open(table) as f:
+            entries = json.load(f)
+        for entry in entries:
+            at = entry["measured_at"]
+            if entry["sketch"].startswith("hydra-"):
+                at["dataset"] = hydra_datasets[(
+                    entry["sketch"], int(entry["sketch_config"]["params"]["cols"]),
+                    at["schema_width"])]
+        with open(table, "w") as f:
+            json.dump(entries, f, indent=2)
     print(f"Done. {table}", file=sys.stderr)
     return 0
 
@@ -380,8 +510,10 @@ def params(config):
 
     means = {}
     for variant, library, config, comparator, data in rows:
-        if library == "exact" or float(params(config[1]).get("heap", TOPK_HEAPS[0])) != \
-                TOPK_HEAPS[0]:
+        # Hydra rows keep their one pass: the optimizer reads Hydra's accuracy
+        # from hydra_saturation.csv.
+        if library == "exact" or variant.startswith("hydra-") or \
+                float(params(config[1]).get("heap", TOPK_HEAPS[0])) != TOPK_HEAPS[0]:
             continue
         errors = [error(run(args.binary, [
             "--variant", variant, "--library", library, *config,
@@ -401,13 +533,264 @@ def params(config):
         json.dump(entries, f, indent=2)
 
 
+def hydra_labels(spec_text):
+    """The label names of a Hydra datagen spec: its column_label before the
+    value column, in schema order."""
+    match = re.search(r"^column_label:\s*\[(.*)\]", spec_text, re.M)
+    return [name.strip() for name in match.group(1).split(",")][:-1]
+
+
+def hydra_spec(spec_text, records, seed):
+    """`spec_text` at `records` rows, every distribution seed offset by
+    1000 · (seed − 1), so seed 1 is the file as written and each study seed
+    draws every column afresh. Plain text edits: the specs are this repo's."""
+    text = re.sub(r"^row_num:.*$", f"row_num: {records}", spec_text, flags=re.M)
+    return re.sub(r"^(\s*seed:\s*)(\d+)\s*$",
+                  lambda m: f"{m.group(1)}{int(m.group(2)) + 1000 * (seed - 1)}", text,
+                  flags=re.M)
+
+
+def hydra_groupings(labels):
+    """Every non-empty subset of `labels`, in schema order, smallest first."""
+    return [list(c) for n in range(1, len(labels) + 1)
+            for c in itertools.combinations(labels, n)]
+
+
+def hydra_config(variant, w):
+    return " ".join(filter(None, [f"rows={HYDRA_R} cols={w}", HYDRA_CELLS[variant]]))
+
+
+def hydra_jobs(sweeps, specs, records=None):
+    """Every run of `sweeps`, in one fixed order (so `--hydra-shard i/n` picks
+    the same runs on every node): (dataset, dtype, variant, config, W,
+    grouping, records, merge shards, seed). `specs` maps a yaml stem to its
+    text; `records`, if given, replaces each sweep's record counts."""
+    jobs = []
+    for sweep in sweeps:
+        datasets, ws, ns, shards, seeds = HYDRA_SWEEPS[sweep]
+        ns = records or ns
+        for stem, dtype, variants in datasets:
+            groupings = hydra_groupings(hydra_labels(specs[stem]))
+            for variant, w, grouping, n, m, seed in itertools.product(
+                    variants, ws, groupings, ns, shards, range(1, seeds + 1)):
+                jobs.append((stem, dtype, variant, hydra_config(variant, w), w, grouping, n,
+                             m, seed))
+    return jobs
+
+
+def hydra_mem_bytes(job):
+    """A run's estimated peak memory: its m shard grids plus the data."""
+    _, _, variant, _, w, _, n, m, _ = job
+    cell = HYDRA_CELL_BYTES.get(variant, HYDRA_UNIVMON_CELL_BYTES)
+    return m * HYDRA_R * w * cell + n * HYDRA_RECORD_BYTES
+
+
+def read_per_group(path):
+    """[(n_q, error or None)] from a --per-group-out CSV; an empty error is a
+    group the comparator could not score."""
+    with open(path, newline="") as f:
+        return [(int(r["n_q"]), float(r["error"]) if r["error"] else None)
+                for r in csv.DictReader(f)]
+
+
+def covered_max(groups, records, tau):
+    """The worst error over the groups holding >= `tau` of the `records`,
+    plus always the largest group (so no coverage is empty): the eval's
+    covers_share rule. Unscored groups count for "largest" but have no
+    error; "" when no covered group was scored."""
+    if not groups:
+        return ""
+    largest = max(range(len(groups)), key=lambda i: groups[i][0])
+    errors = [e for i, (n_q, e) in enumerate(groups)
+              if e is not None and (i == largest or n_q / records >= tau)]
+    return max(errors) if errors else ""
+
+
+def hydra_row(job, labels, record, groups):
+    """One hydra_saturation.csv row from a run's record and per-group file."""
+    stem, _, variant, config, w, grouping, n, m, seed = job
+    acc = record["bench"]["accuracy"]
+    records = int(acc["records"])
+    return {
+        "variant": variant, "statistic": HYDRA_STATISTICS[variant], "config": config,
+        "R": HYDRA_R, "W": w, "dataset": stem, "group_columns": ",".join(grouping),
+        "schema_width": int(acc["schema_width"]), "records": records,
+        "fanned_mass": int(acc["fanned_mass"]), "merge_shards": m,
+        "merge_split": record["bench"].get("merge_split", ""), "seed": seed,
+        **{k: acc.get(k, "") for k in ("err_mean", "err_p50", "err_p90", "err_max")},
+        "groups_scored": int(acc["groups_scored"]),
+        **{f"err_max_cov_{tau}": covered_max(groups, records, tau) for tau in HYDRA_COVERS},
+    }
+
+
+# hydra_saturation.csv's integer columns; the error columns are numbers or
+# "" (no scored group), merge_split is "" for a query run, the rest are text.
+HYDRA_INT_COLUMNS = ["R", "W", "schema_width", "records", "fanned_mass", "merge_shards",
+                     "seed", "groups_scored"]
+
+
+def hydra_row_parses(row):
+    """Whether every column of a hydra_saturation.csv row holds a value of
+    its type (a row cut off mid-write does not)."""
+    try:
+        for c in HYDRA_COLUMNS:
+            if row.get(c) is None:
+                return False
+            if c in HYDRA_INT_COLUMNS:
+                int(row[c])
+            elif c.startswith("err_"):
+                if row[c] != "":
+                    float(row[c])
+            elif c != "merge_split" and not row[c]:
+                return False
+    except ValueError:
+        return False
+    return True
+
+
+def hydra_key(row):
+    """A hydra_saturation.csv row's identity, from a row dict (values as the
+    CSV holds them)."""
+    return tuple(str(row[c]) for c in HYDRA_KEY)
+
+
+class MemoryBudget:
+    """Admits a run while the running runs' estimates fit in `total` bytes;
+    a run larger than the whole budget runs alone."""
+
+    def __init__(self, total):
+        self.total, self.used, self.cond = total, 0, threading.Condition()
+
+    def acquire(self, n):
+        with self.cond:
+            self.cond.wait_for(lambda: self.used == 0 or self.used + n <= self.total)
+            self.used += n
+
+    def release(self, n):
+        with self.cond:
+            self.used -= n
+            self.cond.notify_all()
+
+
+def hydra(args):
+    """--phase hydra: every run of --hydra-sweeps (this node's share with
+    --hydra-shard i/n), one row each in --out/hydra_saturation.csv, rows
+    written as runs finish. --resume keeps the rows already there and runs
+    the rest; a failed run is reported and the phase exits 1 at the end."""
+    os.environ.setdefault("BENCH_WARMUP_SECS", "0")
+    # --jobs runs share the cores; polars' own pool would oversubscribe them.
+    os.environ.setdefault("POLARS_MAX_THREADS", "1")
+    os.makedirs(args.out, exist_ok=True)
+    specs = {}
+    for sweep in args.hydra_sweeps.split(","):
+        for stem, _, _ in HYDRA_SWEEPS[sweep][0]:
+            with open(os.path.join(args.datagen, f"{stem}.yaml")) as f:
+                specs[stem] = f.read()
+    ns = [int(float(n)) for n in args.hydra_ns.split(",")] if args.hydra_ns else None
+    jobs = hydra_jobs(args.hydra_sweeps.split(","), specs, ns)
+    ws = {int(w) for w in args.hydra_ws.split(",")} if args.hydra_ws else None
+    groupings = ({tuple(g.split(",")) for g in args.hydra_groupings.split(";")}
+                 if args.hydra_groupings else None)
+    jobs = [j for j in jobs
+            if (ws is None or j[4] in ws)
+            and (groupings is None or tuple(j[5]) in groupings)
+            and j[8] <= (args.hydra_seeds or j[8])
+            and (not args.hydra_variants or j[2] in args.hydra_variants.split(","))
+            and (not args.hydra_datasets or j[0] in args.hydra_datasets.split(","))]
+    shard, n_shards = (int(x) for x in args.hydra_shard.split("/"))
+    jobs = jobs[shard::n_shards]
+
+    table = os.path.join(args.out, HYDRA_TABLE)
+    kept = []
+    if args.resume and os.path.exists(table):
+        with open(table, newline="") as f:
+            text = f.read()
+        # A row cut off mid-write (no line end, or a field that does not
+        # parse) reruns.
+        kept = [r for r in csv.DictReader(io.StringIO(text[:text.rfind("\n") + 1]))
+                if hydra_row_parses(r)]
+    done = {hydra_key(r) for r in kept}
+    todo = [j for j in jobs if hydra_key({
+        "variant": j[2], "config": j[3], "dataset": j[0], "group_columns": ",".join(j[5]),
+        "records": j[6], "merge_shards": j[7], "seed": j[8]}) not in done]
+    # Largest first, so the long runs do not trail the pool.
+    todo.sort(key=lambda j: (j[6], j[7], j[4]), reverse=True)
+    print(f"{len(jobs)} runs in shard {shard}/{n_shards}, {len(jobs) - len(todo)} kept, "
+          f"{len(todo)} to run", file=sys.stderr)
+
+    spec_dir = os.path.join(args.out, "hydra_specs")
+    group_dir = os.path.join(args.out, "hydra_groups")
+    os.makedirs(spec_dir, exist_ok=True)
+    os.makedirs(group_dir, exist_ok=True)
+    for stem, n, seed in {(j[0], j[6], j[8]) for j in todo}:
+        with open(os.path.join(spec_dir, f"{stem}_n{n}_s{seed}.yaml"), "w") as f:
+            f.write(hydra_spec(specs[stem], n, seed))
+    budget = MemoryBudget(args.hydra_mem_gb * 1e9)
+
+    def measure(i, job):
+        stem, dtype, variant, config, _, grouping, n, m, seed = job
+        labels = hydra_labels(specs[stem])
+        per_group = os.path.join(group_dir, f"{shard}_{i}.csv")
+        operation = ["query"] if m == 1 else [
+            "merge", "--merge-shards", str(m), "--merge-split", "interleaved"]
+        need = hydra_mem_bytes(job)
+        budget.acquire(need)
+        try:
+            record = run(args.binary, [
+                "--variant", variant, "--library", "lib", "--config", config,
+                "--spec", os.path.join(spec_dir, f"{stem}_n{n}_s{seed}.yaml"),
+                "--dtype", dtype,
+                "--group-columns", ",".join(str(labels.index(g)) for g in grouping),
+                "--per-group-out", per_group, "--operations", *operation,
+                "--metrics", "accuracy", "--runs", "1", "--warmup-runs", "0"])
+        finally:
+            budget.release(need)
+        groups = read_per_group(per_group)
+        os.remove(per_group)
+        return record, hydra_row(job, labels, record, groups)
+
+    # The kept rows go to a new file that then replaces the table, so a
+    # crash mid-rewrite cannot lose them; new rows append to it.
+    with open(table + ".tmp", "w", newline="") as f, \
+            open(os.path.join(args.out, "hydra_accuracy.jsonl"), "a") as raw:
+        writer = csv.DictWriter(f, HYDRA_COLUMNS)
+        writer.writeheader()
+        writer.writerows(kept)
+        f.flush()
+        os.fsync(f.fileno())
+        os.replace(table + ".tmp", table)
+        failed = 0
+        with ThreadPoolExecutor(max_workers=args.jobs) as pool:
+            futures = {pool.submit(measure, i, j): j for i, j in enumerate(todo)}
+            for k, future in enumerate(as_completed(futures), 1):
+                job = futures[future]
+                try:
+                    record, row = future.result()
+                except (subprocess.CalledProcessError, OSError, KeyError, ValueError) as e:
+                    failed += 1
+                    detail = getattr(e, "stderr", "") or e
+                    print(f"  FAILED {job}: {str(detail).strip()[-500:]}", file=sys.stderr)
+                    continue
+                raw.write(json.dumps(record) + "\n")
+                writer.writerow(row)
+                f.flush()
+                if k % 50 == 0 or k == len(todo):
+                    print(f"  {k}/{len(todo)} runs", file=sys.stderr)
+    print(f"Done. {table}", file=sys.stderr)
+    if failed:
+        print(f"{failed} run(s) failed; --resume reruns them", file=sys.stderr)
+        return 1
+    return 0
+
+
 def main():
     parser = argparse.ArgumentParser(
         description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter
     )
     parser.add_argument("--binary", default="./target/release/approxbench")
     parser.add_argument("--out", default="out")
-    parser.add_argument("--families", default="frequency,topk,cardinality,quantile")
+    parser.add_argument("--families", default="frequency,topk,cardinality,quantile,hydra",
+                        help="hydra: --phase optimizer-cost's Hydra rows")
     parser.add_argument("--one-config", action="store_true",
                         help="only the first config of each sketch")
     parser.add_argument("--thetas", default="0,0.5,0.8,1.0,1.2,1.5,2.0")
@@ -427,12 +810,14 @@ def main():
     parser.add_argument("--jobs", type=int, default=1,
                         help="parallel accuracy runs (cost runs are always serial)")
     parser.add_argument("--phase",
-                        choices=["accuracy", "cost", "crossover", "all", "optimizer-cost"],
+                        choices=["accuracy", "cost", "crossover", "all", "optimizer-cost",
+                                 "hydra"],
                         default="all",
                         help="cost reads saturation_curve.csv from --out; crossover "
                              "recomputes crossover.csv from the cost phase's JSONL in --out; "
                              "optimizer-cost writes the optimizer's cost table to --out "
-                             "(the grid's configs at one shape, no curves)")
+                             "(the grid's configs at one shape, no curves); hydra writes "
+                             "hydra_saturation.csv (HYDRA_SWEEPS, its own grid)")
     parser.add_argument("--points-from",
                         help="CSV with sketch,config,dist,param,cardinality columns "
                              "(e.g. a filtered saturation.csv): run only those points")
@@ -445,7 +830,28 @@ def main():
     parser.add_argument("--merge-shards-list", default="",
                         help="accuracy: e.g. 1,4,16,64: also score the sketch merged from "
                              "m shards, for each m>1, into saturation_merge_curve.csv")
+    parser.add_argument("--datagen", default=HYDRA_DATAGEN,
+                        help="hydra: where the sweeps' datagen specs are")
+    parser.add_argument("--hydra-sweeps", default="eval,generic",
+                        help="hydra: which of HYDRA_SWEEPS to run")
+    parser.add_argument("--hydra-shard", default="0/1",
+                        help="hydra: i/n runs every n-th run from the i-th (0-based), "
+                             "to split one sweep across n machines")
+    parser.add_argument("--hydra-mem-gb", type=float, default=200,
+                        help="hydra: start a run only while the running runs' estimated "
+                             "memory fits")
+    # Narrow the sweeps (smoke runs, reruns); each keeps only what it names.
+    parser.add_argument("--hydra-ns", help="hydra: these record counts instead")
+    parser.add_argument("--hydra-ws", help="hydra: only these grid widths")
+    parser.add_argument("--hydra-seeds", type=int, help="hydra: only seeds 1..this")
+    parser.add_argument("--hydra-variants", help="hydra: only these variants")
+    parser.add_argument("--hydra-datasets", help="hydra: only these yaml stems")
+    parser.add_argument("--hydra-groupings",
+                        help="hydra: only these groupings, ';'-separated, labels by name "
+                             "joined by ',' in schema order")
     args = parser.parse_args()
+    if args.phase == "hydra":
+        return hydra(args)
 
     families = args.families.split(",")
     thetas = [float(t) for t in args.thetas.split(",")]
diff --git a/scripts/test_study_saturation.py b/scripts/test_study_saturation.py
index c8166ba0..fdc5e226 100644
--- a/scripts/test_study_saturation.py
+++ b/scripts/test_study_saturation.py
@@ -8,6 +8,7 @@
 import csv
 import json
 import os
+import re
 import subprocess
 import sys
 import tempfile
@@ -16,7 +17,9 @@
 
 sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
 
-from study_saturation import checkpoints, n_saturation, n_star, read_curve  # noqa: E402
+from study_saturation import (  # noqa: E402
+    HYDRA_COLUMNS, checkpoints, covered_max, hydra_groupings, hydra_labels, hydra_spec,
+    n_saturation, n_star, read_curve, read_per_group)
 
 SCRIPT = os.path.join(os.path.dirname(os.path.abspath(__file__)), "study_saturation.py")
 
@@ -484,8 +487,14 @@ def test_a_grid_that_already_has_the_cost_shape_adds_nothing(self):
                 v = c[c.index("--variant") + 1]
                 cfg = c[c.index("--config") + 1] if "--config" in c else ""
                 params = {k: float(x) for k, x in (kv.split("=") for kv in cfg.split())}
+                # The reducer reads a Hydra row's width off its scores.
+                at = {}
+                if v.startswith("hydra-"):
+                    labels = open(c[c.index("--spec") + 1]).read().split("column_label: [")[1]
+                    at = {"schema_width": labels.split("]")[0].count(",")}
                 rows.append({"sketch": v, "sketch_config": {"params": params},
-                             "query_accuracy": {metrics[v]: 0.5}, "accuracy_metric": metrics[v]})
+                             "query_accuracy": {metrics[v]: 0.5}, "accuracy_metric": metrics[v],
+                             "measured_at": at})
         json.dump(rows, open(value("--output"), "w"))
         kept = int(os.environ["FAKE_KEPT"])
         sys.stderr.write(f"approxbench atomic-costs: {kept} row(s), 0 skipped\\n")
@@ -498,7 +507,7 @@ def test_a_grid_that_already_has_the_cost_shape_adds_nothing(self):
 
 
 class OptimizerCostTest(unittest.TestCase):
-    def run_phase(self, d, kept):
+    def run_phase(self, d, kept, families="topk,quantile"):
         binary = os.path.join(d, "approxbench")
         with open(binary, "w") as f:
             f.write("#!" + sys.executable + "\n" + FAKE_COST_APPROXBENCH)
@@ -506,7 +515,7 @@ def run_phase(self, d, kept):
         log = os.path.join(d, "log.jsonl")
         result = subprocess.run([
             sys.executable, SCRIPT, "--binary", binary, "--out", os.path.join(d, "out"),
-            "--phase", "optimizer-cost", "--families", "topk,quantile", "--one-config",
+            "--phase", "optimizer-cost", "--families", families, "--one-config",
             "--seeds", "3",
         ], capture_output=True, env={**os.environ, "FAKE_LOG": log, "FAKE_KEPT": str(kept)})
         with open(log) as f:
@@ -558,6 +567,34 @@ def test_sketch_rows_take_the_seed_mean_accuracy(self):
             if row["sketch"].startswith("exact-"):
                 self.assertEqual(row["query_accuracy"][row["accuracy_metric"]], 0.5)
 
+    def test_hydra_rows_on_the_eval_datasets_name_their_dataset(self):
+        with tempfile.TemporaryDirectory() as d:
+            # hydra-hll and hydra-univmon-cardinality on http and flows,
+            # hydra-kll on http_latency, at the first W; plus 5 exact.
+            result, calls = self.run_phase(d, 10, families="hydra")
+            with open(os.path.join(d, "out", "rqe_atomic_costs.json")) as f:
+                rows = json.load(f)
+        self.assertEqual(result.returncode, 0, result.stderr)
+        value = lambda c, flag: c[c.index(flag) + 1]
+        hydra = [c for c in calls if c[0] == "sketchbench" and "--report" in c
+                 and value(c, "--variant").startswith("hydra-")]
+        self.assertEqual(len(hydra), 10)
+        cost = hydra[0]
+        self.assertEqual((value(cost, "--runs"), value(cost, "--warmup-runs"),
+                          value(cost, "--merge-shards")), ("3", "1", "16"))
+        self.assertEqual(value(cost, "--config"), "rows=3 cols=1024")
+        self.assertTrue(value(cost, "--spec").endswith("hydra_http_n1000000.yaml"))
+        self.assertNotIn("--comparator", hydra[1])
+        # No seed runs: Hydra accuracy comes from hydra_saturation.csv.
+        self.assertFalse([c for c in calls if c[0] == "sketchbench" and "--report" not in c])
+        at = {(r["sketch"], r["measured_at"]["schema_width"]): r["measured_at"]["dataset"]
+              for r in rows if r["sketch"].startswith("hydra-")}
+        self.assertEqual(at, {
+            ("hydra-hll", 4): "hydra_http", ("hydra-hll", 3): "hydra_flows",
+            ("hydra-univmon-cardinality", 4): "hydra_http",
+            ("hydra-univmon-cardinality", 3): "hydra_flows",
+            ("hydra-kll", 4): "hydra_http_latency"})
+
     def test_a_skipped_row_fails(self):
         with tempfile.TemporaryDirectory() as d:
             result, _ = self.run_phase(d, 10)
@@ -565,5 +602,227 @@ def test_a_skipped_row_fails(self):
         self.assertIn(b"did not keep all 11 rows", result.stderr)
 
 
+class CoveredMaxTest(unittest.TestCase):
+    # (n_q, error) of 100 records: shares 0.6, 0.3, 0.06, 0.04.
+    GROUPS = [(60, 0.1), (30, 0.2), (6, 0.5), (4, 0.9)]
+
+    def test_groups_at_or_above_the_share_count(self):
+        self.assertEqual(covered_max(self.GROUPS, 100, 0.05), 0.5)
+        self.assertEqual(covered_max(self.GROUPS, 100, 0.04), 0.9)
+        self.assertEqual(covered_max(self.GROUPS, 100, 0.3), 0.2)
+
+    def test_the_largest_group_is_always_covered(self):
+        self.assertEqual(covered_max(self.GROUPS, 100, 0.9), 0.1)
+        # Wherever it sits in the file.
+        self.assertEqual(covered_max(self.GROUPS[::-1], 100, 0.9), 0.1)
+
+    def test_unscored_groups_have_no_error_but_can_be_the_largest(self):
+        groups = [(60, None), (30, 0.2), (6, 0.5)]
+        self.assertEqual(covered_max(groups, 100, 0.05), 0.5)
+        # The largest is unscored, so nothing covered has an error.
+        self.assertEqual(covered_max(groups, 100, 0.9), "")
+        # An unscored group does not hand "largest" to the next one.
+        self.assertEqual(covered_max([(60, None), (30, 0.2)], 100, 0.5), "")
+
+    def test_no_groups(self):
+        self.assertEqual(covered_max([], 100, 0.01), "")
+
+
+SPEC = textwrap.dedent("""\
+    column_num: 3
+    column_label: [region, service, user_id]
+    row_num: 200000
+    column_spec:
+      - data_type: string
+        distribution:
+          kind: uniform
+          seed: 1
+      - data_type: i64
+        distribution: {kind: zipf, skewness: 0.8, population_size: 10, seed: 5}
+""")
+
+
+class HydraSpecTest(unittest.TestCase):
+    def test_rows_and_seeds(self):
+        self.assertEqual(hydra_spec(SPEC, 1000, 1).replace("row_num: 1000", "row_num: 200000"),
+                         SPEC)
+        text = hydra_spec(SPEC, 1000, 3)
+        self.assertIn("row_num: 1000\n", text)
+        self.assertIn("      seed: 2001\n", text)
+        # A flow mapping's seed is not a line of its own; it stays.
+        self.assertIn("seed: 5}", text)
+
+    def test_labels_and_groupings_in_schema_order(self):
+        self.assertEqual(hydra_labels(SPEC), ["region", "service"])
+        self.assertEqual(hydra_groupings(["a", "b", "c"]),
+                         [["a"], ["b"], ["c"], ["a", "b"], ["a", "c"], ["b", "c"],
+                          ["a", "b", "c"]])
+
+    def test_the_eval_specs_match_their_schemas(self):
+        root = os.path.join(os.path.dirname(SCRIPT), "..", "configs", "datagen")
+
+        def read(stem):
+            with open(os.path.join(root, f"{stem}.yaml")) as f:
+                return f.read()
+        http = ["region", "service", "endpoint", "status"]
+        self.assertEqual(hydra_labels(read("hydra_http")), http)
+        self.assertEqual(hydra_labels(read("hydra_http_latency")), http)
+        self.assertEqual(hydra_labels(read("hydra_flows")), ["dst_subnet", "dst_port", "proto"])
+        self.assertEqual(len(hydra_groupings(http)), 15)
+
+        def labels(stem):
+            """Each label column's (cardinality, skew, child_of), as datagen
+            reads them; a uniform column's cardinality is its upper bound."""
+            text = read(stem).split("column_spec:")[1]
+            out = []
+            for column in text.split("\n  - ")[1:]:
+                field = dict(re.findall(r"(\w+): ([\w.]+)\n", column + "\n"))
+                if field["data_type"] != "string":
+                    continue
+                if field["kind"] == "uniform":
+                    out.append((int(float(field["upper_bound"])), 0.0, None))
+                else:
+                    child_of = field.get("child_of")
+                    out.append((int(field["population_size"]), float(field["skewness"]),
+                                int(child_of) if child_of else None))
+                    if child_of:
+                        self.assertEqual(field["fan_out"], field["population_size"])
+            return out
+        # PR #196's SCHEMAS (origin/eval/hydra-workload,
+        # scripts/export_autosketch_eval_table.py): (cardinality, skew,
+        # child_of as a column index) per label, in schema order.
+        http_schema = [(4, 0.5, None), (25, 1.1, None), (25, 1.1, 1), (4, 2.0, None)]
+        self.assertEqual(labels("hydra_http"), http_schema)
+        self.assertEqual(labels("hydra_http_latency"), http_schema)
+        self.assertEqual(labels("hydra_flows"), [(1000, 1.1, None), (1000, 1.2, None),
+                                                 (3, 0.0, None)])
+
+
+# A Hydra run: writes a fixed per-group file for the grouping asked, and a
+# record of N = 100 records over 2 labels. Merged runs score 0.1 worse.
+FAKE_HYDRA_APPROXBENCH = textwrap.dedent("""\
+    import json, os, sys
+    a = sys.argv
+    value = lambda flag: a[a.index(flag) + 1] if flag in a else None
+    with open(os.environ["FAKE_LOG"], "a") as f:
+        f.write(json.dumps(a[1:]) + "\\n")
+    if value("--variant") == "hydra-cms" and value("--group-columns") == "1":
+        sys.exit("boom")
+    m = int(value("--merge-shards") or 1)
+    worse = 0.1 if m > 1 else 0.0
+    groups = [(60, 0.1 + worse), (30, None), (6, 0.5 + worse), (4, 0.9 + worse)]
+    with open(value("--per-group-out"), "w") as f:
+        f.write("group_key,n_q,error\\n")
+        for i, (n, e) in enumerate(groups):
+            f.write(f'"label0:g{i},x",{n},{"" if e is None else e}\\n')
+    errs = sorted(e for _, e in groups if e is not None)
+    bench = {"accuracy": {"err_mean": sum(errs) / 3, "err_p50": errs[1], "err_p90": errs[2],
+                          "err_max": errs[2], "groups_scored": 3.0, "schema_width": 2.0,
+                          "records": 100.0, "fanned_mass": 300.0}}
+    if m > 1:
+        bench.update(merge_shards=m, merge_split=value("--merge-split"))
+    print(json.dumps({"bench": bench}))
+""")
+
+
+class HydraPhaseTest(unittest.TestCase):
+    def run_phase(self, d, *extra, out="out"):
+        binary = os.path.join(d, "approxbench")
+        with open(binary, "w") as f:
+            f.write("#!" + sys.executable + "\n" + FAKE_HYDRA_APPROXBENCH)
+        os.chmod(binary, 0o755)
+        datagen = os.path.join(d, "datagen")
+        os.makedirs(datagen, exist_ok=True)
+        for stem in ("hydra_hier_d2", "hydra_hier_d3", "hydra_hier"):
+            with open(os.path.join(datagen, f"{stem}.yaml"), "w") as f:
+                f.write(SPEC)
+        log = os.path.join(d, "log.jsonl")
+        result = subprocess.run([
+            sys.executable, SCRIPT, "--binary", binary, "--out", os.path.join(d, out),
+            "--phase", "hydra", "--datagen", datagen, "--hydra-sweeps", "generic",
+            "--hydra-datasets", "hydra_hier_d2", "--hydra-variants", "hydra-hll,hydra-cms",
+            "--hydra-ws", "1024", "--hydra-ns", "100", "--hydra-seeds", "1", "--jobs", "2",
+            *extra,
+        ], capture_output=True, env={**os.environ, "FAKE_LOG": log})
+        table = os.path.join(d, out, "hydra_saturation.csv")
+        with open(table, newline="") as f:
+            header = next(csv.reader(f))
+        with open(table, newline="") as f:
+            rows = list(csv.DictReader(f))
+        with open(log) as f:
+            return result, rows, header, [json.loads(line) for line in f]
+
+    def test_one_row_per_run_with_coverage_columns(self):
+        with tempfile.TemporaryDirectory() as d:
+            result, rows, header, calls = self.run_phase(d)
+        # 2 variants x 3 groupings x shards {1, 4}; hydra-cms at {service}
+        # fails (2 runs), is reported, and the phase exits 1.
+        self.assertEqual(result.returncode, 1)
+        self.assertIn(b"2 run(s) failed", result.stderr)
+        self.assertEqual(len(calls), 12)
+        self.assertEqual(header, HYDRA_COLUMNS)
+        self.assertEqual(len(rows), 10)
+        row = {(r["variant"], r["group_columns"], r["merge_shards"]): r for r in rows}
+        query = row[("hydra-hll", "region,service", "1")]
+        self.assertEqual(
+            {k: query[k] for k in ("statistic", "config", "R", "W", "dataset", "schema_width",
+                                   "records", "fanned_mass", "merge_split", "seed",
+                                   "err_max", "groups_scored")},
+            {"statistic": "cardinality", "config": "rows=3 cols=1024", "R": "3",
+             "W": "1024", "dataset": "hydra_hier_d2", "schema_width": "2", "records": "100",
+             "fanned_mass": "300", "merge_split": "", "seed": "1", "err_max": "0.9",
+             "groups_scored": "3"})
+        # Shares 0.6, 0.3 (unscored), 0.06, 0.04.
+        self.assertEqual((query["err_max_cov_0.01"], query["err_max_cov_0.05"]), ("0.9", "0.5"))
+        merged = row[("hydra-hll", "region,service", "4")]
+        self.assertEqual((merged["merge_split"], merged["err_max_cov_0.05"]),
+                         ("interleaved", "0.6"))
+        # Groupings are asked by column index, in schema order.
+        hll = [c for c in calls if c[c.index("--variant") + 1] == "hydra-hll"]
+        self.assertEqual(sorted({c[c.index("--group-columns") + 1] for c in hll}),
+                         ["0", "0,1", "1"])
+        spec = hll[0][hll[0].index("--spec") + 1]
+        self.assertTrue(spec.endswith("hydra_hier_d2_n100_s1.yaml"))
+
+    def test_shards_split_the_runs_and_resume_reruns_only_the_missing(self):
+        with tempfile.TemporaryDirectory() as d:
+            _, first, _, _ = self.run_phase(d, "--hydra-shard", "0/2")
+            os.remove(os.path.join(d, "log.jsonl"))
+            _, other, _, _ = self.run_phase(d, "--hydra-shard", "1/2", out="other")
+            os.remove(os.path.join(d, "log.jsonl"))
+            _, rows, _, calls = self.run_phase(d, "--resume")
+        self.assertEqual(len(first) + len(other), 10)
+        key = lambda r: (r["variant"], r["group_columns"], r["merge_shards"])
+        self.assertFalse({key(r) for r in first} & {key(r) for r in other})
+        # The rerun runs the 6 runs shard 0 lacked; the 2 failures fail again.
+        self.assertEqual(len(calls), 12 - len(first))
+        self.assertEqual(len(rows), 10)
+
+    def test_resume_reruns_a_cut_off_row_and_replaces_the_table(self):
+        with tempfile.TemporaryDirectory() as d:
+            self.run_phase(d)
+            os.remove(os.path.join(d, "log.jsonl"))
+            table = os.path.join(d, "out", "hydra_saturation.csv")
+            with open(table, newline="") as f:
+                lines = f.read().splitlines(keepends=True)
+            # Two rows cut off mid-write: one with every field there but an
+            # empty groups_scored, and the last one with no line end.
+            fields = lines[1].split(",")
+            fields[HYDRA_COLUMNS.index("groups_scored")] = ""
+            lines[1] = ",".join(fields)
+            lines[-1] = lines[-1].rstrip("\r\n")
+            with open(table, "w", newline="") as f:
+                f.writelines(lines)
+            # The rewrite goes to a new file, so the old one survives.
+            os.link(table, table + ".old")
+            _, rows, _, calls = self.run_phase(d, "--resume")
+            with open(table + ".old", newline="") as f:
+                self.assertEqual(f.read(), "".join(lines))
+            self.assertFalse(os.path.exists(table + ".tmp"))
+        # Those 2 rerun, and the 2 failures again.
+        self.assertEqual(len(calls), 4)
+        self.assertEqual(len(rows), 10)
+
+
 if __name__ == "__main__":
     unittest.main()
diff --git a/sketch-bench/src/request.rs b/sketch-bench/src/request.rs
index 67d89371..c3383d84 100644
--- a/sketch-bench/src/request.rs
+++ b/sketch-bench/src/request.rs
@@ -122,12 +122,12 @@ pub struct Requirement {
     pub workers: usize,
     /// How many shards a merge measurement folds
     pub merge_shards: usize,
-    /// How a merge measurement splits the stream into those shards
-    pub merge_split: crate::wrappers::MergeSplit,
     /// The label columns a grouped (`hydra-*`) row asks and scores, ascending
     pub group_columns: Vec<usize>,
     /// Where a grouped row writes each scored group's error, if anywhere
     pub per_group_out: Option<std::path::PathBuf>,
+    /// How a merge measurement splits the stream into those shards
+    pub merge_split: crate::wrappers::MergeSplit,
     /// A named comparator, or `None` for the row's default.
     pub comparator: Option<String>,
     /// Measured runs, and the warm-ups run and discarded before them