tawasuyu/llimphi
Sovereign retained-mode GPU UI framework in Rust (wgpu+vello+taffy+parley), Elm loop, 40+ widgets.
1//! `ecuacion` — **motor de Leyes autorables**: una Ley deja de ser un autómata
2//! cableado en Rust y pasa a ser un **sistema de ecuaciones de actualización sobre
3//! los campos de estado de la celda**, escrito por el usuario y **compilado** una vez
4//! a opcodes para correr en el loop caliente.
5//!
6//! El modelo: cada material lleva un vector chico de **campos** escalares
7//! ([`FieldDef`]: `cantidad`, `energía`, `temperatura`, `edad`…). Una Ley es una lista
8//! de asignaciones `Δcampoₖ/dt = fₖ(campos, vecinos, params)` — un [`Program`]. Cada
9//! tick, el motor integra `campo' = clamp(campo + dt·fₖ, min, max)` sobre toda la
10//! rejilla.
11//!
12//! **Determinista** (regla §1.5 de la sim): se evalúa estilo **Jacobi** — `fₖ` lee
13//! sólo el estado del tick anterior (buffer `cur`) y escribe en `next`, así el
14//! resultado no depende del orden de barrido. Mismas condiciones → misma evolución,
15//! sin azar y sin sesgo de recorrido (a diferencia de Gauss–Seidel).
16//!
17//! Esto **generaliza** la ley `Crecer` (un campo `altura` que sube) y abre reacción‑
18//! difusión (Gray–Scott), difusión de calor, fuego, erosión… con la misma máquina.
19//! NO reemplaza a [`WaterSim`](crate::WaterSim): ese es un autómata que *mueve celdas*
20//! (conserva masa), otra discretización, y sigue como fast‑path nativo aparte.
21//!
22//! El mismo AST [`Expr`] alimenta las **dos vistas de autoría**: la barra de fórmula
23//! (parser [`Expr::parse`] ↔ printer [`Expr::to_source`]) y, a futuro, el grafo de
24//! nodos. Ambas serializan a este `Expr`, así se editan indistintamente.
25
26use serde::{Deserialize, Serialize};
27
28/// Índice de un campo de estado dentro del material (posición en su `Vec<FieldDef>`).
29pub type FieldId = u16;
30/// Índice de un parámetro del material (los sliders que ya expone la UI).
31pub type ParamId = u16;
32
33/// Un **campo de estado** de la celda: un escalar con nombre, valor inicial y rango.
34/// El material declara los suyos; las ecuaciones los leen y escriben por nombre.
35#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
36pub struct FieldDef {
37 pub name: String,
38 pub init: f32,
39 pub min: f32,
40 pub max: f32,
41}
42
43impl FieldDef {
44 pub fn new(name: impl Into<String>, init: f32, min: f32, max: f32) -> Self {
45 Self { name: name.into(), init, min, max }
46 }
47}
48
49/// Dirección de un vecino de 6‑conectividad (caras del cubo).
50#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
51pub enum Dir {
52 Abajo,
53 Arriba,
54 Este,
55 Oeste,
56 Norte,
57 Sur,
58}
59
60impl Dir {
61 /// Desplazamiento `(dx,dy,dz)` de la dirección.
62 #[inline]
63 pub fn delta(self) -> (i32, i32, i32) {
64 match self {
65 Dir::Abajo => (0, -1, 0),
66 Dir::Arriba => (0, 1, 0),
67 Dir::Este => (1, 0, 0),
68 Dir::Oeste => (-1, 0, 0),
69 Dir::Norte => (0, 0, 1),
70 Dir::Sur => (0, 0, -1),
71 }
72 }
73 fn nombre(self) -> &'static str {
74 match self {
75 Dir::Abajo => "abajo",
76 Dir::Arriba => "arriba",
77 Dir::Este => "este",
78 Dir::Oeste => "oeste",
79 Dir::Norte => "norte",
80 Dir::Sur => "sur",
81 }
82 }
83 fn de_nombre(s: &str) -> Option<Dir> {
84 Some(match s {
85 "abajo" => Dir::Abajo,
86 "arriba" => Dir::Arriba,
87 "este" => Dir::Este,
88 "oeste" => Dir::Oeste,
89 "norte" => Dir::Norte,
90 "sur" => Dir::Sur,
91 _ => return None,
92 })
93 }
94}
95
96/// Reducción sobre los 6 vecinos de una celda.
97#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
98pub enum Reduce {
99 Promedio,
100 Suma,
101 Min,
102 Max,
103}
104
105impl Reduce {
106 fn nombre(self) -> &'static str {
107 match self {
108 Reduce::Promedio => "avg",
109 Reduce::Suma => "sum6",
110 Reduce::Min => "min6",
111 Reduce::Max => "max6",
112 }
113 }
114}
115
116/// Operador unario.
117#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
118pub enum UnOp {
119 Neg,
120 Abs,
121 Exp,
122 Sqrt,
123}
124
125/// Operador binario.
126#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
127pub enum BinOp {
128 Add,
129 Sub,
130 Mul,
131 Div,
132 Min,
133 Max,
134 /// `a > b → 1.0`, si no `0.0` (umbral).
135 Gt,
136 /// `a < b → 1.0`, si no `0.0`.
137 Lt,
138}
139
140/// **Árbol de la ecuación**: la expresión que da el lado derecho de un `Δcampo/dt`.
141#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
142pub enum Expr {
143 Const(f32),
144 /// Valor del campo en esta celda.
145 Field(FieldId),
146 /// Parámetro del material (slider).
147 Param(ParamId),
148 /// Paso de tiempo `dt`.
149 Dt,
150 /// Laplaciano discreto 6‑conexo del campo: `Σvecinos − 6·centro`.
151 Laplacian(FieldId),
152 /// Reducción sobre los 6 vecinos.
153 Vecinos(Reduce, FieldId),
154 /// Valor del campo en un vecino concreto.
155 Dir(Dir, FieldId),
156 Un(UnOp, Box<Expr>),
157 Bin(BinOp, Box<Expr>, Box<Expr>),
158 /// `clamp(x, lo, hi)`.
159 Clamp(Box<Expr>, Box<Expr>, Box<Expr>),
160}
161
162/// Una **asignación**: `Δcampo/dt = expr`. El motor integra `campo += dt·expr`.
163#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
164pub struct Assign {
165 pub campo: FieldId,
166 pub expr: Expr,
167}
168
169// =============================================================================
170// Compilación a opcodes (stack machine)
171// =============================================================================
172
173/// Opcode de la máquina de pila que ejecuta una ecuación por celda (compilado una
174/// vez al editar la Ley; corre en el loop caliente sin recorrer el árbol).
175#[derive(Debug, Clone, Copy, PartialEq)]
176enum Op {
177 Const(f32),
178 Field(FieldId),
179 Param(ParamId),
180 Dt,
181 Laplacian(FieldId),
182 Vecinos(Reduce, FieldId),
183 Dir(Dir, FieldId),
184 Un(UnOp),
185 Bin(BinOp),
186 Clamp,
187}
188
189/// Ecuación compilada de un campo: `(campo destino, opcodes en RPN)`.
190#[derive(Debug, Clone)]
191struct Compiled {
192 campo: FieldId,
193 ops: Vec<Op>,
194}
195
196fn compilar(expr: &Expr, out: &mut Vec<Op>) {
197 match expr {
198 Expr::Const(c) => out.push(Op::Const(*c)),
199 Expr::Field(f) => out.push(Op::Field(*f)),
200 Expr::Param(p) => out.push(Op::Param(*p)),
201 Expr::Dt => out.push(Op::Dt),
202 Expr::Laplacian(f) => out.push(Op::Laplacian(*f)),
203 Expr::Vecinos(r, f) => out.push(Op::Vecinos(*r, *f)),
204 Expr::Dir(d, f) => out.push(Op::Dir(*d, *f)),
205 Expr::Un(op, a) => {
206 compilar(a, out);
207 out.push(Op::Un(*op));
208 }
209 Expr::Bin(op, a, b) => {
210 compilar(a, out);
211 compilar(b, out);
212 out.push(Op::Bin(*op));
213 }
214 Expr::Clamp(x, lo, hi) => {
215 compilar(x, out);
216 compilar(lo, out);
217 compilar(hi, out);
218 out.push(Op::Clamp);
219 }
220 }
221}
222
223/// Un **programa de Ley compilado**: las ecuaciones de todos los campos, listas para
224/// correr. Se obtiene con [`Program::compile`] a partir de las [`Assign`] autoradas.
225#[derive(Debug, Clone)]
226pub struct Program {
227 ecuaciones: Vec<Compiled>,
228}
229
230impl Program {
231 /// Compila la lista de asignaciones (una Ley) a opcodes.
232 pub fn compile(asignaciones: &[Assign]) -> Program {
233 let ecuaciones = asignaciones
234 .iter()
235 .map(|a| {
236 let mut ops = Vec::new();
237 compilar(&a.expr, &mut ops);
238 Compiled { campo: a.campo, ops }
239 })
240 .collect();
241 Program { ecuaciones }
242 }
243
244 /// Cantidad de ecuaciones (campos que la Ley modifica).
245 pub fn len(&self) -> usize {
246 self.ecuaciones.len()
247 }
248 pub fn is_empty(&self) -> bool {
249 self.ecuaciones.is_empty()
250 }
251}
252
253// =============================================================================
254// Motor de campo
255// =============================================================================
256
257/// **Motor de campos**: mantiene el estado escalar por celda de todos los campos de
258/// un material y avanza la simulación aplicando un [`Program`] estilo Jacobi.
259///
260/// Layout de memoria **campo‑mayor**: los `ncells` valores de un campo son
261/// contiguos (`buf[f*ncells + i]`), así muestrear vecinos de un campo es cache‑friendly.
262pub struct FieldEngine {
263 dim: [u32; 3],
264 ncells: usize,
265 defs: Vec<FieldDef>,
266 cur: Vec<f32>,
267 next: Vec<f32>,
268}
269
270impl FieldEngine {
271 /// Crea el motor con los campos dados, cada celda inicializada a `FieldDef::init`.
272 pub fn new(dim: [u32; 3], defs: Vec<FieldDef>) -> Self {
273 let ncells = (dim[0] as usize) * (dim[1] as usize) * (dim[2] as usize);
274 let nf = defs.len();
275 let mut cur = vec![0.0f32; nf * ncells];
276 for (f, d) in defs.iter().enumerate() {
277 for v in cur[f * ncells..(f + 1) * ncells].iter_mut() {
278 *v = d.init;
279 }
280 }
281 let next = cur.clone();
282 Self { dim, ncells, defs, cur, next }
283 }
284
285 pub fn dim(&self) -> [u32; 3] {
286 self.dim
287 }
288 pub fn fields(&self) -> &[FieldDef] {
289 &self.defs
290 }
291
292 #[inline]
293 fn idx(&self, x: u32, y: u32, z: u32) -> usize {
294 (x as usize) + (y as usize) * (self.dim[0] as usize)
295 + (z as usize) * (self.dim[0] as usize) * (self.dim[1] as usize)
296 }
297
298 /// Lee el campo `f` en `(x,y,z)`.
299 #[inline]
300 pub fn get(&self, f: FieldId, x: u32, y: u32, z: u32) -> f32 {
301 self.cur[(f as usize) * self.ncells + self.idx(x, y, z)]
302 }
303
304 /// Fija el campo `f` en `(x,y,z)` (clampeado al rango del campo). Útil para
305 /// sembrar condiciones iniciales (una gota, un foco de calor…).
306 pub fn set(&mut self, f: FieldId, x: u32, y: u32, z: u32, v: f32) {
307 let d = &self.defs[f as usize];
308 let vi = v.clamp(d.min, d.max);
309 let i = (f as usize) * self.ncells + self.idx(x, y, z);
310 self.cur[i] = vi;
311 }
312
313 /// Muestra el campo `f` en `(x,y,z)` con **borde Neumann** (fuera de la rejilla =
314 /// el valor de la celda misma → sin flujo artificial en el borde).
315 #[inline]
316 fn sample(&self, f: usize, x: i32, y: i32, z: i32) -> f32 {
317 let xc = x.clamp(0, self.dim[0] as i32 - 1) as u32;
318 let yc = y.clamp(0, self.dim[1] as i32 - 1) as u32;
319 let zc = z.clamp(0, self.dim[2] as i32 - 1) as u32;
320 self.cur[f * self.ncells + self.idx(xc, yc, zc)]
321 }
322
323 /// Evalúa una ecuación compilada en la celda `(x,y,z)` sobre el buffer `cur`.
324 fn eval(&self, ops: &[Op], params: &[f32], dt: f32, x: i32, y: i32, z: i32, stack: &mut Vec<f32>) -> f32 {
325 stack.clear();
326 for op in ops {
327 match *op {
328 Op::Const(c) => stack.push(c),
329 Op::Field(f) => stack.push(self.sample(f as usize, x, y, z)),
330 Op::Param(p) => stack.push(params.get(p as usize).copied().unwrap_or(0.0)),
331 Op::Dt => stack.push(dt),
332 Op::Laplacian(f) => {
333 let f = f as usize;
334 let c = self.sample(f, x, y, z);
335 let s = self.sample(f, x + 1, y, z)
336 + self.sample(f, x - 1, y, z)
337 + self.sample(f, x, y + 1, z)
338 + self.sample(f, x, y - 1, z)
339 + self.sample(f, x, y, z + 1)
340 + self.sample(f, x, y, z - 1);
341 stack.push(s - 6.0 * c);
342 }
343 Op::Vecinos(r, f) => {
344 let f = f as usize;
345 let n = [
346 self.sample(f, x + 1, y, z),
347 self.sample(f, x - 1, y, z),
348 self.sample(f, x, y + 1, z),
349 self.sample(f, x, y - 1, z),
350 self.sample(f, x, y, z + 1),
351 self.sample(f, x, y, z - 1),
352 ];
353 let v = match r {
354 Reduce::Suma => n.iter().sum(),
355 Reduce::Promedio => n.iter().sum::<f32>() / 6.0,
356 Reduce::Min => n.iter().copied().fold(f32::INFINITY, f32::min),
357 Reduce::Max => n.iter().copied().fold(f32::NEG_INFINITY, f32::max),
358 };
359 stack.push(v);
360 }
361 Op::Dir(d, f) => {
362 let (dx, dy, dz) = d.delta();
363 stack.push(self.sample(f as usize, x + dx, y + dy, z + dz));
364 }
365 Op::Un(u) => {
366 let a = stack.pop().unwrap_or(0.0);
367 stack.push(match u {
368 UnOp::Neg => -a,
369 UnOp::Abs => a.abs(),
370 UnOp::Exp => a.exp(),
371 UnOp::Sqrt => a.max(0.0).sqrt(),
372 });
373 }
374 Op::Bin(b) => {
375 let rhs = stack.pop().unwrap_or(0.0);
376 let lhs = stack.pop().unwrap_or(0.0);
377 stack.push(match b {
378 BinOp::Add => lhs + rhs,
379 BinOp::Sub => lhs - rhs,
380 BinOp::Mul => lhs * rhs,
381 BinOp::Div => {
382 if rhs == 0.0 {
383 0.0
384 } else {
385 lhs / rhs
386 }
387 }
388 BinOp::Min => lhs.min(rhs),
389 BinOp::Max => lhs.max(rhs),
390 BinOp::Gt => (lhs > rhs) as i32 as f32,
391 BinOp::Lt => (lhs < rhs) as i32 as f32,
392 });
393 }
394 Op::Clamp => {
395 let hi = stack.pop().unwrap_or(0.0);
396 let lo = stack.pop().unwrap_or(0.0);
397 let xv = stack.pop().unwrap_or(0.0);
398 stack.push(xv.clamp(lo.min(hi), lo.max(hi)));
399 }
400 }
401 }
402 stack.pop().unwrap_or(0.0)
403 }
404
405 /// Avanza **un tick** integrando `campo += dt·f` para cada ecuación, estilo Jacobi
406 /// (todas leen el estado anterior). Determinista: independiente del orden de barrido.
407 pub fn step(&mut self, program: &Program, params: &[f32], dt: f32) {
408 // `next` parte del estado actual: los campos sin ecuación quedan intactos.
409 self.next.copy_from_slice(&self.cur);
410 let (dx, dy, dz) = (self.dim[0] as i32, self.dim[1] as i32, self.dim[2] as i32);
411 let mut stack: Vec<f32> = Vec::with_capacity(16);
412 for ec in &program.ecuaciones {
413 let f = ec.campo as usize;
414 let d = &self.defs[f];
415 let base = f * self.ncells;
416 for z in 0..dz {
417 for y in 0..dy {
418 for x in 0..dx {
419 let i = self.idx(x as u32, y as u32, z as u32);
420 let cur_v = self.cur[base + i];
421 let rate = self.eval(&ec.ops, params, dt, x, y, z, &mut stack);
422 self.next[base + i] = (cur_v + dt * rate).clamp(d.min, d.max);
423 }
424 }
425 }
426 }
427 std::mem::swap(&mut self.cur, &mut self.next);
428 }
429
430 /// Suma de todos los valores de un campo (para tests de conservación).
431 pub fn total(&self, f: FieldId) -> f32 {
432 let base = (f as usize) * self.ncells;
433 self.cur[base..base + self.ncells].iter().sum()
434 }
435}
436
437// =============================================================================
438// Parser texto → Expr y printer Expr → texto (base de las dos vistas)
439// =============================================================================
440
441/// Tabla de símbolos: resuelve nombres de campo y de parámetro a sus ids. La comparte
442/// el editor de la Ley (viene de los `FieldDef` del material y sus params).
443#[derive(Debug, Clone, Default)]
444pub struct Symbols {
445 pub campos: Vec<String>,
446 pub params: Vec<String>,
447}
448
449impl Symbols {
450 fn campo_id(&self, n: &str) -> Option<FieldId> {
451 self.campos.iter().position(|c| c == n).map(|i| i as FieldId)
452 }
453 fn param_id(&self, n: &str) -> Option<ParamId> {
454 self.params.iter().position(|p| p == n).map(|i| i as ParamId)
455 }
456 fn campo_nombre(&self, id: FieldId) -> String {
457 self.campos.get(id as usize).cloned().unwrap_or_else(|| format!("campo{id}"))
458 }
459 fn param_nombre(&self, id: ParamId) -> String {
460 self.params.get(id as usize).cloned().unwrap_or_else(|| format!("param{id}"))
461 }
462}
463
464#[derive(Debug, Clone, PartialEq)]
465enum Tok {
466 Num(f32),
467 Ident(String),
468 LParen,
469 RParen,
470 Comma,
471 Plus,
472 Minus,
473 Star,
474 Slash,
475 Lt,
476 Gt,
477}
478
479fn lex(s: &str) -> Result<Vec<Tok>, String> {
480 let mut toks = Vec::new();
481 let bytes: Vec<char> = s.chars().collect();
482 let mut i = 0;
483 while i < bytes.len() {
484 let c = bytes[i];
485 if c.is_whitespace() {
486 i += 1;
487 continue;
488 }
489 match c {
490 '(' => {
491 toks.push(Tok::LParen);
492 i += 1;
493 }
494 ')' => {
495 toks.push(Tok::RParen);
496 i += 1;
497 }
498 ',' => {
499 toks.push(Tok::Comma);
500 i += 1;
501 }
502 '+' => {
503 toks.push(Tok::Plus);
504 i += 1;
505 }
506 '-' => {
507 toks.push(Tok::Minus);
508 i += 1;
509 }
510 '*' => {
511 toks.push(Tok::Star);
512 i += 1;
513 }
514 '/' => {
515 toks.push(Tok::Slash);
516 i += 1;
517 }
518 '<' => {
519 toks.push(Tok::Lt);
520 i += 1;
521 }
522 '>' => {
523 toks.push(Tok::Gt);
524 i += 1;
525 }
526 _ if c.is_ascii_digit() || c == '.' => {
527 let start = i;
528 while i < bytes.len() && (bytes[i].is_ascii_digit() || bytes[i] == '.') {
529 i += 1;
530 }
531 let txt: String = bytes[start..i].iter().collect();
532 let n: f32 = txt.parse().map_err(|_| format!("número inválido: {txt}"))?;
533 toks.push(Tok::Num(n));
534 }
535 _ if c.is_alphabetic() || c == '_' => {
536 let start = i;
537 while i < bytes.len() && (bytes[i].is_alphanumeric() || bytes[i] == '_') {
538 i += 1;
539 }
540 toks.push(Tok::Ident(bytes[start..i].iter().collect()));
541 }
542 _ => return Err(format!("carácter inesperado: {c}")),
543 }
544 }
545 Ok(toks)
546}
547
548struct Parser<'a> {
549 toks: Vec<Tok>,
550 pos: usize,
551 sym: &'a Symbols,
552}
553
554impl<'a> Parser<'a> {
555 fn peek(&self) -> Option<&Tok> {
556 self.toks.get(self.pos)
557 }
558 fn next(&mut self) -> Option<Tok> {
559 let t = self.toks.get(self.pos).cloned();
560 self.pos += 1;
561 t
562 }
563 fn eat(&mut self, t: &Tok) -> Result<(), String> {
564 if self.peek() == Some(t) {
565 self.pos += 1;
566 Ok(())
567 } else {
568 Err(format!("se esperaba {t:?}, hay {:?}", self.peek()))
569 }
570 }
571
572 // suma/resta (precedencia baja)
573 fn parse_expr(&mut self) -> Result<Expr, String> {
574 let mut lhs = self.parse_term()?;
575 while let Some(t) = self.peek() {
576 let op = match t {
577 Tok::Plus => BinOp::Add,
578 Tok::Minus => BinOp::Sub,
579 Tok::Lt => BinOp::Lt,
580 Tok::Gt => BinOp::Gt,
581 _ => break,
582 };
583 self.pos += 1;
584 let rhs = self.parse_term()?;
585 lhs = Expr::Bin(op, Box::new(lhs), Box::new(rhs));
586 }
587 Ok(lhs)
588 }
589
590 // mult/div (precedencia alta)
591 fn parse_term(&mut self) -> Result<Expr, String> {
592 let mut lhs = self.parse_factor()?;
593 while let Some(t) = self.peek() {
594 let op = match t {
595 Tok::Star => BinOp::Mul,
596 Tok::Slash => BinOp::Div,
597 _ => break,
598 };
599 self.pos += 1;
600 let rhs = self.parse_factor()?;
601 lhs = Expr::Bin(op, Box::new(lhs), Box::new(rhs));
602 }
603 Ok(lhs)
604 }
605
606 fn parse_factor(&mut self) -> Result<Expr, String> {
607 match self.next() {
608 Some(Tok::Minus) => Ok(Expr::Un(UnOp::Neg, Box::new(self.parse_factor()?))),
609 Some(Tok::Num(n)) => Ok(Expr::Const(n)),
610 Some(Tok::LParen) => {
611 let e = self.parse_expr()?;
612 self.eat(&Tok::RParen)?;
613 Ok(e)
614 }
615 Some(Tok::Ident(name)) => {
616 // ¿llamada a función?
617 if self.peek() == Some(&Tok::LParen) {
618 self.pos += 1;
619 let mut args = Vec::new();
620 if self.peek() != Some(&Tok::RParen) {
621 args.push(self.parse_expr()?);
622 while self.peek() == Some(&Tok::Comma) {
623 self.pos += 1;
624 args.push(self.parse_expr()?);
625 }
626 }
627 self.eat(&Tok::RParen)?;
628 self.build_call(&name, args)
629 } else {
630 // constante especial, campo o parámetro
631 if name == "dt" {
632 return Ok(Expr::Dt);
633 }
634 if let Some(id) = self.sym.campo_id(&name) {
635 Ok(Expr::Field(id))
636 } else if let Some(id) = self.sym.param_id(&name) {
637 Ok(Expr::Param(id))
638 } else {
639 Err(format!("símbolo desconocido: {name}"))
640 }
641 }
642 }
643 other => Err(format!("token inesperado: {other:?}")),
644 }
645 }
646
647 /// Resuelve una llamada `nombre(args)` cuyo primer arg suele ser un campo.
648 fn build_call(&self, name: &str, args: Vec<Expr>) -> Result<Expr, String> {
649 // funciones sobre un campo (arg = nombre de campo → id)
650 let campo_arg = |args: &[Expr]| -> Result<FieldId, String> {
651 match args.first() {
652 Some(Expr::Field(f)) => Ok(*f),
653 _ => Err(format!("{name}(...) espera un campo como argumento")),
654 }
655 };
656 match name {
657 "lap" => Ok(Expr::Laplacian(campo_arg(&args)?)),
658 "avg" => Ok(Expr::Vecinos(Reduce::Promedio, campo_arg(&args)?)),
659 "sum6" => Ok(Expr::Vecinos(Reduce::Suma, campo_arg(&args)?)),
660 "min6" => Ok(Expr::Vecinos(Reduce::Min, campo_arg(&args)?)),
661 "max6" => Ok(Expr::Vecinos(Reduce::Max, campo_arg(&args)?)),
662 "abs" => Ok(Expr::Un(UnOp::Abs, Box::new(single(name, args)?))),
663 "exp" => Ok(Expr::Un(UnOp::Exp, Box::new(single(name, args)?))),
664 "sqrt" => Ok(Expr::Un(UnOp::Sqrt, Box::new(single(name, args)?))),
665 "min" => {
666 let (a, b) = pair(name, args)?;
667 Ok(Expr::Bin(BinOp::Min, Box::new(a), Box::new(b)))
668 }
669 "max" => {
670 let (a, b) = pair(name, args)?;
671 Ok(Expr::Bin(BinOp::Max, Box::new(a), Box::new(b)))
672 }
673 "clamp" => {
674 if args.len() != 3 {
675 return Err("clamp(x, lo, hi) espera 3 argumentos".into());
676 }
677 let mut it = args.into_iter();
678 Ok(Expr::Clamp(
679 Box::new(it.next().unwrap()),
680 Box::new(it.next().unwrap()),
681 Box::new(it.next().unwrap()),
682 ))
683 }
684 _ => {
685 if let Some(d) = Dir::de_nombre(name) {
686 Ok(Expr::Dir(d, campo_arg(&args)?))
687 } else {
688 Err(format!("función desconocida: {name}"))
689 }
690 }
691 }
692 }
693}
694
695fn single(name: &str, args: Vec<Expr>) -> Result<Expr, String> {
696 if args.len() != 1 {
697 return Err(format!("{name}(x) espera 1 argumento"));
698 }
699 Ok(args.into_iter().next().unwrap())
700}
701fn pair(name: &str, args: Vec<Expr>) -> Result<(Expr, Expr), String> {
702 if args.len() != 2 {
703 return Err(format!("{name}(a, b) espera 2 argumentos"));
704 }
705 let mut it = args.into_iter();
706 Ok((it.next().unwrap(), it.next().unwrap()))
707}
708
709impl Expr {
710 /// Parsea una fórmula de texto a un [`Expr`], resolviendo campos y params por
711 /// nombre con la tabla `sym`. Ej: `"grav*abajo(agua) + horiz*(avg(agua) - agua)"`.
712 pub fn parse(src: &str, sym: &Symbols) -> Result<Expr, String> {
713 let toks = lex(src)?;
714 let mut p = Parser { toks, pos: 0, sym };
715 let e = p.parse_expr()?;
716 if p.pos != p.toks.len() {
717 return Err(format!("sobra entrada tras la expresión (token {:?})", p.peek()));
718 }
719 Ok(e)
720 }
721
722 /// Imprime el `Expr` como fórmula de texto (inversa de [`parse`](Self::parse)),
723 /// usando los nombres de `sym`. Base de la barra de fórmula y del round‑trip.
724 pub fn to_source(&self, sym: &Symbols) -> String {
725 self.print(sym, 0)
726 }
727
728 /// Evalúa la expresión como **escalar sobre un entorno plano** (sin rejilla): un
729 /// vector de variables (`Field(id) → vars[id]`) y otro de parámetros
730 /// (`Param(id) → params[id]`), más `dt`. Es la evaluación que usan los **agentes**
731 /// (conducta autorable): reusa el mismo AST/parser/printer que las Leyes de rejilla,
732 /// pero las variables son *percepts* del agente, no campos espaciales. Los términos
733 /// espaciales (`lap`/`avg`/`dir…`) no tienen sentido aquí y evalúan a `0`. La semántica
734 /// numérica (división por cero → 0, `clamp` con límites invertidos) coincide con el
735 /// evaluador de opcodes de [`FieldEngine`].
736 pub fn eval_scalar(&self, vars: &[f32], params: &[f32], dt: f32) -> f32 {
737 match self {
738 Expr::Const(c) => *c,
739 Expr::Field(f) => vars.get(*f as usize).copied().unwrap_or(0.0),
740 Expr::Param(p) => params.get(*p as usize).copied().unwrap_or(0.0),
741 Expr::Dt => dt,
742 // Términos espaciales: sin significado para un agente puntual.
743 Expr::Laplacian(_) | Expr::Vecinos(_, _) | Expr::Dir(_, _) => 0.0,
744 Expr::Un(op, a) => {
745 let a = a.eval_scalar(vars, params, dt);
746 match op {
747 UnOp::Neg => -a,
748 UnOp::Abs => a.abs(),
749 UnOp::Exp => a.exp(),
750 UnOp::Sqrt => a.max(0.0).sqrt(),
751 }
752 }
753 Expr::Bin(op, a, b) => {
754 let l = a.eval_scalar(vars, params, dt);
755 let r = b.eval_scalar(vars, params, dt);
756 match op {
757 BinOp::Add => l + r,
758 BinOp::Sub => l - r,
759 BinOp::Mul => l * r,
760 BinOp::Div => {
761 if r == 0.0 {
762 0.0
763 } else {
764 l / r
765 }
766 }
767 BinOp::Min => l.min(r),
768 BinOp::Max => l.max(r),
769 BinOp::Gt => (l > r) as i32 as f32,
770 BinOp::Lt => (l < r) as i32 as f32,
771 }
772 }
773 Expr::Clamp(x, lo, hi) => {
774 let x = x.eval_scalar(vars, params, dt);
775 let lo = lo.eval_scalar(vars, params, dt);
776 let hi = hi.eval_scalar(vars, params, dt);
777 x.clamp(lo.min(hi), lo.max(hi))
778 }
779 }
780 }
781
782 // prec: 0 top, 1 suma/cmp, 2 mult, 3 unario/atómico — para poner paréntesis mínimos
783 fn print(&self, sym: &Symbols, prec: u8) -> String {
784 match self {
785 Expr::Const(c) => fmt_num(*c),
786 Expr::Field(f) => sym.campo_nombre(*f),
787 Expr::Param(p) => sym.param_nombre(*p),
788 Expr::Dt => "dt".to_string(),
789 Expr::Laplacian(f) => format!("lap({})", sym.campo_nombre(*f)),
790 Expr::Vecinos(r, f) => format!("{}({})", r.nombre(), sym.campo_nombre(*f)),
791 Expr::Dir(d, f) => format!("{}({})", d.nombre(), sym.campo_nombre(*f)),
792 Expr::Un(UnOp::Neg, a) => {
793 let s = format!("-{}", a.print(sym, 3));
794 paren_if(prec > 1, s)
795 }
796 Expr::Un(op, a) => {
797 let name = match op {
798 UnOp::Abs => "abs",
799 UnOp::Exp => "exp",
800 UnOp::Sqrt => "sqrt",
801 UnOp::Neg => unreachable!(),
802 };
803 format!("{name}({})", a.print(sym, 0))
804 }
805 Expr::Bin(op, a, b) => {
806 let (sym_op, my_prec) = match op {
807 BinOp::Add => ("+", 1),
808 BinOp::Sub => ("-", 1),
809 BinOp::Lt => ("<", 1),
810 BinOp::Gt => (">", 1),
811 BinOp::Mul => ("*", 2),
812 BinOp::Div => ("/", 2),
813 // min/max/... como funciones
814 BinOp::Min => {
815 return format!("min({}, {})", a.print(sym, 0), b.print(sym, 0));
816 }
817 BinOp::Max => {
818 return format!("max({}, {})", a.print(sym, 0), b.print(sym, 0));
819 }
820 };
821 let s = format!("{} {sym_op} {}", a.print(sym, my_prec), b.print(sym, my_prec + 1));
822 paren_if(prec > my_prec, s)
823 }
824 Expr::Clamp(x, lo, hi) => format!(
825 "clamp({}, {}, {})",
826 x.print(sym, 0),
827 lo.print(sym, 0),
828 hi.print(sym, 0)
829 ),
830 }
831 }
832}
833
834fn fmt_num(c: f32) -> String {
835 if c.fract() == 0.0 && c.abs() < 1e7 {
836 format!("{}", c as i64)
837 } else {
838 let s = format!("{c}");
839 s
840 }
841}
842fn paren_if(cond: bool, s: String) -> String {
843 if cond {
844 format!("({s})")
845 } else {
846 s
847 }
848}
849
850#[cfg(test)]
851mod tests {
852 use super::*;
853
854 fn sym(campos: &[&str], params: &[&str]) -> Symbols {
855 Symbols {
856 campos: campos.iter().map(|s| s.to_string()).collect(),
857 params: params.iter().map(|s| s.to_string()).collect(),
858 }
859 }
860
861 #[test]
862 fn round_trip_parse_print() {
863 let s = sym(&["agua"], &["grav", "horiz"]);
864 let casos = [
865 "grav * abajo(agua) + horiz * (avg(agua) - agua)",
866 "lap(agua) - agua * agua",
867 "clamp(agua + dt, 0, 1)",
868 "min(agua, 1) + max(agua, 0)",
869 ];
870 for c in casos {
871 let e = Expr::parse(c, &s).expect("parsea");
872 let printed = e.to_source(&s);
873 let e2 = Expr::parse(&printed, &s).expect("re-parsea");
874 assert_eq!(e, e2, "round-trip estable para: {c} → {printed}");
875 }
876 }
877
878 #[test]
879 fn eval_scalar_entorno_plano() {
880 // Las mismas fórmulas que usa una conducta autorable, evaluadas sobre un
881 // entorno de percepts (Field) + params (Param), sin rejilla.
882 let s = sym(&["cercania_amenaza"], &["miedo", "gregario"]);
883 // peso de huida por defecto: miedo * 2.5 * cercania_amenaza
884 let e = Expr::parse("miedo * 2.5 * cercania_amenaza", &s).unwrap();
885 // vars[0]=cercania_amenaza=0.8 ; params: miedo=1, gregario=0.35
886 let w = e.eval_scalar(&[0.8], &[1.0, 0.35], 1.0 / 60.0);
887 assert!((w - 2.0).abs() < 1e-4, "peso={w}");
888 // términos espaciales evalúan a 0 (no rompen).
889 let e2 = Expr::parse("lap(cercania_amenaza) + miedo", &s).unwrap();
890 assert!((e2.eval_scalar(&[0.8], &[1.0, 0.35], 0.0) - 1.0).abs() < 1e-4);
891 // clamp y umbral funcionan igual que en la rejilla.
892 let e3 = Expr::parse("clamp(gregario, 0, 1) + (cercania_amenaza > 0.5)", &s).unwrap();
893 assert!((e3.eval_scalar(&[0.8], &[1.0, 0.35], 0.0) - 1.35).abs() < 1e-4);
894 }
895
896 #[test]
897 fn parse_errores_utiles() {
898 let s = sym(&["agua"], &[]);
899 assert!(Expr::parse("desconocido", &s).is_err());
900 assert!(Expr::parse("agua +", &s).is_err());
901 assert!(Expr::parse("lap()", &s).is_err());
902 assert!(Expr::parse("agua ) (", &s).is_err());
903 }
904
905 #[test]
906 fn difusion_de_calor_conserva_y_converge() {
907 // Un campo `t` con ecuación de calor: dt/dt = k*lap(t). Neumann conserva el
908 // total y converge a uniforme.
909 let s = sym(&["t"], &["k"]);
910 let e = Expr::parse("k * lap(t)", &s).unwrap();
911 let prog = Program::compile(&[Assign { campo: 0, expr: e }]);
912 let mut eng = FieldEngine::new([8, 1, 8], vec![FieldDef::new("t", 0.0, 0.0, 1000.0)]);
913 eng.set(0, 4, 0, 4, 100.0); // foco de calor
914 let total0 = eng.total(0);
915 for _ in 0..200 {
916 eng.step(&prog, &[0.1], 1.0);
917 }
918 let total1 = eng.total(0);
919 assert!((total0 - total1).abs() < 1e-2, "conserva calor: {total0} vs {total1}");
920 // convergió: el foco se enfrió y el borde se calentó (menos varianza).
921 let centro = eng.get(0, 4, 0, 4);
922 let borde = eng.get(0, 0, 0, 0);
923 assert!(centro < 50.0, "el foco se difundió (centro={centro})");
924 assert!(borde > 0.1, "el calor llegó al borde (borde={borde})");
925 }
926
927 #[test]
928 fn crecer_como_campo() {
929 // La ley Crecer generalizada: la altura sube a `vel` hasta un tope. Ecuación:
930 // dh/dt = vel * (h < tope). (umbral → 0/1)
931 let s = sym(&["h"], &["vel", "tope"]);
932 let e = Expr::parse("vel * (h < tope)", &s).unwrap();
933 let prog = Program::compile(&[Assign { campo: 0, expr: e }]);
934 let mut eng = FieldEngine::new([1, 1, 1], vec![FieldDef::new("h", 0.0, 0.0, 100.0)]);
935 for _ in 0..20 {
936 eng.step(&prog, &[1.0, 5.0], 1.0);
937 }
938 // sube a razón 1/tick, se detiene al llegar a tope=5.
939 let h = eng.get(0, 0, 0, 0);
940 assert!((h - 5.0).abs() < 1.5, "creció y se detuvo cerca del tope: {h}");
941 }
942
943 #[test]
944 fn reaccion_difusion_gray_scott_hace_patron() {
945 // Gray–Scott: dos campos u,v.
946 // du/dt = Du*lap(u) - u*v*v + F*(1-u)
947 // dv/dt = Dv*lap(v) + u*v*v - (F+k)*v
948 let s = sym(&["u", "v"], &["Du", "Dv", "F", "k"]);
949 let eu = Expr::parse("Du * lap(u) - u * v * v + F * (1 - u)", &s).unwrap();
950 let ev = Expr::parse("Dv * lap(v) + u * v * v - (F + k) * v", &s).unwrap();
951 let prog = Program::compile(&[
952 Assign { campo: 0, expr: eu },
953 Assign { campo: 1, expr: ev },
954 ]);
955 let mut eng = FieldEngine::new(
956 [32, 1, 32],
957 vec![
958 FieldDef::new("u", 1.0, 0.0, 1.0),
959 FieldDef::new("v", 0.0, 0.0, 1.0),
960 ],
961 );
962 // sembrar una mancha de v en el centro
963 for z in 14..18 {
964 for x in 14..18 {
965 eng.set(1, x, 0, z, 0.5);
966 eng.set(0, x, 0, z, 0.25);
967 }
968 }
969 let params = [0.16, 0.08, 0.06, 0.062];
970 for _ in 0..400 {
971 eng.step(&prog, ¶ms, 1.0);
972 }
973 // Hubo estructura: v tiene varianza espacial (no todo cero ni todo igual).
974 let mut vmin = f32::INFINITY;
975 let mut vmax = f32::NEG_INFINITY;
976 for z in 0..32 {
977 for x in 0..32 {
978 let v = eng.get(1, x, 0, z);
979 vmin = vmin.min(v);
980 vmax = vmax.max(v);
981 assert!(v.is_finite(), "estable numéricamente");
982 }
983 }
984 assert!(vmax - vmin > 0.05, "reacción-difusión formó patrón (rango v={})", vmax - vmin);
985 }
986
987 #[test]
988 fn determinismo_independiente_del_orden() {
989 // Dos corridas idénticas dan exactamente lo mismo (Jacobi, sin azar).
990 let s = sym(&["t"], &["k"]);
991 let e = Expr::parse("k * lap(t)", &s).unwrap();
992 let prog = Program::compile(&[Assign { campo: 0, expr: e }]);
993 let run = || {
994 let mut eng = FieldEngine::new([6, 1, 6], vec![FieldDef::new("t", 0.0, 0.0, 100.0)]);
995 eng.set(0, 1, 0, 2, 50.0);
996 eng.set(0, 4, 0, 3, 30.0);
997 for _ in 0..50 {
998 eng.step(&prog, &[0.2], 1.0);
999 }
1000 (0..36).map(|i| eng.cur[i]).collect::<Vec<_>>()
1001 };
1002 assert_eq!(run(), run(), "misma evolución bit a bit");
1003 }
1004}