Skip to main content

gam_problem/
indexed_response.rs

1//! Structural geometry and likelihood measures for separable indexed responses.
2//!
3//! A response cell can be structurally absent (missing, not at risk, or outside
4//! a declared output domain) independently of the numerical weight attached to
5//! an observed cell. Conflating those concepts in a single zero-valued weight
6//! loses the response geometry and makes it impossible to validate event/risk
7//! sets or preserve them in a fitted model. This module keeps them distinct.
8
9use ndarray::{Array1, Array2, ArrayView1, ArrayView2};
10use std::error::Error;
11use std::fmt::{Display, Formatter};
12
13/// Invalid indexed-response geometry or likelihood measure.
14#[derive(Clone, Debug, PartialEq, Eq)]
15pub struct IndexedResponseError {
16    reason: String,
17}
18
19impl IndexedResponseError {
20    fn new(reason: impl Into<String>) -> Self {
21        Self {
22            reason: reason.into(),
23        }
24    }
25
26    /// Human-readable invariant violation.
27    pub fn reason(&self) -> &str {
28        &self.reason
29    }
30}
31
32impl Display for IndexedResponseError {
33    fn fmt(&self, formatter: &mut Formatter<'_>) -> std::fmt::Result {
34        formatter.write_str(&self.reason)
35    }
36}
37
38impl Error for IndexedResponseError {}
39
40/// A sorted, duplicate-free set of `(row, output)` cells in CSR form.
41///
42/// The representation costs `O(n_rows + n_cells)` and supports both sparse
43/// inclusion sets and sparse exclusion sets through [`StructuralCells`].
44#[derive(Clone, Debug, PartialEq, Eq)]
45pub struct IndexedCellSet {
46    n_rows: usize,
47    n_outputs: usize,
48    row_offsets: Vec<usize>,
49    output_indices: Vec<usize>,
50}
51
52impl IndexedCellSet {
53    /// Construct from CSR arrays. Output indices inside every row must be
54    /// strictly increasing; duplicates are rejected rather than silently
55    /// changing the declared response geometry.
56    pub fn new(
57        n_rows: usize,
58        n_outputs: usize,
59        row_offsets: Vec<usize>,
60        output_indices: Vec<usize>,
61    ) -> Result<Self, IndexedResponseError> {
62        if row_offsets.len() != n_rows.saturating_add(1) {
63            return Err(IndexedResponseError::new(format!(
64                "indexed cell row_offsets length {} does not equal n_rows + 1 = {}",
65                row_offsets.len(),
66                n_rows.saturating_add(1)
67            )));
68        }
69        if row_offsets.first().copied() != Some(0) {
70            return Err(IndexedResponseError::new(
71                "indexed cell row_offsets must begin at zero",
72            ));
73        }
74        if row_offsets.last().copied() != Some(output_indices.len()) {
75            return Err(IndexedResponseError::new(format!(
76                "indexed cell final row offset {:?} does not equal cell count {}",
77                row_offsets.last(),
78                output_indices.len()
79            )));
80        }
81        for row in 0..n_rows {
82            let start = row_offsets[row];
83            let end = row_offsets[row + 1];
84            if start > end || end > output_indices.len() {
85                return Err(IndexedResponseError::new(format!(
86                    "indexed cell row {row} has invalid CSR range {start}..{end} for {} cells",
87                    output_indices.len()
88                )));
89            }
90            let outputs = &output_indices[start..end];
91            for (position, &output) in outputs.iter().enumerate() {
92                if output >= n_outputs {
93                    return Err(IndexedResponseError::new(format!(
94                        "indexed cell row {row} output {output} is outside 0..{n_outputs}"
95                    )));
96                }
97                if position > 0 && outputs[position - 1] >= output {
98                    return Err(IndexedResponseError::new(format!(
99                        "indexed cell outputs in row {row} must be strictly increasing; found {} then {output}",
100                        outputs[position - 1]
101                    )));
102                }
103            }
104        }
105        Ok(Self {
106            n_rows,
107            n_outputs,
108            row_offsets,
109            output_indices,
110        })
111    }
112
113    /// Construct from unordered cell coordinates. Coordinates are sorted into
114    /// canonical row-major order; duplicates remain an error.
115    pub fn from_cells(
116        n_rows: usize,
117        n_outputs: usize,
118        mut cells: Vec<(usize, usize)>,
119    ) -> Result<Self, IndexedResponseError> {
120        cells.sort_unstable();
121        if let Some(pair) = cells.windows(2).find(|pair| pair[0] == pair[1]) {
122            return Err(IndexedResponseError::new(format!(
123                "indexed response cell ({}, {}) was declared more than once",
124                pair[0].0, pair[0].1
125            )));
126        }
127        let mut row_offsets = vec![0usize; n_rows.saturating_add(1)];
128        let mut output_indices = Vec::with_capacity(cells.len());
129        for (row, output) in cells {
130            if row >= n_rows {
131                return Err(IndexedResponseError::new(format!(
132                    "indexed cell row {row} is outside 0..{n_rows}"
133                )));
134            }
135            if output >= n_outputs {
136                return Err(IndexedResponseError::new(format!(
137                    "indexed cell row {row} output {output} is outside 0..{n_outputs}"
138                )));
139            }
140            row_offsets[row + 1] += 1;
141            output_indices.push(output);
142        }
143        for row in 0..n_rows {
144            row_offsets[row + 1] += row_offsets[row];
145        }
146        Self::new(n_rows, n_outputs, row_offsets, output_indices)
147    }
148
149    /// Number of rows in the declared response grid.
150    pub fn n_rows(&self) -> usize {
151        self.n_rows
152    }
153
154    /// Number of outputs in the declared response grid.
155    pub fn n_outputs(&self) -> usize {
156        self.n_outputs
157    }
158
159    /// Number of represented cells.
160    pub fn len(&self) -> usize {
161        self.output_indices.len()
162    }
163
164    /// Whether the set contains no cells.
165    pub fn is_empty(&self) -> bool {
166        self.output_indices.is_empty()
167    }
168
169    /// Sorted output indices represented on `row`.
170    pub fn row_outputs(&self, row: usize) -> Option<&[usize]> {
171        if row >= self.n_rows {
172            return None;
173        }
174        Some(&self.output_indices[self.row_offsets[row]..self.row_offsets[row + 1]])
175    }
176
177    /// Whether `(row, output)` belongs to the set.
178    pub fn contains(&self, row: usize, output: usize) -> bool {
179        self.position(row, output).is_some()
180    }
181
182    /// Row-major storage position of `(row, output)`, when present.
183    pub fn position(&self, row: usize, output: usize) -> Option<usize> {
184        if row >= self.n_rows || output >= self.n_outputs {
185            return None;
186        }
187        let start = self.row_offsets[row];
188        let end = self.row_offsets[row + 1];
189        self.output_indices[start..end]
190            .binary_search(&output)
191            .ok()
192            .map(|within_row| start + within_row)
193    }
194
195    fn validate_shape(&self, n_rows: usize, n_outputs: usize) -> Result<(), IndexedResponseError> {
196        if (self.n_rows, self.n_outputs) != (n_rows, n_outputs) {
197            return Err(IndexedResponseError::new(format!(
198                "indexed cell set geometry ({}, {}) does not match response geometry ({n_rows}, {n_outputs})",
199                self.n_rows, self.n_outputs
200            )));
201        }
202        Ok(())
203    }
204}
205
206/// Structural presence of cells in a separable response grid.
207#[derive(Clone, Copy, Debug)]
208pub enum StructuralCells<'a> {
209    /// Every `(row, output)` cell is structurally present.
210    All,
211    /// A dense Boolean activity mask with shape `(N, M)`.
212    Dense(ArrayView2<'a, bool>),
213    /// Only cells in the sparse set are structurally present.
214    Only(&'a IndexedCellSet),
215    /// Every cell except those in the sparse set is structurally present.
216    AllExcept(&'a IndexedCellSet),
217}
218
219/// Numerical likelihood weighting, independent of structural presence.
220#[derive(Clone, Copy, Debug)]
221pub enum LikelihoodWeights<'a> {
222    /// Unit weight for every structurally present cell.
223    Uniform,
224    /// One weight per observation row, shared across outputs.
225    ByRow(ArrayView1<'a, f64>),
226    /// One weight per response cell, shape `(N, M)`.
227    ByCell(ArrayView2<'a, f64>),
228}
229
230/// Owned structural presence for a separable response grid.
231///
232/// This is the lifetime-free counterpart of [`StructuralCells`]. Sparse
233/// inclusion and exclusion sets remain sparse; converting a borrowed measure
234/// for storage never silently allocates an `N × M` mask.
235#[derive(Clone, Debug, PartialEq, Eq)]
236pub enum OwnedStructuralCells {
237    All,
238    Dense(Array2<bool>),
239    Only(IndexedCellSet),
240    AllExcept(IndexedCellSet),
241}
242
243/// Owned numerical likelihood measure, independent of structural presence.
244#[derive(Clone, Debug, PartialEq)]
245pub enum OwnedLikelihoodWeights {
246    Uniform,
247    ByRow(Array1<f64>),
248    ByCell(Array2<f64>),
249}
250
251/// Owned values over a declared `(row, output)` grid.
252///
253/// `ConstantWithOverrides` is the canonical sparse-event representation: zero
254/// is implicit everywhere and only event cells are stored. Row and output
255/// broadcasts retain separable fields such as quadrature exposure or an
256/// endpoint-specific threshold without allocating the Cartesian grid.
257#[derive(Clone, Debug, PartialEq)]
258pub enum OwnedCellValues {
259    Dense(Array2<f64>),
260    ByRow {
261        values: Array1<f64>,
262        n_outputs: usize,
263    },
264    ByOutput {
265        n_rows: usize,
266        values: Array1<f64>,
267    },
268    Constant {
269        n_rows: usize,
270        n_outputs: usize,
271        value: f64,
272    },
273    ConstantWithOverrides {
274        n_rows: usize,
275        n_outputs: usize,
276        default: f64,
277        cells: IndexedCellSet,
278        values: Vec<f64>,
279    },
280}
281
282impl OwnedCellValues {
283    pub fn dense(values: Array2<f64>) -> Self {
284        Self::Dense(values)
285    }
286
287    /// Broadcast one value per row over every output.
288    pub fn by_row(values: Array1<f64>, n_outputs: usize) -> Self {
289        Self::ByRow { values, n_outputs }
290    }
291
292    /// Broadcast one value per output over every row.
293    pub fn by_output(n_rows: usize, values: Array1<f64>) -> Self {
294        Self::ByOutput { n_rows, values }
295    }
296
297    pub fn constant(n_rows: usize, n_outputs: usize, value: f64) -> Self {
298        Self::Constant {
299            n_rows,
300            n_outputs,
301            value,
302        }
303    }
304
305    pub fn constant_with_overrides(
306        n_rows: usize,
307        n_outputs: usize,
308        default: f64,
309        mut overrides: Vec<(usize, usize, f64)>,
310    ) -> Result<Self, IndexedResponseError> {
311        overrides.sort_unstable_by_key(|&(row, output, _)| (row, output));
312        if let Some(pair) = overrides
313            .windows(2)
314            .find(|pair| (pair[0].0, pair[0].1) == (pair[1].0, pair[1].1))
315        {
316            return Err(IndexedResponseError::new(format!(
317                "indexed value cell ({}, {}) was overridden more than once",
318                pair[0].0, pair[0].1,
319            )));
320        }
321        let cells = IndexedCellSet::from_cells(
322            n_rows,
323            n_outputs,
324            overrides
325                .iter()
326                .map(|&(row, output, _)| (row, output))
327                .collect(),
328        )?;
329        let values = overrides.into_iter().map(|(_, _, value)| value).collect();
330        Ok(Self::ConstantWithOverrides {
331            n_rows,
332            n_outputs,
333            default,
334            cells,
335            values,
336        })
337    }
338
339    pub fn n_rows(&self) -> usize {
340        match self {
341            Self::Dense(values) => values.nrows(),
342            Self::ByRow { values, .. } => values.len(),
343            Self::ByOutput { n_rows, .. } => *n_rows,
344            Self::Constant { n_rows, .. } | Self::ConstantWithOverrides { n_rows, .. } => *n_rows,
345        }
346    }
347
348    pub fn n_outputs(&self) -> usize {
349        match self {
350            Self::Dense(values) => values.ncols(),
351            Self::ByRow { n_outputs, .. } => *n_outputs,
352            Self::ByOutput { values, .. } => values.len(),
353            Self::Constant { n_outputs, .. } | Self::ConstantWithOverrides { n_outputs, .. } => {
354                *n_outputs
355            }
356        }
357    }
358
359    pub fn value(&self, row: usize, output: usize) -> Option<f64> {
360        if row >= self.n_rows() || output >= self.n_outputs() {
361            return None;
362        }
363        Some(match self {
364            Self::Dense(values) => values[[row, output]],
365            Self::ByRow { values, .. } => values[row],
366            Self::ByOutput { values, .. } => values[output],
367            Self::Constant { value, .. } => *value,
368            Self::ConstantWithOverrides {
369                default,
370                cells,
371                values,
372                ..
373            } => cells
374                .position(row, output)
375                .map(|position| values[position])
376                .unwrap_or(*default),
377        })
378    }
379}
380
381/// Lifetime-free structural geometry and likelihood measure.
382///
383/// The realized `(n_rows, n_outputs)` shape is part of the value. This makes a
384/// stored measure self-validating and prevents reusing a sparse activity set or
385/// weight vector against a different response grid.
386#[derive(Clone, Debug, PartialEq)]
387pub struct OwnedSeparableCellMeasure {
388    n_rows: usize,
389    n_outputs: usize,
390    structural: OwnedStructuralCells,
391    likelihood_weights: OwnedLikelihoodWeights,
392}
393
394impl OwnedSeparableCellMeasure {
395    /// Construct and validate an owned response measure.
396    pub fn new(
397        n_rows: usize,
398        n_outputs: usize,
399        structural: OwnedStructuralCells,
400        likelihood_weights: OwnedLikelihoodWeights,
401    ) -> Result<Self, IndexedResponseError> {
402        let measure = Self {
403            n_rows,
404            n_outputs,
405            structural,
406            likelihood_weights,
407        };
408        measure.as_borrowed().validate(n_rows, n_outputs)?;
409        Ok(measure)
410    }
411
412    /// All cells present with unit numerical weight.
413    pub fn uniform(n_rows: usize, n_outputs: usize) -> Self {
414        Self {
415            n_rows,
416            n_outputs,
417            structural: OwnedStructuralCells::All,
418            likelihood_weights: OwnedLikelihoodWeights::Uniform,
419        }
420    }
421
422    pub fn n_rows(&self) -> usize {
423        self.n_rows
424    }
425
426    pub fn n_outputs(&self) -> usize {
427        self.n_outputs
428    }
429
430    /// Borrow this owned value through the zero-copy evaluation API.
431    pub fn as_borrowed(&self) -> SeparableCellMeasure<'_> {
432        let structural = match &self.structural {
433            OwnedStructuralCells::All => StructuralCells::All,
434            OwnedStructuralCells::Dense(mask) => StructuralCells::Dense(mask.view()),
435            OwnedStructuralCells::Only(cells) => StructuralCells::Only(cells),
436            OwnedStructuralCells::AllExcept(cells) => StructuralCells::AllExcept(cells),
437        };
438        let likelihood_weights = match &self.likelihood_weights {
439            OwnedLikelihoodWeights::Uniform => LikelihoodWeights::Uniform,
440            OwnedLikelihoodWeights::ByRow(weights) => LikelihoodWeights::ByRow(weights.view()),
441            OwnedLikelihoodWeights::ByCell(weights) => LikelihoodWeights::ByCell(weights.view()),
442        };
443        SeparableCellMeasure::new(structural, likelihood_weights)
444    }
445
446    pub fn is_active(&self, row: usize, output: usize) -> bool {
447        self.as_borrowed().is_active(row, output)
448    }
449
450    pub fn active_weight(&self, row: usize, output: usize) -> Option<f64> {
451        self.as_borrowed().active_weight(row, output)
452    }
453
454    /// Visit every structurally active cell in deterministic row-major order.
455    /// Sparse inclusion geometry runs in `O(active cells)`; zero numerical
456    /// weights remain visible to the visitor because they do not erase model
457    /// structure.
458    pub fn try_for_each_active<E>(
459        &self,
460        mut visitor: impl FnMut(usize, usize, f64) -> Result<(), E>,
461    ) -> Result<(), E> {
462        let weight = |row: usize, output: usize| match &self.likelihood_weights {
463            OwnedLikelihoodWeights::Uniform => 1.0,
464            OwnedLikelihoodWeights::ByRow(weights) => weights[row],
465            OwnedLikelihoodWeights::ByCell(weights) => weights[[row, output]],
466        };
467        match &self.structural {
468            OwnedStructuralCells::All => {
469                for row in 0..self.n_rows {
470                    for output in 0..self.n_outputs {
471                        visitor(row, output, weight(row, output))?;
472                    }
473                }
474            }
475            OwnedStructuralCells::Dense(active) => {
476                for ((row, output), &is_active) in active.indexed_iter() {
477                    if is_active {
478                        visitor(row, output, weight(row, output))?;
479                    }
480                }
481            }
482            OwnedStructuralCells::Only(cells) => {
483                for row in 0..self.n_rows {
484                    for &output in cells
485                        .row_outputs(row)
486                        .expect("owned sparse cell geometry was validated at construction")
487                    {
488                        visitor(row, output, weight(row, output))?;
489                    }
490                }
491            }
492            OwnedStructuralCells::AllExcept(excluded) => {
493                for row in 0..self.n_rows {
494                    for output in 0..self.n_outputs {
495                        if !excluded.contains(row, output) {
496                            visitor(row, output, weight(row, output))?;
497                        }
498                    }
499                }
500            }
501        }
502        Ok(())
503    }
504}
505
506/// Structural activity plus numerical likelihood weights for a separable
507/// response. An inactive cell has no likelihood contribution; an active cell
508/// with weight zero remains an observed member of the response geometry.
509#[derive(Clone, Copy, Debug)]
510pub struct SeparableCellMeasure<'a> {
511    /// Which response cells exist in the likelihood.
512    pub structural: StructuralCells<'a>,
513    /// Numerical measure applied to existing cells.
514    pub likelihood_weights: LikelihoodWeights<'a>,
515}
516
517impl<'a> SeparableCellMeasure<'a> {
518    /// All cells present with unit numerical weight.
519    pub const fn uniform() -> Self {
520        Self {
521            structural: StructuralCells::All,
522            likelihood_weights: LikelihoodWeights::Uniform,
523        }
524    }
525
526    /// All cells present with a shared weight per row.
527    pub const fn row_weighted(weights: ArrayView1<'a, f64>) -> Self {
528        Self {
529            structural: StructuralCells::All,
530            likelihood_weights: LikelihoodWeights::ByRow(weights),
531        }
532    }
533
534    /// Construct from an explicit structural layout and numerical measure.
535    pub const fn new(
536        structural: StructuralCells<'a>,
537        likelihood_weights: LikelihoodWeights<'a>,
538    ) -> Self {
539        Self {
540            structural,
541            likelihood_weights,
542        }
543    }
544
545    /// Validate every shape, weight, and sparse-set index against `(N, M)`.
546    pub fn validate(&self, n_rows: usize, n_outputs: usize) -> Result<(), IndexedResponseError> {
547        match self.structural {
548            StructuralCells::All => {}
549            StructuralCells::Dense(mask) => {
550                if mask.dim() != (n_rows, n_outputs) {
551                    return Err(IndexedResponseError::new(format!(
552                        "structural cell mask shape {:?} does not match ({n_rows}, {n_outputs})",
553                        mask.dim()
554                    )));
555                }
556            }
557            StructuralCells::Only(cells) | StructuralCells::AllExcept(cells) => {
558                cells.validate_shape(n_rows, n_outputs)?;
559            }
560        }
561        match self.likelihood_weights {
562            LikelihoodWeights::Uniform => {}
563            LikelihoodWeights::ByRow(weights) => {
564                if weights.len() != n_rows {
565                    return Err(IndexedResponseError::new(format!(
566                        "row likelihood weights length {} does not match N={n_rows}",
567                        weights.len()
568                    )));
569                }
570                for (row, &weight) in weights.iter().enumerate() {
571                    validate_weight(weight, format!("row likelihood weight[{row}]"))?;
572                }
573            }
574            LikelihoodWeights::ByCell(weights) => {
575                if weights.dim() != (n_rows, n_outputs) {
576                    return Err(IndexedResponseError::new(format!(
577                        "cell likelihood weights shape {:?} does not match ({n_rows}, {n_outputs})",
578                        weights.dim()
579                    )));
580                }
581                for ((row, output), &weight) in weights.indexed_iter() {
582                    validate_weight(weight, format!("cell likelihood weight[{row},{output}]"))?;
583                }
584            }
585        }
586        Ok(())
587    }
588
589    /// Whether `(row, output)` is structurally present.
590    pub fn is_active(&self, row: usize, output: usize) -> bool {
591        match self.structural {
592            StructuralCells::All => true,
593            StructuralCells::Dense(mask) => mask[[row, output]],
594            StructuralCells::Only(cells) => cells.contains(row, output),
595            StructuralCells::AllExcept(cells) => !cells.contains(row, output),
596        }
597    }
598
599    /// Numerical weight of an active cell, or `None` when the cell is
600    /// structurally absent. Returning `Some(0.0)` preserves the distinction
601    /// between a present zero-weight observation and an absent cell.
602    pub fn active_weight(&self, row: usize, output: usize) -> Option<f64> {
603        if !self.is_active(row, output) {
604            return None;
605        }
606        Some(match self.likelihood_weights {
607            LikelihoodWeights::Uniform => 1.0,
608            LikelihoodWeights::ByRow(weights) => weights[row],
609            LikelihoodWeights::ByCell(weights) => weights[[row, output]],
610        })
611    }
612
613    /// Copy this borrowed view into a lifetime-free measure while preserving
614    /// sparse structural representations.
615    pub fn to_owned(
616        &self,
617        n_rows: usize,
618        n_outputs: usize,
619    ) -> Result<OwnedSeparableCellMeasure, IndexedResponseError> {
620        self.validate(n_rows, n_outputs)?;
621        let structural = match self.structural {
622            StructuralCells::All => OwnedStructuralCells::All,
623            StructuralCells::Dense(mask) => OwnedStructuralCells::Dense(mask.to_owned()),
624            StructuralCells::Only(cells) => OwnedStructuralCells::Only(cells.clone()),
625            StructuralCells::AllExcept(cells) => OwnedStructuralCells::AllExcept(cells.clone()),
626        };
627        let likelihood_weights = match self.likelihood_weights {
628            LikelihoodWeights::Uniform => OwnedLikelihoodWeights::Uniform,
629            LikelihoodWeights::ByRow(weights) => OwnedLikelihoodWeights::ByRow(weights.to_owned()),
630            LikelihoodWeights::ByCell(weights) => {
631                OwnedLikelihoodWeights::ByCell(weights.to_owned())
632            }
633        };
634        OwnedSeparableCellMeasure::new(
635            n_rows,
636            n_outputs,
637            structural,
638            likelihood_weights,
639        )
640    }
641}
642
643fn validate_weight(weight: f64, context: String) -> Result<(), IndexedResponseError> {
644    if !(weight.is_finite() && weight >= 0.0) {
645        return Err(IndexedResponseError::new(format!(
646            "{context} must be finite and non-negative (got {weight})"
647        )));
648    }
649    Ok(())
650}
651
652#[cfg(test)]
653mod tests {
654    use super::*;
655
656    #[test]
657    fn sparse_inclusion_and_exclusion_preserve_structural_geometry() {
658        let cells = IndexedCellSet::from_cells(3, 4, vec![(2, 3), (0, 1), (2, 0)])
659            .expect("valid sparse cell set");
660        assert_eq!(cells.row_outputs(0), Some(&[1][..]));
661        assert_eq!(cells.row_outputs(1), Some(&[][..]));
662        assert_eq!(cells.row_outputs(2), Some(&[0, 3][..]));
663
664        let only =
665            SeparableCellMeasure::new(StructuralCells::Only(&cells), LikelihoodWeights::Uniform);
666        only.validate(3, 4).expect("matching inclusion geometry");
667        assert_eq!(only.active_weight(2, 3), Some(1.0));
668        assert_eq!(only.active_weight(2, 2), None);
669
670        let except = SeparableCellMeasure::new(
671            StructuralCells::AllExcept(&cells),
672            LikelihoodWeights::Uniform,
673        );
674        except.validate(3, 4).expect("matching exclusion geometry");
675        assert_eq!(except.active_weight(2, 3), None);
676        assert_eq!(except.active_weight(2, 2), Some(1.0));
677    }
678
679    #[test]
680    fn structural_absence_is_distinct_from_zero_likelihood_weight() {
681        let excluded = IndexedCellSet::from_cells(1, 2, vec![(0, 0)]).expect("valid exclusion set");
682        let weights = ndarray::array![0.0];
683        let measure = SeparableCellMeasure::new(
684            StructuralCells::AllExcept(&excluded),
685            LikelihoodWeights::ByRow(weights.view()),
686        );
687        measure.validate(1, 2).expect("valid measure");
688        assert_eq!(measure.active_weight(0, 0), None);
689        assert_eq!(measure.active_weight(0, 1), Some(0.0));
690    }
691
692    #[test]
693    fn owned_measure_round_trip_preserves_sparse_geometry_and_cell_weights() {
694        let active = IndexedCellSet::from_cells(2, 3, vec![(0, 2), (1, 0)])
695            .expect("valid sparse activity set");
696        let weights = ndarray::array![[7.0, 8.0, 0.0], [2.5, 9.0, 10.0]];
697        let borrowed = SeparableCellMeasure::new(
698            StructuralCells::Only(&active),
699            LikelihoodWeights::ByCell(weights.view()),
700        );
701        let owned = borrowed.to_owned(2, 3).expect("owned response measure");
702
703        assert_eq!((owned.n_rows(), owned.n_outputs()), (2, 3));
704        assert_eq!(owned.active_weight(0, 2), Some(0.0));
705        assert_eq!(owned.active_weight(1, 0), Some(2.5));
706        assert_eq!(owned.active_weight(0, 0), None);
707        assert_eq!(owned.active_weight(1, 2), None);
708
709        let wrong_shape = borrowed
710            .to_owned(3, 3)
711            .expect_err("sparse geometry cannot be relabeled with another shape");
712        assert!(wrong_shape.reason().contains("does not match response geometry"));
713    }
714
715    #[test]
716    fn duplicate_cells_and_malformed_measures_are_rejected() {
717        let duplicate = IndexedCellSet::from_cells(2, 2, vec![(0, 1), (0, 1)])
718            .expect_err("duplicate structural declarations must fail");
719        assert!(duplicate.reason().contains("more than once"));
720
721        let bad_weights = ndarray::array![[1.0, -1.0]];
722        let measure = SeparableCellMeasure::new(
723            StructuralCells::All,
724            LikelihoodWeights::ByCell(bad_weights.view()),
725        );
726        let error = measure
727            .validate(1, 2)
728            .expect_err("negative likelihood weight must fail");
729        assert!(error.reason().contains("non-negative"));
730    }
731
732    #[test]
733    fn constant_values_with_sparse_overrides_preserve_row_major_identity() {
734        let values = OwnedCellValues::constant_with_overrides(
735            3,
736            4,
737            0.0,
738            vec![(2, 3, 9.0), (0, 1, 5.0), (2, 0, 7.0)],
739        )
740        .expect("valid sparse value field");
741
742        // Coordinates and values are canonicalized together.
743        assert_eq!(values.value(0, 1), Some(5.0));
744        assert_eq!(values.value(2, 0), Some(7.0));
745        assert_eq!(values.value(2, 3), Some(9.0));
746        assert_eq!(values.value(1, 2), Some(0.0));
747        assert_eq!(values.value(3, 0), None);
748    }
749
750    #[test]
751    fn row_and_output_broadcasts_preserve_declared_grid_geometry() {
752        let by_row = OwnedCellValues::by_row(ndarray::array![0.25, 1.5], 3);
753        assert_eq!((by_row.n_rows(), by_row.n_outputs()), (2, 3));
754        assert_eq!(by_row.value(0, 0), Some(0.25));
755        assert_eq!(by_row.value(0, 2), Some(0.25));
756        assert_eq!(by_row.value(1, 1), Some(1.5));
757        assert_eq!(by_row.value(2, 0), None);
758
759        let by_output = OwnedCellValues::by_output(2, ndarray::array![3.0, 5.0, 7.0]);
760        assert_eq!((by_output.n_rows(), by_output.n_outputs()), (2, 3));
761        assert_eq!(by_output.value(0, 1), Some(5.0));
762        assert_eq!(by_output.value(1, 1), Some(5.0));
763        assert_eq!(by_output.value(0, 3), None);
764    }
765
766    #[test]
767    fn sparse_value_overrides_reject_duplicate_coordinates() {
768        let error = OwnedCellValues::constant_with_overrides(
769            2,
770            3,
771            0.0,
772            vec![(1, 2, 4.0), (1, 2, 9.0)],
773        )
774        .expect_err("duplicate override coordinates must fail");
775        assert!(error.reason().contains("overridden more than once"));
776    }
777
778    #[test]
779    fn sparse_activity_visitor_keeps_zero_mass_cells_and_row_major_order() {
780        let active = IndexedCellSet::from_cells(3, 4, vec![(2, 3), (0, 1), (2, 0)])
781            .expect("valid active cells");
782        let weights = ndarray::array![
783            [1.0, 0.0, 3.0, 4.0],
784            [5.0, 6.0, 7.0, 8.0],
785            [9.0, 10.0, 11.0, 12.0]
786        ];
787        let measure = OwnedSeparableCellMeasure::new(
788            3,
789            4,
790            OwnedStructuralCells::Only(active),
791            OwnedLikelihoodWeights::ByCell(weights),
792        )
793        .expect("valid sparse measure");
794        let mut visited = Vec::new();
795        measure
796            .try_for_each_active::<std::convert::Infallible>(|row, output, weight| {
797                visited.push((row, output, weight));
798                Ok(())
799            })
800            .expect("infallible visit");
801
802        assert_eq!(visited, vec![(0, 1, 0.0), (2, 0, 9.0), (2, 3, 12.0)]);
803    }
804}