1#[derive(Debug, Clone)]
10pub struct GeoForm {
11 pub data: Vec<f64>, pub alpha: f64, pub damping: f64, }
15
16impl GeoForm {
17 pub fn from_leonov(len: usize, seed: u64) -> Self {
19 use rand::{Rng, SeedableRng};
20 use rand_pcg::Pcg64;
21
22 let mut rng = Pcg64::seed_from_u64(seed);
23 let mut data = vec![1.0, 1.0];
24
25 while data.len() < len {
26 let n = data.len();
27 let a = rng.random_range(-2.0..2.0);
28 let b = rng.random_range(-2.0..2.0);
29 let c = rng.random_range(-2.0..2.0);
30 let d = rng.random_range(1..=n.min(10));
31 let val = (a * data[n - 1] + b * data[n - 2] + c * data[n - d] / (1.0 + n as f64)).abs();
32 data.push(val);
33 }
34
35 Self {
36 data,
37 alpha: 0.1,
38 damping: 0.01,
39 }
40 }
41
42 pub fn adapt(&mut self, wave: &[f64], dt: f64) {
45 for (f, &psi) in self.data.iter_mut().zip(wave.iter()) {
46 *f += self.alpha * psi * psi * dt; *f -= self.damping * (*f - 1.0); }
49 }
50}
51
52#[derive(Debug, Clone)]
54pub struct Wave {
55 pub psi_prev: Vec<f64>, pub psi: Vec<f64>, pub alpha: f64, }
59
60impl Wave {
61 pub fn new(len: usize, impulse_pos: usize) -> Self {
63 let mut psi = vec![0.0; len];
64 psi[impulse_pos] = 1.0; Self {
67 psi_prev: psi.clone(),
68 psi,
69 alpha: 0.5,
70 }
71 }
72
73 pub fn step(&mut self, form: &GeoForm) {
75 let len = self.psi.len();
76 let mut next = vec![0.0; len];
77
78 for i in 1..len - 1 {
79 let laplace = (self.psi[i + 1] - 2.0 * self.psi[i] + self.psi[i - 1]) / form.data[i].max(1e-6);
80 next[i] = 2.0 * self.psi[i] - self.psi_prev[i] + self.alpha * laplace;
81 }
82
83 self.psi_prev = self.psi.clone();
84 self.psi = next;
85 }
86
87 pub fn energy(&self) -> Vec<f64> {
89 self.psi.iter().map(|x| x * x).collect()
90 }
91}
92
93#[cfg(test)]
94mod tests {
95 use super::*;
96
97 #[test]
99 fn test_adaptive_wave_memory() {
100 let mut form = GeoForm::from_leonov(64, 42);
101 let mut wave = Wave::new(64, 32);
102
103 for _ in 0..50 {
104 wave.step(&form);
105 let energy = wave.energy();
106 form.adapt(&energy, 0.1);
107 }
108
109 assert!(form.data.iter().any(|&f| f > 1.01));
111 }
112}