sc_neurocore_engine/
connectome.rs1use rand::{RngExt, SeedableRng};
12use rand_xoshiro::Xoshiro256PlusPlus;
13
14pub fn watts_strogatz(n: usize, k: usize, p_rewire: f64, seed: u64) -> Vec<u8> {
18 let mut rng = Xoshiro256PlusPlus::seed_from_u64(seed);
19 let mut adj = vec![0u8; n * n];
20
21 if k >= n {
22 for i in 0..n {
23 for j in 0..n {
24 if i != j {
25 adj[i * n + j] = 1;
26 }
27 }
28 }
29 return adj;
30 }
31
32 for i in 0..n {
34 for d in 1..=k / 2 {
35 let j = (i + d) % n;
36 adj[i * n + j] = 1;
37 adj[j * n + i] = 1;
38 }
39 }
40
41 for i in 0..n {
43 for d in 1..=k / 2 {
44 let j = (i + d) % n;
45 if rng.random::<f64>() < p_rewire {
46 adj[i * n + j] = 0;
47 let mut new_j = i;
48 while new_j == i || adj[i * n + new_j] == 1 {
49 new_j = (rng.random::<f64>() * n as f64) as usize % n;
50 }
51 adj[i * n + new_j] = 1;
52 }
53 }
54 }
55
56 adj
57}
58
59#[allow(clippy::needless_range_loop)]
60pub fn barabasi_albert(n: usize, seed: u64) -> Vec<u8> {
64 let mut rng = Xoshiro256PlusPlus::seed_from_u64(seed);
65 let mut adj = vec![0u8; n * n];
66 let mut degree = vec![0u64; n];
67
68 if n < 2 {
69 return adj;
70 }
71
72 adj[1] = 1;
74 adj[n] = 1;
75 degree[0] = 1;
76 degree[1] = 1;
77
78 for i in 2..n {
79 let total_degree: u64 = degree[..i].iter().sum();
80 if total_degree == 0 {
81 adj[i * n] = 1;
82 degree[i] += 1;
83 degree[0] += 1;
84 continue;
85 }
86
87 let r = rng.random::<f64>() * total_degree as f64;
88 let mut cumulative = 0.0;
89 let mut target = 0;
90 for j in 0..i {
91 cumulative += degree[j] as f64;
92 if cumulative > r {
93 target = j;
94 break;
95 }
96 }
97
98 adj[i * n + target] = 1;
99 degree[i] += 1;
100 degree[target] += 1;
101 }
102
103 adj
104}
105
106#[cfg(test)]
107mod tests {
108 use super::*;
109
110 #[test]
111 fn ws_no_self_loops() {
112 let adj = watts_strogatz(20, 4, 0.3, 42);
113 for i in 0..20 {
114 assert_eq!(adj[i * 20 + i], 0);
115 }
116 }
117
118 #[test]
119 fn ws_connected() {
120 let adj = watts_strogatz(10, 4, 0.0, 0);
121 for i in 0..10 {
122 let row_sum: u8 = (0..10).map(|j| adj[i * 10 + j]).sum();
123 assert!(row_sum > 0, "node {i} has no connections");
124 }
125 }
126
127 #[test]
128 fn barabasi_albert_grows_correctly() {
129 let adj = barabasi_albert(10, 99);
130 for i in 0..10 {
131 assert_eq!(adj[i * 10 + i], 0);
132 }
133 for i in 2..10 {
135 let row_sum: u8 = (0..10).map(|j| adj[i * 10 + j]).sum();
136 assert!(row_sum > 0, "node {i} has no outgoing edges");
137 }
138 }
139
140 #[test]
141 fn deterministic() {
142 let a = watts_strogatz(10, 4, 0.3, 123);
143 let b = watts_strogatz(10, 4, 0.3, 123);
144 assert_eq!(a, b);
145 }
146}