tawasuyu/brahman

1//! Difusión y entropía de los campos de la grilla.
2//!
3//! Ecuación de fluidos discreta: cada celda intercambia una fracción de
4//! su valor con sus 4 vecinas, y luego pierde una fracción al ambiente
5//! (entropía). Difunden los 3 campos dinámicos — materia, psique,
6//! poder. `oro` (materia sólida) y `degradacion` (cicatriz permanente)
7//! no difunden.
8
9use dominium_core::{Grid, SimParams};
10
11/// Difunde una sola capa: `new[c] = c + rate·(media_vecinos − c)`, y luego
12/// aplica la entropía. Usa un buffer de lectura separado (la difusión
13/// debe leer el estado viejo).
14fn diffuse_layer(layer: &mut [f32], width: usize, height: usize, rate: f32, entropy: f32) {
15    let old = layer.to_vec();
16    for y in 0..height {
17        for x in 0..width {
18            let c = y * width + x;
19            let mut sum = 0.0f32;
20            let mut count = 0.0f32;
21            // 4-vecindad (von Neumann), bordes sin wrap.
22            if x > 0 {
23                sum += old[c - 1];
24                count += 1.0;
25            }
26            if x + 1 < width {
27                sum += old[c + 1];
28                count += 1.0;
29            }
30            if y > 0 {
31                sum += old[c - width];
32                count += 1.0;
33            }
34            if y + 1 < height {
35                sum += old[c + width];
36                count += 1.0;
37            }
38            let neighbor_avg = if count > 0.0 { sum / count } else { old[c] };
39            let diffused = old[c] + rate * (neighbor_avg - old[c]);
40            layer[c] = diffused * (1.0 - entropy);
41        }
42    }
43}
44
45/// Aplica un paso de difusión + entropía a los 3 campos dinámicos.
46pub fn diffuse(grid: &mut Grid, p: &SimParams) {
47    let (w, h) = (grid.width, grid.height);
48    let (rate, ent) = (p.diffusion_rate, p.entropy_rate);
49    diffuse_layer(&mut grid.materia, w, h, rate, ent);
50    diffuse_layer(&mut grid.psique, w, h, rate, ent);
51    diffuse_layer(&mut grid.poder, w, h, rate, ent);
52}
53
54#[cfg(test)]
55mod tests {
56    use super::*;
57
58    #[test]
59    fn diffusion_spreads_a_spike_to_neighbors() {
60        let mut g = Grid::new(5, 5);
61        let center = g.idx(2, 2);
62        g.materia[center] = 100.0;
63        let p = SimParams::default();
64        diffuse(&mut g, &p);
65        // El pico bajó; las vecinas subieron desde 0.
66        assert!(g.materia[center] < 100.0);
67        assert!(g.materia[g.idx(1, 2)] > 0.0);
68        assert!(g.materia[g.idx(3, 2)] > 0.0);
69    }
70
71    #[test]
72    fn entropy_decays_a_uniform_field() {
73        let mut g = Grid::new(4, 4);
74        for v in g.psique.iter_mut() {
75            *v = 10.0;
76        }
77        let p = SimParams::default();
78        diffuse(&mut g, &p);
79        // Campo uniforme: la difusión no cambia nada, pero la entropía sí.
80        for &v in &g.psique {
81            assert!(v < 10.0 && v > 9.0);
82        }
83    }
84
85    #[test]
86    fn diffusion_conserves_mass_minus_entropy() {
87        let mut g = Grid::new(6, 6);
88        let c = g.idx(3, 3);
89        g.materia[c] = 60.0;
90        let total_before: f32 = g.materia.iter().sum();
91        let mut p = SimParams::default();
92        p.entropy_rate = 0.0; // sin pérdida → masa conservada
93        diffuse(&mut g, &p);
94        let total_after: f32 = g.materia.iter().sum();
95        assert!((total_before - total_after).abs() < 1e-2);
96    }
97}