entropy_conservation/
flow.rs1use std::collections::HashMap;
5
6
7pub type ModuleId = String;
9
10#[derive(Debug, Clone, serde::Serialize, serde::Deserialize)]
12pub struct EntropyFlow {
13 pub source: ModuleId,
14 pub sink: ModuleId,
15 pub rate: f64,
16}
17
18#[derive(Debug, Clone, serde::Serialize, serde::Deserialize)]
20pub struct FlowNetwork {
21 pub modules: Vec<ModuleId>,
23 pub flows: Vec<EntropyFlow>,
25}
26
27impl FlowNetwork {
28 pub fn new(flows: Vec<EntropyFlow>) -> Self {
30 let mut module_set = std::collections::HashSet::new();
31 for f in &flows {
32 module_set.insert(f.source.clone());
33 module_set.insert(f.sink.clone());
34 }
35 let mut modules: Vec<_> = module_set.into_iter().collect();
36 modules.sort();
37 Self { modules, flows }
38 }
39
40 pub fn net_flow(&self) -> HashMap<ModuleId, f64> {
44 let mut net: HashMap<ModuleId, f64> = HashMap::new();
45 for m in &self.modules {
46 net.insert(m.clone(), 0.0);
47 }
48 for f in &self.flows {
49 *net.get_mut(&f.source).unwrap() -= f.rate;
50 *net.get_mut(&f.sink).unwrap() += f.rate;
51 }
52 net
53 }
54
55 pub fn check_conservation(&self) -> (bool, f64) {
59 let net = self.net_flow();
60 let max_violation = net.values().map(|v| v.abs()).fold(0.0_f64, f64::max);
61 (max_violation < 1e-10, max_violation)
62 }
63
64 pub fn total_flow(&self) -> f64 {
66 self.flows.iter().map(|f| f.rate).sum()
67 }
68
69 pub fn adjacency_matrix(&self) -> (usize, ndarray::Array2<f64>) {
71 let n = self.modules.len();
72 let mut mat = ndarray::Array2::<f64>::zeros((n, n));
73 let idx: HashMap<&ModuleId, usize> = self.modules.iter().enumerate().map(|(i, m)| (m, i)).collect();
74
75 for f in &self.flows {
76 if let (Some(&i), Some(&j)) = (idx.get(&f.source), idx.get(&f.sink)) {
77 mat[[i, j]] += f.rate;
78 }
79 }
80 (n, mat)
81 }
82
83 pub fn laplacian(&self) -> ndarray::Array2<f64> {
85 let (n, a) = self.adjacency_matrix();
86 let mut d = ndarray::Array2::<f64>::zeros((n, n));
87 for i in 0..n {
88 let mut row_sum = 0.0;
89 for j in 0..n {
90 row_sum += a[[i, j]];
91 }
92 d[[i, i]] = row_sum;
93 }
94 d - a
95 }
96}
97
98#[cfg(test)]
99mod tests {
100 use super::*;
101
102 #[test]
103 fn kirchhoff_conserved_cycle() {
104 let net = FlowNetwork::new(vec![
106 EntropyFlow { source: "A".into(), sink: "B".into(), rate: 1.0 },
107 EntropyFlow { source: "B".into(), sink: "C".into(), rate: 1.0 },
108 EntropyFlow { source: "C".into(), sink: "A".into(), rate: 1.0 },
109 ]);
110 let (ok, violation) = net.check_conservation();
111 assert!(ok, "cycle should be conserved, max violation = {violation}");
112 }
113
114 #[test]
115 fn kirchhoff_source_sink_violates() {
116 let net = FlowNetwork::new(vec![
118 EntropyFlow { source: "A".into(), sink: "B".into(), rate: 1.0 },
119 ]);
120 let (ok, _) = net.check_conservation();
121 assert!(!ok, "unbalanced flow should violate conservation");
122 }
123
124 #[test]
125 fn kirchhoff_multi_flow_conserved() {
126 let net = FlowNetwork::new(vec![
128 EntropyFlow { source: "A".into(), sink: "B".into(), rate: 2.0 },
129 EntropyFlow { source: "A".into(), sink: "C".into(), rate: 1.0 },
130 EntropyFlow { source: "B".into(), sink: "C".into(), rate: 2.0 },
131 EntropyFlow { source: "C".into(), sink: "A".into(), rate: 3.0 },
132 ]);
133 let (ok, _) = net.check_conservation();
134 assert!(ok, "balanced multi-flow should be conserved");
135
136 let net_flow = net.net_flow();
137 for (m, v) in &net_flow {
138 assert!(v.abs() < 1e-10, "node {m}: net flow = {v}");
139 }
140 }
141
142 #[test]
143 fn laplacian_rows_sum_zero() {
144 let net = FlowNetwork::new(vec![
145 EntropyFlow { source: "A".into(), sink: "B".into(), rate: 1.0 },
146 EntropyFlow { source: "B".into(), sink: "C".into(), rate: 1.0 },
147 EntropyFlow { source: "C".into(), sink: "A".into(), rate: 1.0 },
148 ]);
149 let l = net.laplacian();
150 let n = l.nrows();
151 for i in 0..n {
152 let row_sum: f64 = (0..n).map(|j| l[[i, j]]).sum();
153 assert!((row_sum).abs() < 1e-12, "Laplacian row {i} sums to {row_sum}, expected 0");
154 }
155 }
156
157 #[test]
158 fn adjacency_matrix_symmetry() {
159 let net = FlowNetwork::new(vec![
160 EntropyFlow { source: "A".into(), sink: "B".into(), rate: 1.0 },
161 EntropyFlow { source: "B".into(), sink: "A".into(), rate: 1.0 },
162 ]);
163 let (_, a) = net.adjacency_matrix();
164 assert!((a[[0, 1]] - a[[1, 0]]).abs() < 1e-12);
165 }
166}