Skip to main content

cas_expr/
lib.rs

1//! cas 表达式层(P0-M1):arena + hash-consing、规范形构造、确定性全序。
2//!
3//! # 核心不变量
4//!
5//! - **规范形即身份**:构造经 intern 表去重,同一 Context 内结构相等 ⇔
6//!   索引相等([`Expr::raw_id`] 相同),相等性 O(1)。
7//! - **构造即规范**:`Add`/`Mul` 在构造器内完成扁平化、全序排序、精确数值
8//!   折叠、同类项/同底幂合并;任何时刻取出的 `Expr` 都是规范形。
9//! - **确定性**:全序见 [`order`] 模块文档(`ORDER_VERSION = 1`);全库哈希
10//!   用固定键(禁 `RandomState`),哈希迭代不进入输出。
11//! - **浮点边界**:Float 仅作字面量(值恒非负,负号由系数 -1 携带),
12//!   不参与任何算术折叠;`-0.0` 归一为 `0.0`;NaN 拒绝。
13//!
14//! arena 只增不减(P0 不回收,P4 评估压缩);`Context` 按问题作用域创建。
15//! 详见 `cas/DESIGN.md` D1/D6。
16
17use std::cell::{Ref, RefCell};
18use std::cmp::Ordering;
19use std::collections::HashMap;
20use std::fmt;
21use std::rc::Rc;
22
23mod assume;
24mod calculus;
25mod canonical;
26mod eval;
27mod factor;
28mod hash;
29mod node;
30mod order;
31mod poly_bridge;
32#[cfg(feature = "test-gen")]
33pub mod test_gen;
34mod transform;
35
36use hash::FxBuild;
37pub(crate) use node::Node;
38
39pub use crate::assume::{Assumptions, Predicate, Trinary};
40pub use cas_domain::{Integer, Rational};
41
42/// arena 与全部查表的宿主。类型名公开是为 `IntoExpr` 等签名可提及;
43/// 字段全部 crate 私有,外部无法构造或窥视。
44pub struct Inner {
45    pub(crate) nodes: Vec<Node>,
46    pub(crate) hashes: Vec<u64>,
47    pub(crate) intern: HashMap<u64, Vec<u32>, FxBuild>,
48    pub(crate) args_slab: Vec<u32>,
49    pub(crate) syms: HashMap<Box<str>, u32, FxBuild>,
50    pub(crate) sym_names: Vec<Box<str>>,
51    pub(crate) fn_heads: HashMap<Box<str>, u32, FxBuild>,
52    pub(crate) fn_names: Vec<Box<str>>,
53    pub(crate) sym_assumptions: HashMap<u32, assume::Assumptions>,
54}
55
56impl Inner {
57    fn new() -> Self {
58        Inner {
59            nodes: vec![Node::Int(Integer::zero())], // id 0 占位,永不引用
60            hashes: vec![0],
61            intern: HashMap::with_hasher(FxBuild::default()),
62            args_slab: Vec::new(),
63            syms: HashMap::with_hasher(FxBuild::default()),
64            sym_names: Vec::new(),
65            fn_heads: HashMap::with_hasher(FxBuild::default()),
66            fn_names: Vec::new(),
67            sym_assumptions: HashMap::new(),
68        }
69    }
70}
71
72/// 表达式构造与求值的宿主。按问题作用域创建;`Clone` 得到共享同一 arena
73/// 的另一个入口(不复制数据)。
74#[derive(Clone)]
75pub struct Context {
76    inner: Rc<RefCell<Inner>>,
77}
78
79impl Context {
80    pub fn new() -> Self {
81        Context {
82            inner: Rc::new(RefCell::new(Inner::new())),
83        }
84    }
85
86    fn wrap(&self, id: u32) -> Expr {
87        Expr {
88            ctx: Rc::clone(&self.inner),
89            id,
90        }
91    }
92
93    fn with(&self, f: impl FnOnce(&mut Inner) -> u32) -> Expr {
94        let id = f(&mut self.inner.borrow_mut());
95        self.wrap(id)
96    }
97
98    /// 声明符号。名字须匹配 `[A-Za-z_][A-Za-z0-9_]*`(违反属编程错误,直接 panic)。
99    pub fn sym(&self, name: &str) -> Expr {
100        validate_name(name, "符号");
101        self.with(|c| c.sym_node(name))
102    }
103
104    pub fn int(&self, v: i64) -> Expr {
105        self.with(|c| c.lit_int(v))
106    }
107
108    pub fn integer(&self, v: &Integer) -> Expr {
109        self.with(|c| c.lit_integer(v))
110    }
111
112    /// 带假设声明符号(D5):闭包补全(even⇒integer⇒…)+ 冲突检测;
113    /// 同名重设不同假设报错(返回 Err 由调用方决定——MATLAB assume
114    /// 的覆盖语义是有意不采用的,见 DESIGN D5)。
115    pub fn sym_with(&self, name: &str, preds: &[assume::Predicate]) -> Result<Expr, String> {
116        validate_name(name, "符号");
117        let a = assume::Assumptions::union(preds);
118        let e = self.sym(name);
119        let sid = self.inner.borrow().syms.get(name).copied();
120        if let Some(sid) = sid {
121            self.inner
122                .borrow_mut()
123                .sym_assumptions
124                .entry(sid)
125                .and_modify(|old| {
126                    assert!(
127                        *old == a,
128                        "符号 {name} 的假设只能一次性设定(当前 {:?},重设 {:?})",
129                        old,
130                        a
131                    );
132                })
133                .or_insert(a);
134        }
135        Ok(e)
136    }
137
138    /// 三值假设查询(D5):True / False / Unknown。
139    pub fn query(&self, e: &Expr, p: assume::Predicate) -> assume::Trinary {
140        let inner = self.inner.borrow();
141        inner.query_at(e.id, p, 0)
142    }
143
144    /// 浮点字面量。负值返回 `-1 * |v|` 的规范形(Float 节点恒非负);
145    /// NaN 属编程错误。
146    pub fn float(&self, v: f64) -> Expr {
147        assert!(!v.is_nan(), "Float 字面量禁止 NaN");
148        if v < 0.0 {
149            self.with(|c| {
150                let f = c.lit_float(-v);
151                let neg = c.lit_int(-1);
152                c.make_mul(&[neg, f])
153            })
154        } else {
155            self.with(|c| c.lit_float(v))
156        }
157    }
158
159    /// 有理数字面量;分母为零返回 `None`。
160    pub fn rational(&self, n: i64, d: i64) -> Option<Expr> {
161        let r = Rational::from_ints(&Integer::from_i64(n), &Integer::from_i64(d))?;
162        Some(self.with(|c| c.lit_rational(&r)))
163    }
164
165    /// n 元加法(运算符重载的批量化入口)。
166    pub fn add(&self, args: &[Expr]) -> Expr {
167        let ids: Vec<u32> = args.iter().map(|e| e.id).collect();
168        self.with(|c| c.make_add(&ids))
169    }
170
171    /// n 元乘法。
172    pub fn mul(&self, args: &[Expr]) -> Expr {
173        let ids: Vec<u32> = args.iter().map(|e| e.id).collect();
174        self.with(|c| c.make_mul(&ids))
175    }
176
177    pub fn pow(&self, base: &Expr, exp: &Expr) -> Expr {
178        self.with(|c| c.make_pow(base.id, exp.id))
179    }
180
181    /// 函数应用。参数保持语义次序(不排序);head 命名规则同符号。
182    pub fn call(&self, head: &str, args: &[Expr]) -> Expr {
183        validate_name(head, "函数");
184        let ids: Vec<u32> = args.iter().map(|e| e.id).collect();
185        self.with(|c| c.fn_node(head, &ids))
186    }
187
188    /// 表达式全序(同一 Context 内;跨 Context 用 `Expr::eq` 的结构比较)。
189    pub fn cmp_expr(&self, a: &Expr, b: &Expr) -> Ordering {
190        let inner = self.inner.borrow();
191        inner.cmp_ids(a.id, b.id)
192    }
193
194    pub fn eq_expr(&self, a: &Expr, b: &Expr) -> bool {
195        self.cmp_expr(a, b) == Ordering::Equal
196    }
197
198    /// arena 规模(节点数;诊断与后续基准用)。
199    pub fn node_count(&self) -> usize {
200        self.inner.borrow().nodes.len()
201    }
202
203    /// 只读检视。持有 Inspector 期间不可构造新表达式(构造方独占借用 arena)。
204    pub fn inspect<R>(&self, f: impl FnOnce(&Inspector<'_>) -> R) -> R {
205        f(&Inspector {
206            inner: self.inner.borrow(),
207        })
208    }
209
210    // ── L1 显式变换与求值 ──────────────────────────────────────
211
212    /// 展开(L1):分配乘积 over 和式、展开非负整数幂(带规模守卫)、
213    /// 函数参数内部展开。每步经规范形构造器,输出为规范形。
214    pub fn expand(&self, e: &Expr) -> Expr {
215        self.with(|c| c.expand_at(e.id, 0))
216    }
217
218    /// 有理函数约化(L1):分子分母的多项式公因子(gcd)约去。
219    /// 非多项式因子原样保留;无可约化时返回与输入同一节点。
220    pub fn cancel(&self, e: &Expr) -> Expr {
221        self.with(|c| c.cancel_at(e.id))
222    }
223
224    /// 多项式因式分解(M4/M5):`cont · Π 因子^重数`,重建为规范形。
225    /// 1/2 变元完全分解,≥3 变元部分(内容+平方自由);非多项式输入返回原节点。
226    pub fn factor(&self, e: &Expr) -> Expr {
227        self.with(|c| c.factor_at(e.id))
228    }
229
230    /// 符号求导(P1):和/积/链式/幂(整数、有理、一般指数)与初等
231    /// 函数表(sin cos tan exp log sqrt,abs→sign)。输出规范形。
232    pub fn diff(&self, e: &Expr, x: &Expr) -> Expr {
233        match x_name_of(self, x) {
234            Some(name) => {
235                let sid = self.inner.borrow().syms.get(name.as_str()).copied();
236                match sid {
237                    Some(s) => self.with(|c| c.diff_at(e.id, s, 0)),
238                    None => self.with(|c| c.lit_int(0)),
239                }
240            }
241            None => self.with(|c| c.lit_int(0)),
242        }
243    }
244
245    /// Taylor 级数(P1):x=a 处的多项式部分到 (x−a)^order(含)。
246    pub fn taylor(&self, e: &Expr, x: &Expr, at: i64, order: u32) -> Expr {
247        match x_name_of(self, x) {
248            Some(name) => {
249                let sid = self.inner.borrow().syms.get(name.as_str()).copied();
250                match sid {
251                    Some(s) => self.with(|c| c.taylor_at(e.id, s, at, order)),
252                    None => e.clone(),
253                }
254            }
255            None => e.clone(),
256        }
257    }
258
259    /// 定向化简(P1,L2):常量折叠 + 无条件恒等式(sin²+cos²→1),
260    /// 自底向上单遍,确定性输出。
261    pub fn simplify(&self, e: &Expr) -> Expr {
262        self.with(|c| c.simplify_at(e.id, 0))
263    }
264
265    /// 代换:按符号名替换子表达式(替换值须属同一 Context)。
266    /// 重建经规范形构造器,`subst(e, x→x)` 与 `e` 同节点。
267    pub fn subst(&self, e: &Expr, map: &[(&str, Expr)]) -> Expr {
268        let m = {
269            let inner = self.inner.borrow();
270            map.iter()
271                .filter_map(|(name, val)| inner.syms.get(*name).map(|&s| (s, val.id)))
272                .collect::<std::collections::HashMap<u32, u32>>()
273        };
274        self.with(|c| c.subst_at(e.id, &m, 0))
275    }
276
277    /// ℚ 上精确求值。Float 字面量与函数节点返回 `None`(D2 域边界)。
278    pub fn eval_rational(&self, e: &Expr, vals: &[(&str, Rational)]) -> Option<Rational> {
279        let inner = self.inner.borrow();
280        let m = vals
281            .iter()
282            .filter_map(|(n, v)| inner.syms.get(*n).map(|&s| (s, v.clone())))
283            .collect::<std::collections::HashMap<u32, Rational>>();
284        eval::eval_rational_at(&inner, e.id, &m, 0)
285    }
286
287    /// f64 数值求值(差分测试的语义判据通道)。支持初等函数头。
288    pub fn eval_float(&self, e: &Expr, vals: &[(&str, f64)]) -> Option<f64> {
289        let inner = self.inner.borrow();
290        let m = vals
291            .iter()
292            .filter_map(|(n, v)| inner.syms.get(*n).map(|&s| (s, *v)))
293            .collect::<std::collections::HashMap<u32, f64>>();
294        eval::eval_float_at(&inner, e.id, &m, 0)
295    }
296
297    /// 把 `e` 视作 `vars` 的多项式,返回 `(指数向量, 系数)` 列表,指数按 `vars` 顺序。
298    ///
299    /// - 依赖 `e` 已 [`Context::expand`](未展开也能跑,但同幂次可能分散多项)。
300    /// - 同指数的项会合并;系数保持符号形式(不数值化)。
301    /// - `vars` 之外的自由符号被视为系数的一部分。
302    /// - 变量的非多项式用法返回 [`CasErrorKind::NonPolynomial`],负指数返回
303    ///   [`CasErrorKind::NegativeExponent`],变量名不存在返回 [`CasErrorKind::UnknownSymbol`]。
304    pub fn monomial_coeffs(
305        &self,
306        e: &Expr,
307        vars: &[&str],
308    ) -> Result<Vec<(Vec<u32>, Expr)>, CasError> {
309        let var_ids: Vec<u32> = {
310            let inner = self.inner.borrow();
311            let mut ids = Vec::with_capacity(vars.len());
312            for name in vars {
313                match inner.syms.get(*name) {
314                    Some(&s) => ids.push(s),
315                    None => {
316                        return Err(CasError::new(
317                            CasErrorKind::UnknownSymbol,
318                            format!("变量 {name} 不在该 Context 中"),
319                        ));
320                    }
321                }
322            }
323            ids
324        };
325        let monos = {
326            let inner = self.inner.borrow();
327            transform::monomials_at(&inner, e.id, &var_ids, 0)?
328        };
329
330        // 按指数分组:组内每个单项式的系数 = 其因子之积;同类项之间**相加**。
331        let mut groups: HashMap<Vec<u32>, Vec<Vec<u32>>> = HashMap::new();
332        let mut order: Vec<Vec<u32>> = Vec::new();
333        for m in monos {
334            match groups.get_mut(&m.exps) {
335                Some(list) => list.push(m.coef),
336                None => {
337                    order.push(m.exps.clone());
338                    groups.insert(m.exps, vec![m.coef]);
339                }
340            }
341        }
342        let mut out = Vec::with_capacity(order.len());
343        for exps in order {
344            let terms = groups.remove(&exps).unwrap_or_default();
345            let coef = self.with(|c| {
346                let ids: Vec<u32> = terms
347                    .iter()
348                    .map(|factors| {
349                        if factors.is_empty() {
350                            c.lit_int(1)
351                        } else {
352                            c.make_mul(factors)
353                        }
354                    })
355                    .collect();
356                if ids.len() == 1 {
357                    ids[0]
358                } else {
359                    c.make_add(&ids)
360                }
361            });
362            out.push((exps, coef));
363        }
364        Ok(out)
365    }
366}
367
368/// 从(预期为 Sym 的)Expr 提取符号名(非 Sym 返回 None)。
369fn x_name_of(ctx: &Context, x: &Expr) -> Option<String> {
370    ctx.inspect(|i| match i.kind(x.raw_id()) {
371        Kind::Sym(name) => Some(name.to_string()),
372        _ => None,
373    })
374}
375
376impl Default for Context {
377    fn default() -> Self {
378        Self::new()
379    }
380}
381
382fn validate_name(name: &str, what: &str) {
383    let mut chars = name.chars();
384    let valid = matches!(chars.next(), Some(c) if c.is_ascii_alphabetic() || c == '_')
385        && chars.all(|c| c.is_ascii_alphanumeric() || c == '_');
386    assert!(
387        valid && !name.is_empty(),
388        "{what}名非法: {name:?},须匹配 [A-Za-z_][A-Za-z0-9_]*"
389    );
390}
391
392/// 表达式句柄:arena 索引 + 共享上下文引用。
393///
394/// `Clone` 廉价(引用计数 + 4 字节索引)。同一 Context 内结构相等 ⇔ 索引相等;
395/// 跨 Context 的 `==` 走结构比较。不实现 `Hash`/`Ord`:句柄的序只在所属
396/// Context 内有意义,用 [`Context::cmp_expr`]。
397#[derive(Clone)]
398pub struct Expr {
399    ctx: Rc<RefCell<Inner>>,
400    id: u32,
401}
402
403impl Expr {
404    /// arena 索引。仅在同一 Context 内有意义(相等判据/诊断用)。
405    pub fn raw_id(&self) -> u32 {
406        self.id
407    }
408
409    fn bin(self, other: Expr, f: impl FnOnce(&mut Inner, u32, u32) -> u32) -> Expr {
410        let ctx = Rc::clone(&self.ctx);
411        let id = f(&mut ctx.borrow_mut(), self.id, other.id);
412        Expr { ctx, id }
413    }
414
415    /// 幂。指数可传 `Expr`/`&Expr`/整数字面量/浮点字面量。
416    pub fn pow<E: IntoExpr>(self, exp: E) -> Expr {
417        let ctx = Rc::clone(&self.ctx);
418        let mut inner = ctx.borrow_mut();
419        let eid = exp.into_id(&mut inner);
420        let id = inner.make_pow(self.id, eid);
421        drop(inner);
422        Expr { ctx, id }
423    }
424}
425
426/// 四种(自有/引用 × 自有/引用)组合的二元运算符转发。
427macro_rules! forward_binop {
428    ($tr:ident, $method:ident, $inner:expr) => {
429        impl std::ops::$tr<Expr> for Expr {
430            type Output = Expr;
431            fn $method(self, rhs: Expr) -> Expr {
432                self.bin(rhs, $inner)
433            }
434        }
435        impl std::ops::$tr<&Expr> for Expr {
436            type Output = Expr;
437            fn $method(self, rhs: &Expr) -> Expr {
438                self.bin(rhs.clone(), $inner)
439            }
440        }
441        impl std::ops::$tr<Expr> for &Expr {
442            type Output = Expr;
443            fn $method(self, rhs: Expr) -> Expr {
444                self.clone().bin(rhs, $inner)
445            }
446        }
447        impl std::ops::$tr<&Expr> for &Expr {
448            type Output = Expr;
449            fn $method(self, rhs: &Expr) -> Expr {
450                self.clone().bin(rhs.clone(), $inner)
451            }
452        }
453    };
454}
455
456forward_binop!(Add, add, |c, a, b| c.make_add(&[a, b]));
457forward_binop!(Sub, sub, |c, a, b| {
458    let neg1 = c.lit_int(-1);
459    let neg = c.make_mul(&[neg1, b]);
460    c.make_add(&[a, neg])
461});
462forward_binop!(Mul, mul, |c, a, b| c.make_mul(&[a, b]));
463forward_binop!(Div, div, |c, a, b| {
464    let neg1 = c.lit_int(-1);
465    let inv = c.make_pow(b, neg1);
466    c.make_mul(&[a, inv])
467});
468
469impl std::ops::Neg for Expr {
470    type Output = Expr;
471    fn neg(self) -> Expr {
472        let ctx = Rc::clone(&self.ctx);
473        let id = {
474            let mut c = ctx.borrow_mut();
475            let neg1 = c.lit_int(-1);
476            c.make_mul(&[neg1, self.id])
477        };
478        Expr { ctx, id }
479    }
480}
481
482impl std::ops::Neg for &Expr {
483    type Output = Expr;
484    fn neg(self) -> Expr {
485        (*self).clone().neg()
486    }
487}
488
489impl PartialEq for Expr {
490    fn eq(&self, other: &Self) -> bool {
491        if Rc::ptr_eq(&self.ctx, &other.ctx) {
492            return self.id == other.id;
493        }
494        let x = self.ctx.borrow();
495        let y = other.ctx.borrow();
496        order::deep_eq(&x, self.id, &y, other.id)
497    }
498}
499
500impl Eq for Expr {}
501
502impl fmt::Debug for Expr {
503    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
504        f.write_str(
505            &Inspector {
506                inner: self.ctx.borrow(),
507            }
508            .dump(self.id),
509        )
510    }
511}
512
513/// [`Expr::pow`] 可接受的指数类型。
514pub trait IntoExpr {
515    fn into_id(self, c: &mut Inner) -> u32;
516}
517
518impl IntoExpr for Expr {
519    fn into_id(self, _c: &mut Inner) -> u32 {
520        self.id
521    }
522}
523
524impl IntoExpr for &Expr {
525    fn into_id(self, _c: &mut Inner) -> u32 {
526        self.id
527    }
528}
529
530impl IntoExpr for i64 {
531    fn into_id(self, c: &mut Inner) -> u32 {
532        c.lit_int(self)
533    }
534}
535
536impl IntoExpr for i32 {
537    fn into_id(self, c: &mut Inner) -> u32 {
538        c.lit_int(self as i64)
539    }
540}
541
542impl IntoExpr for u32 {
543    fn into_id(self, c: &mut Inner) -> u32 {
544        c.lit_int(self as i64)
545    }
546}
547
548impl IntoExpr for f64 {
549    fn into_id(self, c: &mut Inner) -> u32 {
550        assert!(self >= 0.0, "pow 的浮点指数须非负(负浮点请用 Expr 形态)");
551        c.lit_float(self)
552    }
553}
554
555/// 按幂次提取系数时的错误类别。
556#[derive(Clone, Copy, Debug, PartialEq, Eq)]
557pub enum CasErrorKind {
558    /// 出现变量的非多项式用法(如 `sqrt(x)`、变量作分母)。
559    NonPolynomial,
560    /// 变量出现负指数(Laurent 情形,当前不支持)。
561    NegativeExponent,
562    /// 传入的变量名在该 Context 中不存在。
563    UnknownSymbol,
564}
565
566/// `monomial_coeffs` 的错误。
567#[derive(Clone, Debug)]
568pub struct CasError {
569    pub kind: CasErrorKind,
570    pub msg: String,
571}
572
573impl CasError {
574    pub fn new(kind: CasErrorKind, msg: impl Into<String>) -> Self {
575        CasError {
576            kind,
577            msg: msg.into(),
578        }
579    }
580}
581
582impl std::fmt::Display for CasError {
583    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
584        write!(f, "{:?}: {}", self.kind, self.msg)
585    }
586}
587
588impl std::error::Error for CasError {}
589
590/// 只读检视器:在 [`Context::inspect`] 闭包内遍历表达式结构。
591pub struct Inspector<'a> {
592    inner: Ref<'a, Inner>,
593}
594
595/// 节点的只读视图。复合节点给子节点索引,调用方继续经 Inspector 递归。
596#[derive(Debug)]
597pub enum Kind<'a> {
598    Int(&'a Integer),
599    Rat(&'a Rational),
600    Float(f64),
601    Sym(&'a str),
602    Fn { head: &'a str, args: &'a [u32] },
603    Pow { base: u32, exp: u32 },
604    Mul(&'a [u32]),
605    Add(&'a [u32]),
606}
607
608impl Inspector<'_> {
609    pub fn kind(&self, id: u32) -> Kind<'_> {
610        match &self.inner.nodes[id as usize] {
611            Node::Int(v) => Kind::Int(v),
612            Node::Rat(v) => Kind::Rat(v),
613            Node::Float { bits, .. } => Kind::Float(f64::from_bits(*bits)),
614            Node::Sym(s) => Kind::Sym(&self.inner.sym_names[*s as usize]),
615            Node::Fn { head, args } => Kind::Fn {
616                head: &self.inner.fn_names[*head as usize],
617                args: self.inner.node_args(*args),
618            },
619            Node::Pow { base, exp } => Kind::Pow {
620                base: *base,
621                exp: *exp,
622            },
623            Node::Mul { args } => Kind::Mul(self.inner.node_args(*args)),
624            Node::Add { args } => Kind::Add(self.inner.node_args(*args)),
625        }
626    }
627
628    pub fn cmp(&self, a: u32, b: u32) -> Ordering {
629        self.inner.cmp_ids(a, b)
630    }
631
632    /// 诊断转储:`add[mul[int(3), pow[sym(x), int(2)]], ...]` 风格。
633    pub fn dump(&self, id: u32) -> String {
634        let mut s = String::new();
635        self.dump_at(id, &mut s, 0);
636        s
637    }
638
639    fn dump_at(&self, id: u32, out: &mut String, depth: u32) {
640        assert!(depth <= 10_000, "表达式嵌套过深");
641        let sub = |c: &Self, cid: u32, o: &mut String, d: u32| c.dump_at(cid, o, d);
642        match self.kind(id) {
643            Kind::Int(v) => out.push_str(&format!("int({v})")),
644            Kind::Rat(r) => out.push_str(&format!("rat({r})")),
645            Kind::Float(v) => out.push_str(&format!("flt({v})")),
646            Kind::Sym(n) => out.push_str(&format!("sym({n})")),
647            Kind::Fn { head, args } => {
648                out.push_str(&format!("fn({head}, ["));
649                for (i, &a) in args.iter().enumerate() {
650                    if i > 0 {
651                        out.push_str(", ");
652                    }
653                    sub(self, a, out, depth + 1);
654                }
655                out.push_str("])");
656            }
657            Kind::Pow { base, exp } => {
658                out.push_str("pow(");
659                sub(self, base, out, depth + 1);
660                out.push_str(", ");
661                sub(self, exp, out, depth + 1);
662                out.push(')');
663            }
664            Kind::Mul(args) => {
665                out.push_str("mul[");
666                for (i, &a) in args.iter().enumerate() {
667                    if i > 0 {
668                        out.push_str(", ");
669                    }
670                    sub(self, a, out, depth + 1);
671                }
672                out.push(']');
673            }
674            Kind::Add(args) => {
675                out.push_str("add[");
676                for (i, &a) in args.iter().enumerate() {
677                    if i > 0 {
678                        out.push_str(", ");
679                    }
680                    sub(self, a, out, depth + 1);
681                }
682                out.push(']');
683            }
684        }
685    }
686}
687
688/// 批量声明符号:`let (x, y, z) = sym!(&ctx, x, y, z);`(单个名字直接返回 Expr)。
689#[macro_export]
690macro_rules! sym {
691    ($ctx:expr, $name:ident) => {
692        $ctx.sym(stringify!($name))
693    };
694    ($ctx:expr, $($name:ident),+ $(,)?) => {
695        ($($ctx.sym(stringify!($name))),+)
696    };
697}
698
699#[cfg(test)]
700mod tests {
701    use super::*;
702
703    #[test]
704    fn 规范形_同类合并与数值折叠() {
705        let ctx = Context::new();
706        let x = ctx.sym("x");
707        let y = ctx.sym("y");
708
709        // x + x 与 2*x 是同一节点
710        let a = x.clone() + x.clone();
711        let b = ctx.int(2) * x.clone();
712        assert!(a == b && a.raw_id() == b.raw_id());
713
714        // (a+b)+c 扁平化,c+(b+a) 同一规范形
715        let c = ctx.sym("c");
716        let p = (x.clone() + y.clone()) + c.clone();
717        let q = c + (y + x);
718        assert_eq!(p.raw_id(), q.raw_id());
719
720        // 数值折叠
721        let n = ctx.int(2) + ctx.int(3);
722        assert!(
723            ctx.inspect(
724                |i| matches!(i.kind(n.raw_id()), Kind::Int(v) if *v == Integer::from_i64(5))
725            )
726        );
727    }
728
729    #[test]
730    fn 规范形_同底幂合并() {
731        let ctx = Context::new();
732        let x = ctx.sym("x");
733        let a = ctx.sym("a");
734        let b = ctx.sym("b");
735
736        // x * x == x^2
737        let p = x.clone() * x.clone();
738        let q = x.clone().pow(2);
739        assert_eq!(p.raw_id(), q.raw_id());
740
741        // x^2 * x^3 == x^5
742        let p = x.clone().pow(2) * x.clone().pow(3);
743        let q = x.clone().pow(5);
744        assert_eq!(p.raw_id(), q.raw_id());
745
746        // x^a * x^b == x^(a+b)
747        let p = x.clone().pow(&a) * x.clone().pow(&b);
748        let q = x.clone().pow(&a + &b);
749        assert_eq!(p.raw_id(), q.raw_id());
750
751        // (x^a)^2 == x^(2*a);(x^a)^b 保持嵌套(非整数外幂)
752        let p = x.clone().pow(&a).pow(2);
753        let q = x.clone().pow(&ctx.int(2) * &a);
754        assert_eq!(p.raw_id(), q.raw_id());
755        let nested = x.clone().pow(&a).pow(&b);
756        assert!(ctx.inspect(|i| matches!(i.kind(nested.raw_id()), Kind::Pow { .. })));
757    }
758
759    #[test]
760    fn 规范形_数值幂折叠与守卫() {
761        let ctx = Context::new();
762        let p = ctx.int(2).pow(100);
763        assert!(ctx.inspect(|i| matches!(i.kind(p.raw_id()), Kind::Int(v)
764            if *v == Integer::parse("1267650600228229401496703205376").unwrap())));
765
766        // 负整数幂 → 有理数
767        let p = ctx.int(2).pow(-2);
768        assert!(ctx.inspect(|i| matches!(i.kind(p.raw_id()), Kind::Rat(_))));
769
770        // 0^0 = 1,0^3 = 0,0^-1 无定义:保留 Pow 节点
771        assert_eq!(ctx.int(0).pow(0).raw_id(), ctx.int(1).raw_id());
772        assert_eq!(ctx.int(0).pow(3).raw_id(), ctx.int(0).raw_id());
773        let zero_neg_pow = ctx.int(0).pow(-1);
774        assert!(ctx.inspect(|i| matches!(i.kind(zero_neg_pow.raw_id()), Kind::Pow { .. })));
775
776        // 规模守卫:超大指数不折叠
777        let huge = ctx.int(2).pow(1 << 30);
778        assert!(ctx.inspect(|i| matches!(i.kind(huge.raw_id()), Kind::Pow { .. })));
779    }
780
781    #[test]
782    fn 规范形_除法与负浮点() {
783        let ctx = Context::new();
784        let x = ctx.sym("x");
785
786        // x / 2 == 1/2 * x
787        let p = x.clone() / ctx.int(2);
788        let q = ctx.mul(&[ctx.rational(1, 2).unwrap(), x.clone()]);
789        assert_eq!(p.raw_id(), q.raw_id());
790
791        // 负浮点 = -1 * |v|
792        let p = ctx.float(-2.5);
793        let q = -ctx.float(2.5);
794        assert_eq!(p.raw_id(), q.raw_id());
795
796        // x/x == 1(sympy 同约定的退化点)
797        assert_eq!((x.clone() / x.clone()).raw_id(), ctx.int(1).raw_id());
798    }
799
800    #[test]
801    fn 全序_数值与符号() {
802        let ctx = Context::new();
803        let x = ctx.sym("x");
804        let y = ctx.sym("y");
805        let two = ctx.int(2);
806        let half = ctx.rational(1, 2).unwrap();
807        let three = ctx.int(3);
808        let f1 = ctx.float(1.5);
809
810        // 精确数按值:1/2 < 2 < 3;Float 排在精确数之后
811        assert_eq!(ctx.cmp_expr(&half, &two), Ordering::Less);
812        assert_eq!(ctx.cmp_expr(&two, &three), Ordering::Less);
813        assert_eq!(ctx.cmp_expr(&three, &f1), Ordering::Less);
814
815        // 符号按名字节序,与创建先后无关
816        let big = ctx.sym("zz");
817        assert_eq!(ctx.cmp_expr(&x, &y), Ordering::Less);
818        assert_eq!(ctx.cmp_expr(&y, &big), Ordering::Less);
819
820        // 数值在 Mul 中居首
821        let m = ctx.int(3) * x.clone() * y.clone();
822        let first = ctx.inspect(|i| match i.kind(m.raw_id()) {
823            Kind::Mul(args) => args[0],
824            k => panic!("期望 Mul: {k:?}"),
825        });
826        assert!(ctx.cmp_expr(&ctx.wrap(first), &x) == Ordering::Less);
827    }
828
829    #[test]
830    fn 全序_复合节点字典序() {
831        let ctx = Context::new();
832        let x = ctx.sym("x");
833        let y = ctx.sym("y");
834
835        let x2 = x.clone().pow(2);
836        let xy = x.clone() * y.clone();
837        let sum = x.clone() + y.clone();
838
839        // Sym < Pow < Mul < Add
840        assert_eq!(ctx.cmp_expr(&x, &x2), Ordering::Less);
841        assert_eq!(ctx.cmp_expr(&x2, &xy), Ordering::Less);
842        assert_eq!(ctx.cmp_expr(&xy, &sum), Ordering::Less);
843    }
844
845    // ── L1 变换与求值 ─────────────────────────────────────────
846
847    #[test]
848    fn 展开_黄金快照() {
849        // 断言用与手工构造表达式的节点恒等(比字符串更强)
850        let ctx = Context::new();
851        let x = ctx.sym("x");
852        let y = ctx.sym("y");
853
854        let e = (x.clone() + y.clone()) * (x.clone() - y.clone());
855        let g = ctx.expand(&e);
856        let expect = x.clone().pow(2) - y.clone().pow(2);
857        assert_eq!(g.raw_id(), expect.raw_id());
858
859        let e = (x.clone() + y.clone()).pow(2);
860        let g = ctx.expand(&e);
861        let expect = x.clone().pow(2) + y.clone().pow(2) + ctx.int(2) * x.clone() * y.clone();
862        assert_eq!(g.raw_id(), expect.raw_id());
863
864        // 立方:交叉项系数正确(含 3 与 1/3 的往返核对)
865        let e = (x.clone() + y.clone()).pow(3);
866        let g = ctx.expand(&e);
867        let expect = x.clone().pow(3)
868            + y.clone().pow(3)
869            + ctx.int(3) * x.clone() * y.clone().pow(2)
870            + ctx.int(3) * x.clone().pow(2) * y.clone();
871        assert_eq!(g.raw_id(), expect.raw_id());
872    }
873
874    #[test]
875    fn 展开_语义一致与幂等() {
876        let ctx = Context::new();
877        let x = ctx.sym("x");
878        let y = ctx.sym("y");
879        let z = ctx.sym("z");
880
881        let e = ((x.clone() + y.clone()).pow(3) - z.clone() * (x.clone() + y.clone()))
882            * (x.clone() - z.clone());
883        let g = ctx.expand(&e);
884        // 幂等
885        let g2 = ctx.expand(&g);
886        assert_eq!(g.raw_id(), g2.raw_id());
887        // 与原式在随机点精确同值
888        for (xv, yv, zv) in [(2, -3, 5), (-1, 4, 7), (0, 9, -11)] {
889            let pts = [
890                ("x", Rational::from_integer(&Integer::from_i64(xv))),
891                ("y", Rational::from_integer(&Integer::from_i64(yv))),
892                ("z", Rational::from_integer(&Integer::from_i64(zv))),
893            ];
894            let a = ctx.eval_rational(&e, &pts);
895            let b = ctx.eval_rational(&g, &pts);
896            assert_eq!(a, b, "展开改变语义");
897        }
898
899        // (规模守卫在展开_四元二十次幂与 huge 用例覆盖)
900        let huge = (x.clone() + y.clone()).pow(20_000);
901        assert!(ctx.inspect(|i| matches!(i.kind(huge.raw_id()), Kind::Pow { .. })));
902    }
903
904    #[test]
905    fn 展开_四元二十次幂() {
906        // 设计基准负载:项数 C(23,3) = 1771
907        let ctx = Context::new();
908        let (x, y, z, w) = crate::sym!(&ctx, x, y, z, w);
909        let base = ctx.int(1) + x + y + z + w;
910        let g = ctx.expand(&base.pow(20));
911        let n = ctx.inspect(|i| match i.kind(g.raw_id()) {
912            Kind::Add(args) => args.len(),
913            _ => 0,
914        });
915        assert_eq!(n, 10_626); // C(24,4)
916    }
917
918    #[test]
919    fn 代换() {
920        let ctx = Context::new();
921        let x = ctx.sym("x");
922        let y = ctx.sym("y");
923
924        let e = x.clone() + y.clone();
925        let g = ctx.subst(&e, &[("x", x.clone().pow(2))]);
926        let expect = x.clone().pow(2) + y.clone();
927        assert_eq!(g.raw_id(), expect.raw_id());
928
929        // 恒等代换回原节点(intern 去重)
930        let same = ctx.subst(&e, &[("x", x.clone())]);
931        assert_eq!(same.raw_id(), e.raw_id());
932
933        // 对换
934        let p = x.clone() * y.clone();
935        let q = ctx.subst(&p, &[("x", y.clone()), ("y", x.clone())]);
936        assert_eq!(q.raw_id(), p.raw_id()); // x*y 交换后仍同节点
937    }
938
939    #[test]
940    fn 展开快慢路径同构() {
941        // M3 快路径(poly 域)与朴素 DAG 分配必须产出同一节点
942        let ctx = Context::new();
943        let (x, y, z, w) = crate::sym!(&ctx, x, y, z, w);
944        let syms = [x.clone(), y.clone(), z.clone(), w.clone()];
945
946        let mut xs = 777u64;
947        let mut nxt = move || {
948            xs ^= xs << 13;
949            xs ^= xs >> 7;
950            xs ^= xs << 17;
951            xs
952        };
953        fn gen_poly(
954            ctx: &Context,
955            syms: &[Expr],
956            nxt: &mut impl FnMut() -> u64,
957            depth: u32,
958        ) -> Expr {
959            if depth == 0 || nxt() % 3 == 0 {
960                match nxt() % 3 {
961                    0 => ctx.int((nxt() % 19) as i64 - 9),
962                    1 => ctx
963                        .rational((nxt() % 15) as i64 - 7, (nxt() % 8) as i64 + 2)
964                        .unwrap(),
965                    _ => syms[(nxt() % syms.len() as u64) as usize].clone(),
966                }
967            } else {
968                match nxt() % 3 {
969                    0 => gen_poly(ctx, syms, nxt, depth - 1) + gen_poly(ctx, syms, nxt, depth - 1),
970                    1 => gen_poly(ctx, syms, nxt, depth - 1) * gen_poly(ctx, syms, nxt, depth - 1),
971                    _ => gen_poly(ctx, syms, nxt, depth - 1).pow((nxt() % 5) as i64),
972                }
973            }
974        }
975        let slow_of = |ctx: &Context, e: &Expr| -> u32 {
976            let mut inner = ctx.inner.borrow_mut();
977            inner.expand_impl(e.raw_id(), 0, false)
978        };
979
980        for _ in 0..60 {
981            let e = gen_poly(&ctx, &syms, &mut nxt, 4);
982            let fast = ctx.expand(&e);
983            let slow = slow_of(&ctx, &e);
984            assert_eq!(fast.raw_id(), slow, "快慢路径不同构: {e:?}");
985        }
986
987        // 混合表达式(函数/浮点使快路径部分让位)仍须同构
988        let e = ctx.call("sin", std::slice::from_ref(&x)) * (y.clone() + z.clone()).pow(5)
989            + ctx.float(1.5) * (x.clone() + w.clone()).pow(3);
990        assert_eq!(ctx.expand(&e).raw_id(), slow_of(&ctx, &e));
991    }
992
993    #[test]
994    fn 约化() {
995        let ctx = Context::new();
996        let x = ctx.sym("x");
997        let y = ctx.sym("y");
998        let z = ctx.sym("z");
999
1000        // (x^2 - y^2)/(x - y) → x + y(节点恒等)
1001        let e = (x.clone().pow(2) - y.clone().pow(2)) / (x.clone() - y.clone());
1002        let expect = x.clone() + y.clone();
1003        assert_eq!(ctx.cancel(&e).raw_id(), expect.raw_id());
1004
1005        // x*y/(x*z) → y/z(分母为复合 Mul 的情形)
1006        let e = (x.clone() * y.clone()) / (x.clone() * z.clone());
1007        let expect = y.clone() / z.clone();
1008        assert_eq!(ctx.cancel(&e).raw_id(), expect.raw_id());
1009
1010        // (x*y + x*z)/x → y + z
1011        let e = (x.clone() * y.clone() + x.clone() * z.clone()) / x.clone();
1012        let expect = y.clone() + z.clone();
1013        assert_eq!(ctx.cancel(&e).raw_id(), expect.raw_id());
1014
1015        // 互素:返回原节点
1016        let e = (x.clone() + y.clone()) / x.clone();
1017        assert_eq!(ctx.cancel(&e).raw_id(), e.raw_id());
1018
1019        // 语义保持:约化前后在随机点精确同值
1020        let e = (x.clone().pow(2) * y.clone() - y.clone().pow(3))
1021            / (x.clone() * y.clone() + y.clone().pow(2));
1022        let g = ctx.cancel(&e);
1023        for (xv, yv) in [(2, 3), (-4, 5), (7, -2)] {
1024            let pts = [
1025                ("x", Rational::from_integer(&Integer::from_i64(xv))),
1026                ("y", Rational::from_integer(&Integer::from_i64(yv))),
1027            ];
1028            assert_eq!(ctx.eval_rational(&e, &pts), ctx.eval_rational(&g, &pts));
1029        }
1030
1031        // 浮点因子保留:1.5*(x^2 - 1)/(x - 1) → 1.5*(x + 1)
1032        let e = ctx.float(1.5) * ((x.clone().pow(2) - ctx.int(1)) / (x.clone() - ctx.int(1)));
1033        let g = ctx.cancel(&e);
1034        let v = ctx.eval_float(&g, &[("x", 3.0)]).unwrap();
1035        assert!((v - 1.5 * 4.0).abs() < 1e-12);
1036    }
1037
1038    #[test]
1039    fn 因式分解() {
1040        let ctx = Context::new();
1041        let x = ctx.sym("x");
1042
1043        // x^2 − 1 → (x − 1)(x + 1)(节点恒等)
1044        let e = x.clone().pow(2) - ctx.int(1);
1045        let g = ctx.factor(&e);
1046        let expect = (x.clone() - ctx.int(1)) * (x.clone() + ctx.int(1));
1047        assert_eq!(g.raw_id(), expect.raw_id());
1048
1049        // (x^2 − 1)^2 → (x − 1)^2 (x + 1)^2
1050        let e = x.clone().pow(2) - ctx.int(1);
1051        let e = e.pow(2);
1052        let g = ctx.factor(&e);
1053        let expect = (x.clone() - ctx.int(1)).pow(2) * (x.clone() + ctx.int(1)).pow(2);
1054        assert_eq!(g.raw_id(), expect.raw_id());
1055
1056        // x^4 + 4 在 ℚ 上可约(Sophie Germain)
1057        let e = x.clone().pow(4) + ctx.int(4);
1058        let g = ctx.factor(&e);
1059        let expect = (x.clone().pow(2) - ctx.int(2) * x.clone() + ctx.int(2))
1060            * (x.clone().pow(2) + ctx.int(2) * x.clone() + ctx.int(2));
1061        assert_eq!(g.raw_id(), expect.raw_id());
1062
1063        // x^4 + 1 不可约:重建后应与原式同节点
1064        let e = x.clone().pow(4) + ctx.int(1);
1065        let g = ctx.factor(&e);
1066        let expect = x.clone().pow(4) + ctx.int(1);
1067        assert_eq!(g.raw_id(), expect.raw_id());
1068
1069        // 双变元:x^2 − y^2 → 两个一次因子(节点形态随 pp 归一约定,
1070        // 断言用语义恒等 + 因子可观测形态)
1071        let y = ctx.sym("y");
1072        let e = x.clone().pow(2) - y.clone().pow(2);
1073        let g = ctx.factor(&e);
1074        let s = ctx.inspect(|i| i.dump(g.raw_id()));
1075        assert_eq!(
1076            s,
1077            "mul[int(-1), add[sym(x), sym(y)], add[sym(y), mul[int(-1), sym(x)]]]"
1078        );
1079        for (xv, yv) in [(3, 2), (-5, 7)] {
1080            let pts = [
1081                ("x", Rational::from_integer(&Integer::from_i64(xv))),
1082                ("y", Rational::from_integer(&Integer::from_i64(yv))),
1083            ];
1084            assert_eq!(ctx.eval_rational(&e, &pts), ctx.eval_rational(&g, &pts));
1085        }
1086    }
1087
1088    #[test]
1089    fn 精确与数值求值() {
1090        let ctx = Context::new();
1091        let x = ctx.sym("x");
1092
1093        let e = ctx.rational(3, 2).unwrap() * x.clone() + ctx.int(2);
1094        let v = ctx.eval_rational(&e, &[("x", Rational::from_integer(&Integer::from_i64(4)))]);
1095        assert_eq!(v.map(|r| r.to_string()), Some("8".to_string()));
1096
1097        // 负幂
1098        let e = x.clone().pow(-2);
1099        let v = ctx.eval_rational(&e, &[("x", Rational::from_integer(&Integer::from_i64(2)))]);
1100        assert_eq!(v.map(|r| r.to_string()), Some("1/4".to_string()));
1101
1102        // Float 与函数不进精确求值
1103        assert_eq!(ctx.eval_rational(&ctx.float(1.5), &[]), None);
1104        let s = ctx.call("sin", std::slice::from_ref(&x));
1105        assert_eq!(ctx.eval_rational(&s, &[("x", Rational::zero())]), None);
1106
1107        // 数值求值:初等函数 + 与展开一致性
1108        let sv = ctx.eval_float(&s, &[("x", 0.5)]);
1109        assert!((sv.unwrap() - 0.5f64.sin()).abs() < 1e-15);
1110        let e = (x.clone() + ctx.int(1)).pow(10);
1111        let g = ctx.expand(&e);
1112        let pts = [("x", 1.7f64)];
1113        let a = ctx.eval_float(&e, &pts).unwrap();
1114        let b = ctx.eval_float(&g, &pts).unwrap();
1115        assert!((a - b).abs() < 1e-9 * a.abs().max(1.0));
1116
1117        // log 非正数 → None
1118        let lg = ctx.call("log", std::slice::from_ref(&x));
1119        assert_eq!(ctx.eval_float(&lg, &[("x", -1.0)]), None);
1120    }
1121}