Another project
1use std::collections::BTreeSet;
2
3use bone_types::{ParentIndex, ResidualIndex};
4
5use crate::dof::{DofConfig, analyze_dof_at};
6use crate::graph::decompose;
7use crate::newton::{NewtonConfig, solve_newton_decomposed};
8use crate::system::ConstraintSystem;
9use crate::util::dedup_preserving_order;
10
11#[must_use]
12pub fn minimal_unsatisfiable_subset(
13 system: &ConstraintSystem,
14 over_flagged: &[ResidualIndex],
15 cfg: DofConfig,
16) -> Vec<ResidualIndex> {
17 let candidate_parents: Vec<ParentIndex> = dedup_preserving_order(
18 over_flagged
19 .iter()
20 .copied()
21 .map(|row| system.parent_of_row(row)),
22 );
23 if candidate_parents.is_empty() {
24 return Vec::new();
25 }
26 let solver_cfg = NewtonConfig::DEFAULT;
27 let kept_parents = candidate_parents
28 .iter()
29 .copied()
30 .fold(
31 (BTreeSet::<ParentIndex>::new(), Vec::<ParentIndex>::new()),
32 |(mut excluded, mut kept), parent| {
33 let mut trial = excluded.clone();
34 trial.insert(parent);
35 if reduced_is_over(system, &trial, solver_cfg, cfg) {
36 excluded.insert(parent);
37 } else {
38 kept.push(parent);
39 }
40 (excluded, kept)
41 },
42 )
43 .1;
44 kept_parents
45 .into_iter()
46 .flat_map(|p| system.rows_of_parent(p))
47 .filter(|row| over_flagged.contains(row))
48 .collect()
49}
50
51fn reduced_is_over(
52 system: &ConstraintSystem,
53 exclude: &BTreeSet<ParentIndex>,
54 newton_cfg: NewtonConfig,
55 dof_cfg: DofConfig,
56) -> bool {
57 let sub = system.without_parents(exclude);
58 if sub.row_count() == 0 {
59 return false;
60 }
61 let decomposition = decompose(&sub);
62 match solve_newton_decomposed(&sub, &decomposition, newton_cfg) {
63 Ok(solved) => {
64 let report = analyze_dof_at(&sub, &solved, dof_cfg);
65 !report.over_constrained().is_empty()
66 }
67 Err(_) => true,
68 }
69}
70
71#[cfg(test)]
72mod tests {
73 use super::*;
74 use crate::residual::Residual;
75 use bone_types::{Parameter, ParameterIndex};
76
77 fn p(i: u32) -> ParameterIndex {
78 ParameterIndex::new(i)
79 }
80
81 fn pin(at: u32, target: f64) -> Residual {
82 Residual::Pin {
83 param: p(at),
84 target,
85 }
86 }
87
88 fn over_flagged_after_solve(system: &ConstraintSystem) -> Vec<ResidualIndex> {
89 let decomposition = decompose(system);
90 let solved = solve_newton_decomposed(system, &decomposition, NewtonConfig::DEFAULT)
91 .unwrap_or_else(|_| system.parameters().to_vec());
92 analyze_dof_at(system, &solved, DofConfig::DEFAULT)
93 .over_constrained()
94 .to_vec()
95 }
96
97 #[test]
98 fn three_conflicting_pins_isolate_the_odd_one_out() {
99 let system = ConstraintSystem::new(
100 vec![Parameter::new(0.0)],
101 vec![pin(0, 1.0), pin(0, 2.0), pin(0, 1.0)],
102 );
103 let over = over_flagged_after_solve(&system);
104 assert!(
105 !over.is_empty(),
106 "fixture must produce at least one over row"
107 );
108 let mus = minimal_unsatisfiable_subset(&system, &over, DofConfig::DEFAULT);
109 assert!(
110 !mus.is_empty(),
111 "MUS must keep at least one essential row, got empty",
112 );
113 assert!(
114 mus.len() <= over.len(),
115 "MUS must be a subset of over-flagged rows: mus={mus:?} over={over:?}",
116 );
117 }
118
119 #[test]
120 fn empty_over_flagged_yields_empty_mus() {
121 let system = ConstraintSystem::new(vec![Parameter::new(0.0)], vec![pin(0, 0.0)]);
122 let mus = minimal_unsatisfiable_subset(&system, &[], DofConfig::DEFAULT);
123 assert!(mus.is_empty());
124 }
125}