use rustc_hash::{FxHashMap, FxHashSet};
use smallvec::SmallVec;
use num_traits::One;
use num_traits::Signed;
use num_traits::Zero;
use crate::base::arena::{Arena, FN_CHI, FN_ERFI, FN_SHI};
use crate::base::node::{ExprId, ExprNode, SymbolId};
use crate::base::numeric::Q;
use crate::base::stage::stage;
use crate::poly::Poly;
pub(crate) fn integrate(arena: &mut Arena, expr: ExprId, var: ExprId) -> ExprId {
let call = IntegrateCall::enter(arena);
if call.depth > MAX_INTEGRATE_NESTING {
tracing::debug!("integrate: nesting limit reached");
return arena.intern(ExprNode::Integral(expr, var));
}
let _jumps = crate::transforms::trig_integ::JumpMemo::enter();
stage!(arena, "integrate", expr, integrate_impl(arena, expr, var))
}
const MAX_INTEGRATE_NESTING: usize = 48;
struct IntegrateCall {
outermost: bool,
depth: usize,
}
impl IntegrateCall {
fn enter(arena: &Arena) -> Self {
let depth = NODE_BUDGET_DEPTH.with(|d| {
let depth = d.get();
d.set(depth + 1);
depth
});
let outermost = depth == 0;
if outermost {
NODE_BUDGET_USED.with(|u| u.set(0));
NODE_BUDGET_CAP.with(|c| c.set(usize::MAX));
ARENA_BUDGET_START.with(|s| s.set(arena.node_count()));
INTEGRATE_MEMO.with(|m| m.borrow_mut().clear());
VERDICT_MEMO.with(|m| m.borrow_mut().clear());
CLEARED_MEMO.with(|m| m.borrow_mut().clear());
VARIES_MEMO.with(|m| m.borrow_mut().clear());
}
IntegrateCall { outermost, depth }
}
}
impl Drop for IntegrateCall {
fn drop(&mut self) {
let _ = NODE_BUDGET_DEPTH.try_with(|d| d.set(d.get().saturating_sub(1)));
if self.outermost {
let _ = INTEGRATE_MEMO.try_with(|m| {
if let Ok(mut m) = m.try_borrow_mut() {
m.clear();
}
});
let _ = VERDICT_MEMO.try_with(|m| {
if let Ok(mut m) = m.try_borrow_mut() {
m.clear();
}
});
let _ = CLEARED_MEMO.try_with(|m| {
if let Ok(mut m) = m.try_borrow_mut() {
m.clear();
}
});
let _ = VARIES_MEMO.try_with(|m| {
if let Ok(mut m) = m.try_borrow_mut() {
m.clear();
}
});
}
}
}
const INTEGRATE_NODE_BUDGET: usize = 20_000;
const INTEGRATE_ARENA_BUDGET: usize = 2_000_000;
thread_local! {
static NODE_BUDGET_DEPTH: std::cell::Cell<usize> = const { std::cell::Cell::new(0) };
static NODE_BUDGET_USED: std::cell::Cell<usize> = const { std::cell::Cell::new(0) };
static ARENA_BUDGET_START: std::cell::Cell<usize> = const { std::cell::Cell::new(0) };
static NODE_BUDGET_CAP: std::cell::Cell<usize> = const { std::cell::Cell::new(usize::MAX) };
}
struct StepCap(usize);
impl StepCap {
fn enter(steps: usize) -> Self {
let used = NODE_BUDGET_USED.with(std::cell::Cell::get);
let old = NODE_BUDGET_CAP.with(|c| {
let old = c.get();
c.set(old.min(used.saturating_add(steps)));
old
});
StepCap(old)
}
fn active() -> bool {
NODE_BUDGET_CAP.with(std::cell::Cell::get) != usize::MAX
}
}
impl Drop for StepCap {
fn drop(&mut self) {
let _ = NODE_BUDGET_CAP.try_with(|c| c.set(self.0));
}
}
fn arena_budget_exhausted(arena: &Arena) -> bool {
let start = ARENA_BUDGET_START.with(std::cell::Cell::get);
let exhausted = arena.node_count().saturating_sub(start) > INTEGRATE_ARENA_BUDGET;
if exhausted {
tracing::debug!("integrate: arena budget exhausted");
}
exhausted
}
fn node_budget_exhausted() -> bool {
NODE_BUDGET_USED.with(|u| {
let used = u.get() + 1;
u.set(used);
used > INTEGRATE_NODE_BUDGET || used > NODE_BUDGET_CAP.with(std::cell::Cell::get)
})
}
fn integrate_impl(arena: &mut Arena, expr: ExprId, var: ExprId) -> ExprId {
let var_sym = match arena.node(var) {
ExprNode::Symbol(sid) => *sid,
_ => {
return arena.intern(ExprNode::Integral(expr, var));
}
};
let expr = crate::transforms::eval::eval(arena, expr);
if crate::base::walk::post_order_ids(arena, expr)
.iter()
.any(|&id| {
matches!(arena.node(id), ExprNode::RootSum(..)) && contains_var(arena, id, var_sym)
})
{
return arena.intern(ExprNode::Integral(expr, var));
}
if let Some((constants, f)) = split_constant_factors(arena, expr, var_sym) {
let big_f = integrate_stages(arena, f, var, var_sym);
if matches!(arena.node(big_f), ExprNode::Integral(..)) {
return arena.intern(ExprNode::Integral(expr, var));
}
return wrap_with_constants(arena, big_f, &constants);
}
integrate_stages(arena, expr, var, var_sym)
}
fn integrate_stages(arena: &mut Arena, expr: ExprId, var: ExprId, var_sym: SymbolId) -> ExprId {
let result = integrate_stages_once(arena, expr, var, var_sym);
let result = if let ExprNode::Integral(_, _) = arena.node(result)
&& let Some(reduced) = crate::calculus::risch::reduced_fraction_terms(arena, expr, var)
{
let r = stage!(
arena,
"reduced_fraction_retry",
reduced,
integrate_stages_once(arena, reduced, var, var_sym)
);
if crate::base::walk::has_unevaluated(arena, r) {
result
} else {
accept_or_unevaluated(arena, expr, r, var, var_sym)
}
} else {
result
};
let result = if crate::base::walk::has_unevaluated(arena, result) {
integrate_terms_separately(arena, expr, var, var_sym).unwrap_or(result)
} else {
result
};
stage!(
arena,
"piecewise_wrap",
result,
try_piecewise_wrap(arena, result, expr, var, var_sym)
)
}
fn integrate_terms_separately(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let sum = match arena.node(expr) {
ExprNode::Add(_) => expr,
ExprNode::Mul(children)
if children
.iter()
.any(|&c| matches!(arena.node(c), ExprNode::Add(_))) =>
{
crate::transforms::expand::expand(arena, expr)
}
_ => return None,
};
let ExprNode::Add(terms) = arena.node(sum).clone() else {
return None;
};
if terms.len() < 2 || terms.len() > MAX_SEPARATE_TERMS {
return None;
}
let mut parts = Vec::with_capacity(terms.len());
for &t in &terms {
if !contains_var(arena, t, var_sym) {
parts.push(arena.mul(&[t, var]));
continue;
}
let (constants, f) = match split_constant_factors(arena, t, var_sym) {
Some((c, f)) => (c, f),
None => (SmallVec::new(), t),
};
let big_f = integrate_stages(arena, f, var, var_sym);
if crate::base::walk::has_unevaluated(arena, big_f) {
return None;
}
parts.push(wrap_with_constants(arena, big_f, &constants));
}
tracing::debug!("integrate: the terms of a sum, each through every stage");
let sum = arena.add(&parts);
let r = accept_or_unevaluated(arena, expr, sum, var, var_sym);
(!crate::base::walk::has_unevaluated(arena, r)).then_some(r)
}
const MAX_SEPARATE_TERMS: usize = 8;
fn integrate_stages_once(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> ExprId {
let result = integrate_node(arena, expr, var, var_sym, 20);
let result = accept_or_unevaluated(arena, expr, result, var, var_sym);
let result = if let ExprNode::Integral(_, _) = arena.node(result)
&& let Some(r) = stage!(
arena,
"substitutions",
expr,
try_substitution_strategies(arena, expr, var, var_sym)
) {
accept_or_unevaluated(arena, expr, r, var, var_sym)
} else {
result
};
if let ExprNode::Integral(_, _) = arena.node(result) {
match stage!(
arena,
"risch_tower",
expr,
crate::calculus::risch::try_risch_tower(arena, expr, var)
) {
crate::calculus::risch::TowerResult::Elementary(id) => {
accept_or_unevaluated(arena, expr, id, var, var_sym)
}
crate::calculus::risch::TowerResult::NonElementary => {
result
}
crate::calculus::risch::TowerResult::NotApplicable => {
if let Some(heurisch_result) = stage!(
arena,
"heurisch",
expr,
crate::transforms::heurisch::heurisch_integrate(arena, expr, var, var_sym)
) {
accept_or_unevaluated(arena, expr, heurisch_result, var, var_sym)
} else {
result
}
}
}
} else {
result
}
}
fn accept_or_unevaluated(
arena: &mut Arena,
f: ExprId,
candidate: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> ExprId {
let unevaluated = arena.intern(ExprNode::Integral(f, var));
if candidate == unevaluated {
return candidate;
}
let candidate = analytic_log_of_complex_arguments(arena, candidate, var_sym);
if candidate_rejected(arena, f, candidate, var, var_sym) {
tracing::debug!("integrate: rejecting an unverified closed form");
return unevaluated;
}
candidate
}
fn candidate_rejected(
arena: &mut Arena,
f: ExprId,
big_f: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> bool {
if introduces_root_of(arena, f, big_f) {
tracing::debug!("integrate: rejecting a closed form over implicit roots");
return true;
}
let big_f = analytic_log_of_complex_arguments(arena, big_f, var_sym);
antiderivative_rejected(arena, f, big_f, var, var_sym)
}
fn introduces_root_of(arena: &Arena, f: ExprId, big_f: ExprId) -> bool {
let root_ofs = |e: ExprId| -> Vec<ExprId> {
crate::base::walk::post_order_ids(arena, e)
.into_iter()
.filter(|&id| matches!(arena.node(id), ExprNode::RootOf(..)))
.collect()
};
let new = root_ofs(big_f);
!new.is_empty() && {
let old: FxHashSet<ExprId> = root_ofs(f).into_iter().collect();
new.iter().any(|id| !old.contains(id))
}
}
fn analytic_log_of_complex_arguments(arena: &mut Arena, e: ExprId, var_sym: SymbolId) -> ExprId {
let mut out = e;
let mut reals = crate::base::assumptions::AssumptionCache::new();
for id in crate::base::walk::post_order_ids(arena, e) {
if let ExprNode::Ln(inner) = arena.node(id).clone()
&& let ExprNode::Abs(u) = arena.node(inner).clone()
&& contains_var(arena, u, var_sym)
&& !log_argument_is_real(arena, u, var_sym, &mut reals)
{
let plain = arena.ln(u);
out = arena.subs_structural(out, id, plain);
}
}
out
}
fn log_argument_is_real(
arena: &mut Arena,
u: ExprId,
var_sym: SymbolId,
reals: &mut crate::base::assumptions::AssumptionCache,
) -> bool {
let mut stack = vec![u];
let mut seen: FxHashSet<ExprId> = FxHashSet::default();
while let Some(id) = stack.pop() {
if !seen.insert(id) {
continue;
}
if contains_var(arena, id, var_sym) {
arena.node(id).for_each_child(|c| stack.push(c));
} else if !constant_is_real(arena, id, reals) {
return false;
}
}
true
}
fn constant_is_real(
arena: &mut Arena,
c: ExprId,
reals: &mut crate::base::assumptions::AssumptionCache,
) -> bool {
crate::transforms::realness::constant_realness_as_declared(arena, c, FTC_CHECK_DIGITS, reals)
== Some(true)
}
fn is_by_parts_candidate(arena: &Arena, expr: ExprId, var: ExprId, var_sym: SymbolId) -> bool {
if is_polynomial_in(arena, expr, var, var_sym) {
return true;
}
match arena.node(expr) {
ExprNode::Ln(inner)
| ExprNode::Asin(inner)
| ExprNode::Acos(inner)
| ExprNode::Atan(inner) => contains_var(arena, *inner, var_sym),
_ => false,
}
}
fn liate_rank(arena: &Arena, expr: ExprId, var: ExprId, var_sym: SymbolId) -> u8 {
match arena.node(expr) {
ExprNode::Ln(_) => 1,
ExprNode::Asin(_) | ExprNode::Acos(_) | ExprNode::Atan(_) => 2,
_ if is_polynomial_in(arena, expr, var, var_sym) => 3,
ExprNode::Sin(_) | ExprNode::Cos(_) | ExprNode::Tan(_) => 4,
ExprNode::Exp(_) | ExprNode::Sinh(_) | ExprNode::Cosh(_) | ExprNode::Tanh(_) => 5,
_ => 6,
}
}
fn is_var_squared(arena: &Arena, expr: ExprId, var: ExprId) -> bool {
if let ExprNode::Pow(b, e) = arena.node(expr)
&& *b == var
&& let Some(n) = arena.as_num(*e)
{
return *n == num_rational::Ratio::from_integer(2.into());
}
false
}
fn is_neg_var_squared(arena: &Arena, expr: ExprId, var: ExprId) -> bool {
match arena.node(expr).clone() {
ExprNode::Neg(inner) => is_var_squared(arena, inner, var),
ExprNode::Mul(ref children) if children.len() == 2 => {
let neg_one_val = num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
let has_neg_one = children
.iter()
.any(|&c| arena.as_num(c).is_some_and(|n| *n == neg_one_val));
let has_var_sq = children.iter().any(|&c| is_var_squared(arena, c, var));
has_neg_one && has_var_sq
}
_ => false,
}
}
fn try_standard_form_integral(
arena: &mut Arena,
_expr: ExprId,
base: ExprId,
exp: ExprId,
var: ExprId,
_var_sym: SymbolId,
) -> Option<ExprId> {
let exp_val = arena.as_num(exp)?.clone();
let neg_one = num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
let neg_half = num_rational::Ratio::<num_bigint::BigInt>::new((-1).into(), 2.into());
let is_neg_one = exp_val == neg_one;
let is_neg_half = exp_val == neg_half;
if !is_neg_one && !is_neg_half {
return None;
}
let children = match arena.node(base).clone() {
ExprNode::Add(c) if c.len() == 2 => c,
_ => return None,
};
let mut const_val: Option<Q> = None;
let mut has_pos_x2 = false;
let mut has_neg_x2 = false;
for &child in children.iter() {
if let Some(n) = arena.as_num(child) {
const_val = Some(n.clone());
} else if is_var_squared(arena, child, var) {
has_pos_x2 = true;
} else if is_neg_var_squared(arena, child, var) {
has_neg_x2 = true;
} else {
return None;
}
}
let c_val = const_val?;
if !has_pos_x2 && !has_neg_x2 {
return None;
}
let a_squared = c_val.abs();
if a_squared.is_zero() {
return None;
}
let a_sq_is_one = a_squared == num_rational::Ratio::<num_bigint::BigInt>::one();
let x_over_a = if a_sq_is_one {
var
} else {
let a_sq_id = arena.num_ratio(a_squared.clone());
let nh = arena.rational(-1, 2);
let a_inv = arena.pow(a_sq_id, nh); arena.mul(&[var, a_inv])
};
let one_over_a = if a_sq_is_one {
None
} else {
let a_sq_id = arena.num_ratio(a_squared.clone());
let nh = arena.rational(-1, 2);
Some(arena.pow(a_sq_id, nh))
};
if has_pos_x2 && c_val.is_positive() {
if is_neg_one {
let atan_val = arena.atan(x_over_a);
return Some(match one_over_a {
Some(inv_a) => arena.mul(&[inv_a, atan_val]),
None => atan_val,
});
}
if is_neg_half {
let asinh_val = arena.asinh(x_over_a);
return Some(asinh_val);
}
}
if has_neg_x2 && c_val.is_positive() {
if is_neg_half {
let asin_val = arena.asin(x_over_a);
return Some(asin_val);
}
if is_neg_one {
let atanh_val = arena.atanh(x_over_a);
return Some(match one_over_a {
Some(inv_a) => arena.mul(&[inv_a, atanh_val]),
None => atanh_val,
});
}
}
if has_pos_x2 && c_val.is_negative() && is_neg_half {
return Some(ln_abs_x_plus_sqrt(arena, var, base));
}
None
}
fn ln_abs_x_plus_sqrt(arena: &mut Arena, v: ExprId, q: ExprId) -> ExprId {
let half = arena.rational(1, 2);
let sqrt_q = arena.pow(q, half);
let sum = arena.add(&[v, sqrt_q]);
let abs_sum = arena.abs(sum);
arena.ln(abs_sum)
}
fn try_complete_square_integral(
arena: &mut Arena,
base: ExprId,
exp: ExprId,
var: ExprId,
_var_sym: SymbolId,
) -> Option<ExprId> {
let exp_r = arena.as_num(exp)?.clone();
let neg_one = num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
let neg_half = num_rational::Ratio::<num_bigint::BigInt>::new((-1).into(), 2.into());
let is_neg_one = exp_r == neg_one;
let is_neg_half = exp_r == neg_half;
if !is_neg_one && !is_neg_half {
return None;
}
let poly_opt = crate::poly::polybridge::expr_to_poly(arena, base, var);
if poly_opt.is_none() && is_neg_one {
tracing::debug!(
"try_complete_square: expr_to_poly failed, trying symbolic quadratic coefficients"
);
if let Some((c_id, d_id, e_id)) = symbolic_quadratic_coeffs(arena, base, var, _var_sym) {
tracing::debug!(
"try_complete_square: symbolic quadratic coefficients extracted, completing the square"
);
let b_sym = arena.div(d_id, c_id);
let c_norm = arena.div(e_id, c_id);
let two = arena.int(2);
let half_b = arena.div(b_sym, two);
let half_b_sq = arena.mul(&[half_b, half_b]);
let disc = arena.sub(c_norm, half_b_sq);
let disc = crate::transforms::eval::eval(arena, disc);
let disc_sign = crate::poly::algebraic::sign_checked(arena, disc);
tracing::debug!(?disc_sign, "try_complete_square: discriminant sign");
if disc_sign != Some(1) {
tracing::debug!(
"try_complete_square: discriminant non-positive, atan form not applicable"
);
return None;
}
let half = arena.rational(1, 2);
let sqrt_disc = arena.pow(disc, half);
let shifted = arena.add(&[var, half_b]);
let ratio = arena.div(shifted, sqrt_disc);
let atan_val = arena.atan(ratio);
let sqrt_disc2 = arena.pow(disc, half);
let a_sqrt_d = arena.mul(&[c_id, sqrt_disc2]);
tracing::debug!(
"try_complete_square: symbolic completing-the-square succeeded → atan form"
);
return Some(arena.div(atan_val, a_sqrt_d));
}
}
let poly = poly_opt?;
if poly.degree()? != 2 {
return None;
}
let a_coeff = poly.coeff(2);
let b_coeff = poly.coeff(1);
let c_coeff = poly.coeff(0);
if a_coeff.is_zero() {
return None;
}
if is_neg_one {
let b = &b_coeff / &a_coeff;
let c = &c_coeff / &a_coeff;
if b.is_zero() {
return None;
}
let half_b = &b / &num_rational::Ratio::from_integer(2.into());
let d = &c - &(&half_b * &half_b);
if d.is_zero() || d.is_negative() {
return None; }
let half_b_id = {
let nid = arena.intern_num(half_b.clone());
arena.intern(crate::base::node::ExprNode::Num(nid))
};
let shifted = arena.add(&[var, half_b_id]);
let d_id = {
let nid = arena.intern_num(d.clone());
arena.intern(crate::base::node::ExprNode::Num(nid))
};
let half = arena.rational(1, 2);
let sqrt_d = arena.pow(d_id, half);
let ratio = arena.div(shifted, sqrt_d);
let atan_result = arena.atan(ratio);
let a_id = {
let nid = arena.intern_num(a_coeff);
arena.intern(crate::base::node::ExprNode::Num(nid))
};
let sqrt_d2 = arena.pow(d_id, half);
let a_sqrt_d = arena.mul(&[a_id, sqrt_d2]);
return Some(arena.div(atan_result, a_sqrt_d));
}
let two_r = num_rational::Ratio::<num_bigint::BigInt>::from_integer(2.into());
let four_r = num_rational::Ratio::<num_bigint::BigInt>::from_integer(4.into());
let shift = &b_coeff / &(&a_coeff * &two_r); let d = &c_coeff - &(&b_coeff * &b_coeff / &(&a_coeff * &four_r));
let u_expr = if shift.is_zero() {
var
} else {
let shift_id = arena.num_ratio(shift.clone());
arena.add(&[var, shift_id])
};
let half = arena.rational(1, 2);
if a_coeff.is_positive() {
let a_id = arena.num_ratio(a_coeff.clone());
let sqrt_a = arena.pow(a_id, half); let inv_sqrt_a = {
let neg_half_e = arena.rational(-1, 2);
arena.pow(a_id, neg_half_e)
};
if d.is_positive() {
let d_id = arena.num_ratio(d.clone());
let sqrt_d = arena.pow(d_id, half);
let u_sqrt_a = arena.mul(&[u_expr, sqrt_a]);
let arg = arena.div(u_sqrt_a, sqrt_d);
let asinh_val = arena.asinh(arg);
return Some(arena.mul(&[inv_sqrt_a, asinh_val]));
} else if d.is_negative() {
let u_sqrt_a = arena.mul(&[u_expr, sqrt_a]);
let ln_val = ln_abs_x_plus_sqrt(arena, u_sqrt_a, base);
return Some(arena.mul(&[inv_sqrt_a, ln_val]));
} else {
let abs_u = arena.abs(u_expr);
let ln_u = arena.ln(abs_u);
let sign_u = arena.sign(u_expr);
return Some(arena.mul(&[inv_sqrt_a, sign_u, ln_u]));
}
} else if a_coeff.is_negative() && d.is_positive() {
let abs_a = a_coeff.abs();
let abs_a_id = arena.num_ratio(abs_a.clone());
let sqrt_abs_a = arena.pow(abs_a_id, half);
let inv_sqrt_abs_a = {
let neg_half_e = arena.rational(-1, 2);
arena.pow(abs_a_id, neg_half_e)
};
let d_id = arena.num_ratio(d.clone());
let sqrt_d = arena.pow(d_id, half);
let u_sqrt_abs_a = arena.mul(&[u_expr, sqrt_abs_a]);
let arg = arena.div(u_sqrt_abs_a, sqrt_d);
let asin_val = arena.asin(arg);
return Some(arena.mul(&[inv_sqrt_abs_a, asin_val]));
}
None
}
fn try_trig_recip_product(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
if dependent.len() != 2 {
return None;
}
let neg_two = num_rational::Ratio::<num_bigint::BigInt>::from_integer((-2).into());
for (i, j) in [(0usize, 1usize), (1, 0)] {
let node_i = arena.node(dependent[i]).clone();
let node_j = arena.node(dependent[j]).clone();
if let ExprNode::Sin(inner_sin) = node_i
&& let ExprNode::Pow(base_j, exp_j) = node_j
&& let ExprNode::Cos(inner_cos) = arena.node(base_j).clone()
&& inner_sin == inner_cos
&& let Some(e) = arena.as_num(exp_j)
&& *e == neg_two
{
let neg_one_e = arena.int(-1);
if inner_sin == var {
return Some(arena.pow(base_j, neg_one_e));
} else if let Some((a_expr, _)) =
symbolic_linear_coeff_of(arena, inner_sin, var, var_sym)
{
let recip = arena.pow(base_j, neg_one_e);
return Some(arena.div(recip, a_expr));
}
}
if let ExprNode::Cos(inner_cos) = node_i
&& let ExprNode::Pow(base_j, exp_j) = node_j
&& let ExprNode::Sin(inner_sin) = arena.node(base_j).clone()
&& inner_cos == inner_sin
&& let Some(e) = arena.as_num(exp_j)
&& *e == neg_two
{
let neg_one_e = arena.int(-1);
if inner_cos == var {
let recip = arena.pow(base_j, neg_one_e);
return Some(arena.neg(recip));
} else if let Some((a_expr, _)) =
symbolic_linear_coeff_of(arena, inner_cos, var, var_sym)
{
let recip = arena.pow(base_j, neg_one_e);
let neg_recip = arena.neg(recip);
return Some(arena.div(neg_recip, a_expr));
}
}
}
None
}
fn try_x_over_sqrt_quadratic(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
if dependent.len() != 2 {
return None;
}
let pow_idx = if dependent[0] == var {
1
} else if dependent[1] == var {
0
} else {
return None;
};
let (base, exp_id) = match arena.node(dependent[pow_idx]).clone() {
ExprNode::Pow(b, e) => (b, e),
_ => return None,
};
let exp_r = arena.as_num(exp_id)?.clone();
let neg_half = num_rational::Ratio::<num_bigint::BigInt>::new((-1).into(), 2.into());
if exp_r != neg_half {
return None;
}
let poly = crate::poly::polybridge::expr_to_poly(arena, base, var)?;
if poly.degree()? != 2 {
return None;
}
let a_coeff = poly.coeff(2);
let b_coeff = poly.coeff(1);
if a_coeff.is_zero() {
return None;
}
let half = arena.rational(1, 2);
let sqrt_r = arena.pow(base, half);
let a_id = arena.num_ratio(a_coeff.clone());
let first_term = arena.div(sqrt_r, a_id);
if b_coeff.is_zero() {
return Some(first_term);
}
let i_0 = integrate_node(
arena,
dependent[pow_idx],
var,
var_sym,
depth.saturating_sub(1),
);
if matches!(arena.node(i_0), ExprNode::Integral(_, _)) {
return None;
}
let two_r = num_rational::Ratio::<num_bigint::BigInt>::from_integer(2.into());
let b_over_2a = &b_coeff / &(&a_coeff * &two_r);
let b_over_2a_id = arena.num_ratio(b_over_2a.clone());
let second_term = arena.mul(&[b_over_2a_id, i_0]);
Some(arena.sub(first_term, second_term))
}
fn try_trig_sub_sqrt_integral(
arena: &mut Arena,
base: ExprId,
exp: ExprId,
var: ExprId,
_var_sym: SymbolId,
) -> Option<ExprId> {
let exp_val = arena.as_num(exp)?.clone();
let pos_half = num_rational::Ratio::<num_bigint::BigInt>::new(1.into(), 2.into());
if exp_val != pos_half {
return None;
}
let children = match arena.node(base).clone() {
ExprNode::Add(c) if c.len() == 2 => c,
_ => return None,
};
let mut const_val: Option<Q> = None;
let mut has_pos_x2 = false;
let mut has_neg_x2 = false;
for &child in children.iter() {
if let Some(n) = arena.as_num(child) {
const_val = Some(n.clone());
} else if is_var_squared(arena, child, var) {
has_pos_x2 = true;
} else if is_neg_var_squared(arena, child, var) {
has_neg_x2 = true;
} else {
return None;
}
}
let c_val = const_val?;
if !has_pos_x2 && !has_neg_x2 {
return None;
}
let a_squared = c_val.abs();
if a_squared.is_zero() {
return None;
}
let half = arena.rational(1, 2);
let a_sq_is_one = a_squared == num_rational::Ratio::<num_bigint::BigInt>::one();
let sqrt_base = arena.pow(base, half);
let a_sq_expr = if a_sq_is_one {
arena.one
} else {
arena.num_ratio(a_squared.clone())
};
let x_over_a = if a_sq_is_one {
var
} else {
let neg_half = arena.rational(-1, 2);
let a_inv = arena.pow(a_sq_expr, neg_half);
arena.mul(&[var, a_inv])
};
if has_neg_x2 && c_val.is_positive() {
let x_sqrt = arena.mul(&[var, sqrt_base]);
let asin_term = arena.asin(x_over_a);
let a_sq_asin = if a_sq_is_one {
asin_term
} else {
arena.mul(&[a_sq_expr, asin_term])
};
let sum = arena.add(&[x_sqrt, a_sq_asin]);
return Some(arena.mul(&[half, sum]));
}
if has_pos_x2 && c_val.is_positive() {
let x_sqrt = arena.mul(&[var, sqrt_base]);
let asinh_term = arena.asinh(x_over_a);
let a_sq_asinh = if a_sq_is_one {
asinh_term
} else {
arena.mul(&[a_sq_expr, asinh_term])
};
let sum = arena.add(&[x_sqrt, a_sq_asinh]);
return Some(arena.mul(&[half, sum]));
}
if has_pos_x2 && c_val.is_negative() {
let x_sqrt = arena.mul(&[var, sqrt_base]);
let ln_term = ln_abs_x_plus_sqrt(arena, var, base);
let a_sq_ln = if a_sq_is_one {
ln_term
} else {
arena.mul(&[a_sq_expr, ln_term])
};
let diff = arena.sub(x_sqrt, a_sq_ln);
return Some(arena.mul(&[half, diff]));
}
None
}
fn try_linear_over_quadratic(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
if dependent.len() != 2 {
return None;
}
let (linear_idx, pow_idx) = {
let mut li = None;
let mut pi = None;
for (i, &d) in dependent.iter().enumerate() {
if let ExprNode::Pow(_, _) = arena.node(d) {
if pi.is_none() {
pi = Some(i);
}
} else if li.is_none() {
li = Some(i);
}
}
(li?, pi?)
};
let linear = dependent[linear_idx];
let (pow_base, pow_exp) = match arena.node(dependent[pow_idx]).clone() {
ExprNode::Pow(b, e) => (b, e),
_ => return None,
};
let exp_val = arena.as_num(pow_exp)?.clone();
let neg_one = num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
if exp_val != neg_one {
return None;
}
let lin_poly = crate::poly::polybridge::expr_to_poly(arena, linear, var);
let quad_poly = crate::poly::polybridge::expr_to_poly(arena, pow_base, var);
if let (Some(lp), Some(qp)) = (&lin_poly, &quad_poly) {
tracing::trace!("try_linear_over_quadratic: rational coefficient path");
if lp.degree() == Some(1) && qp.degree() == Some(2) {
let a_coeff = lp.coeff(1);
let b_coeff = lp.coeff(0);
let c_coeff = qp.coeff(2);
let d_coeff = qp.coeff(1);
use num_traits::Zero;
if !c_coeff.is_zero() && !a_coeff.is_zero() {
let two_r = num_rational::Ratio::<num_bigint::BigInt>::from_integer(2.into());
let a_over_2c = &a_coeff / &(&c_coeff * &two_r);
let remainder = &b_coeff - &(&a_coeff * &d_coeff / &(&c_coeff * &two_r));
let mut terms: Vec<ExprId> = Vec::new();
if !a_over_2c.is_zero() {
let coeff_id = arena.num_ratio(a_over_2c.clone());
let abs_quad = arena.abs(pow_base);
let ln_quad = arena.ln(abs_quad);
terms.push(arena.mul(&[coeff_id, ln_quad]));
}
if !remainder.is_zero() {
let inv_quad = arena.pow(pow_base, pow_exp); let inv_integral =
integrate_node(arena, inv_quad, var, var_sym, depth.saturating_sub(1));
if matches!(arena.node(inv_integral), ExprNode::Integral(_, _)) {
return None;
}
let rem_id = arena.num_ratio(remainder.clone());
terms.push(arena.mul(&[rem_id, inv_integral]));
}
return match terms.len() {
0 => Some(arena.zero),
1 => Some(terms[0]),
_ => Some(arena.add(&terms)),
};
}
}
return None;
}
tracing::debug!(
"try_linear_over_quadratic: expr_to_poly failed, trying symbolic coefficient path"
);
let (a_id, b_id) = symbolic_linear_coeff_of(arena, linear, var, var_sym)?;
let (c_id, d_id, _e_id) = symbolic_quadratic_coeffs(arena, pow_base, var, var_sym)?;
tracing::debug!(
"try_linear_over_quadratic: symbolic coefficients extracted for linear/quadratic decomposition"
);
let two = arena.int(2);
let two_c = arena.mul(&[two, c_id]);
let log_coeff = arena.div(a_id, two_c);
let log_coeff = crate::transforms::eval::eval(arena, log_coeff);
let a_d = arena.mul(&[a_id, d_id]);
let a_d_over_2c = arena.div(a_d, two_c);
let remainder_expr = arena.sub(b_id, a_d_over_2c);
let remainder_expr = crate::transforms::eval::eval(arena, remainder_expr);
let mut terms: Vec<ExprId> = Vec::new();
let known_zero = |arena: &mut Arena, c: ExprId| {
arena.is_zero_structural(c)
|| crate::transforms::evalf::eval_const_f64(arena, c).is_some_and(|v| v.abs() <= 1e-14)
};
if !known_zero(arena, log_coeff) {
let abs_quad = arena.abs(pow_base);
let ln_quad = arena.ln(abs_quad);
terms.push(arena.mul(&[log_coeff, ln_quad]));
}
if !known_zero(arena, remainder_expr) {
let inv_quad = arena.pow(pow_base, pow_exp); let inv_integral = integrate_node(arena, inv_quad, var, var_sym, depth.saturating_sub(1));
if crate::base::walk::has_unevaluated(arena, inv_integral) {
return None;
}
terms.push(arena.mul(&[remainder_expr, inv_integral]));
}
match terms.len() {
0 => Some(arena.zero),
1 => Some(terms[0]),
_ => Some(arena.add(&terms)),
}
}
fn is_rational_trig(arena: &Arena, expr: ExprId, var: ExprId, var_sym: SymbolId) -> bool {
if !contains_var(arena, expr, var_sym) {
return true; }
match arena.node(expr).clone() {
ExprNode::Sin(inner) if inner == var => true,
ExprNode::Cos(inner) if inner == var => true,
ExprNode::Add(children) => children
.iter()
.all(|&c| is_rational_trig(arena, c, var, var_sym)),
ExprNode::Mul(children) => children
.iter()
.all(|&c| is_rational_trig(arena, c, var, var_sym)),
ExprNode::Neg(inner) => is_rational_trig(arena, inner, var, var_sym),
ExprNode::Pow(base, exp) => {
if !contains_var(arena, exp, var_sym)
&& let Some(e) = arena.as_num(exp)
&& e.is_integer()
{
return is_rational_trig(arena, base, var, var_sym);
}
false
}
_ => false,
}
}
fn try_weierstrass_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
if depth < 3 {
return None;
}
if !is_rational_trig(arena, expr, var, var_sym) {
return None;
}
let _jumps = crate::transforms::trig_integ::RealParameterJumps::enable(false);
if has_only_real_parameters(arena, expr, var_sym)
&& let Some(result) = try_bioche_tan_substitution(arena, expr, var, var_sym, depth, false)
{
return Some(result);
}
tracing::debug!("trying Weierstrass substitution");
let t = arena.symbol("__wt");
let ExprNode::Symbol(t_sym) = *arena.node(t) else {
return None;
};
let two = arena.int(2);
let one = arena.one;
let t_sq = arena.pow(t, two);
let one_plus_t_sq = arena.add(&[one, t_sq]);
let sin_var = arena.sin(var);
let cos_var = arena.cos(var);
let two_t = arena.mul(&[two, t]);
let sin_sub = arena.div(two_t, one_plus_t_sq);
let one_minus_t_sq = arena.sub(one, t_sq);
let cos_sub = arena.div(one_minus_t_sq, one_plus_t_sq);
let dx_factor = arena.div(two, one_plus_t_sq);
let mut sub_expr = arena.subs_structural(expr, sin_var, sin_sub);
sub_expr = arena.subs_structural(sub_expr, cos_var, cos_sub);
let integrand_t = arena.mul(&[sub_expr, dx_factor]);
let integrand_t = crate::transforms::eval::eval(arena, integrand_t);
let integrand_t = crate::transforms::expand::expand(arena, integrand_t);
let integrand_t = crate::transforms::eval::eval(arena, integrand_t);
let integrand_t = arena.cancel_expr(integrand_t, t);
let integrand_t = crate::transforms::eval::eval(arena, integrand_t);
let integrand_t = clear_nested_fractions(arena, integrand_t, t);
tracing::debug!(integrand_t = %arena.display(integrand_t), "weierstrass: integrand in t");
let integral_t = integrate_node(arena, integrand_t, t, t_sym, depth.saturating_sub(2));
if crate::base::walk::has_unevaluated(arena, integral_t) {
return None;
}
if candidate_rejected(arena, integrand_t, integral_t, t, t_sym) {
tracing::debug!("weierstrass: closed form in t failed verification");
return None;
}
let half = arena.rational(1, 2);
let half_var = arena.mul(&[half, var]);
match crate::transforms::trig_integ::continuous_through_tan(
arena, integral_t, t, t_sym, half_var,
) {
Some(result) => Some(result),
None if undecided_jump_kept(arena, integrand_t, integral_t, t, t_sym) => {
tracing::debug!("weierstrass: the jump at x = π is not decided; uncorrected form");
let tan_half = arena.tan(half_var);
let back = arena.subs_structural(integral_t, t, tan_half);
Some(crate::transforms::trig_integ::atan_tan_terms_to_argument(
arena, back, half_var,
))
}
None => {
tracing::debug!("weierstrass: the jump of G(tan(x/2)) at x = π is not decided");
None
}
}
}
fn has_only_real_parameters(arena: &Arena, expr: ExprId, var_sym: SymbolId) -> bool {
let params: Vec<SymbolId> = crate::base::walk::free_symbols(arena, expr)
.into_iter()
.filter_map(|s| match *arena.node(s) {
ExprNode::Symbol(sid) if sid != var_sym => Some(sid),
_ => None,
})
.collect();
!params.is_empty()
&& params
.iter()
.all(|&sid| crate::transforms::realness::symbol_declared_real(arena, sid))
}
fn try_bioche_tan_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
depth: usize,
keep_undecided: bool,
) -> Option<ExprId> {
let t = arena.symbol("__bt");
let c = arena.symbol("__bc");
let ExprNode::Symbol(t_sym) = *arena.node(t) else {
return None;
};
let sin_var = arena.sin(var);
let cos_var = arena.cos(var);
let t_c = arena.mul(&[t, c]);
let substituted = arena.subs_structural(expr, sin_var, t_c);
let substituted = arena.subs_structural(substituted, cos_var, c);
if contains_var(arena, substituted, var_sym) {
return None;
}
let in_c = clear_nested_fractions(arena, substituted, c);
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, in_c);
let num_c = crate::transforms::solve::symbolic_poly_coeffs(arena, num, c)?;
let den_c = crate::transforms::solve::symbolic_poly_coeffs(arena, den, c)?;
let vanishes = |arena: &mut Arena, cs: &[ExprId], parity: usize| {
cs.iter()
.skip(parity)
.step_by(2)
.all(|&k| crate::poly::algebraic::is_zero_checked(arena, k) == Some(true))
};
let shift = if vanishes(arena, &num_c, 1) && vanishes(arena, &den_c, 1) {
0
} else if vanishes(arena, &num_c, 0) && vanishes(arena, &den_c, 0) {
1
} else {
return None;
};
let one = arena.one;
let two = arena.int(2);
let t_sq = arena.pow(t, two);
let one_plus_t_sq = arena.add(&[one, t_sq]);
let in_t = |arena: &mut Arena, cs: &[ExprId]| -> ExprId {
let mut terms: Vec<ExprId> = Vec::new();
for (j, &k) in cs.iter().enumerate().skip(shift).step_by(2) {
let half_degree = i64::try_from((j - shift) / 2).unwrap_or(0);
let e = arena.int(-half_degree);
let c_pow = arena.pow(one_plus_t_sq, e);
terms.push(arena.mul(&[k, c_pow]));
}
arena.add(&terms)
};
let num_t = in_t(arena, &num_c);
let den_t = in_t(arena, &den_c);
let dt = arena.mul(&[den_t, one_plus_t_sq]);
let integrand_t = arena.div(num_t, dt);
let integrand_t = crate::transforms::eval::eval(arena, integrand_t);
let integrand_t = clear_nested_fractions(arena, integrand_t, t);
if contains_var(arena, integrand_t, var_sym)
|| crate::base::walk::contains(arena, integrand_t, c)
{
return None;
}
tracing::debug!(integrand_t = %arena.display(integrand_t), "bioche: integrand in tan x");
let integral_t = integrate_node(arena, integrand_t, t, t_sym, depth.saturating_sub(2));
if crate::base::walk::has_unevaluated(arena, integral_t)
|| candidate_rejected(arena, integrand_t, integral_t, t, t_sym)
{
return None;
}
let result = match crate::transforms::trig_integ::continuous_through_tan(
arena, integral_t, t, t_sym, var,
) {
Some(result) => result,
None if keep_undecided && undecided_jump_kept(arena, integrand_t, integral_t, t, t_sym) => {
let tan_var = arena.tan(var);
let back = arena.subs_structural(integral_t, t, tan_var);
crate::transforms::trig_integ::atan_tan_terms_to_argument(arena, back, var)
}
None => return None,
};
if candidate_rejected(arena, expr, result, var, var_sym) {
tracing::debug!("bioche: the antiderivative in x failed verification");
return None;
}
Some(result)
}
const MAX_ANGLE_MULTIPLE: i64 = 12;
fn multiple_angle(arena: &mut Arena, k: i64, s: ExprId, c: ExprId) -> (ExprId, ExprId) {
let two = arena.int(2);
let two_c = arena.mul(&[two, c]);
let minus_one = arena.int(-1);
let (mut t_prev, mut t_cur) = (arena.one, c);
let (mut u_prev, mut u_cur) = (arena.zero, arena.one);
for _ in 1..k {
let t_lead = arena.mul(&[two_c, t_cur]);
let t_back = arena.mul(&[minus_one, t_prev]);
let t_next = arena.add(&[t_lead, t_back]);
let t_next = crate::transforms::expand::expand(arena, t_next);
let u_lead = arena.mul(&[two_c, u_cur]);
let u_back = arena.mul(&[minus_one, u_prev]);
let u_next = arena.add(&[u_lead, u_back]);
let u_next = crate::transforms::expand::expand(arena, u_next);
(t_prev, t_cur) = (t_cur, t_next);
(u_prev, u_cur) = (u_cur, u_next);
}
let sin_k = arena.mul(&[s, u_cur]);
let sin_k = crate::transforms::expand::expand(arena, sin_k);
(sin_k, t_cur)
}
fn parity_in(arena: &mut Arena, p: ExprId, v: ExprId) -> Option<usize> {
let coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, v)?;
let mut nonzero = [false, false];
for (j, &k) in coeffs.iter().enumerate() {
let zero = arena.is_zero_structural(k)
|| (crate::base::walk::free_symbols(arena, k).is_empty()
&& crate::poly::algebraic::is_zero_checked(arena, k) == Some(true));
if !zero {
nonzero[j % 2] = true;
}
}
match nonzero {
[true, false] => Some(0),
[false, true] => Some(1),
_ => None,
}
}
fn parameter_content(arena: &mut Arena, p: ExprId, vars: &[ExprId]) -> Option<(ExprId, ExprId)> {
let mut monomials: Vec<(ExprId, ExprId)> = vec![(p, arena.one)];
for &v in vars {
let mut next: Vec<(ExprId, ExprId)> = Vec::new();
for (q, mono) in monomials {
let coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, q, v)?;
for (j, k) in coeffs.into_iter().enumerate() {
if arena.is_zero_structural(k) {
continue;
}
let e = arena.int(i64::try_from(j).ok()?);
let v_j = arena.pow(v, e);
next.push((k, arena.mul(&[mono, v_j])));
}
}
monomials = next;
}
let &(lead, _) = monomials.first()?;
if crate::base::walk::free_symbols(arena, lead).is_empty() {
return None;
}
let mut terms: Vec<ExprId> = Vec::with_capacity(monomials.len());
for (k, mono) in monomials {
let ratio = arena.div(k, lead);
let ratio = crate::transforms::eval::eval(arena, ratio);
if !crate::base::walk::free_symbols(arena, ratio).is_empty() {
return None;
}
terms.push(arena.mul(&[ratio, mono]));
}
Some((lead, arena.add(&terms)))
}
fn odd_part_in_t(
arena: &mut Arena,
fraction: (ExprId, ExprId),
parities: (usize, usize),
w: ExprId,
o: ExprId,
t: ExprId,
) -> Option<ExprId> {
let one = arena.one;
let two = arena.int(2);
let t_sq = arena.pow(t, two);
let minus_t_sq = arena.neg(t_sq);
let w_sq = arena.add(&[one, minus_t_sq]);
let reduce = |arena: &mut Arena, p: ExprId, parity: usize| -> Option<ExprId> {
let coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, w)?;
let mut terms: Vec<ExprId> = Vec::new();
for (j, &k) in coeffs.iter().enumerate().skip(parity).step_by(2) {
let half = arena.int(i64::try_from((j - parity) / 2).ok()?);
let k_t = arena.subs_structural(k, o, t);
let power = arena.pow(w_sq, half);
terms.push(arena.mul(&[k_t, power]));
}
Some(arena.add(&terms))
};
let n = reduce(arena, fraction.0, parities.0)?;
let d = reduce(arena, fraction.1, parities.1)?;
let d = if parities.0 == 1 {
d
} else {
arena.mul(&[w_sq, d])
};
Some(arena.div(n, d))
}
fn is_rational_in_all(arena: &Arena, e: ExprId, vars: &[ExprId]) -> bool {
vars.iter().all(|&v| match *arena.node(v) {
ExprNode::Symbol(sid) => is_rational_in(arena, e, sid),
_ => false,
})
}
fn try_trig_rational_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
use num_integer::Integer;
let mut trig: Vec<(ExprId, ExprId, ExprId)> = Vec::new();
for id in crate::base::walk::post_order_ids(arena, expr) {
let arg = match *arena.node(id) {
ExprNode::Sin(a) | ExprNode::Cos(a) | ExprNode::Tan(a) => a,
_ => continue,
};
if !contains_var(arena, arg, var_sym) {
continue;
}
let (alpha, beta) = symbolic_linear_coeff_of(arena, arg, var, var_sym)?;
trig.push((id, alpha, beta));
}
let &(_, alpha0, beta0) = trig.first()?;
let mut ratios: Vec<Q> = Vec::with_capacity(trig.len());
for &(_, alpha, _) in &trig {
let r = arena.div(alpha, alpha0);
let r = crate::transforms::eval::eval(arena, r);
ratios.push(arena.as_num(r)?.clone());
}
let mut num_gcd = num_bigint::BigInt::zero();
let mut den_lcm = num_bigint::BigInt::one();
for r in &ratios {
num_gcd = num_gcd.gcd(r.numer());
den_lcm = den_lcm.lcm(r.denom());
}
if num_gcd.is_zero() {
return None;
}
let g = Q::new(num_gcd, den_lcm);
let u = arena.symbol("__tru");
let s = arena.symbol("__trs");
let c = arena.symbol("__trc");
let ExprNode::Symbol(u_sym) = *arena.node(u) else {
return None;
};
let mut sub = expr;
for (&(node, _, beta), r) in trig.iter().zip(&ratios) {
let k = r / &g;
if !k.is_integer() {
return None;
}
let k = i64::try_from(k.to_integer()).ok()?;
if k == 0 || k.abs() > MAX_ANGLE_MULTIPLE {
return None;
}
let r_id = arena.num_ratio(r.clone());
let r_beta0 = arena.mul(&[r_id, beta0]);
let gamma = arena.sub(beta, r_beta0);
let gamma = crate::transforms::eval::eval(arena, gamma);
let (sin_k, cos_k) = multiple_angle(arena, k.abs(), s, c);
let sin_k = if k < 0 { arena.neg(sin_k) } else { sin_k };
let (sin_l, cos_l) = if arena.is_zero_structural(gamma) {
(sin_k, cos_k)
} else {
let sin_g = arena.sin(gamma);
let cos_g = arena.cos(gamma);
let a = arena.mul(&[sin_k, cos_g]);
let b = arena.mul(&[cos_k, sin_g]);
let sin_l = arena.add(&[a, b]);
let a = arena.mul(&[cos_k, cos_g]);
let b = arena.mul(&[sin_k, sin_g]);
let cos_l = arena.sub(a, b);
(sin_l, cos_l)
};
let repl = match arena.node(node) {
ExprNode::Sin(_) => sin_l,
ExprNode::Cos(_) => cos_l,
_ => arena.div(sin_l, cos_l),
};
sub = arena.subs_structural(sub, node, repl);
}
let rate = arena.num_ratio(g.clone());
let rate = arena.mul(&[rate, alpha0]);
let shift = arena.num_ratio(g);
let shift = arena.mul(&[shift, beta0]);
let shift = crate::transforms::eval::eval(arena, shift);
if contains_var(arena, sub, var_sym) {
return trig_polynomial_substitution(arena, expr, sub, (s, c), (rate, shift), var, var_sym);
}
let undeclared = |arena: &Arena, e: ExprId| {
crate::base::walk::free_symbols(arena, e)
.into_iter()
.any(|v| match *arena.node(v) {
ExprNode::Symbol(sid) => {
!crate::transforms::realness::symbol_declared_real(arena, sid)
}
_ => true,
})
};
let plain_logs = undeclared(arena, rate) || undeclared(arena, shift);
let _jumps = crate::transforms::trig_integ::RealParameterJumps::enable(plain_logs);
let _strict = crate::transforms::trig_integ::StrictJumps::enable();
let (kappa, result_u) = if is_rational_in_all(arena, sub, &[s, c]) {
tracing::debug!(
"integrate: a rational function of the trigonometric functions of a linear argument"
);
let r = clear_nested_fractions(arena, sub, s);
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, r);
let num = crate::transforms::expand::expand(arena, num);
let den = crate::transforms::expand::expand(arena, den);
let (num_k, num) = parameter_content(arena, num, &[s, c]).unwrap_or((arena.one, num));
let (den_k, den) = parameter_content(arena, den, &[s, c]).unwrap_or((arena.one, den));
let kappa = arena.div(num_k, den_k);
(
kappa,
trig_rational_in_u(arena, (num, den), (s, c), u, u_sym)?,
)
} else {
let result = match odd_algebraic_in_u(arena, sub, (s, c), u) {
Some(r) => r,
None => match tan_algebraic_in_u(arena, sub, (s, c), u) {
Some(r) => r,
None => half_angle_radical_in_u(arena, sub, (s, c), u)?,
},
};
(arena.one, result)
};
let l = arena.mul(&[rate, var]);
let l = arena.add(&[l, shift]);
let back = arena.subs_structural(result_u, u, l);
let minus_one = arena.int(-1);
let inv_rate = arena.pow(rate, minus_one);
let back = arena.mul(&[kappa, inv_rate, back]);
let back = crate::transforms::eval::eval(arena, back);
if candidate_rejected(arena, expr, back, var, var_sym) {
tracing::debug!(
"integrate: rejecting an unverified trigonometric-substitution closed form"
);
return None;
}
Some(back)
}
fn parameter_terms(arena: &mut Arena, p: ExprId, vars: &[ExprId]) -> Vec<(ExprId, ExprId)> {
let p = crate::transforms::expand::expand(arena, p);
let terms: Vec<ExprId> = match arena.node(p) {
ExprNode::Add(children) => children.to_vec(),
_ => vec![p],
};
let mut groups: Vec<(ExprId, Vec<ExprId>)> = Vec::new();
for term in terms {
let factors: Vec<ExprId> = match arena.node(term) {
ExprNode::Mul(children) => children.to_vec(),
_ => vec![term],
};
let mut key: Vec<ExprId> = Vec::new();
let mut part: Vec<ExprId> = Vec::new();
for f in factors {
let dependent = vars.iter().any(|&v| contains_var_id(arena, f, v));
if !dependent && !crate::base::walk::free_symbols(arena, f).is_empty() {
key.push(f);
} else {
part.push(f);
}
}
let key = arena.mul(&key);
let part = arena.mul(&part);
match groups.iter_mut().find(|(k, _)| *k == key) {
Some((_, parts)) => parts.push(part),
None => groups.push((key, vec![part])),
}
}
groups
.into_iter()
.map(|(k, parts)| (k, arena.add(&parts)))
.collect()
}
fn trig_polynomial_substitution(
arena: &mut Arena,
expr: ExprId,
sub: ExprId,
sc: (ExprId, ExprId),
linear: (ExprId, ExprId),
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let (s, c) = sc;
let (rate, shift) = linear;
if !is_rational_in_all(arena, sub, &[s, c])
|| crate::transforms::solve::symbolic_poly_coeffs(arena, sub, var).is_none()
{
return None;
}
let u = arena.symbol("__trp");
let minus_shift = arena.neg(shift);
let u_minus = arena.add(&[u, minus_shift]);
let x_of_u = arena.div(u_minus, rate);
let in_u = arena.subs_structural(sub, var, x_of_u);
let sin_u = arena.sin(u);
let cos_u = arena.cos(u);
let in_u = arena.subs_structural(in_u, s, sin_u);
let in_u = arena.subs_structural(in_u, c, cos_u);
let in_u = arena.div(in_u, rate);
let in_u = crate::transforms::eval::eval(arena, in_u);
let same = arena.subs_structural(in_u, u, var);
if same == expr {
return None;
}
tracing::debug!(integrand_u = %arena.display(in_u), "integrate: polynomial times trigonometric functions of a linear argument");
let res = integrate_nested_generic(arena, in_u, u)?;
let l = arena.mul(&[rate, var]);
let l = arena.add(&[l, shift]);
let back = arena.subs_structural(res, u, l);
let back = crate::transforms::eval::eval(arena, back);
if candidate_rejected(arena, expr, back, var, var_sym) {
return None;
}
Some(back)
}
fn trig_rational_in_u(
arena: &mut Arena,
fraction: (ExprId, ExprId),
sc: (ExprId, ExprId),
u: ExprId,
u_sym: SymbolId,
) -> Option<ExprId> {
let (num, den) = fraction;
let (s, c) = sc;
let parametric = |arena: &Arena, e: ExprId| {
crate::base::walk::free_symbols(arena, e)
.into_iter()
.any(|v| v != s && v != c)
};
let groups = if !parametric(arena, den) && parametric(arena, num) {
parameter_terms(arena, num, &[s, c])
} else {
Vec::new()
};
if groups.len() < 2 {
return trig_rational_single(arena, fraction, sc, u, u_sym);
}
let mut sum: Vec<ExprId> = Vec::with_capacity(groups.len());
for (key, part) in groups {
let r = trig_rational_single(arena, (part, den), sc, u, u_sym)?;
sum.push(arena.mul(&[key, r]));
}
Some(arena.add(&sum))
}
fn trig_rational_single(
arena: &mut Arena,
fraction: (ExprId, ExprId),
sc: (ExprId, ExprId),
u: ExprId,
u_sym: SymbolId,
) -> Option<ExprId> {
let (num, den) = fraction;
let (s, c) = sc;
let ratio = arena.div(num, den);
if !contains_var_id(arena, ratio, s) && !contains_var_id(arena, ratio, c) {
return Some(arena.mul(&[ratio, u]));
}
let t = arena.symbol("__trt");
let ExprNode::Symbol(t_sym) = *arena.node(t) else {
return None;
};
let sin_u = arena.sin(u);
let cos_u = arena.cos(u);
for (w, o, back, sign) in [(s, c, cos_u, -1), (c, s, sin_u, 1)] {
let (Some(pn), Some(pd)) = (parity_in(arena, num, w), parity_in(arena, den, w)) else {
continue;
};
if pn == pd {
continue;
}
let Some(in_t) = odd_part_in_t(arena, (num, den), (pn, pd), w, o, t) else {
continue;
};
let sign = arena.int(sign);
let in_t = arena.mul(&[sign, in_t]);
let in_t = clear_nested_fractions(arena, in_t, t);
tracing::debug!(integrand_t = %arena.display(in_t), "bioche: integrand in sin u or cos u");
let g = integrate_node(arena, in_t, t, t_sym, 18);
if crate::base::walk::has_unevaluated(arena, g)
|| candidate_rejected(arena, in_t, g, t, t_sym)
{
continue;
}
return Some(arena.subs_structural(g, t, back));
}
let r_s = arena.subs_structural(ratio, s, sin_u);
let r_u = arena.subs_structural(r_s, c, cos_u);
let r_u = crate::transforms::eval::eval(arena, r_u);
let lam = arena.symbol("__trl");
let joint = |arena: &mut Arena, p: ExprId| -> Option<usize> {
let ls = arena.mul(&[lam, s]);
let lc = arena.mul(&[lam, c]);
let p = arena.subs_structural(p, s, ls);
let p = arena.subs_structural(p, c, lc);
parity_in(arena, p, lam)
};
if let (Some(pn), Some(pd)) = (joint(arena, num), joint(arena, den))
&& pn == pd
&& let Some(result) = try_bioche_tan_substitution(arena, r_u, u, u_sym, 20, true)
{
return Some(result);
}
try_weierstrass_substitution(arena, r_u, u, u_sym, 20)
}
fn try_gaussian_rational(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let i = arena.i_unit();
if depth < 2
|| !crate::base::walk::contains(arena, expr, i)
|| !is_rational_in(arena, expr, var_sym)
{
return None;
}
let iota = arena.symbol("__gi");
let e = arena.subs_structural(expr, i, iota);
if crate::base::walk::contains(arena, e, i) {
return None;
}
let e = clear_nested_fractions(arena, e, var);
let (n, d) = crate::poly::polybridge::as_numer_denom(arena, e);
let split = |arena: &mut Arena, p: ExprId| -> Option<(ExprId, ExprId)> {
let coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, iota)?;
let mut re: Vec<ExprId> = Vec::new();
let mut im: Vec<ExprId> = Vec::new();
for (k, &c) in coeffs.iter().enumerate() {
let c = if (k / 2) % 2 == 0 { c } else { arena.neg(c) };
if k % 2 == 0 {
re.push(c);
} else {
im.push(c);
}
}
let re = arena.add(&re);
let im = arena.add(&im);
Some((
crate::transforms::eval::eval(arena, re),
crate::transforms::eval::eval(arena, im),
))
};
let (n0, n1) = split(arena, n)?;
let (d0, d1) = split(arena, d)?;
let (p, q, den) = if arena.is_zero_structural(d1) {
(n0, n1, d0)
} else {
let a = arena.mul(&[n0, d0]);
let b = arena.mul(&[n1, d1]);
let p = arena.add(&[a, b]);
let a = arena.mul(&[n1, d0]);
let b = arena.mul(&[n0, d1]);
let q = arena.sub(a, b);
let two = arena.int(2);
let a = arena.pow(d0, two);
let b = arena.pow(d1, two);
(p, q, arena.add(&[a, b]))
};
let mut parts: Vec<ExprId> = Vec::with_capacity(2);
for (num, factor) in [(p, arena.one), (q, i)] {
let num = crate::transforms::expand::expand(arena, num);
if arena.is_zero_structural(num) {
continue;
}
let f = arena.div(num, den);
let f = clear_nested_fractions(arena, f, var);
if crate::base::walk::contains(arena, f, iota) {
return None;
}
let big_f = integrate_node(arena, f, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, big_f) {
return None;
}
let terms: Vec<ExprId> = match arena.node(big_f) {
ExprNode::Add(children) => children.to_vec(),
_ => vec![big_f],
};
for term in terms {
parts.push(arena.mul(&[factor, term]));
}
}
tracing::debug!(
"integrate: a rational function with Gaussian coefficients, by real and imaginary parts"
);
let result = arena.add(&parts);
if candidate_rejected(arena, expr, result, var, var_sym) {
return None;
}
Some(result)
}
fn half_angle_radical_in_u(
arena: &mut Arena,
e: ExprId,
sc: (ExprId, ExprId),
u: ExprId,
) -> Option<ExprId> {
let (s, c) = sc;
let mut base_kind: Option<(ExprId, ExprId)> = None;
let mut radicals: Vec<(ExprId, ExprId, i64, ExprId)> = Vec::new();
for id in crate::base::walk::post_order_ids(arena, e) {
let ExprNode::Pow(base, ex) = *arena.node(id) else {
continue;
};
let Some(r) = arena.as_num(ex).cloned() else {
continue;
};
if r.is_integer() || !(contains_var_id(arena, base, s) || contains_var_id(arena, base, c)) {
continue;
}
if *r.denom() != num_bigint::BigInt::from(2) {
return None;
}
let k = i64::try_from(r.numer()).ok()?;
let (w, other) = if contains_var_id(arena, base, c) {
(c, s)
} else {
(s, c)
};
if contains_var_id(arena, base, other) {
return None;
}
let coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, base, w)?;
let [kappa, b1] = coeffs.as_slice() else {
return None;
};
let eps = arena.div(*b1, *kappa);
let eps = crate::transforms::eval::eval(arena, eps);
let eps_one = arena.as_num(eps).filter(|q| q.abs().is_one())?.clone();
let eps = arena.num_ratio(eps_one);
match base_kind {
None => base_kind = Some((w, eps)),
Some(kind) if kind == (w, eps) => {}
Some(_) => return None,
}
radicals.push((id, *kappa, k, w));
}
let (w, eps) = base_kind?;
let eps_positive = arena.as_num(eps).is_some_and(|q| q.is_positive());
let rad = arena.symbol("__hrad");
let mut sub = e;
for &(id, kappa, k, _) in &radicals {
let half_k = arena.rational(k, 2);
let kappa_pow = arena.pow(kappa, half_k);
let k_id = arena.int(k);
let rad_pow = arena.pow(rad, k_id);
let repl = arena.mul(&[kappa_pow, rad_pow]);
sub = arena.subs_structural(sub, id, repl);
}
let ExprNode::Symbol(rad_sym) = *arena.node(rad) else {
return None;
};
if !is_rational_in_all(arena, sub, &[s, c, rad]) || !is_rational_in(arena, sub, rad_sym) {
return None;
}
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, sub);
let one = arena.one;
let eps_w = arena.mul(&[eps, w]);
let rad_sq = arena.add(&[one, eps_w]);
let reduce_rad = |arena: &mut Arena, p: ExprId| -> Option<(ExprId, ExprId)> {
let cs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, rad)?;
let mut even: Vec<ExprId> = Vec::new();
let mut odd: Vec<ExprId> = Vec::new();
for (j, &k) in cs.iter().enumerate() {
let e = arena.int(i64::try_from(j / 2).ok()?);
let power = arena.pow(rad_sq, e);
let term = arena.mul(&[k, power]);
if j % 2 == 0 {
even.push(term);
} else {
odd.push(term);
}
}
let even = arena.add(&even);
let odd = arena.add(&odd);
Some((
crate::transforms::expand::expand(arena, even),
crate::transforms::expand::expand(arena, odd),
))
};
let (n_even, n_odd) = reduce_rad(arena, num)?;
let (d_even, d_odd) = reduce_rad(arena, den)?;
let e1 = if arena.is_zero_structural(n_even) && arena.is_zero_structural(d_odd) {
arena.div(n_odd, d_even)
} else if arena.is_zero_structural(n_odd) && arena.is_zero_structural(d_even) {
let d = arena.mul(&[rad_sq, d_odd]);
arena.div(n_even, d)
} else {
return None;
};
let big_s = arena.symbol("__hs");
let big_c = arena.symbol("__hc");
let two = arena.int(2);
let two_sc = arena.mul(&[two, big_s, big_c]);
let s_sq = arena.pow(big_s, two);
let c_sq = arena.pow(big_c, two);
let c2_minus_s2 = arena.sub(c_sq, s_sq);
let neg = |arena: &mut Arena, z: ExprId| arena.neg(z);
let w_is_sin = w == s;
let (sin_u, cos_u, phi) = match (w_is_sin, eps_positive) {
(true, true) => (neg(arena, c2_minus_s2), two_sc, arena.rational(1, 4)),
(true, false) => (c2_minus_s2, neg(arena, two_sc), arena.rational(-1, 4)),
(false, true) => (
neg(arena, two_sc),
neg(arena, c2_minus_s2),
arena.rational(1, 2),
),
(false, false) => (two_sc, c2_minus_s2, arena.zero),
};
let pi = arena.pi();
let phi = arena.mul(&[phi, pi]);
let g = arena.subs_structural(e1, s, sin_u);
let g = arena.subs_structural(g, c, cos_u);
let sqrt2 = arena.sqrt(two);
let g = arena.mul(&[two, sqrt2, big_s, g]);
let g = clear_nested_fractions(arena, g, big_s);
let (g_num, g_den) = crate::poly::polybridge::as_numer_denom(arena, g);
let reduce_c = |arena: &mut Arena, p: ExprId| -> Option<(ExprId, ExprId)> {
let cs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, big_c)?;
let one = arena.one;
let s_sq = arena.pow(big_s, two);
let minus = arena.neg(s_sq);
let one_minus = arena.add(&[one, minus]);
let mut even: Vec<ExprId> = Vec::new();
let mut odd: Vec<ExprId> = Vec::new();
for (j, &k) in cs.iter().enumerate() {
let e = arena.int(i64::try_from(j / 2).ok()?);
let power = arena.pow(one_minus, e);
let term = arena.mul(&[k, power]);
if j % 2 == 0 {
even.push(term);
} else {
odd.push(term);
}
}
let even = arena.add(&even);
let odd = arena.add(&odd);
Some((
crate::transforms::expand::expand(arena, even),
crate::transforms::expand::expand(arena, odd),
))
};
let (a_part, b_part) = reduce_c(arena, g_num)?;
let (den_s, den_c) = reduce_c(arena, g_den)?;
if !arena.is_zero_structural(den_c) {
return None;
}
let den_coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, den_s, big_s)?;
let nonzero: Vec<usize> = (0..den_coeffs.len())
.filter(|&j| !arena.is_zero_structural(den_coeffs[j]))
.collect();
let [m] = nonzero.as_slice() else {
return None;
};
let mut m = *m;
let lowest = |arena: &mut Arena, p: ExprId| -> Option<usize> {
if arena.is_zero_structural(p) {
return Some(usize::MAX);
}
let cs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, big_s)?;
cs.iter().position(|&k| !arena.is_zero_structural(k))
};
let k = lowest(arena, a_part)?.min(lowest(arena, b_part)?).min(m);
let (a_part, b_part, den_s) = if k > 0 {
let shift = |arena: &mut Arena, p: ExprId| -> Option<ExprId> {
if arena.is_zero_structural(p) {
return Some(p);
}
let cs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, big_s)?;
let mut terms: Vec<ExprId> = Vec::new();
for (j, &c) in cs.iter().enumerate().skip(k) {
let e = arena.int(i64::try_from(j - k).ok()?);
let power = arena.pow(big_s, e);
terms.push(arena.mul(&[c, power]));
}
Some(arena.add(&terms))
};
let a = shift(arena, a_part)?;
let b = shift(arena, b_part)?;
let d = shift(arena, den_s)?;
m -= k;
(a, b, d)
} else {
(a_part, b_part, den_s)
};
let v = arena.symbol("__hv");
let ExprNode::Symbol(v_sym) = *arena.node(v) else {
return None;
};
let b_times_c = arena.mul(&[big_c, b_part]);
let half = arena.rational(1, 2);
let half_u = arena.mul(&[half, u]);
let v_of_u = arena.add(&[half_u, phi]);
if m > 0 {
let n = arena.add(&[a_part, b_times_c]);
let big_g = trig_rational_in_u(arena, (n, den_s), (big_s, big_c), v, v_sym)?;
let back = arena.subs_structural(big_g, v, v_of_u);
let sin_v = arena.sin(v_of_u);
let sqrt2 = arena.sqrt(two);
let w_u = if w_is_sin { arena.sin(u) } else { arena.cos(u) };
let eps_w_u = arena.mul(&[eps, w_u]);
let one_eps_w = arena.add(&[one, eps_w_u]);
let rad_u = arena.sqrt(one_eps_w);
let minus_one = arena.int(-1);
let inv_rad = arena.pow(rad_u, minus_one);
tracing::debug!("integrate: a radical of 1 ± sin u (cos u) through the half angle");
let big_f = arena.mul(&[sqrt2, sin_v, inv_rad, back]);
return Some(crate::transforms::eval::eval(arena, big_f));
}
let odd_or_zero = |arena: &mut Arena, p: ExprId, parity: usize| {
arena.is_zero_structural(p) || parity_in(arena, p, big_s) == Some(parity)
};
if !odd_or_zero(arena, a_part, 1) || !odd_or_zero(arena, b_part, 0) {
return None;
}
let sin_u_node = arena.sin(u);
let cos_u_node = arena.cos(u);
let w_u = if w_is_sin { sin_u_node } else { cos_u_node };
let eps_w_u = arena.mul(&[eps, w_u]);
let one_eps_w = arena.add(&[one, eps_w_u]);
let s2_u = arena.mul(&[half, one_eps_w]);
let minus_s2 = arena.neg(s2_u);
let c2_u = arena.add(&[one, minus_s2]);
let sin_2v = match (w_is_sin, eps_positive) {
(true, true) => cos_u_node,
(true, false) => arena.neg(cos_u_node),
(false, true) => arena.neg(sin_u_node),
(false, false) => sin_u_node,
};
let mut over_s: Vec<ExprId> = Vec::new();
let mut over_c: Vec<ExprId> = Vec::new();
let mut at_zero: Vec<ExprId> = Vec::new();
let sin_v = arena.sin(v);
let cos_v = arena.cos(v);
for (part, fv, sym, sq, out) in [
(b_times_c, sin_v, big_s, s2_u, &mut over_s),
(a_part, cos_v, big_c, c2_u, &mut over_c),
] {
if arena.is_zero_structural(part) {
continue;
}
let big_g = trig_rational_in_u(arena, (part, den_s), (big_s, big_c), v, v_sym)?;
let p = arena.subs_structural(big_g, fv, sym);
let p = crate::transforms::eval::eval(arena, p);
if contains_var_id(arena, p, v) || parity_in(arena, p, sym) != Some(1) {
return None;
}
let coeffs = crate::transforms::solve::symbolic_poly_coeffs(arena, p, sym)?;
for (j, &k) in coeffs.iter().enumerate().skip(1).step_by(2) {
let e = arena.int(i64::try_from(j / 2).ok()?);
let power = arena.pow(sq, e);
out.push(arena.mul(&[k, power]));
if sym == big_c {
at_zero.push(k);
}
}
}
let rad_u = arena.sqrt(one_eps_w);
let inv_sqrt2 = {
let r = arena.sqrt(two);
let minus_one = arena.int(-1);
arena.pow(r, minus_one)
};
let p_s = arena.add(&over_s);
let p_c = arena.add(&over_c);
let ratio = arena.div(sin_2v, one_eps_w);
let c_term = arena.mul(&[ratio, p_c]);
let bracket = arena.add(&[p_s, c_term]);
let main = arena.mul(&[rad_u, inv_sqrt2, bracket]);
let g0 = arena.add(&at_zero);
let g0 = crate::transforms::eval::eval(arena, g0);
let minus_two = arena.int(-2);
let minus_one = arena.int(-1);
let inv_pi = arena.pow(pi, minus_one);
let v_over_pi = arena.mul(&[v_of_u, inv_pi]);
let fl = arena.floor(v_over_pi);
let correction = arena.mul(&[minus_two, g0, fl]);
tracing::debug!("integrate: a radical of 1 ± sin u (cos u) through the half angle");
let big_f = arena.add(&[main, correction]);
Some(crate::transforms::eval::eval(arena, big_f))
}
fn odd_algebraic_in_u(
arena: &mut Arena,
e: ExprId,
sc: (ExprId, ExprId),
u: ExprId,
) -> Option<ExprId> {
let (s, c) = sc;
let t = arena.symbol("__toa");
let radical_bases: Vec<ExprId> = crate::base::walk::post_order_ids(arena, e)
.into_iter()
.filter_map(|id| match *arena.node(id) {
ExprNode::Pow(base, ex) if !arena.as_num(ex).is_some_and(|r| r.is_integer()) => {
Some(base)
}
_ => None,
})
.collect();
let sin_u = arena.sin(u);
let cos_u = arena.cos(u);
for (w, o, back, sign) in [(c, s, sin_u, 1), (s, c, cos_u, -1)] {
let ExprNode::Symbol(w_sym) = *arena.node(w) else {
continue;
};
if radical_bases.iter().any(|&b| contains_var_id(arena, b, w))
|| !is_rational_in(arena, e, w_sym)
{
continue;
}
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, e);
let (Some(pn), Some(pd)) = (parity_in(arena, num, w), parity_in(arena, den, w)) else {
continue;
};
if pn == pd {
continue;
}
let Some(in_t) = odd_part_in_t(arena, (num, den), (pn, pd), w, o, t) else {
continue;
};
let sign = arena.int(sign);
let in_t = arena.mul(&[sign, in_t]);
let in_t = crate::transforms::eval::eval(arena, in_t);
if contains_var_id(arena, in_t, s) || contains_var_id(arena, in_t, c) {
continue;
}
tracing::debug!(integrand_t = %arena.display(in_t), "integrate: an algebraic function through t = sin u or cos u");
let Some(g) = integrate_nested_generic(arena, in_t, t) else {
continue;
};
return Some(arena.subs_structural(g, t, back));
}
None
}
fn tan_algebraic_in_u(
arena: &mut Arena,
e: ExprId,
sc: (ExprId, ExprId),
u: ExprId,
) -> Option<ExprId> {
let (s, c) = sc;
let ExprNode::Symbol(c_sym) = *arena.node(c) else {
return None;
};
let t = arena.symbol("__tat");
let ExprNode::Symbol(t_sym) = *arena.node(t) else {
return None;
};
let tc = arena.mul(&[t, c]);
let e = arena.subs_structural(e, s, tc);
let e = crate::transforms::eval::eval(arena, e);
if contains_var_id(arena, e, s) || !is_rational_in(arena, e, c_sym) {
return None;
}
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, e);
let num_c = crate::transforms::solve::symbolic_poly_coeffs(arena, num, c)?;
let den_c = crate::transforms::solve::symbolic_poly_coeffs(arena, den, c)?;
let parity = |arena: &mut Arena, cs: &[ExprId]| -> Option<usize> {
let mut nonzero = [false, false];
for (j, &k) in cs.iter().enumerate() {
if !arena.is_zero_structural(k) {
nonzero[j % 2] = true;
}
}
match nonzero {
[true, false] => Some(0),
[false, true] => Some(1),
_ => None,
}
};
let shift = parity(arena, &num_c)?;
if parity(arena, &den_c)? != shift {
return None;
}
tracing::debug!("integrate: an algebraic function of tan u, through t = tan u");
let one = arena.one;
let two = arena.int(2);
let t_sq = arena.pow(t, two);
let one_plus_t_sq = arena.add(&[one, t_sq]);
let in_t = |arena: &mut Arena, cs: &[ExprId]| -> Option<ExprId> {
let mut terms: Vec<ExprId> = Vec::new();
for (j, &k) in cs.iter().enumerate().skip(shift).step_by(2) {
let e = arena.int(-i64::try_from((j - shift) / 2).ok()?);
let c_pow = arena.pow(one_plus_t_sq, e);
terms.push(arena.mul(&[k, c_pow]));
}
Some(arena.add(&terms))
};
let num_t = in_t(arena, &num_c)?;
let den_t = in_t(arena, &den_c)?;
let dt = arena.mul(&[den_t, one_plus_t_sq]);
let integrand_t = arena.div(num_t, dt);
let integrand_t = crate::transforms::eval::eval(arena, integrand_t);
if contains_var_id(arena, integrand_t, c) {
return None;
}
tracing::debug!(integrand_t = %arena.display(integrand_t), "tan u: algebraic integrand in t");
let g = integrate_nested_generic(arena, integrand_t, t)?;
if let Some(r) = crate::transforms::trig_integ::continuous_through_tan(arena, g, t, t_sym, u) {
return Some(r);
}
if tan_pole_continuous(arena, integrand_t, t) != Some(false) {
tracing::debug!("tan u: the jump at the poles of tan u is not decided");
return None;
}
let tan_u = arena.tan(u);
let back = arena.subs_structural(g, t, tan_u);
Some(crate::transforms::trig_integ::atan_tan_terms_to_argument(
arena, back, u,
))
}
fn tan_pole_continuous(arena: &mut Arena, h: ExprId, t: ExprId) -> Option<bool> {
let one = arena.one;
let two = arena.int(2);
let t_sq = arena.pow(t, two);
let one_plus_t_sq = arena.add(&[one, t_sq]);
let mut k = arena.mul(&[h, one_plus_t_sq]);
let mut params: Vec<(String, ExprId)> = Vec::new();
for s in crate::base::walk::free_symbols(arena, k) {
if s == t {
continue;
}
let ExprNode::Symbol(sid) = *arena.node(s) else {
return None;
};
params.push((arena.symbol_name(sid).to_owned(), s));
}
if params.len() > FTC_PARAMETER_VALUES.len() {
return None;
}
params.sort();
for (j, &(_, s)) in params.iter().enumerate() {
let (p, d) = FTC_PARAMETER_VALUES[j];
let value = arena.rational(p, d);
k = crate::transforms::subs::subs(arena, k, s, value);
}
let k = crate::transforms::eval::eval(arena, k);
let mut values = Vec::with_capacity(2);
for point in [arena.infinity(), arena.neg_infinity()] {
let v = arena.limit_expr(k, t, point).ok()?;
let v = crate::transforms::eval::eval(arena, v);
if crate::base::walk::has_unevaluated(arena, v) || contains_var_id(arena, v, t) {
return None;
}
if contains_non_finite(arena, v) {
return Some(false);
}
let Ok(z) = crate::transforms::evalf::evalf_complex64(arena, v) else {
return None;
};
if !z.re.is_finite() || !z.im.is_finite() {
return Some(false);
}
values.push(z);
}
let gap = (values[0] - values[1]).norm();
Some(gap <= 1e-9 * values[0].norm().max(1.0))
}
fn contains_var_id(arena: &Arena, e: ExprId, v: ExprId) -> bool {
match *arena.node(v) {
ExprNode::Symbol(sid) => contains_var(arena, e, sid),
_ => true,
}
}
fn undecided_jump_kept(
arena: &mut Arena,
h: ExprId,
g: ExprId,
t: ExprId,
t_sym: SymbolId,
) -> bool {
jump_kept_undecided(arena, g, t_sym)
&& (!crate::transforms::trig_integ::strict_jumps_active()
|| tan_pole_continuous(arena, h, t) == Some(false))
}
fn jump_kept_undecided(arena: &Arena, g: ExprId, t_sym: SymbolId) -> bool {
crate::base::walk::free_symbols(arena, g)
.iter()
.any(|&s| !matches!(arena.node(s), ExprNode::Symbol(sid) if *sid == t_sym))
|| crate::base::walk::post_order_ids(arena, g)
.iter()
.any(|&id| matches!(arena.node(id), ExprNode::RootSum(..)))
}
pub(crate) fn clear_nested_fractions(arena: &mut Arena, expr: ExprId, var: ExprId) -> ExprId {
let inside = NODE_BUDGET_DEPTH.with(std::cell::Cell::get) > 0;
let key = (std::ptr::from_ref::<Arena>(arena) as usize, expr, var);
if inside && let Some(found) = CLEARED_MEMO.with(|m| m.borrow().get(&key).copied()) {
return found;
}
let cleared = clear_nested_fractions_uncached(arena, expr, var);
if inside {
CLEARED_MEMO.with(|m| m.borrow_mut().insert(key, cleared));
}
cleared
}
thread_local! {
static CLEARED_MEMO: std::cell::RefCell<FxHashMap<(usize, ExprId, ExprId), ExprId>> =
std::cell::RefCell::new(FxHashMap::default());
}
fn clear_nested_fractions_uncached(arena: &mut Arena, expr: ExprId, var: ExprId) -> ExprId {
let mut current = expr;
for _ in 0..16 {
let order = crate::base::walk::post_order_ids(arena, current);
let mut changed = false;
for id in order {
if id == current {
continue;
}
if let ExprNode::Add(_) = arena.node(id) {
let t = crate::poly::polybridge::together(arena, id);
if t != id {
current = arena.subs_structural(current, id, t);
changed = true;
break;
}
}
}
if !changed {
break;
}
}
let together = crate::poly::polybridge::together(arena, current);
let (n, d) = crate::poly::polybridge::as_numer_denom(arena, together);
let n = crate::transforms::expand::expand(arena, n);
let d = crate::transforms::expand::expand(arena, d);
let n = crate::transforms::eval::eval(arena, n);
let d = crate::transforms::eval::eval(arena, d);
let ratio = arena.div(n, d);
let cancelled = arena.cancel_expr(ratio, var);
crate::transforms::eval::eval(arena, cancelled)
}
fn try_cyclic_ibp(
arena: &mut Arena,
factors: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
if factors.len() != 2 || depth < 2 {
return None;
}
let (u_idx, dv_idx) = {
let r0 = liate_rank(arena, factors[0], var, var_sym);
let r1 = liate_rank(arena, factors[1], var, var_sym);
if r0 <= r1 { (0, 1) } else { (1, 0) }
};
let u1 = factors[u_idx];
let dv1 = factors[dv_idx];
if !contains_var(arena, u1, var_sym) || !contains_var(arena, dv1, var_sym) {
return None;
}
let v1 = integrate_node(arena, dv1, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, v1) {
tracing::trace!("try_cyclic_ibp: v1 has unevaluated nodes, bailing");
return None;
}
let du1 = crate::transforms::diff::diff(arena, u1, var);
let boundary1 = arena.mul(&[u1, v1]);
let v2 = integrate_node(arena, v1, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, v2) {
tracing::trace!("try_cyclic_ibp: v2 has unevaluated nodes, bailing");
return None;
}
let du2 = crate::transforms::diff::diff(arena, du1, var); let boundary2 = arena.mul(&[du1, v2]);
let remaining2 = arena.mul(&[v2, du2]);
let original = arena.mul(&[factors[0], factors[1]]);
let sum = arena.add(&[remaining2, original]);
if sum == arena.zero {
tracing::debug!("cyclic IBP detected (c = -1)");
let numerator = arena.sub(boundary1, boundary2);
let two = arena.int(2);
return Some(arena.div(numerator, two));
}
let diff_check = arena.sub(remaining2, original);
if diff_check == arena.zero {
return None;
}
let ratio = arena.div(remaining2, original);
let cancelled = arena.cancel_expr(ratio, var);
if !contains_var(arena, cancelled, var_sym) && cancelled != arena.one {
tracing::debug!("cyclic IBP detected (general c)");
let numerator = arena.sub(boundary1, boundary2);
let one = arena.one;
let one_minus_c = arena.sub(one, cancelled);
return Some(arena.div(numerator, one_minus_c));
}
None
}
#[derive(Clone, Copy)]
enum MemoOutcome {
Failed { depth: usize },
Solved(ExprId),
}
#[derive(Clone, Copy, PartialEq, Eq, Hash)]
struct MemoKey {
arena: usize,
expr: ExprId,
var: ExprId,
in_risch: bool,
subst_depth: u8,
}
thread_local! {
static INTEGRATE_MEMO: std::cell::RefCell<FxHashMap<MemoKey, MemoOutcome>> =
std::cell::RefCell::new(FxHashMap::default());
}
fn integrate_node(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> ExprId {
let key = MemoKey {
arena: std::ptr::from_ref::<Arena>(arena) as usize,
expr,
var,
in_risch: crate::calculus::risch::rational_guard_active(),
subst_depth: SUBST_DEPTH.with(std::cell::Cell::get),
};
match INTEGRATE_MEMO.with(|m| m.borrow().get(&key).copied()) {
Some(MemoOutcome::Solved(result)) => return result,
Some(MemoOutcome::Failed { depth: failed }) if depth <= failed => {
return arena.intern(ExprNode::Integral(expr, var));
}
_ => {}
}
let result = if depth > 0
&& !matches!(arena.node(expr), ExprNode::Add(_))
&& let Some(r) = try_log_substitution(arena, expr, var, var_sym)
{
r
} else {
integrate_node_uncached(arena, expr, var, var_sym, depth)
};
let result = analytic_log_of_complex_arguments(arena, result, var_sym);
let outcome = if matches!(arena.node(result), ExprNode::Integral(e, v) if *e == expr && *v == var)
{
(!StepCap::active()).then_some(MemoOutcome::Failed { depth })
} else if !crate::base::walk::has_unevaluated(arena, result) {
Some(MemoOutcome::Solved(result))
} else {
None
};
if let Some(outcome) = outcome {
INTEGRATE_MEMO.with(|m| {
let mut m = m.borrow_mut();
let keep_old = matches!(
(m.get(&key), outcome),
(Some(MemoOutcome::Failed { depth: old }), MemoOutcome::Failed { depth: new }) if *old >= new
);
if !keep_old {
m.insert(key, outcome);
}
});
}
result
}
fn integrate_node_uncached(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> ExprId {
tracing::trace!(depth = depth, "integrate_node entered");
if depth == 0 || node_budget_exhausted() || arena_budget_exhausted(arena) {
return arena.intern(ExprNode::Integral(expr, var));
}
if let Some((constants, f)) = split_constant_factors(arena, expr, var_sym) {
let d = if split_constant_factors(arena, f, var_sym).is_some() {
depth - 1
} else {
depth
};
let big_f = integrate_node(arena, f, var, var_sym, d);
if matches!(arena.node(big_f), ExprNode::Integral(..)) {
return arena.intern(ExprNode::Integral(expr, var));
}
return wrap_with_constants(arena, big_f, &constants);
}
if let Some(result) =
crate::transforms::trig_integ::try_trig_power_integral(arena, expr, var, var_sym)
{
if !matches!(arena.node(result), ExprNode::Integral(_, _)) {
return result;
}
}
if let ExprNode::Pow(trig_base, trig_exp) = arena.node(expr).clone()
&& let Some(n_val) = arena.as_num(trig_exp)
{
let two_r = num_rational::Ratio::<num_bigint::BigInt>::from_integer(2.into());
if *n_val == two_r {
if let ExprNode::Tan(inner) = arena.node(trig_base).clone() {
let cos_inner = arena.cos(inner);
let neg_two = arena.int(-2);
let sec_sq = arena.pow(cos_inner, neg_two);
let rewritten = arena.sub(sec_sq, arena.one);
return integrate_node(arena, rewritten, var, var_sym, depth - 1);
}
if let ExprNode::Tanh(inner) = arena.node(trig_base).clone() {
let cosh_inner = arena.cosh(inner);
let neg_two = arena.int(-2);
let sech_sq = arena.pow(cosh_inner, neg_two);
let rewritten = arena.sub(arena.one, sech_sq);
return integrate_node(arena, rewritten, var, var_sym, depth - 1);
}
}
}
let risch_rejected = match stage!(
arena,
"risch_rational",
expr,
crate::calculus::risch::try_risch_rational(arena, expr, var)
) {
Some(result) if !matches!(arena.node(result), ExprNode::Integral(_, _)) => {
if !candidate_rejected(arena, expr, result, var, var_sym) {
return result;
}
tracing::debug!("integrate: rejecting unverified rational-function closed form");
true
}
_ => false,
};
if let Some(result) = stage!(
arena,
"biquadratic",
expr,
try_biquadratic_rational(arena, expr, var, var_sym)
) {
return result;
}
if risch_rejected {
return arena.intern(ExprNode::Integral(expr, var));
}
if let Some(result) = try_poly_over_symbolic_linear(arena, expr, var, var_sym) {
return result;
}
if let Some(result) = stage!(
arena,
"risch_rational_params",
expr,
crate::calculus::risch::try_risch_rational_params(arena, expr, var)
) && !candidate_rejected(arena, expr, result, var, var_sym)
{
return result;
}
if let Some(result) = stage!(
arena,
"gaussian_rational",
expr,
try_gaussian_rational(arena, expr, var, var_sym, depth)
) {
return result;
}
if let Some(result) = try_poly_over_symbolic_quadratic(arena, expr, var, var_sym, depth) {
return result;
}
let node = arena.node(expr).clone();
match node {
_ if !matches!(node, ExprNode::Add(_)) && !contains_var(arena, expr, var_sym) => {
arena.mul(&[expr, var])
}
ExprNode::Num(_)
| ExprNode::Pi
| ExprNode::E
| ExprNode::ImaginaryUnit
| ExprNode::Infinity
| ExprNode::NegInfinity
| ExprNode::ComplexInfinity
| ExprNode::NaN => {
arena.mul(&[expr, var])
}
ExprNode::Symbol(sid) => {
if sid == var_sym {
let two = arena.int(2);
let x_sq = arena.pow(var, two);
let half = arena.rational(1, 2);
arena.mul(&[half, x_sq])
} else {
arena.mul(&[expr, var])
}
}
ExprNode::Add(ref children) => {
let integrals: SmallVec<[ExprId; 6]> = children
.iter()
.map(|&child| integrate_node(arena, child, var, var_sym, depth - 1))
.collect();
arena.add(&integrals)
}
ExprNode::Mul(ref children) => {
let (mut constants, dependent) = partition_factors(arena, children, var_sym);
if dependent.is_empty() {
return arena.mul(&[expr, var]);
}
if let Some(expanded) = expand_polynomial_product(arena, &dependent, var_sym) {
let expanded = wrap_with_constants(arena, expanded, &constants);
let result = integrate_node(arena, expanded, var, var_sym, depth - 1);
if !crate::base::walk::has_unevaluated(arena, result) {
return result;
}
}
if let Some(merged) = merge_exponential_factors(arena, &dependent, var, var_sym) {
let merged = wrap_with_constants(arena, merged, &constants);
let result = integrate_node(arena, merged, var, var_sym, depth - 1);
if !crate::base::walk::has_unevaluated(arena, result) {
return result;
}
}
if !constants.is_empty() && dependent.len() == 1 {
let inner_integral = integrate_node(arena, dependent[0], var, var_sym, depth - 1);
if let ExprNode::Integral(_, _) = arena.node(inner_integral) {
return arena.intern(ExprNode::Integral(expr, var));
}
constants.push(inner_integral);
return arena.mul(&constants);
}
if dependent.len() == 2
&& let Some(result) = try_tan_sec_patterns(arena, &dependent, var, var_sym)
{
return wrap_with_constants(arena, result, &constants);
}
if dependent.len() == 2
&& let Some(result) = try_trig_recip_product(arena, &dependent, var, var_sym)
{
if constants.is_empty() {
return result;
} else {
let mut all = constants.clone();
all.push(result);
return arena.mul(&all);
}
}
if dependent.len() == 2
&& let Some(result) =
try_x_over_sqrt_quadratic(arena, &dependent, var, var_sym, depth)
{
if constants.is_empty() {
return result;
} else {
let mut all = constants.clone();
all.push(result);
return arena.mul(&all);
}
}
if dependent.len() == 2
&& let Some(result) =
try_linear_over_quadratic(arena, &dependent, var, var_sym, depth)
{
if constants.is_empty() {
return result;
} else {
let mut all = constants.clone();
all.push(result);
return arena.mul(&all);
}
}
{
let all_children: SmallVec<[ExprId; 6]> = children.clone();
for (i, &child) in all_children.iter().enumerate() {
if let ExprNode::DiracDelta(delta_arg) = arena.node(child).clone() {
tracing::debug!(
"integrate: detected DiracDelta factor in Mul, attempting sifting property"
);
let other_factors: SmallVec<[ExprId; 4]> = all_children
.iter()
.enumerate()
.filter(|&(j, _)| j != i)
.map(|(_, &c)| c)
.collect();
let f_expr = if other_factors.len() == 1 {
other_factors[0]
} else if other_factors.is_empty() {
arena.one
} else {
arena.mul(&other_factors)
};
if delta_arg == var {
let f_at_0 =
crate::transforms::subs::subs(arena, f_expr, var, arena.zero);
let f_at_0_eval = crate::transforms::eval::eval(arena, f_at_0);
let heaviside = arena.intern(ExprNode::Heaviside(var));
return arena.mul(&[f_at_0_eval, heaviside]);
}
let solutions = crate::transforms::solve::solve(arena, delta_arg, var);
if solutions.len() == 1 {
let root = solutions[0].value;
let f_at_root = crate::transforms::subs::subs(arena, f_expr, var, root);
let f_at_root_eval = crate::transforms::eval::eval(arena, f_at_root);
let heaviside = arena.intern(ExprNode::Heaviside(delta_arg));
return arena.mul(&[f_at_root_eval, heaviside]);
}
break;
}
}
}
if let Some(result) = try_abs_sign_product(arena, &dependent, var, var_sym, depth) {
return wrap_with_constants(arena, result, &constants);
}
if dependent.len() == 2
&& let Some(result) = stage!(
arena,
"by_parts",
expr,
try_by_parts_pair(arena, &dependent, var, var_sym, depth)
)
{
return wrap_with_constants(arena, result, &constants);
}
if dependent.len() == 2
&& let Some(result) = stage!(
arena,
"cyclic_ibp",
expr,
try_cyclic_ibp(arena, &dependent, var, var_sym, depth)
)
{
if constants.is_empty() {
return result;
} else {
let mut all = constants.clone();
all.push(result);
return arena.mul(&all);
}
}
if let Some(result) = try_poly_times_half_power(arena, &dependent, var, var_sym, depth)
{
return wrap_with_constants(arena, result, &constants);
}
if let Some(result) =
try_reciprocal_sqrt_substitution(arena, &dependent, var, var_sym, depth)
{
return wrap_with_constants(arena, result, &constants);
}
if dependent.len() == 3
&& let Some(result) =
try_by_parts_poly_times_pair(arena, &dependent, var, var_sym, depth)
{
return wrap_with_constants(arena, result, &constants);
}
if dependent.len() >= 3
&& let Some(result) =
try_merge_polynomial_factors(arena, &dependent, var, var_sym, depth)
{
return wrap_with_constants(arena, result, &constants);
}
if dependent.len() >= 2
&& let Some(result) = stage!(
arena,
"trig_product_to_sum",
expr,
try_trig_product_to_sum(arena, expr, &dependent, var, var_sym, depth)
)
{
return result;
}
{
let (_numer, denom) = crate::poly::polybridge::as_numer_denom(arena, expr);
if denom != arena.one {
let decomposed = stage!(
arena,
"apart",
expr,
crate::transforms::apart::apart(arena, expr, var)
);
if decomposed != expr && !introduces_root_of(arena, expr, decomposed) {
let result = integrate_node(arena, decomposed, var, var_sym, depth - 1);
if !matches!(arena.node(result), ExprNode::Integral(_, _))
&& !candidate_rejected(arena, expr, result, var, var_sym)
{
return result;
}
}
}
}
if let Some(result) = stage!(
arena,
"u_substitution",
expr,
try_u_substitution(arena, &dependent, var, var_sym, depth - 1)
) {
if constants.is_empty() {
return result;
} else {
let mut all = constants.clone();
all.push(result);
return arena.mul(&all);
}
}
if dependent.len() >= 2
&& let Some(linear) = crate::transforms::trig_integ::linearize_trig_powers(
arena, &dependent, var, var_sym,
)
{
let result = integrate_node(arena, linear, var, var_sym, depth - 1);
if !crate::base::walk::has_unevaluated(arena, result) {
return wrap_with_constants(arena, result, &constants);
}
}
if let Some(result) = stage!(
arena,
"distribute_sum_factor",
expr,
try_distribute_sum_factor(arena, &dependent, var, var_sym, depth)
) {
return wrap_with_constants(arena, result, &constants);
}
let dependent_product = if dependent.len() == 1 {
dependent[0]
} else {
arena.mul(&dependent)
};
if let Some(result) = stage!(
arena,
"weierstrass",
dependent_product,
try_weierstrass_substitution(arena, dependent_product, var, var_sym, depth)
) {
return wrap_with_constants(arena, result, &constants);
}
if let Some(result) =
try_special_function_integral(arena, &dependent, &constants, var, var_sym)
{
return result;
}
tracing::debug!("integration: no strategy succeeded, returning unevaluated");
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Neg(inner) => {
let inner_int = integrate_node(arena, inner, var, var_sym, depth - 1);
arena.neg(inner_int)
}
ExprNode::Pow(base, exp) => {
if let Some(flattened) = flatten_nested_pow(arena, expr, var_sym)
.or_else(|| exp_power_as_exp(arena, expr, var_sym))
&& flattened != expr
{
return integrate_node(arena, flattened, var, var_sym, depth - 1);
}
if let ExprNode::Abs(g) = arena.node(base).clone()
&& contains_var(arena, g, var_sym)
&& real_on_reals(arena, g, var_sym)
&& let Some(n_val) = arena.as_num(exp).cloned()
&& n_val.is_integer()
&& n_val >= num_rational::Ratio::from_integer(2.into())
{
let n_int = n_val.to_integer();
let n_is_even = (&n_int % num_bigint::BigInt::from(2)).is_zero();
let rewritten = if n_is_even {
arena.pow(g, exp)
} else {
let n_minus_1 = arena.num_ratio(&n_val - &Q::one());
let g_pow = arena.pow(g, n_minus_1);
arena.mul(&[base, g_pow])
};
let result = integrate_node(arena, rewritten, var, var_sym, depth - 1);
if !crate::base::walk::has_unevaluated(arena, result) {
return result;
}
}
let base_is_var = base == var;
let exp_has_var = contains_var(arena, exp, var_sym);
let base_has_var = contains_var(arena, base, var_sym);
if base_is_var && !exp_has_var {
if let Some(n) = arena.as_num(exp) {
let n = n.clone();
if n == num_rational::Ratio::from_integer((-1).into()) {
let abs_x = arena.abs(var);
return arena.ln(abs_x);
}
let one = num_rational::Ratio::<num_bigint::BigInt>::one();
let n_plus_1 = &n + &one;
let n_plus_1_id = {
let nid = arena.intern_num(n_plus_1.clone());
arena.intern(ExprNode::Num(nid))
};
let x_pow = arena.pow(var, n_plus_1_id);
let recip = {
let inv = one / n_plus_1;
let nid = arena.intern_num(inv);
arena.intern(ExprNode::Num(nid))
};
return arena.mul(&[recip, x_pow]);
} else {
let one_id = arena.one;
let n_plus_1 = arena.add(&[exp, one_id]);
let x_pow = arena.pow(var, n_plus_1);
return arena.div(x_pow, n_plus_1);
}
}
if !base_has_var && !exp_has_var {
return arena.mul(&[expr, var]);
}
if !base_has_var && exp_has_var && base != arena.e_const() {
let ln_c = arena.ln(base);
let ln_c = crate::transforms::eval::eval(arena, ln_c);
let new_exp = arena.mul(&[exp, ln_c]);
let rewritten = arena.exp(new_exp);
let result = integrate_node(arena, rewritten, var, var_sym, depth - 1);
if !crate::base::walk::has_unevaluated(arena, result) {
return result;
}
return arena.intern(ExprNode::Integral(expr, var));
}
if let ExprNode::Tan(inner) = arena.node(base).clone()
&& let Some(n_val) = arena.as_num(exp).cloned()
&& n_val.is_integer()
&& n_val >= num_rational::Ratio::from_integer(3.into())
&& n_val <= num_rational::Ratio::from_integer(12.into())
{
let n_minus_2 =
arena.num_ratio(&n_val - &num_rational::Ratio::from_integer(2.into()));
let tan_pow = arena.pow(base, n_minus_2);
let cos_inner = arena.cos(inner);
let neg_two = arena.int(-2);
let sec_sq = arena.pow(cos_inner, neg_two);
let term1 = arena.mul(&[tan_pow, sec_sq]);
let rewritten = arena.sub(term1, tan_pow);
let result = integrate_node(arena, rewritten, var, var_sym, depth - 1);
if !crate::base::walk::has_unevaluated(arena, result) {
return result;
}
}
if let ExprNode::Cosh(inner) = arena.node(base).clone()
&& let Some(e_val) = arena.as_num(exp)
{
let neg_two_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer((-2).into());
if *e_val == neg_two_r {
if inner == var {
return arena.tanh(var);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym)
{
let tanh_inner = arena.tanh(inner);
return arena.div(tanh_inner, a_expr);
}
}
}
if let ExprNode::Sinh(inner) = arena.node(base).clone()
&& let Some(e_val) = arena.as_num(exp)
{
let neg_two_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer((-2).into());
if *e_val == neg_two_r {
if inner == var {
let cosh_v = arena.cosh(var);
let sinh_v = arena.sinh(var);
let neg1 = arena.int(-1);
let sinh_inv = arena.pow(sinh_v, neg1);
let coth_v = arena.mul(&[cosh_v, sinh_inv]);
return arena.neg(coth_v);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym)
{
let cosh_i = arena.cosh(inner);
let sinh_i = arena.sinh(inner);
let neg1 = arena.int(-1);
let sinh_inv = arena.pow(sinh_i, neg1);
let coth_i = arena.mul(&[cosh_i, sinh_inv]);
let neg_coth = arena.neg(coth_i);
return arena.div(neg_coth, a_expr);
}
}
}
if let ExprNode::Ln(inner) = arena.node(base).clone()
&& inner == var
&& !exp_has_var
&& let Some(n_val) = arena.as_num(exp)
{
let n_val = n_val.clone();
if n_val.is_integer() && n_val.is_positive() {
let n_i64: i64 = n_val.to_integer().try_into().unwrap_or(0);
if n_i64 >= 2 {
let x_ln_n = arena.mul(&[var, expr]);
let n_id = arena.num_ratio(n_val.clone());
let n_minus_1 = {
let v = &n_val - &num_rational::Ratio::<num_bigint::BigInt>::one();
arena.num_ratio(v.clone())
};
let ln_x = arena.ln(var);
let ln_nm1 = if n_i64 == 2 {
ln_x
} else {
arena.pow(ln_x, n_minus_1)
};
let sub_int = integrate_node(arena, ln_nm1, var, var_sym, depth - 1);
let n_times_sub = arena.mul(&[n_id, sub_int]);
return arena.sub(x_ln_n, n_times_sub);
}
}
}
if !exp_has_var
&& base_has_var
&& let Some((a_expr, _b_expr)) = symbolic_linear_coeff_of(arena, base, var, var_sym)
&& let Some(n) = arena.as_num(exp)
{
let n = n.clone();
let neg_one = num_rational::Ratio::from_integer((-1).into());
if n != neg_one {
let one = num_rational::Ratio::<num_bigint::BigInt>::one();
let n_plus_1 = &n + &one;
let n_plus_1_id = arena.num_ratio(n_plus_1.clone());
let base_pow = arena.pow(base, n_plus_1_id);
let denom = arena.mul(&[a_expr, n_plus_1_id]);
return arena.div(base_pow, denom);
} else {
let abs_base = arena.abs(base);
let ln_base = arena.ln(abs_base);
return arena.div(ln_base, a_expr);
}
}
if base_has_var
&& !exp_has_var
&& let Some(result) =
try_standard_form_integral(arena, expr, base, exp, var, var_sym)
{
return result;
}
if base_has_var
&& !exp_has_var
&& let Some(exp_val) = arena.as_num(exp)
{
let neg_one_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
if *exp_val == neg_one_r
&& let ExprNode::Add(ref ac) = arena.node(base).clone()
&& ac.len() == 2
{
let (mut x2_found, mut k_id) = (false, None);
for &ch in ac.iter() {
if is_var_squared(arena, ch, var) {
x2_found = true;
} else if !contains_var(arena, ch, var_sym) {
k_id = Some(ch);
}
}
if x2_found
&& let Some(k) = k_id
&& arena.as_num(k).is_none()
{
let half = arena.rational(1, 2);
let sqrt_k = arena.pow(k, half);
let x_over_sk = arena.div(var, sqrt_k);
let atan_val = arena.atan(x_over_sk);
let neg_half = arena.rational(-1, 2);
let inv_sk = arena.pow(k, neg_half);
return arena.mul(&[inv_sk, atan_val]);
}
}
}
if base_has_var
&& !exp_has_var
&& let Some(result) = try_complete_square_integral(arena, base, exp, var, var_sym)
{
return result;
}
if base_has_var
&& !exp_has_var
&& let Some(result) = try_trig_sub_sqrt_integral(arena, base, exp, var, var_sym)
{
return result;
}
{
let (_numer, denom) = crate::poly::polybridge::as_numer_denom(arena, expr);
if denom != arena.one {
let decomposed = stage!(
arena,
"apart",
expr,
crate::transforms::apart::apart(arena, expr, var)
);
if decomposed != expr && !introduces_root_of(arena, expr, decomposed) {
let result = integrate_node(arena, decomposed, var, var_sym, depth - 1);
if !matches!(arena.node(result), ExprNode::Integral(_, _))
&& !candidate_rejected(arena, expr, result, var, var_sym)
{
return result;
}
}
}
}
if let ExprNode::Add(_) = arena.node(base)
&& let Some(n) = arena.as_num(exp)
&& n.is_integer()
&& n.is_positive()
{
let n_i64: i64 = n.to_integer().try_into().unwrap_or(0);
if (2..=10).contains(&n_i64) {
let expanded = crate::transforms::expand::expand(arena, expr);
if expanded != expr {
let result = integrate_node(arena, expanded, var, var_sym, depth - 1);
if !matches!(arena.node(result), ExprNode::Integral(_, _)) {
return result;
}
}
}
}
if let Some(result) = stage!(
arena,
"weierstrass",
expr,
try_weierstrass_substitution(arena, expr, var, var_sym, depth)
) {
return result;
}
if let ExprNode::Ln(inner) = arena.node(base).clone()
&& inner == var
&& let Some(n_val) = arena.as_num(exp)
{
let neg_one_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
if *n_val == neg_one_r {
return arena.li(var);
}
}
if let Some(n_val) = arena.as_num(exp) {
let neg_one_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
if *n_val == neg_one_r
&& base_has_var
&& let Some((a_expr, _b_expr)) =
symbolic_linear_coeff_of(arena, base, var, var_sym)
{
let abs_base = arena.abs(base);
let ln_abs = arena.ln(abs_base);
return arena.div(ln_abs, a_expr);
}
}
if let ExprNode::Mul(ref children) = arena.node(base).clone()
&& let Some(e_val) = arena.as_num(exp)
&& e_val.is_negative()
&& e_val.is_integer()
{
let factors: SmallVec<[ExprId; 6]> = children
.iter()
.map(|&child| arena.pow(child, exp))
.collect();
let distributed = arena.mul(&factors);
if distributed != expr {
let result = integrate_node(arena, distributed, var, var_sym, depth - 1);
if !matches!(arena.node(result), ExprNode::Integral(_, _)) {
return result;
}
}
}
tracing::debug!("integration: no strategy succeeded, returning unevaluated");
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Sin(inner) => {
if inner == var {
let cos_x = arena.cos(var);
return arena.neg(cos_x);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let cos_inner = arena.cos(inner);
let neg_cos = arena.neg(cos_inner);
return arena.div(neg_cos, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Cos(inner) => {
if inner == var {
return arena.sin(var);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let sin_inner = arena.sin(inner);
return arena.div(sin_inner, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Tan(inner) => {
if inner == var {
let cos_x = arena.cos(var);
let abs_cos = arena.abs(cos_x);
let ln_abs_cos = arena.ln(abs_cos);
return arena.neg(ln_abs_cos);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let cos_inner = arena.cos(inner);
let abs_cos = arena.abs(cos_inner);
let ln_abs_cos = arena.ln(abs_cos);
let neg_ln = arena.neg(ln_abs_cos);
return arena.div(neg_ln, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Exp(inner) => {
if inner == var {
return arena.exp(var);
}
if let Some(poly) = crate::poly::polybridge::expr_to_poly(arena, inner, var)
&& poly.degree() == Some(2)
{
let coeffs = poly.coeffs(); let a_coeff = &coeffs[2];
let b_coeff = &coeffs[1];
let c_coeff = &coeffs[0];
if a_coeff.is_negative() {
let neg_a = -a_coeff.clone();
let neg_a_expr = arena.num_ratio(neg_a.clone());
let sqrt_neg_a = arena.sqrt(neg_a_expr);
let two = arena.int(2);
let pi_id = arena.pi;
let sqrt_pi = arena.sqrt(pi_id);
let two_sqrt_neg_a = arena.mul(&[two, sqrt_neg_a]);
let front = arena.div(sqrt_pi, two_sqrt_neg_a);
let b_expr = arena.num_ratio(b_coeff.clone());
let a_expr = arena.num_ratio(a_coeff.clone());
let c_expr = arena.num_ratio(c_coeff.clone());
let b_sq = arena.mul(&[b_expr, b_expr]);
let four = arena.int(4);
let four_a = arena.mul(&[four, a_expr]);
let b_sq_over_4a = arena.div(b_sq, four_a);
let exp_arg = arena.sub(c_expr, b_sq_over_4a);
let exp_arg = arena.eval_expr(exp_arg);
let exp_factor = if arena.is_zero_structural(exp_arg) {
arena.int(1)
} else {
arena.exp(exp_arg)
};
let neg_two_a = {
let two_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer(2.into());
let val = -two_r * a_coeff;
arena.num_ratio(val.clone())
};
let neg_2ax = arena.mul(&[neg_two_a, var]);
let erf_numer = arena.sub(neg_2ax, b_expr);
let sqrt_neg_a2 = arena.sqrt(neg_a_expr);
let erf_denom = arena.mul(&[two, sqrt_neg_a2]);
let erf_arg = arena.div(erf_numer, erf_denom);
let erf_term = arena.erf(erf_arg);
return arena.mul(&[front, exp_factor, erf_term]);
}
if a_coeff.is_positive() {
let a_expr = arena.num_ratio(a_coeff.clone());
let sqrt_a = arena.sqrt(a_expr);
let two = arena.int(2);
let pi_id = arena.pi;
let sqrt_pi = arena.sqrt(pi_id);
let two_sqrt_a = arena.mul(&[two, sqrt_a]);
let front = arena.div(sqrt_pi, two_sqrt_a);
let b_expr = arena.num_ratio(b_coeff.clone());
let c_expr = arena.num_ratio(c_coeff.clone());
let b_sq = arena.mul(&[b_expr, b_expr]);
let four = arena.int(4);
let four_a = arena.mul(&[four, a_expr]);
let b_sq_over_4a = arena.div(b_sq, four_a);
let exp_arg = arena.sub(c_expr, b_sq_over_4a);
let exp_arg = arena.eval_expr(exp_arg);
let exp_factor = if arena.is_zero_structural(exp_arg) {
arena.int(1)
} else {
arena.exp(exp_arg)
};
let two_a = {
let two_r =
num_rational::Ratio::<num_bigint::BigInt>::from_integer(2.into());
let val = two_r * a_coeff;
arena.num_ratio(val.clone())
};
let two_ax = arena.mul(&[two_a, var]);
let numer = arena.add(&[two_ax, b_expr]);
let sqrt_a2 = arena.sqrt(a_expr);
let denom = arena.mul(&[two, sqrt_a2]);
let erfi_arg = arena.div(numer, denom);
let erfi_term = special_apply(arena, FN_ERFI, erfi_arg);
return arena.mul(&[front, exp_factor, erfi_term]);
}
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let exp_inner = arena.exp(inner);
return arena.div(exp_inner, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Ln(inner) => {
if inner == var {
let ln_var = arena.ln(var);
let x_ln_x = arena.mul(&[var, ln_var]);
return arena.sub(x_ln_x, var);
}
if let ExprNode::Ln(ln_inner) = arena.node(inner).clone()
&& ln_inner == var
{
let ln_x = arena.ln(var);
let ln_ln_x = arena.ln(ln_x);
let x_ln_ln_x = arena.mul(&[var, ln_ln_x]);
let li_x = arena.li(var);
return arena.sub(x_ln_ln_x, li_x);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let ln_inner = arena.ln(inner);
let inner_times_ln = arena.mul(&[inner, ln_inner]);
let diff = arena.sub(inner_times_ln, inner);
return arena.div(diff, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Sinh(inner) => {
if inner == var {
return arena.intern(ExprNode::Cosh(var));
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let cosh_inner = arena.cosh(inner);
return arena.div(cosh_inner, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Cosh(inner) => {
if inner == var {
return arena.intern(ExprNode::Sinh(var));
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let sinh_inner = arena.sinh(inner);
return arena.div(sinh_inner, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Tanh(inner) => {
if inner == var {
let cosh_x = arena.cosh(var);
return arena.ln(cosh_x);
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let cosh_inner = arena.cosh(inner);
let ln_cosh = arena.ln(cosh_inner);
return arena.div(ln_cosh, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Asin(inner) | ExprNode::Acos(inner) | ExprNode::Atan(inner) => {
let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) else {
return arena.intern(ExprNode::Integral(expr, var));
};
let g = inner;
let g_f = arena.mul(&[g, expr]);
let one = arena.one;
let two = arena.int(2);
let g2 = arena.pow(g, two);
let half = arena.rational(1, 2);
let numer = match node {
ExprNode::Asin(_) => {
let one_minus_g2 = arena.sub(one, g2);
let sqrt_term = arena.pow(one_minus_g2, half);
arena.add(&[g_f, sqrt_term])
}
ExprNode::Acos(_) => {
let one_minus_g2 = arena.sub(one, g2);
let sqrt_term = arena.pow(one_minus_g2, half);
arena.sub(g_f, sqrt_term)
}
_ => {
let one_plus_g2 = arena.add(&[one, g2]);
let ln_term = arena.ln(one_plus_g2);
let half_ln = arena.mul(&[half, ln_term]);
arena.sub(g_f, half_ln)
}
};
if a_expr == one {
return numer;
}
arena.div(numer, a_expr)
}
ExprNode::Erf(inner) | ExprNode::Erfc(inner) => {
let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) else {
return arena.intern(ExprNode::Integral(expr, var));
};
let g = inner;
let g_f = arena.mul(&[g, expr]);
let two = arena.int(2);
let g2 = arena.pow(g, two);
let neg_g2 = arena.neg(g2);
let e = arena.exp(neg_g2);
let pi = arena.pi();
let sqrt_pi = arena.sqrt(pi);
let gauss = arena.div(e, sqrt_pi);
let numer = if matches!(node, ExprNode::Erf(_)) {
arena.add(&[g_f, gauss])
} else {
arena.sub(g_f, gauss)
};
if a_expr == arena.one {
return numer;
}
arena.div(numer, a_expr)
}
ExprNode::Abs(_) | ExprNode::Sign(_) => {
if let Some(r) = try_abs_sign_product(arena, &[expr], var, var_sym, depth) {
return r;
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Piecewise(ref pairs) => {
if let Some(r) = integrate_piecewise(arena, pairs, var, var_sym, depth) {
return r;
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::DiracDelta(inner) => {
if inner == var {
return arena.intern(ExprNode::Heaviside(var));
}
if let Some((a_expr, _)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let h = arena.intern(ExprNode::Heaviside(inner));
let abs_a = arena.abs(a_expr);
return arena.div(h, abs_a);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Heaviside(inner) => {
if inner == var {
return arena.mul(&[var, expr]);
}
if let Some((a_expr, _b_expr)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let h = arena.intern(ExprNode::Heaviside(inner));
let product = arena.mul(&[inner, h]);
return arena.div(product, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Asinh(inner) => {
if inner == var {
tracing::debug!("integrate: matched asinh(x) direct");
let asinh_var = arena.asinh(var);
let x_asinh = arena.mul(&[var, asinh_var]);
let two = arena.int(2);
let x2 = arena.pow(var, two);
let one = arena.one;
let x2_plus_1 = arena.add(&[x2, one]);
let half = arena.rational(1, 2);
let sqrt_term = arena.pow(x2_plus_1, half);
return arena.sub(x_asinh, sqrt_term);
}
if let Some((a_expr, _b_expr)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
tracing::debug!("integrate: matched asinh(ax+b) linear");
let asinh_g = arena.asinh(inner);
let g_asinh = arena.mul(&[inner, asinh_g]);
let two = arena.int(2);
let g2 = arena.pow(inner, two);
let one = arena.one;
let g2_plus_1 = arena.add(&[g2, one]);
let half = arena.rational(1, 2);
let sqrt_term = arena.pow(g2_plus_1, half);
let numer = arena.sub(g_asinh, sqrt_term);
return arena.div(numer, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Acosh(inner) => {
let a_expr = if inner == var {
Some(arena.one)
} else {
symbolic_linear_coeff_of(arena, inner, var, var_sym).map(|(a, _)| a)
};
if let Some(a_expr) = a_expr {
tracing::debug!("integrate: matched acosh(ax+b)");
let acosh_g = arena.acosh(inner);
let g_acosh = arena.mul(&[inner, acosh_g]);
let one = arena.one;
let g_minus_1 = arena.sub(inner, one);
let g_plus_1 = arena.add(&[inner, one]);
let half = arena.rational(1, 2);
let root_minus = arena.pow(g_minus_1, half);
let root_plus = arena.pow(g_plus_1, half);
let roots = arena.mul(&[root_minus, root_plus]);
let numer = arena.sub(g_acosh, roots);
return if a_expr == arena.one {
numer
} else {
arena.div(numer, a_expr)
};
}
arena.intern(ExprNode::Integral(expr, var))
}
ExprNode::Atanh(inner) => {
if inner == var {
tracing::debug!("integrate: matched atanh(x) direct");
let atanh_var = arena.atanh(var);
let x_atanh = arena.mul(&[var, atanh_var]);
let two = arena.int(2);
let x2 = arena.pow(var, two);
let one = arena.one;
let one_minus_x2 = arena.sub(one, x2);
let half = arena.rational(1, 2);
let ln_term = arena.ln(one_minus_x2);
let half_ln = arena.mul(&[half, ln_term]);
return arena.add(&[x_atanh, half_ln]);
}
if let Some((a_expr, _b_expr)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
tracing::debug!("integrate: matched atanh(ax+b) linear");
let atanh_g = arena.atanh(inner);
let g_atanh = arena.mul(&[inner, atanh_g]);
let two = arena.int(2);
let g2 = arena.pow(inner, two);
let one = arena.one;
let one_minus_g2 = arena.sub(one, g2);
let half = arena.rational(1, 2);
let ln_term = arena.ln(one_minus_g2);
let half_ln = arena.mul(&[half, ln_term]);
let numer = arena.add(&[g_atanh, half_ln]);
return arena.div(numer, a_expr);
}
arena.intern(ExprNode::Integral(expr, var))
}
_ => {
tracing::debug!("integration: no strategy succeeded, returning unevaluated");
arena.intern(ExprNode::Integral(expr, var))
}
}
}
fn contains_var(arena: &Arena, expr: ExprId, var: SymbolId) -> bool {
crate::base::walk::has_free_symbol(arena, expr, var)
}
fn is_polynomial_in(arena: &Arena, expr: ExprId, var: ExprId, var_sym: SymbolId) -> bool {
let mut stack = vec![expr];
while let Some(e) = stack.pop() {
if e == var || !contains_var(arena, e, var_sym) {
continue;
}
match arena.node(e) {
ExprNode::Pow(base, exp) => match arena.as_num(*exp) {
Some(r)
if r.is_integer()
&& !r.is_negative()
&& (*base == var || *r <= Q::from_integer(MAX_POLYNOMIAL_POWER.into())) =>
{
stack.push(*base);
}
_ => return false,
},
ExprNode::Mul(children) | ExprNode::Add(children) => {
stack.extend(children.iter().copied())
}
ExprNode::Neg(inner) => stack.push(*inner),
ExprNode::Num(_) => {}
_ => return false,
}
}
true
}
const MAX_POLYNOMIAL_POWER: i64 = 16;
#[allow(dead_code)]
fn linear_coeff_of(arena: &Arena, expr: ExprId, _var: ExprId, _var_sym: SymbolId) -> Option<Q> {
let poly = crate::poly::polybridge::expr_to_poly(arena, expr, _var)?;
if poly.degree()? != 1 {
return None;
}
let a = poly.coeff(1);
if a.is_zero() {
return None;
}
Some(a)
}
fn symbolic_linear_coeff_of(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<(ExprId, ExprId)> {
if let Some(a) = linear_coeff_of(arena, expr, var, var_sym) {
let a_id = arena.num_ratio(a.clone());
if let Some(poly) = crate::poly::polybridge::expr_to_poly(arena, expr, var) {
let b = poly.coeff(0);
let b_id = arena.num_ratio(b.clone());
return Some((a_id, b_id));
}
return Some((a_id, arena.zero));
}
if expr == var {
return Some((arena.one, arena.zero));
}
if let ExprNode::Neg(inner) = arena.node(expr).clone() {
if let Some((coeff, constant)) = symbolic_linear_coeff_of(arena, inner, var, var_sym) {
let neg_coeff = arena.neg(coeff);
let neg_const = arena.neg(constant);
return Some((neg_coeff, neg_const));
}
return None;
}
if let ExprNode::Mul(ref children) = arena.node(expr).clone() {
let mut has_var = false;
let mut other_factors: SmallVec<[ExprId; 4]> = SmallVec::new();
let mut var_count = 0u32;
let mut nontrivial = false;
for &child in children {
if child == var {
var_count += 1;
if var_count > 1 {
return None; }
has_var = true;
} else if contains_var(arena, child, var_sym) {
nontrivial = true; } else {
other_factors.push(child);
}
}
if has_var && var_count == 1 && !nontrivial {
let coeff = match other_factors.len() {
0 => arena.one,
1 => other_factors[0],
_ => arena.mul(&other_factors),
};
return Some((coeff, arena.zero));
}
let dependent: SmallVec<[ExprId; 2]> = children
.iter()
.copied()
.filter(|&c| contains_var(arena, c, var_sym))
.collect();
if let [sum] = dependent[..]
&& matches!(arena.node(sum), ExprNode::Add(_))
&& let Some((a, b)) = linear_sum_parts(arena, sum, var, var_sym)
{
let c: SmallVec<[ExprId; 4]> = children.iter().copied().filter(|&f| f != sum).collect();
let c = arena.mul(&c);
let a = arena.mul(&[c, a]);
let b = arena.mul(&[c, b]);
return Some((a, b));
}
}
if matches!(arena.node(expr), ExprNode::Add(_)) {
return linear_sum_parts(arena, expr, var, var_sym);
}
None
}
fn linear_sum_parts(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<(ExprId, ExprId)> {
{
let ExprNode::Add(ref children) = arena.node(expr).clone() else {
return None;
};
let mut var_terms: SmallVec<[ExprId; 4]> = SmallVec::new();
let mut const_terms: SmallVec<[ExprId; 4]> = SmallVec::new();
for &child in children {
if contains_var(arena, child, var_sym) {
var_terms.push(child);
} else {
const_terms.push(child);
}
}
if var_terms.is_empty() {
return None; }
let mut coeffs: SmallVec<[ExprId; 4]> = SmallVec::new();
for &term in &var_terms {
if term == var {
coeffs.push(arena.one);
continue;
}
let ExprNode::Mul(ref mul_children) = arena.node(term).clone() else {
return None; };
let mut vc = 0u32;
let mut other: SmallVec<[ExprId; 4]> = SmallVec::new();
for &mc in mul_children {
if mc == var {
vc += 1;
} else if contains_var(arena, mc, var_sym) {
return None;
} else {
other.push(mc);
}
}
if vc != 1 {
return None;
}
coeffs.push(match other.len() {
0 => arena.one,
1 => other[0],
_ => arena.mul(&other),
});
}
let coeff = match coeffs.len() {
1 => coeffs[0],
_ => arena.add(&coeffs),
};
if coeff == arena.zero {
return None;
}
if contains_var(arena, coeff, var_sym) {
return None;
}
let constant = match const_terms.len() {
0 => arena.zero,
1 => const_terms[0],
_ => arena.add(&const_terms),
};
if contains_var(arena, constant, var_sym) {
return None;
}
Some((coeff, constant))
}
}
fn symbolic_quadratic_coeffs(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<(ExprId, ExprId, ExprId)> {
tracing::trace!("symbolic_quadratic_coeffs: attempting to extract quadratic coefficients");
if let Some(poly) = crate::poly::polybridge::expr_to_poly(arena, expr, var) {
if poly.degree()? == 2 {
let c = arena.num_ratio(poly.coeff(2).clone());
let d = arena.num_ratio(poly.coeff(1).clone());
let e = arena.num_ratio(poly.coeff(0).clone());
return Some((c, d, e));
}
return None;
}
tracing::trace!("symbolic_quadratic_coeffs: rational path failed, trying symbolic extraction");
let children = match arena.node(expr).clone() {
ExprNode::Add(c) => c.to_vec(),
_ => vec![expr],
};
let mut x2_terms: SmallVec<[ExprId; 4]> = SmallVec::new(); let mut x1_terms: SmallVec<[ExprId; 4]> = SmallVec::new(); let mut x0_terms: SmallVec<[ExprId; 4]> = SmallVec::new();
let two = arena.int(2);
for &child in &children {
if !contains_var(arena, child, var_sym) {
x0_terms.push(child);
continue;
}
if child == arena.pow(var, two) {
x2_terms.push(arena.one);
continue;
}
if let ExprNode::Pow(base, exp) = arena.node(child).clone()
&& base == var
{
if let Some(e) = arena.as_num(exp)
&& *e == num_rational::Ratio::from_integer(2.into())
{
x2_terms.push(arena.one);
continue;
}
return None;
}
if let ExprNode::Mul(ref mul_children) = arena.node(child).clone() {
let mul_children = mul_children.clone();
let mut has_var_sq = false;
let mut var_count = 0u32;
let mut other_factors: SmallVec<[ExprId; 4]> = SmallVec::new();
for &mc in &mul_children {
if mc == var {
var_count += 1;
if var_count > 2 {
return None; }
} else if let ExprNode::Pow(base, exp) = arena.node(mc).clone() {
if base == var {
{
let e = arena.as_num(exp)?;
if *e == num_rational::Ratio::from_integer(2.into()) {
has_var_sq = true;
} else if e.is_integer()
&& *e > num_rational::Ratio::from_integer(2.into())
{
return None; } else {
return None;
}
}
} else if contains_var(arena, mc, var_sym) {
return None; } else {
other_factors.push(mc);
}
} else if contains_var(arena, mc, var_sym) {
return None; } else {
other_factors.push(mc);
}
}
let scalar = match other_factors.len() {
0 => arena.one,
1 => other_factors[0],
_ => arena.mul(&other_factors),
};
if has_var_sq || var_count == 2 {
x2_terms.push(scalar);
} else if var_count == 1 {
x1_terms.push(scalar);
} else {
x0_terms.push(child);
}
continue;
}
if child == var {
x1_terms.push(arena.one);
continue;
}
if let ExprNode::Neg(inner) = arena.node(child).clone() {
if let Some((c, d, e)) = symbolic_quadratic_coeffs(arena, inner, var, var_sym) {
x2_terms.push(arena.neg(c));
x1_terms.push(arena.neg(d));
x0_terms.push(arena.neg(e));
continue;
}
return None;
}
return None;
}
if x2_terms.is_empty() {
tracing::trace!("symbolic_quadratic_coeffs: no x² terms found, not quadratic");
return None;
}
tracing::trace!(
n_x2_terms = x2_terms.len(),
n_x1_terms = x1_terms.len(),
n_x0_terms = x0_terms.len(),
"symbolic_quadratic_coeffs: classified terms"
);
let c_expr = match x2_terms.len() {
1 => x2_terms[0],
_ => arena.add(&x2_terms),
};
let d_expr = match x1_terms.len() {
0 => arena.zero,
1 => x1_terms[0],
_ => arena.add(&x1_terms),
};
let e_expr = match x0_terms.len() {
0 => arena.zero,
1 => x0_terms[0],
_ => arena.add(&x0_terms),
};
if contains_var(arena, c_expr, var_sym) {
return None;
}
Some((c_expr, d_expr, e_expr))
}
fn try_u_substitution(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
for (i, &factor) in dependent.iter().enumerate() {
let candidates = u_sub_candidates(arena, factor, var_sym);
for u_expr in candidates {
if u_expr == var {
continue;
}
let du = crate::transforms::diff::diff(arena, u_expr, var);
if du == arena.zero {
continue;
}
let remaining_expr = remaining_product(arena, dependent, i);
let quotient = arena.div(remaining_expr, du);
let placeholder = arena.symbol("__usub_u");
let mut through_u: Option<ExprId> = None;
let coeff = if !contains_var(arena, quotient, var_sym) {
quotient
} else {
let cancelled = arena.cancel_expr(quotient, var);
if !contains_var(arena, cancelled, var_sym) {
cancelled
} else {
let varies = var_outside_trig_arguments(arena, cancelled, var_sym)
|| varies_numerically(arena, cancelled, var);
let trig = if varies {
cancelled
} else {
let trig = arena.trigsimp_expr(cancelled);
crate::transforms::eval::eval(arena, trig)
};
if !contains_var(arena, trig, var_sym) {
trig
} else {
if !contains_var(arena, du, var_sym) {
continue;
}
let q_in_u = arena.subs_structural(quotient, u_expr, placeholder);
if contains_var(arena, q_in_u, var_sym) {
continue;
}
through_u = Some(quotient);
arena.one
}
}
};
if contains_non_finite(arena, coeff) {
continue;
}
let factor = match through_u {
Some(q) => arena.mul(&[factor, q]),
None => factor,
};
let factor_in_u = arena.subs_structural(factor, u_expr, placeholder);
if contains_var(arena, factor_in_u, var_sym) {
continue;
}
let g_of_var = arena.subs_structural(factor, u_expr, var);
let g_integrated = integrate_node(arena, g_of_var, var, var_sym, depth);
if matches!(arena.node(g_integrated), ExprNode::Integral(_, _)) {
continue;
}
let antideriv = arena.subs_structural(g_integrated, var, u_expr);
let antideriv = match continuous_through_tan_of_u(
arena,
g_integrated,
antideriv,
u_expr,
var,
var_sym,
) {
Some(continuous) => continuous,
None if jump_kept_undecided(arena, g_integrated, var_sym) => {
match tan_argument_of(arena, u_expr, var_sym) {
Some(w) => crate::transforms::trig_integ::atan_tan_terms_to_argument_in(
arena,
antideriv,
w,
Some(var_sym),
),
None => antideriv,
}
}
None => continue,
};
tracing::debug!("u-substitution succeeded");
return Some(arena.mul(&[coeff, antideriv]));
}
}
None
}
fn continuous_through_tan_of_u(
arena: &mut Arena,
g: ExprId,
antideriv: ExprId,
u_expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let Some(w) = tan_argument_of(arena, u_expr, var_sym) else {
return Some(antideriv);
};
let tan_w = arena.tan(w);
let s = arena.symbol("__tan_s");
let ExprNode::Symbol(s_sym) = *arena.node(s) else {
return None;
};
let u_of_s = arena.subs_structural(u_expr, tan_w, s);
if contains_var(arena, u_of_s, var_sym) {
return None;
}
let h = arena.subs_structural(g, var, u_of_s);
let (antideriv, jump) = match crate::transforms::trig_integ::infinity_jump(arena, h, s, s_sym) {
Some(jump) => (antideriv, jump),
None => {
let all = crate::base::walk::free_symbols(arena, w)
.into_iter()
.any(|v| match *arena.node(v) {
ExprNode::Symbol(sid) => {
sid != var_sym
&& !crate::transforms::realness::symbol_declared_real(arena, sid)
}
_ => true,
});
let g = crate::transforms::trig_integ::plain_logs(arena, g, var_sym, all);
let h = arena.subs_structural(g, var, u_of_s);
let jump =
crate::transforms::trig_integ::infinity_jump_real_parameters(arena, h, s, s_sym)?;
(arena.subs_structural(g, var, u_expr), jump)
}
};
Some(crate::transforms::trig_integ::add_tan_floor_in(
arena,
antideriv,
w,
jump,
Some(var_sym),
))
}
fn tan_argument_of(arena: &Arena, u: ExprId, var_sym: SymbolId) -> Option<ExprId> {
crate::base::walk::post_order_ids(arena, u)
.into_iter()
.find_map(|id| match *arena.node(id) {
ExprNode::Tan(w) if contains_var(arena, w, var_sym) => Some(w),
_ => None,
})
}
fn var_outside_trig_arguments(arena: &Arena, q: ExprId, var_sym: SymbolId) -> bool {
let mut stack = vec![q];
let mut seen: FxHashSet<ExprId> = FxHashSet::default();
while let Some(id) = stack.pop() {
if !seen.insert(id) {
continue;
}
match arena.node(id) {
ExprNode::Sin(_)
| ExprNode::Cos(_)
| ExprNode::Tan(_)
| ExprNode::Sinh(_)
| ExprNode::Cosh(_)
| ExprNode::Tanh(_) => {}
ExprNode::Symbol(s) if *s == var_sym => return true,
node => node.for_each_child(|c| stack.push(c)),
}
}
false
}
const VARIES_POINTS: [(i64, i64); 2] = [(3, 7), (5, 11)];
fn varies_numerically(arena: &mut Arena, q: ExprId, var: ExprId) -> bool {
let inside = NODE_BUDGET_DEPTH.with(std::cell::Cell::get) > 0;
let key = (std::ptr::from_ref::<Arena>(arena) as usize, q, var);
if inside && let Some(found) = VARIES_MEMO.with(|m| m.borrow().get(&key).copied()) {
return found;
}
let varies = varies_numerically_uncached(arena, q, var);
if inside {
VARIES_MEMO.with(|m| m.borrow_mut().insert(key, varies));
}
varies
}
thread_local! {
static VARIES_MEMO: std::cell::RefCell<FxHashMap<(usize, ExprId, ExprId), bool>> =
std::cell::RefCell::new(FxHashMap::default());
}
fn varies_numerically_uncached(arena: &mut Arena, q: ExprId, var: ExprId) -> bool {
let mut params: Vec<(String, ExprId)> = Vec::new();
for s in crate::base::walk::free_symbols(arena, q) {
if s == var {
continue;
}
let ExprNode::Symbol(sid) = *arena.node(s) else {
return false;
};
params.push((arena.symbol_name(sid).to_owned(), s));
}
if params.len() > FTC_PARAMETER_VALUES.len() {
return false;
}
params.sort();
let mut bound = q;
for (k, &(_, s)) in params.iter().enumerate() {
let (p, d) = FTC_PARAMETER_VALUES[k];
let value = arena.rational(p, d);
bound = crate::transforms::subs::subs(arena, bound, s, value);
}
let mut values = Vec::with_capacity(VARIES_POINTS.len());
for &(p, d) in &VARIES_POINTS {
let point = arena.rational(p, d);
let at = crate::transforms::subs::subs(arena, bound, var, point);
let at = crate::transforms::eval::eval(arena, at);
if crate::base::walk::has_unevaluated(arena, at) {
return false;
}
let Ok(v) = crate::transforms::evalf::evalf_complex(arena, at, 15) else {
return false;
};
values.push(v);
}
let (Some(scale), Some(gap)) = (
crate::transforms::evalf::abs_to_f64(&values[0]),
crate::transforms::evalf::distance_to_f64(&values[0], &values[1]),
) else {
return false;
};
scale.is_finite() && gap.is_finite() && gap > 1e-8 * scale.max(1.0)
}
fn u_sub_candidates(arena: &Arena, factor: ExprId, var_sym: SymbolId) -> SmallVec<[ExprId; 4]> {
let mut out: SmallVec<[ExprId; 4]> = SmallVec::new();
match arena.node(factor).clone() {
ExprNode::Sin(inner)
| ExprNode::Cos(inner)
| ExprNode::Tan(inner)
| ExprNode::Exp(inner)
| ExprNode::Ln(inner)
| ExprNode::Sinh(inner)
| ExprNode::Cosh(inner)
| ExprNode::Tanh(inner)
| ExprNode::Asin(inner)
| ExprNode::Acos(inner)
| ExprNode::Atan(inner)
| ExprNode::Asinh(inner)
| ExprNode::Acosh(inner)
| ExprNode::Atanh(inner)
| ExprNode::Abs(inner)
if contains_var(arena, inner, var_sym) =>
{
out.push(inner);
out.push(factor);
}
ExprNode::Pow(base, _exp) if contains_var(arena, base, var_sym) => {
out.push(base);
inner_function_candidates(arena, base, var_sym, &mut out);
}
ExprNode::Add(_) => {
inner_function_candidates(arena, factor, var_sym, &mut out);
}
_ => {}
}
out
}
fn inner_function_candidates(
arena: &Arena,
expr: ExprId,
var_sym: SymbolId,
out: &mut SmallVec<[ExprId; 4]>,
) {
for id in crate::base::walk::post_order_ids(arena, expr) {
if out.len() >= 4 {
return;
}
if id == expr || out.contains(&id) {
continue;
}
let is_candidate = match arena.node(id) {
ExprNode::Sin(inner)
| ExprNode::Cos(inner)
| ExprNode::Tan(inner)
| ExprNode::Exp(inner)
| ExprNode::Ln(inner)
| ExprNode::Sinh(inner)
| ExprNode::Cosh(inner)
| ExprNode::Tanh(inner) => contains_var(arena, *inner, var_sym),
_ => false,
};
if is_candidate {
out.push(id);
}
}
}
fn remaining_product(arena: &mut Arena, children: &[ExprId], skip: usize) -> ExprId {
let parts: SmallVec<[ExprId; 4]> = children
.iter()
.enumerate()
.filter(|&(j, _)| j != skip)
.map(|(_, &c)| c)
.collect();
match parts.len() {
0 => arena.one,
1 => parts[0],
_ => arena.mul(&parts),
}
}
thread_local! {
static SUBST_DEPTH: std::cell::Cell<u8> = const { std::cell::Cell::new(0) };
}
const MAX_SUBST_DEPTH: u8 = 3;
thread_local! {
static GENERIC_ONLY: std::cell::Cell<bool> = const { std::cell::Cell::new(false) };
}
fn integrate_nested_generic(arena: &mut Arena, expr: ExprId, var: ExprId) -> Option<ExprId> {
struct Restore(bool);
impl Drop for Restore {
fn drop(&mut self) {
let _ = GENERIC_ONLY.try_with(|g| g.set(self.0));
}
}
let restore = Restore(GENERIC_ONLY.with(|g| g.replace(true)));
let result = integrate_nested(arena, expr, var);
drop(restore);
result
}
fn integrate_nested(arena: &mut Arena, expr: ExprId, var: ExprId) -> Option<ExprId> {
struct Restore(u8);
impl Drop for Restore {
fn drop(&mut self) {
let _ = SUBST_DEPTH.try_with(|d| d.set(self.0));
}
}
let depth = SUBST_DEPTH.with(|d| d.get());
if depth >= MAX_SUBST_DEPTH {
return None;
}
SUBST_DEPTH.with(|d| d.set(depth + 1));
let restore = Restore(depth);
let result = integrate(arena, expr, var);
drop(restore);
if crate::base::walk::has_unevaluated(arena, result) {
None
} else {
Some(result)
}
}
fn try_substitution_strategies(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
if let Some(r) = stage!(
arena,
"trig_rational",
expr,
try_trig_rational_substitution(arena, expr, var, var_sym)
) {
return Some(r);
}
if let Some(r) = try_exp_rational_substitution(arena, expr, var, var_sym) {
return Some(r);
}
if let Some(r) = try_exp_linear_substitution(arena, expr, var, var_sym) {
return Some(r);
}
if let Some(r) = try_radical_substitution(arena, expr, var, var_sym) {
return Some(r);
}
if let Some(r) = try_even_power_radicals(arena, expr, var, var_sym) {
return Some(r);
}
if let Some(r) = try_mobius_radical_substitution(arena, expr, var, var_sym) {
return Some(r);
}
if let Some(r) = try_exp_atanh_rewrite(arena, expr, var, var_sym) {
return Some(r);
}
if let Some(r) = try_hyperbolic_to_exp(arena, expr, var, var_sym) {
return Some(r);
}
None
}
fn monomial_parts(
arena: &mut Arena,
g: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<(ExprId, ExprId)> {
let power = |arena: &Arena, f: ExprId| -> Option<ExprId> {
if f == var {
return Some(arena.one);
}
match arena.node(f) {
ExprNode::Pow(b, e)
if *b == var
&& !contains_var(arena, *e, var_sym)
&& !arena.is_zero_structural(*e) =>
{
Some(*e)
}
_ => None,
}
};
if let Some(n) = power(arena, g) {
return Some((arena.one, n));
}
let ExprNode::Mul(factors) = arena.node(g).clone() else {
return None;
};
let mut n: Option<ExprId> = None;
let mut rest: SmallVec<[ExprId; 4]> = SmallVec::new();
for f in factors {
if let Some(k) = power(arena, f) {
if n.replace(k).is_some() {
return None;
}
} else if contains_var(arena, f, var_sym) {
return None;
} else {
rest.push(f);
}
}
let c = match rest.len() {
0 => arena.one,
1 => rest[0],
_ => arena.mul(&rest),
};
Some((c, n?))
}
fn try_log_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let order = crate::base::walk::post_order_ids(arena, expr);
let mut log: Option<(ExprId, ExprId)> = None;
for &id in &order {
if let ExprNode::Ln(g) = *arena.node(id)
&& contains_var(arena, g, var_sym)
{
match log {
Some((l, _)) if l != id => return None,
_ => log = Some((id, g)),
}
}
}
let (ln_node, g) = log?;
let inside_function = order.iter().any(|&id| {
matches!(arena.node(id),
ExprNode::Sin(a) | ExprNode::Cos(a) | ExprNode::Tan(a) | ExprNode::Sinh(a)
| ExprNode::Cosh(a) | ExprNode::Tanh(a) | ExprNode::Exp(a)
if crate::base::walk::contains(arena, *a, ln_node))
});
if !inside_function {
return None;
}
let w = arena.symbol("__lw");
let ExprNode::Symbol(w_sym) = *arena.node(w) else {
return None;
};
let Some((_, n)) = monomial_parts(arena, g, var, var_sym) else {
return log_of_linear_substitution(arena, expr, ln_node, g, w, var, var_sym);
};
let x_power = |arena: &mut Arena, m: ExprId| -> ExprId {
let mw = arena.mul(&[m, w]);
let mwn = arena.div(mw, n);
arena.exp(mwn)
};
let mut sub = arena.subs_structural(expr, ln_node, w);
for id in crate::base::walk::post_order_ids(arena, sub) {
if let ExprNode::Pow(b, m) = *arena.node(id)
&& b == var
&& !contains_var(arena, m, var_sym)
{
let r = x_power(arena, m);
sub = arena.subs_structural(sub, id, r);
}
}
let one = arena.one;
let x_w = x_power(arena, one);
sub = arena.subs_structural(sub, var, x_w);
if contains_var(arena, sub, var_sym) {
return None;
}
let integrand = arena.mul(&[sub, x_w]);
let integrand = arena.div(integrand, n);
let integrand = crate::transforms::eval::eval(arena, integrand);
let res = integrate_nested(arena, integrand, w)?;
let mut back = res;
for id in crate::base::walk::post_order_ids(arena, res) {
if let ExprNode::Exp(arg) = *arena.node(id)
&& crate::base::walk::has_free_symbol(arena, arg, w_sym)
{
let (q, r) = symbolic_linear_coeff_of(arena, arg, w, w_sym)?;
let qn = arena.mul(&[q, n]);
let x_qn = arena.pow(var, qn);
let e_r = arena.exp(r);
let repl = arena.mul(&[e_r, x_qn]);
back = arena.subs_structural(back, id, repl);
}
}
let back = crate::transforms::eval::eval(arena, back);
let back = arena.subs_structural(back, w, ln_node);
let back = crate::transforms::eval::eval(arena, back);
if candidate_rejected(arena, expr, back, var, var_sym) {
tracing::debug!("integrate: rejecting unverified logarithmic-substitution closed form");
return None;
}
Some(back)
}
fn log_of_linear_substitution(
arena: &mut Arena,
expr: ExprId,
ln_node: ExprId,
g: ExprId,
w: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let ExprNode::Symbol(w_sym) = *arena.node(w) else {
return None;
};
let (alpha, beta) = symbolic_linear_coeff_of(arena, g, var, var_sym)?;
if arena.is_zero_structural(beta) {
return None;
}
let sub = arena.subs_structural(expr, ln_node, w);
if contains_var(arena, sub, var_sym) {
return None;
}
let e_w = arena.exp(w);
let integrand = arena.mul(&[sub, e_w]);
let integrand = arena.div(integrand, alpha);
let integrand = crate::transforms::eval::eval(arena, integrand);
let res = integrate_nested(arena, integrand, w)?;
let mut back = res;
for id in crate::base::walk::post_order_ids(arena, res) {
if let ExprNode::Exp(arg) = *arena.node(id)
&& crate::base::walk::has_free_symbol(arena, arg, w_sym)
{
let (q, r) = symbolic_linear_coeff_of(arena, arg, w, w_sym)?;
let g_q = arena.pow(g, q);
let e_r = arena.exp(r);
let repl = arena.mul(&[e_r, g_q]);
back = arena.subs_structural(back, id, repl);
}
}
let back = crate::transforms::eval::eval(arena, back);
let back = arena.subs_structural(back, w, ln_node);
let back = crate::transforms::eval::eval(arena, back);
if candidate_rejected(arena, expr, back, var, var_sym) {
tracing::debug!("integrate: rejecting unverified logarithmic-substitution closed form");
return None;
}
Some(back)
}
fn collect_exp_multiples(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<Vec<(ExprId, Q)>> {
let order = crate::base::walk::post_order_ids(arena, expr);
let mut out = Vec::new();
for id in order {
if let ExprNode::Exp(arg) = arena.node(id).clone()
&& contains_var(arena, arg, var_sym)
{
let (alpha, beta) = symbolic_linear_coeff_of(arena, arg, var, var_sym)?;
if !arena.is_zero_structural(beta) {
return None;
}
let k = arena.as_num(alpha)?.clone();
out.push((id, k));
}
}
Some(out)
}
fn try_exp_rational_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
use num_integer::Integer;
let exps = collect_exp_multiples(arena, expr, var, var_sym)?;
if exps.is_empty() {
return None;
}
let mut num_gcd = num_bigint::BigInt::zero();
let mut den_lcm = num_bigint::BigInt::one();
for (_, k) in &exps {
num_gcd = num_gcd.gcd(k.numer());
den_lcm = den_lcm.lcm(k.denom());
}
if num_gcd.is_zero() {
return None;
}
let a = num_rational::Ratio::new(num_gcd, den_lcm);
let u = arena.symbol("__eu");
let u_sym = match arena.node(u) {
ExprNode::Symbol(s) => *s,
_ => return None,
};
let mut sub = expr;
for (node, k) in &exps {
let power = k / &a; let power_id = arena.num_ratio(power.clone());
let u_pow = arena.pow(u, power_id);
sub = arena.subs_structural(sub, *node, u_pow);
}
if contains_var(arena, sub, var_sym) {
return None;
}
let a_id = arena.num_ratio(a.clone());
let a_u = arena.mul(&[a_id, u]);
let integrand_u = arena.div(sub, a_u);
let integrand_u = clear_nested_fractions(arena, integrand_u, u);
let (n, d) = crate::poly::polybridge::as_numer_denom(arena, integrand_u);
if crate::poly::polybridge::expr_to_poly(arena, n, u).is_none()
|| crate::poly::polybridge::expr_to_poly(arena, d, u).is_none()
{
return None;
}
let res_u = integrate_node(arena, integrand_u, u, u_sym, 20);
if crate::base::walk::has_unevaluated(arena, res_u) {
return None;
}
let ax = arena.mul(&[a_id, var]);
let e_ax = arena.exp(ax);
let result = arena.subs_structural(res_u, u, e_ax);
Some(crate::transforms::eval::eval(arena, result))
}
fn try_exp_linear_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
use num_integer::Integer;
let mut exps: Vec<(ExprId, ExprId, ExprId)> = Vec::new();
for id in crate::base::walk::post_order_ids(arena, expr) {
if let ExprNode::Exp(arg) = *arena.node(id)
&& contains_var(arena, arg, var_sym)
{
let (alpha, beta) = symbolic_linear_coeff_of(arena, arg, var, var_sym)?;
exps.push((id, alpha, beta));
}
}
let &(_, alpha0, beta0) = exps.first()?;
let mut ratios: Vec<Q> = Vec::with_capacity(exps.len());
for &(_, alpha, _) in &exps {
let r = arena.div(alpha, alpha0);
let r = crate::transforms::eval::eval(arena, r);
ratios.push(arena.as_num(r)?.clone());
}
let mut num_gcd = num_bigint::BigInt::zero();
let mut den_lcm = num_bigint::BigInt::one();
for r in &ratios {
num_gcd = num_gcd.gcd(r.numer());
den_lcm = den_lcm.lcm(r.denom());
}
if num_gcd.is_zero() {
return None;
}
let g = Q::new(num_gcd, den_lcm);
let u = arena.symbol("__elu");
let ExprNode::Symbol(u_sym) = *arena.node(u) else {
return None;
};
let mut constants: Vec<(ExprId, ExprId)> = Vec::new();
let mut sub = expr;
for (&(node, _, beta), r) in exps.iter().zip(&ratios) {
let r_id = arena.num_ratio(r.clone());
let r_beta0 = arena.mul(&[r_id, beta0]);
let gamma = arena.sub(beta, r_beta0);
let gamma = crate::transforms::eval::eval(arena, gamma);
let k = arena.num_ratio(r / &g);
let u_k = arena.pow(u, k);
let repl = if arena.is_zero_structural(gamma) {
u_k
} else {
let e_gamma = arena.exp(gamma);
let kappa = match constants.iter().find(|&&(e, _)| e == e_gamma) {
Some(&(_, s)) => s,
None => {
let s = arena.symbol(&format!("__elk{}", constants.len()));
constants.push((e_gamma, s));
s
}
};
arena.mul(&[kappa, u_k])
};
sub = arena.subs_structural(sub, node, repl);
}
if contains_var(arena, sub, var_sym) {
return None;
}
let g_id = arena.num_ratio(g.clone());
let du = arena.mul(&[g_id, alpha0, u]);
let integrand = arena.div(sub, du);
let integrand = crate::transforms::eval::eval(arena, integrand);
if !is_rational_in(arena, integrand, u_sym) {
return None;
}
let res = integrate_nested_generic(arena, integrand, u)?;
let res = analytic_logs(arena, res, u_sym)?;
let l = arena.mul(&[alpha0, var]);
let l = arena.add(&[l, beta0]);
let gl = arena.mul(&[g_id, l]);
let e_gl = arena.exp(gl);
let mut back = arena.subs_structural(res, u, e_gl);
for &(e_gamma, kappa) in &constants {
back = arena.subs_structural(back, kappa, e_gamma);
}
let back = crate::transforms::eval::eval(arena, back);
if candidate_rejected(arena, expr, back, var, var_sym) {
tracing::debug!("integrate: rejecting unverified exponential-substitution closed form");
return None;
}
Some(back)
}
fn try_even_power_radicals(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let two = num_bigint::BigInt::from(2);
let mut out = expr;
let mut constant: Vec<ExprId> = Vec::new();
let mut odd = false;
for id in crate::base::walk::post_order_ids(arena, expr) {
let ExprNode::Pow(base, ex) = *arena.node(id) else {
continue;
};
let Some(r) = arena.as_num(ex).cloned() else {
continue;
};
if r.is_integer() || !contains_var(arena, base, var_sym) {
continue;
}
let factors: Vec<ExprId> = match arena.node(base) {
ExprNode::Mul(children) => children.to_vec(),
_ => vec![base],
};
let mut power: Option<Q> = None;
let mut kappa: Vec<ExprId> = Vec::new();
for f in factors {
if !contains_var(arena, f, var_sym) {
kappa.push(f);
continue;
}
match *arena.node(f) {
ExprNode::Pow(b, e) if b == var && power.is_none() => {
power = arena.as_num(e).cloned();
}
_ => {
power = None;
break;
}
}
}
let Some(n) = power else {
continue;
};
if !n.is_integer() || n.is_negative() || !(n.numer() % &two).is_zero() {
continue;
}
let m = &n * &r;
if !m.is_integer() {
continue;
}
odd |= !(m.numer() % &two).is_zero();
let m_id = arena.num_ratio(m.clone());
let x_m = arena.pow(var, m_id);
let minus_m = arena.num_ratio(-m);
let x_minus_m = arena.pow(var, minus_m);
constant.push(arena.mul(&[id, x_minus_m]));
out = arena.subs_structural(out, id, x_m);
}
if constant.is_empty() {
return None;
}
let out = crate::transforms::eval::eval(arena, out);
tracing::debug!(rewritten = %arena.display(out), "integrate: radicals of even powers of the variable");
let mut g = integrate_nested(arena, out, var)?;
if odd {
let at_zero = arena.subs_structural(g, var, arena.zero);
let at_zero = crate::transforms::eval::eval(arena, at_zero);
if contains_non_finite(arena, at_zero) || contains_var(arena, at_zero, var_sym) {
return None;
}
let minus = arena.neg(at_zero);
g = arena.add(&[g, minus]);
}
constant.push(g);
let res = arena.mul(&constant);
if candidate_rejected(arena, expr, res, var, var_sym) {
return None;
}
Some(res)
}
fn try_radical_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
use num_integer::Integer;
let order = crate::base::walk::post_order_ids(arena, expr);
let mut radicals: Vec<(ExprId, Q)> = Vec::new();
let mut q_lcm = num_bigint::BigInt::one();
for id in order {
if let ExprNode::Pow(base, e) = arena.node(id).clone()
&& base == var
&& let Some(r) = arena.as_num(e).cloned()
&& !r.is_integer()
{
q_lcm = q_lcm.lcm(r.denom());
radicals.push((id, r));
}
}
if radicals.is_empty() {
return None;
}
let q: i64 = q_lcm.to_string().parse().ok()?;
if !(2..=6).contains(&q) {
return None;
}
let s = arena.symbol("__rs");
let mut sub = expr;
for (node, r) in &radicals {
let k = r * num_rational::Ratio::from_integer(num_bigint::BigInt::from(q));
let k_id = arena.num_ratio(k.clone());
let s_pow = arena.pow(s, k_id);
sub = arena.subs_structural(sub, *node, s_pow);
}
let q_id = arena.int(q);
let s_q = arena.pow(s, q_id);
sub = arena.subs_structural(sub, var, s_q);
if contains_var(arena, sub, var_sym) {
return None;
}
let s_sym = match arena.node(s) {
ExprNode::Symbol(sid) => *sid,
_ => return None,
};
let real_only = crate::base::walk::post_order_ids(arena, sub)
.into_iter()
.any(|id| match arena.node(id) {
ExprNode::Abs(g) | ExprNode::Sign(g) | ExprNode::Heaviside(g) => {
contains_var(arena, *g, s_sym)
}
_ => false,
});
if real_only {
return None;
}
let qm1 = arena.int(q - 1);
let s_qm1 = arena.pow(s, qm1);
let integrand_s = arena.mul(&[q_id, s_qm1, sub]);
let integrand_s = crate::transforms::eval::eval(arena, integrand_s);
let res_s = integrate_nested(arena, integrand_s, s)?;
let res_s = analytic_logs(arena, res_s, s_sym)?;
let inv_q = arena.rational(1, q);
let root = arena.pow(var, inv_q);
let result = arena.subs_structural(res_s, s, root);
let result = crate::transforms::eval::eval(arena, result);
if candidate_rejected(arena, expr, result, var, var_sym) {
tracing::debug!("integrate: rejecting unverified radical-substitution closed form");
return None;
}
Some(result)
}
fn analytic_logs(arena: &mut Arena, e: ExprId, s_sym: SymbolId) -> Option<ExprId> {
let order = crate::base::walk::post_order_ids(arena, e);
let mut out = e;
for id in order {
if let ExprNode::Ln(inner) = arena.node(id).clone()
&& let ExprNode::Abs(g) = arena.node(inner).clone()
&& contains_var(arena, g, s_sym)
{
let plain = arena.ln(g);
out = arena.subs_structural(out, id, plain);
}
}
let leftover = crate::base::walk::post_order_ids(arena, out)
.into_iter()
.any(|id| match arena.node(id) {
ExprNode::Abs(g) | ExprNode::Sign(g) | ExprNode::Heaviside(g) => {
contains_var(arena, *g, s_sym)
}
_ => false,
});
(!leftover).then_some(out)
}
const MAX_MOBIUS_ROOT: i64 = 6;
struct Mobius {
a: ExprId,
b: ExprId,
c: ExprId,
d: ExprId,
}
fn mobius_parts(arena: &mut Arena, m: ExprId, var: ExprId, var_sym: SymbolId) -> Option<Mobius> {
let m = clear_nested_fractions(arena, m, var);
let (numer, denom) = crate::poly::polybridge::as_numer_denom(arena, m);
let mut linear = |e: ExprId| -> Option<(ExprId, ExprId)> {
if contains_var(arena, e, var_sym) {
symbolic_linear_coeff_of(arena, e, var, var_sym)
} else {
Some((arena.zero, e))
}
};
let (a, b) = linear(numer)?;
let (c, d) = linear(denom)?;
if contains_var(arena, a, var_sym)
|| contains_var(arena, b, var_sym)
|| contains_var(arena, c, var_sym)
|| contains_var(arena, d, var_sym)
{
return None;
}
Some(Mobius { a, b, c, d })
}
fn radical_powers(
arena: &Arena,
expr: ExprId,
var_sym: SymbolId,
) -> Option<Vec<(ExprId, ExprId, Q)>> {
let mut out = Vec::new();
for id in crate::base::walk::post_order_ids(arena, expr) {
if let ExprNode::Pow(base, e) = *arena.node(id)
&& contains_var(arena, base, var_sym)
{
let r = arena.as_num(e)?;
if !r.is_integer() {
out.push((id, base, r.clone()));
}
} else if let ExprNode::Pow(_, e) = *arena.node(id)
&& contains_var(arena, e, var_sym)
{
return None;
}
}
Some(out)
}
const FACTORED_FORM_DEGREE: i64 = 10;
fn expanded_degree_in(arena: &Arena, p: ExprId, v: ExprId) -> i64 {
let mut deg = 0;
for id in crate::base::walk::post_order_ids(arena, p) {
if id == v {
deg = deg.max(1);
} else if let ExprNode::Pow(b, k) = *arena.node(id)
&& b == v
&& let Some(k) = arena
.as_num(k)
.and_then(|k| i64::try_from(k.to_integer()).ok())
{
deg = deg.max(k);
}
}
deg
}
fn rational_in_factored_form(arena: &mut Arena, e: ExprId, v: ExprId) -> ExprId {
let e = crate::transforms::eval::eval(arena, e);
let cleared = clear_nested_fractions(arena, e, v);
let has_parameters = crate::base::walk::free_symbols(arena, cleared)
.into_iter()
.any(|s| s != v);
let (_, d) = crate::poly::polybridge::as_numer_denom(arena, cleared);
if !has_parameters || expanded_degree_in(arena, d, v) <= FACTORED_FORM_DEGREE {
return cleared;
}
let (n, d) = crate::poly::polybridge::fraction_parts(arena, e);
let n = crate::transforms::expand::expand(arena, n);
let n = crate::transforms::eval::eval(arena, n);
let factors: SmallVec<[ExprId; 6]> = match arena.node(d) {
ExprNode::Mul(children) => children.clone(),
_ => SmallVec::from_slice(&[d]),
};
let mut expanded: SmallVec<[ExprId; 6]> = SmallVec::new();
for f in factors {
let (base, k) = match *arena.node(f) {
ExprNode::Pow(b, k) => (b, Some(k)),
_ => (f, None),
};
let base = crate::transforms::expand::expand(arena, base);
let base = crate::transforms::eval::eval(arena, base);
expanded.push(match k {
Some(k) => arena.pow(base, k),
None => base,
});
}
let d = arena.mul(&expanded);
arena.div(n, d)
}
fn is_rational_in(arena: &Arena, e: ExprId, v: SymbolId) -> bool {
let mut depends: FxHashMap<ExprId, bool> = FxHashMap::default();
for id in crate::base::walk::post_order_ids(arena, e) {
let node = arena.node(id);
let mut below = matches!(node, ExprNode::Symbol(s) if *s == v);
node.for_each_child(|c| below |= depends.get(&c).copied().unwrap_or(false));
if below {
let ok = match node {
ExprNode::Symbol(_) | ExprNode::Add(_) | ExprNode::Mul(_) | ExprNode::Neg(_) => {
true
}
ExprNode::Pow(_, ex) => arena.as_num(*ex).is_some_and(|r| r.is_integer()),
_ => false,
};
if !ok {
return false;
}
}
depends.insert(id, below);
}
true
}
fn pair_radicals(
arena: &mut Arena,
expr: ExprId,
b1: ExprId,
b2: ExprId,
t: ExprId,
q: i64,
) -> Option<ExprId> {
let radical_of = |arena: &Arena, id: ExprId, base: ExprId| -> Option<(ExprId, Q)> {
match *arena.node(id) {
ExprNode::Pow(b, e) if b == base => {
let r = arena.as_num(e)?;
(!r.is_integer()).then(|| (id, r.clone()))
}
_ => None,
}
};
let q_ratio = Q::from_integer(q.into());
let mut cur = expr;
for _ in 0..64 {
let mut hit = None;
for id in crate::base::walk::post_order_ids(arena, cur) {
if let ExprNode::Mul(children) = arena.node(id) {
let p1 = children.iter().find_map(|&c| radical_of(arena, c, b1));
let p2 = children.iter().find_map(|&c| radical_of(arena, c, b2));
if let (Some(p1), Some(p2)) = (p1, p2) {
hit = Some((id, children.clone(), p1, p2));
break;
}
}
}
let Some((id, children, (n1, e1), (n2, e2))) = hit else {
break;
};
let sum = &e1 + &e2;
if !sum.is_integer() {
return None;
}
let mut rest: SmallVec<[ExprId; 6]> = children
.iter()
.copied()
.filter(|&c| c != n1 && c != n2)
.collect();
let k = arena.num_ratio(&e1 * &q_ratio);
rest.push(arena.pow(t, k));
let s = arena.num_ratio(sum);
rest.push(arena.pow(b2, s));
let replacement = arena.mul(&rest);
cur = arena.subs_structural(cur, id, replacement);
}
let leftover = crate::base::walk::post_order_ids(arena, cur)
.into_iter()
.any(|id| radical_of(arena, id, b1).is_some() || radical_of(arena, id, b2).is_some());
(!leftover).then_some(cur)
}
fn try_mobius_radical_substitution(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let back = match power_substituted(arena, expr, var, var_sym) {
Some((f_w, w, n)) => {
let ExprNode::Symbol(w_sym) = *arena.node(w) else {
return None;
};
let g = mobius_core(arena, f_w, w, w_sym, true)?;
unmap_power(arena, g, w, var, n)
}
None => match mobius_core(arena, expr, var, var_sym, false) {
Some(g) => g,
None => mobius_core_split(arena, expr, var, var_sym)?,
},
};
let back = crate::transforms::eval::eval(arena, back);
if candidate_rejected(arena, expr, back, var, var_sym) {
tracing::debug!("integrate: rejecting unverified Möbius-radical closed form");
return None;
}
Some(back)
}
fn mobius_core_split(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let radicals = radical_powers(arena, expr, var_sym)?;
if radicals.is_empty() {
return None;
}
let mut replaced: Vec<(ExprId, ExprId)> = Vec::new();
let mut constants: Vec<(ExprId, ExprId)> = Vec::new();
let mut double_roots: Vec<(ExprId, ExprId)> = Vec::new();
let two = arena.int(2);
for (node, base, e) in &radicals {
if replaced.iter().any(|&(n, _)| n == *node) {
continue;
}
let mut factors: SmallVec<[(ExprId, i64); 2]> = SmallVec::new();
let mut double: Option<ExprId> = None;
if let Some((alpha, beta)) = symbolic_linear_coeff_of(arena, *base, var, var_sym) {
let r = arena.div(beta, alpha);
let r = arena.neg(r);
let r = crate::transforms::eval::eval(arena, r);
factors.push((arena.sub(var, r), 1));
} else if let Some((a2, a1, a0)) = symbolic_quadratic_coeffs(arena, *base, var, var_sym) {
let a1_sq = arena.pow(a1, two);
let four = arena.int(4);
let prod = arena.mul(&[four, a2, a0]);
let disc = arena.sub(a1_sq, prod);
let disc = crate::transforms::expand::expand(arena, disc);
let disc = crate::transforms::eval::eval(arena, disc);
let s = crate::calculus::risch::polynomial_sqrt(arena, disc)?;
let two_a2 = arena.mul(&[two, a2]);
let neg_a1 = arena.neg(a1);
if arena.is_zero_structural(s) {
let r = arena.div(neg_a1, two_a2);
let r = crate::transforms::eval::eval(arena, r);
factors.push((arena.sub(var, r), 2));
double = Some(r);
} else {
for sign_s in [arena.neg(s), s] {
let num = arena.add(&[neg_a1, sign_s]);
let r = arena.div(num, two_a2);
let r = crate::transforms::eval::eval(arena, r);
factors.push((arena.sub(var, r), 1));
}
}
} else {
return None;
}
let e_id = arena.num_ratio(e.clone());
let mut k_parts: SmallVec<[ExprId; 3]> = SmallVec::new();
k_parts.push(arena.pow(*base, e_id));
let mut new_parts: SmallVec<[ExprId; 3]> = SmallVec::new();
for &(l, mult) in &factors {
let em = e * Q::from_integer(mult.into());
let em_id = arena.num_ratio(em.clone());
new_parts.push(arena.pow(l, em_id));
let neg_em = arena.num_ratio(-em);
k_parts.push(arena.pow(l, neg_em));
}
let k = arena.mul(&k_parts);
let k = crate::transforms::eval::eval(arena, k);
if k != arena.one {
let kappa = arena.symbol(&format!("__mk{}", constants.len()));
constants.push((kappa, k));
if let Some(r) = double {
double_roots.push((kappa, r));
}
new_parts.push(kappa);
} else if double.is_some() {
return None;
}
replaced.push((*node, arena.mul(&new_parts)));
}
if constants.is_empty() {
return None;
}
let rewritten = arena.subs_map_structural(expr, &replaced);
let rewritten = crate::transforms::eval::eval(arena, rewritten);
let mut g = mobius_core(arena, rewritten, var, var_sym, false)?;
for &(kappa, r) in &double_roots {
let g1 = crate::transforms::diff::diff(arena, g, kappa);
let at = arena.subs_structural(g1, var, r);
let at = crate::transforms::eval::eval(arena, at);
if contains_var(arena, at, var_sym)
|| contains_non_finite(arena, at)
|| crate::base::walk::has_unevaluated(arena, at)
{
return None;
}
let shift = arena.mul(&[kappa, at]);
g = arena.sub(g, shift);
}
for &(kappa, k) in &constants {
g = arena.subs_structural(g, kappa, k);
}
Some(g)
}
fn power_substituted(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<(ExprId, ExprId, i64)> {
use num_integer::Integer;
let monomial_exponent = |arena: &Arena, f: ExprId| -> Option<Q> {
if f == var {
return Some(Q::one());
}
match *arena.node(f) {
ExprNode::Pow(b, e) if b == var => arena.as_num(e).cloned(),
_ => None,
}
};
let (m, rest): (Q, SmallVec<[ExprId; 6]>) = match arena.node(expr) {
ExprNode::Mul(children) => {
let children = children.clone();
let mut m = Q::zero();
let mut rest = SmallVec::new();
for c in children {
match monomial_exponent(arena, c) {
Some(k) => m += k,
None => rest.push(c),
}
}
(m, rest)
}
_ => (Q::zero(), SmallVec::from_slice(&[expr])),
};
let mut g = num_bigint::BigInt::zero();
let mut found = false;
for &r in &rest {
for id in crate::base::walk::post_order_ids(arena, r) {
let node = arena.node(id);
if let ExprNode::Pow(b, e) = *node
&& b == var
{
let k = arena.as_num(e)?;
if !k.is_integer() {
return None;
}
g = g.gcd(&k.to_integer());
found = true;
continue;
}
let mut bare = false;
node.for_each_child(|c| bare |= c == var);
if bare {
return None;
}
}
}
let n = i64::try_from(&g).ok()?;
if !found || !(2..=12).contains(&n) {
return None;
}
let w = arena.symbol("__mw");
let n_q = Q::from_integer(n.into());
let mut pairs: Vec<(ExprId, ExprId)> = Vec::new();
for &r in &rest {
for id in crate::base::walk::post_order_ids(arena, r) {
if let ExprNode::Pow(b, e) = *arena.node(id)
&& b == var
&& let Some(k) = arena.as_num(e).cloned()
&& !pairs.iter().any(|&(old, _)| old == id)
{
let k_w = arena.num_ratio(k / &n_q);
let w_k = arena.pow(w, k_w);
pairs.push((id, w_k));
}
}
}
let mut factors: SmallVec<[ExprId; 6]> = SmallVec::new();
for &r in &rest {
factors.push(arena.subs_map_structural(r, &pairs));
}
let e_w = arena.num_ratio((&m + Q::one() - &n_q) / &n_q);
factors.push(arena.pow(w, e_w));
factors.push(arena.rational(1, n));
let f_w = arena.mul(&factors);
if contains_var(arena, f_w, var_sym) {
return None;
}
Some((f_w, w, n))
}
fn unmap_power(arena: &mut Arena, e: ExprId, w: ExprId, var: ExprId, n: i64) -> ExprId {
let n_q = Q::from_integer(n.into());
let mut pairs: Vec<(ExprId, ExprId)> = Vec::new();
for id in crate::base::walk::post_order_ids(arena, e) {
if let ExprNode::Pow(b, r) = *arena.node(id)
&& b == w
&& let Some(r) = arena.as_num(r).cloned()
&& !pairs.iter().any(|&(old, _)| old == id)
{
let nr = arena.num_ratio(r * &n_q);
let x_nr = arena.pow(var, nr);
pairs.push((id, x_nr));
}
}
let e = arena.subs_map_structural(e, &pairs);
let n_id = arena.int(n);
let x_n = arena.pow(var, n_id);
arena.subs_structural(e, w, x_n)
}
fn mobius_core(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
w_first: bool,
) -> Option<ExprId> {
use num_integer::Integer;
let radicals = radical_powers(arena, expr, var_sym)?;
let mut bases: SmallVec<[ExprId; 2]> = SmallVec::new();
let mut q_lcm = num_bigint::BigInt::one();
for (_, b, e) in &radicals {
if !bases.contains(b) {
bases.push(*b);
}
q_lcm = q_lcm.lcm(e.denom());
}
if bases.is_empty() || bases.len() > 2 {
return None;
}
let q = i64::try_from(&q_lcm).ok()?;
if !(2..=MAX_MOBIUS_ROOT).contains(&q) {
return None;
}
let t = arena.symbol("__mt");
let ExprNode::Symbol(t_sym) = *arena.node(t) else {
return None;
};
let q_ratio = Q::from_integer(q.into());
let inv_q = arena.rational(1, q);
let (ratio, sub, t_of_x) = if let [b] = bases[..] {
let mut sub = expr;
for (node, _, e) in &radicals {
let k = arena.num_ratio(e * &q_ratio);
let tk = arena.pow(t, k);
sub = arena.subs_structural(sub, *node, tk);
}
(b, sub, arena.pow(b, inv_q))
} else {
let (b1, b2) = if w_first && bases[1] == var {
(bases[1], bases[0])
} else {
(bases[0], bases[1])
};
let sub = pair_radicals(arena, expr, b1, b2, t, q)?;
let neg_inv_q = arena.rational(-1, q);
let r1 = arena.pow(b1, inv_q);
let r2 = arena.pow(b2, neg_inv_q);
(arena.div(b1, b2), sub, arena.mul(&[r1, r2]))
};
let Mobius { a, b, c, d } = mobius_parts(arena, ratio, var, var_sym)?;
let ad = arena.mul(&[a, d]);
let bc = arena.mul(&[b, c]);
let det = arena.sub(ad, bc);
let det = crate::transforms::eval::eval(arena, det);
if arena.is_zero_structural(det) {
return None;
}
let q_id = arena.int(q);
let t_q = arena.pow(t, q_id);
let d_tq = arena.mul(&[d, t_q]);
let x_num = arena.sub(d_tq, b);
let c_tq = arena.mul(&[c, t_q]);
let x_den = arena.sub(a, c_tq);
let x_of_t = arena.div(x_num, x_den);
let qm1 = arena.int(q - 1);
let t_qm1 = arena.pow(t, qm1);
let two = arena.int(2);
let den_sq = arena.pow(x_den, two);
let dx_num = arena.mul(&[q_id, t_qm1, det]);
let dx_dt = arena.div(dx_num, den_sq);
let sub = arena.subs_structural(sub, var, x_of_t);
if contains_var(arena, sub, var_sym) {
return None;
}
let integrand = arena.mul(&[sub, dx_dt]);
let integrand = rational_in_factored_form(arena, integrand, t);
if !is_rational_in(arena, integrand, t_sym) {
return None;
}
tracing::debug!(integrand_t = %arena.display(integrand), "integrate: Möbius radical substitution");
let res = integrate_nested_generic(arena, integrand, t)?;
let res = analytic_logs(arena, res, t_sym)?;
Some(arena.subs_structural(res, t, t_of_x))
}
fn try_exp_atanh_rewrite(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let mut rewritten = expr;
for id in crate::base::walk::post_order_ids(arena, expr) {
let ExprNode::Exp(arg) = *arena.node(id) else {
continue;
};
if !contains_var(arena, arg, var_sym) {
continue;
}
let (k, u) = match arena.node(arg).clone() {
ExprNode::Atanh(u) => (Q::one(), u),
ExprNode::Mul(fs) if fs.len() == 2 => {
let (k, u) = match (arena.as_num(fs[0]), arena.node(fs[1])) {
(Some(k), ExprNode::Atanh(u)) => (k.clone(), *u),
_ => continue,
};
(k, u)
}
_ => continue,
};
let half_k = arena.num_ratio(&k / Q::from_integer(2.into()));
let neg_half_k = arena.num_ratio(-&k / Q::from_integer(2.into()));
let one = arena.one;
let one_plus = arena.add(&[one, u]);
let one_minus = arena.sub(one, u);
let p1 = arena.pow(one_plus, half_k);
let p2 = arena.pow(one_minus, neg_half_k);
let repl = arena.mul(&[p1, p2]);
rewritten = arena.subs_structural(rewritten, id, repl);
}
if rewritten == expr {
return None;
}
let rewritten = crate::transforms::eval::eval(arena, rewritten);
let res = integrate_nested(arena, rewritten, var)?;
if candidate_rejected(arena, expr, res, var, var_sym) {
tracing::debug!("integrate: rejecting unverified e^(k·atanh u) closed form");
return None;
}
Some(res)
}
fn try_hyperbolic_to_exp(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let order = crate::base::walk::post_order_ids(arena, expr);
let has_hyp = order.iter().any(|&id| {
matches!(
arena.node(id),
ExprNode::Sinh(_) | ExprNode::Cosh(_) | ExprNode::Tanh(_)
) && contains_var(arena, id, var_sym)
});
if !has_hyp {
return None;
}
let rewritten = arena.rewrite_as_exp_expr(expr);
if rewritten == expr {
return None;
}
let rewritten = crate::transforms::eval::eval(arena, rewritten);
try_exp_rational_substitution(arena, rewritten, var, var_sym).or_else(|| {
let r = try_exp_linear_substitution(arena, rewritten, var, var_sym)?;
(!candidate_rejected(arena, expr, r, var, var_sym)).then_some(r)
})
}
fn split_poly_and_half_power(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
) -> Option<(crate::poly::Poly, crate::poly::Poly, ExprId, Q)> {
let mut half: Option<(ExprId, Q)> = None;
let mut poly = crate::poly::Poly::from_int(1);
for &d in dependent {
if let ExprNode::Pow(base, e) = arena.node(d).clone()
&& let Some(r) = arena.as_num(e).cloned()
&& *r.denom() == num_bigint::BigInt::from(2)
{
if half.is_some() {
return None;
}
half = Some((base, r));
continue;
}
let p = crate::poly::polybridge::expr_to_poly(arena, d, var)?;
poly = poly.mul(&p);
}
let (q_expr, k) = half?;
let q = crate::poly::polybridge::expr_to_poly(arena, q_expr, var)?;
let qd = q.degree()?;
if !(1..=2).contains(&qd) {
return None;
}
Some((poly, q, q_expr, k))
}
fn try_poly_times_half_power(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let (p, q, q_expr, k) = split_poly_and_half_power(arena, dependent, var)?;
let two = Q::from_integer(2.into());
let k2 = &k * &two; let k2: i64 = k2.to_integer().to_string().parse().ok()?;
if !matches!(k2, -1 | 1 | 3) {
return None;
}
let mut p_tilde = p;
for _ in 0..((k2 + 1) / 2) {
p_tilde = p_tilde.mul(&q);
}
let pd = p_tilde.degree()?;
let qd = q.degree()?;
if pd == 0 && qd == 2 {
return None; }
let m: i64 = if qd == 2 { pd as i64 - 1 } else { pd as i64 };
let with_c = qd == 2;
let n_unknowns = (m + 1).max(0) as usize + usize::from(with_c);
let n_eq = pd + 1;
if n_unknowns != n_eq {
return None;
}
let q_prime = q.derivative();
let mut rows: Vec<Vec<Q>> = vec![vec![Q::zero(); n_unknowns + 1]; n_eq];
for i in 0..=(m.max(-1)) {
if i < 0 {
break;
}
let iu = i as usize;
let mut mono = vec![Q::zero(); iu + 1];
mono[iu] = Q::one();
let r_i = crate::poly::Poly::from_coeffs(mono);
let term = r_i
.derivative()
.mul(&q)
.scale(&two)
.add(&r_i.mul(&q_prime))
.scale(&Q::new(1.into(), 2.into()));
for (j, row) in rows.iter_mut().enumerate() {
row[iu] = term.coeff(j);
}
}
if with_c {
rows[0][n_unknowns - 1] = Q::one();
}
for (j, row) in rows.iter_mut().enumerate() {
row[n_unknowns] = p_tilde.coeff(j);
}
let sol = solve_linear_system(rows, n_unknowns)?;
let r_coeffs: Vec<Q> = sol[..(m.max(-1) + 1) as usize].to_vec();
let r_poly = crate::poly::Poly::from_coeffs(r_coeffs);
let r_expr = crate::poly::polybridge::poly_to_expr(arena, &r_poly, var);
let half = arena.rational(1, 2);
let sqrt_q = arena.pow(q_expr, half);
let mut result = arena.mul(&[r_expr, sqrt_q]);
if with_c {
let c = sol[n_unknowns - 1].clone();
if !c.is_zero() {
let neg_half = arena.rational(-1, 2);
let inv_sqrt = arena.pow(q_expr, neg_half);
let base_int = integrate_node(arena, inv_sqrt, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, base_int) {
return None;
}
let c_id = arena.num_ratio(c.clone());
let c_term = arena.mul(&[c_id, base_int]);
result = arena.add(&[result, c_term]);
}
}
Some(result)
}
fn solve_linear_system(mut rows: Vec<Vec<Q>>, n: usize) -> Option<Vec<Q>> {
let m = rows.len();
let mut pivot_row = 0;
let mut pivot_cols: Vec<usize> = Vec::new();
for col in 0..n {
let Some(p) = (pivot_row..m).find(|&r| !rows[r][col].is_zero()) else {
continue;
};
rows.swap(pivot_row, p);
let inv = Q::one() / rows[pivot_row][col].clone();
for v in rows[pivot_row].iter_mut() {
*v = &*v * &inv;
}
let pivot = rows[pivot_row].clone();
for (r, row) in rows.iter_mut().enumerate() {
if r != pivot_row && !row[col].is_zero() {
let f = row[col].clone();
for (c, cell) in row.iter_mut().enumerate() {
let sub = &pivot[c] * &f;
*cell = &*cell - ⊂
}
}
}
pivot_cols.push(col);
pivot_row += 1;
if pivot_row == m {
break;
}
}
for row in rows.iter().skip(pivot_row) {
if !row[n].is_zero() {
return None;
}
}
let mut sol = vec![Q::zero(); n];
for (r, &col) in pivot_cols.iter().enumerate() {
sol[col] = rows[r][n].clone();
}
Some(sol)
}
fn try_reciprocal_sqrt_substitution(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
if dependent.len() != 2 {
return None;
}
let mut n_neg: Option<i64> = None;
let mut half: Option<(ExprId, Q)> = None;
for &d in dependent {
if let ExprNode::Pow(base, e) = arena.node(d).clone()
&& let Some(r) = arena.as_num(e).cloned()
{
if base == var && r.is_integer() && r.is_negative() {
n_neg = Some(-r.to_integer().to_string().parse::<i64>().ok()?);
continue;
}
if *r.denom() == num_bigint::BigInt::from(2) && contains_var(arena, base, var_sym) {
half = Some((base, r));
continue;
}
}
return None;
}
let n = n_neg?;
let (q_expr, k) = half?;
if k != Q::new((-1).into(), 2.into()) || !(1..=4).contains(&n) {
return None;
}
let q = crate::poly::polybridge::expr_to_poly(arena, q_expr, var)?;
if q.degree()? != 2 {
return None;
}
let t = arena.symbol("__rt");
let t_sym = match arena.node(t) {
ExprNode::Symbol(s) => *s,
_ => return None,
};
let coeffs = q.coeffs().to_vec(); let reversed = crate::poly::Poly::from_coeffs(coeffs.iter().rev().cloned().collect());
let q_rev = crate::poly::polybridge::poly_to_expr(arena, &reversed, t);
let neg_half = arena.rational(-1, 2);
let q_rev_pow = arena.pow(q_rev, neg_half);
let t_pow_id = arena.int(n - 1);
let t_pow = arena.pow(t, t_pow_id);
let integrand_t = arena.mul(&[t_pow, q_rev_pow]);
let neg_integrand = arena.neg(integrand_t);
let res_t = integrate_node(arena, neg_integrand, t, t_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, res_t) {
return None;
}
let inv_x = arena.pow(var, arena.neg_one);
let g_of_x = arena.subs_structural(res_t, t, inv_x);
let sgn = arena.sign(var);
let result = arena.mul(&[sgn, g_of_x]);
Some(crate::transforms::eval::eval(arena, result))
}
fn try_tan_sec_patterns(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
if dependent.len() != 2 {
return None;
}
let mut tan_part: Option<(ExprId, Q)> = None;
let mut sec_part: Option<(ExprId, Q)> = None;
for &d in dependent {
let (base, e) = if let ExprNode::Pow(b, e) = arena.node(d).clone() {
(b, arena.as_num(e)?.clone())
} else {
(d, Q::one())
};
match arena.node(base).clone() {
ExprNode::Tan(g) if e.is_integer() && e.is_positive() => tan_part = Some((g, e)),
ExprNode::Cos(g) if e.is_integer() && e.is_negative() => sec_part = Some((g, -e)),
_ => return None,
}
}
let (g, m) = tan_part?;
let (g2, n) = sec_part?;
if g != g2 {
return None;
}
let (a_expr, _) = symbolic_linear_coeff_of(arena, g, var, var_sym)?;
let two = Q::from_integer(2.into());
if n == two {
let m1 = &m + &Q::one();
let m1_id = arena.num_ratio(m1.clone());
let tan_g = arena.tan(g);
let tp = arena.pow(tan_g, m1_id);
let denom = arena.mul(&[m1_id, a_expr]);
return Some(arena.div(tp, denom));
}
if m == Q::one() {
let cos_g = arena.cos(g);
let neg_n = arena.num_ratio(-n.clone());
let sec_n = arena.pow(cos_g, neg_n);
let n_id = arena.num_ratio(n.clone());
let denom = arena.mul(&[n_id, a_expr]);
return Some(arena.div(sec_n, denom));
}
None
}
fn try_by_parts_poly_times_pair(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let idx = dependent
.iter()
.position(|&d| is_polynomial_in(arena, d, var, var_sym))?;
let u = dependent[idx];
let dv = remaining_product(arena, dependent, idx);
let v = integrate_node(arena, dv, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, v) {
return None;
}
let du = crate::transforms::diff::diff(arena, u, var);
let v_du = arena.mul(&[v, du]);
let u_dv = arena.mul(&[u, dv]);
match by_parts_cycle(arena, u, v, u_dv, v_du, var_sym) {
Cycle::Solved(result) => return Some(result),
Cycle::Hopeless => return None,
Cycle::No => {}
}
let v_du = crate::transforms::expand::expand(arena, v_du);
let rest = integrate_node(arena, v_du, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, rest) {
return None;
}
let uv = arena.mul(&[u, v]);
Some(arena.sub(uv, rest))
}
enum Cycle {
No,
Solved(ExprId),
Hopeless,
}
fn by_parts_cycle(
arena: &mut Arena,
u: ExprId,
v: ExprId,
u_dv: ExprId,
v_du: ExprId,
var_sym: SymbolId,
) -> Cycle {
let ratio = arena.div(v_du, u_dv);
let c = crate::transforms::eval::eval(arena, ratio);
if contains_var(arena, c, var_sym)
|| contains_non_finite(arena, c)
|| arena.is_zero_structural(c)
{
return Cycle::No;
}
let one = arena.one;
let one_plus_c = arena.add(&[one, c]);
let one_plus_c = crate::transforms::eval::eval(arena, one_plus_c);
if arena.is_zero_structural(one_plus_c) {
return Cycle::Hopeless;
}
tracing::debug!("by parts: the remaining integral is a multiple of the integrand");
let uv = arena.mul(&[u, v]);
Cycle::Solved(arena.div(uv, one_plus_c))
}
fn try_by_parts_pair(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let orderings = {
let r0 = liate_rank(arena, dependent[0], var, var_sym);
let r1 = liate_rank(arena, dependent[1], var, var_sym);
tracing::debug!(u_rank = r0, dv_rank = r1, "by-parts LIATE ordering");
if r0 <= r1 {
[(0usize, 1usize), (1, 0)]
} else {
[(1, 0), (0, 1)]
}
};
for (u_idx, dv_idx) in orderings {
let u = dependent[u_idx];
let dv = dependent[dv_idx];
if !is_by_parts_candidate(arena, u, var, var_sym) {
continue;
}
let v = integrate_node(arena, dv, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, v) {
tracing::trace!("by-parts: v = ∫dv has unevaluated nodes, skipping this ordering");
continue; }
let du = crate::transforms::diff::diff(arena, u, var);
let v_du = arena.mul(&[v, du]);
let u_dv = arena.mul(&[u, dv]);
match by_parts_cycle(arena, u, v, u_dv, v_du, var_sym) {
Cycle::Solved(result) => return Some(result),
Cycle::Hopeless => continue,
Cycle::No => {}
}
let integral_v_du = integrate_node(arena, v_du, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, integral_v_du) {
tracing::trace!("by-parts: ∫v·du has unevaluated nodes, skipping this ordering");
continue; }
tracing::debug!("integration by parts succeeded");
let u_v = arena.mul(&[u, v]);
return Some(arena.sub(u_v, integral_v_du));
}
None
}
fn try_merge_polynomial_factors(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let (polys, mut rest): (Vec<ExprId>, Vec<ExprId>) = dependent
.iter()
.partition(|&&d| is_polynomial_in(arena, d, var, var_sym));
if polys.len() < 2 || rest.is_empty() {
return None;
}
let p = arena.mul(&polys);
let p = crate::transforms::expand::expand(arena, p);
rest.push(p);
let merged = arena.mul(&rest);
let r = integrate_node(arena, merged, var, var_sym, depth - 1);
(!crate::base::walk::has_unevaluated(arena, r)).then_some(r)
}
fn try_distribute_sum_factor(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
if dependent.len() < 2 || depth < 2 {
return None;
}
let radical = dependent.iter().any(|&d| {
crate::base::walk::post_order_ids(arena, d)
.iter()
.any(|&id| {
matches!(arena.node(id), ExprNode::Pow(b, e)
if !arena.as_num(*e).is_some_and(|q| q.is_integer())
&& contains_var(arena, *b, var_sym))
})
});
if radical {
return None;
}
let product = arena.mul(dependent);
let (numer, denom) = crate::poly::polybridge::as_numer_denom(arena, product);
if crate::poly::polybridge::expr_to_poly(arena, numer, var).is_some()
&& crate::poly::polybridge::expr_to_poly(arena, denom, var).is_some()
{
return None;
}
let is_sum = |arena: &Arena, d: ExprId| matches!(arena.node(d), ExprNode::Add(_));
let is_power_of_sum = |arena: &Arena, d: ExprId| {
matches!(arena.node(d), ExprNode::Pow(b, n)
if is_sum(arena, *b)
&& !is_polynomial_in(arena, *b, var, var_sym)
&& arena.as_num(*n).is_some_and(|q| {
q.is_integer() && *q >= Q::from_integer(2.into()) && *q <= Q::from_integer(4.into())
}))
};
let idx = dependent
.iter()
.position(|&d| is_sum(arena, d) && !is_polynomial_in(arena, d, var, var_sym))
.or_else(|| dependent.iter().position(|&d| is_sum(arena, d)))
.or_else(|| dependent.iter().position(|&d| is_power_of_sum(arena, d)))?;
let sum = if is_sum(arena, dependent[idx]) {
dependent[idx]
} else {
crate::transforms::expand::expand(arena, dependent[idx])
};
let ExprNode::Add(terms) = arena.node(sum).clone() else {
return None;
};
let rest = remaining_product(arena, dependent, idx);
let mut parts: SmallVec<[ExprId; 6]> = SmallVec::new();
for &term in terms.iter() {
let piece = arena.mul(&[term, rest]);
let big_f = integrate_node(arena, piece, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, big_f) {
return None;
}
parts.push(big_f);
}
Some(arena.add(&parts))
}
fn try_trig_product_to_sum(
arena: &mut Arena,
expr: ExprId,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let trig_count = dependent
.iter()
.filter(|&&d| matches!(arena.node(d), ExprNode::Sin(_) | ExprNode::Cos(_)))
.count();
if trig_count < 2 {
return None;
}
let combined = arena.trig_combine_expr(expr);
let combined = crate::transforms::eval::eval(arena, combined);
if combined == expr {
return None;
}
let combined = crate::transforms::expand::expand(arena, combined);
let result = integrate_node(arena, combined, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, result) {
return None;
}
Some(result)
}
fn real_roots(arena: &mut Arena, g: ExprId, var: ExprId) -> Option<Vec<ExprId>> {
let var_sym = match arena.node(var) {
ExprNode::Symbol(s) => *s,
_ => return None,
};
let poly = crate::poly::polybridge::expr_to_poly(arena, g, var);
if poly.is_none()
&& crate::base::walk::post_order_ids(arena, g)
.into_iter()
.any(|id| match arena.node(id) {
ExprNode::Sin(a) | ExprNode::Cos(a) => contains_var(arena, *a, var_sym),
_ => false,
})
{
return None;
}
let sols = crate::transforms::solve::solve(arena, g, var);
if sols.is_empty() && poly.is_none() {
return None;
}
let mut roots = Vec::new();
for s in sols {
if !crate::base::walk::free_symbols(arena, s.value).is_empty() {
return None; }
match crate::transforms::evalf::eval_const_f64(arena, s.value) {
Some(v) if v.is_finite() => roots.push(s.value),
Some(_) => return None,
None => {
let value = crate::transforms::eval::eval(arena, s.value);
match crate::transforms::evalf::evalf_complex64(arena, value) {
Ok(z) if z.im != 0.0 && z.re.is_finite() && z.im.is_finite() => {} _ => return None,
}
}
}
}
if let Some(p) = &poly
&& p.degree().unwrap_or(0) > 0
&& crate::poly::sturm::SturmChain::new(p).count_real_roots() != roots.len()
{
return None;
}
Some(roots)
}
fn try_poly_over_symbolic_linear(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, expr);
if !contains_var(arena, den, var_sym) {
return None;
}
let (a, b) = symbolic_linear_coeff_of(arena, den, var, var_sym)?;
if contains_var(arena, a, var_sym)
|| contains_var(arena, b, var_sym)
|| arena.is_zero_structural(a)
|| arena.as_num(a).is_some() && arena.as_num(b).is_some()
{
return None;
}
let p = crate::poly::polybridge::expr_to_poly(arena, num, var)?;
let coeffs: Vec<ExprId> = p
.coeffs()
.iter()
.map(|c| arena.num_ratio(c.clone()))
.collect();
let n = coeffs.len().checked_sub(1)?;
let neg_b = arena.neg(b);
let x0 = arena.div(neg_b, a);
let mut q = vec![arena.zero; n];
let mut carry = coeffs[n];
for k in (1..=n).rev() {
q[k - 1] = carry;
let prod = arena.mul(&[x0, carry]);
carry = arena.add(&[coeffs[k - 1], prod]);
}
let remainder = crate::transforms::eval::eval(arena, carry);
let mut terms = Vec::with_capacity(n + 1);
for (k, &qk) in q.iter().enumerate() {
let e = arena.int(k as i64 + 1);
let xp = arena.pow(var, e);
let inv = arena.rational(1, k as i64 + 1);
terms.push(arena.mul(&[inv, qk, xp]));
}
if !arena.is_zero_structural(remainder) {
let abs_den = arena.abs(den);
let log = arena.ln(abs_den);
terms.push(arena.mul(&[remainder, log]));
}
let sum = arena.add(&terms);
let inv_a = arena.pow(a, arena.neg_one());
let result = arena.mul(&[inv_a, sum]);
let result = crate::transforms::eval::eval(arena, result);
(!candidate_rejected(arena, expr, result, var, var_sym)).then_some(result)
}
fn try_poly_over_symbolic_quadratic(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
const MAX_NUMERATOR_DEGREE: usize = 16;
let (num, den) = crate::poly::polybridge::as_numer_denom(arena, expr);
if !contains_var(arena, den, var_sym) {
return None;
}
if crate::poly::polybridge::expr_to_poly(arena, den, var).is_some()
&& crate::poly::polybridge::expr_to_poly(arena, num, var).is_some()
{
return None;
}
let dense = |arena: &mut Arena, e: ExprId| -> Option<Vec<ExprId>> {
let terms = crate::poly::polybridge::symbolic_multipoly_terms(arena, e, &[var])?;
let degree = terms.iter().map(|(exps, _)| exps[0] as usize).max()?;
if degree > MAX_NUMERATOR_DEGREE {
return None;
}
let mut coeffs = vec![arena.zero; degree + 1];
for (exps, c) in terms {
if contains_var(arena, c, var_sym) {
return None;
}
coeffs[exps[0] as usize] = c;
}
Some(coeffs)
};
let d = dense(arena, den)?;
if d.len() != 3 || arena.is_zero_structural(d[2]) {
return None;
}
let mut rem = dense(arena, num)?;
if rem.len() < 3 {
return None;
}
let mut quotient = vec![arena.zero; rem.len() - 2];
for k in (2..rem.len()).rev() {
let qk = arena.div(rem[k], d[2]);
let qk = crate::transforms::eval::eval(arena, qk);
quotient[k - 2] = qk;
for (j, &dj) in d.iter().enumerate().take(2) {
let prod = arena.mul(&[qk, dj]);
let diff = arena.sub(rem[k - 2 + j], prod);
rem[k - 2 + j] = crate::transforms::eval::eval(arena, diff);
}
rem[k] = arena.zero;
}
let mut terms: Vec<ExprId> = Vec::with_capacity(quotient.len() + 1);
for (j, &qj) in quotient.iter().enumerate() {
let e = arena.int(i64::try_from(j).ok()? + 1);
let xp = arena.pow(var, e);
let inv = arena.rational(1, i64::try_from(j).ok()? + 1);
terms.push(arena.mul(&[inv, qj, xp]));
}
let r1x = arena.mul(&[rem[1], var]);
let linear = arena.add(&[r1x, rem[0]]);
let linear = crate::transforms::eval::eval(arena, linear);
if !arena.is_zero_structural(linear) {
let inv_den = arena.pow(den, arena.neg_one());
let piece = arena.mul(&[linear, inv_den]);
let piece_integral = integrate_node(arena, piece, var, var_sym, depth.saturating_sub(1));
if crate::base::walk::has_unevaluated(arena, piece_integral) {
return None;
}
terms.push(piece_integral);
}
let result = arena.add(&terms);
let result = crate::transforms::eval::eval(arena, result);
(!candidate_rejected(arena, expr, result, var, var_sym)).then_some(result)
}
fn contains_non_finite(arena: &Arena, e: ExprId) -> bool {
crate::base::walk::post_order_ids(arena, e)
.into_iter()
.any(|id| {
matches!(
arena.node(id),
ExprNode::NaN
| ExprNode::ComplexInfinity
| ExprNode::Infinity
| ExprNode::NegInfinity
)
})
}
fn real_on_reals(arena: &Arena, e: ExprId, var_sym: SymbolId) -> bool {
crate::base::walk::post_order_ids(arena, e)
.into_iter()
.all(|id| match arena.node(id) {
ExprNode::Symbol(s) => {
*s == var_sym || crate::transforms::realness::symbol_declared_real(arena, *s)
}
ExprNode::Num(_)
| ExprNode::Pi
| ExprNode::E
| ExprNode::EulerGamma
| ExprNode::Catalan
| ExprNode::GoldenRatio
| ExprNode::Add(_)
| ExprNode::Mul(_)
| ExprNode::Sin(_)
| ExprNode::Cos(_)
| ExprNode::Exp(_)
| ExprNode::Atan(_)
| ExprNode::Sinh(_)
| ExprNode::Cosh(_)
| ExprNode::Tanh(_)
| ExprNode::Abs(_) => true,
ExprNode::Pow(base, exp) => arena.as_num(*exp).is_some_and(|r| {
r.is_integer() && (!r.is_negative() || !contains_var(arena, *base, var_sym))
}),
_ => false,
})
}
fn try_abs_sign_product(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let idx = dependent.iter().position(|&d| {
matches!(
arena.node(d),
ExprNode::Abs(_) | ExprNode::Sign(_) | ExprNode::Heaviside(_)
) && contains_var(arena, d, var_sym)
})?;
let node = arena.node(dependent[idx]).clone();
let g = match node {
ExprNode::Abs(g) | ExprNode::Sign(g) | ExprNode::Heaviside(g) => g,
_ => return None,
};
if !real_on_reals(arena, g, var_sym) {
return None;
}
let rest = remaining_product(arena, dependent, idx);
let roots = real_roots(arena, g, var)?;
if roots.len() > 1 {
return None;
}
let root = roots.first().copied();
let smooth = match node {
ExprNode::Abs(_) => arena.mul(&[rest, g]),
_ => rest,
};
let smooth_int = integrate_node(arena, smooth, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, smooth_int) {
return None;
}
match root {
Some(r) => {
let at_r = crate::transforms::subs::subs(arena, smooth_int, var, r);
let at_r = crate::transforms::eval::eval(arena, at_r);
if contains_non_finite(arena, at_r) {
return None;
}
let shifted = arena.sub(smooth_int, at_r);
let factor = match node {
ExprNode::Heaviside(_) => arena.heaviside(g),
_ => arena.sign(g),
};
Some(arena.mul(&[factor, shifted]))
}
None => {
let sample = crate::transforms::subs::subs(arena, g, var, arena.zero);
let sample = crate::transforms::eval::eval(arena, sample);
let sgn = crate::transforms::evalf::eval_const_f64(arena, sample)?;
if sgn == 0.0 || !sgn.is_finite() {
return None;
}
let positive = sgn > 0.0;
match node {
ExprNode::Abs(_) | ExprNode::Sign(_) => {
if positive {
Some(smooth_int)
} else {
Some(arena.neg(smooth_int))
}
}
_ => {
if positive {
Some(smooth_int)
} else {
Some(arena.zero)
}
}
}
}
}
}
fn integrate_piecewise(
arena: &mut Arena,
pairs: &[(ExprId, ExprId)],
var: ExprId,
var_sym: SymbolId,
depth: usize,
) -> Option<ExprId> {
let mut antis: Vec<ExprId> = Vec::with_capacity(pairs.len());
for &(val, _) in pairs {
let f = integrate_node(arena, val, var, var_sym, depth - 1);
if crate::base::walk::has_unevaluated(arena, f) {
return None;
}
antis.push(f);
}
let mut breakpoints: Vec<ExprId> = Vec::new();
let mut chain = true;
for &(_, cond) in &pairs[..pairs.len().saturating_sub(1)] {
match arena.node(cond).clone() {
ExprNode::Gt(c, v) | ExprNode::Ge(c, v)
if v == var && !contains_var(arena, c, var_sym) =>
{
breakpoints.push(c);
}
_ => {
chain = false;
break;
}
}
}
let mut out: Vec<(ExprId, ExprId)> = Vec::with_capacity(pairs.len());
if chain && pairs.len() >= 2 {
let mut k = arena.zero;
out.push((antis[0], pairs[0].1));
for i in 1..pairs.len() {
let c = breakpoints[i - 1];
let prev_at_c = crate::transforms::subs::subs(arena, antis[i - 1], var, c);
let cur_at_c = crate::transforms::subs::subs(arena, antis[i], var, c);
let diff = arena.sub(prev_at_c, cur_at_c);
k = arena.add(&[k, diff]);
k = crate::transforms::eval::eval(arena, k);
let branch = arena.add(&[antis[i], k]);
out.push((branch, pairs[i].1));
}
} else {
for (i, &(_, cond)) in pairs.iter().enumerate() {
out.push((antis[i], cond));
}
}
Some(arena.piecewise(&out))
}
fn try_piecewise_wrap(
arena: &mut Arena,
result: ExprId,
original_integrand: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> ExprId {
if matches!(arena.node(result), ExprNode::Integral(_, _)) {
return result;
}
if GENERIC_ONLY.with(std::cell::Cell::get) {
return result;
}
let mut generic = result;
let mut handled: Vec<(ExprId, ExprId)> = Vec::new();
while let ExprNode::Piecewise(pairs) = arena.node(generic).clone()
&& pairs.len() == 2
&& pairs[1].1 == arena.bool_true
&& let ExprNode::Ne(p, v) = *arena.node(pairs[0].1)
{
handled.push((p, v));
handled.push((v, p));
generic = pairs[0].0;
}
let denoms = collect_denominators(arena, generic);
if denoms.is_empty() {
return result;
}
let mut wrapped = result;
let mut handled_by_denom: Vec<(ExprId, ExprId, ExprId)> = Vec::new();
let mut failed: Vec<(ExprId, ExprId)> = Vec::new();
for denom in &denoms {
let denom_syms = crate::base::walk::free_symbols(arena, *denom);
for sym_expr in &denom_syms {
if let ExprNode::Symbol(sid) = arena.node(*sym_expr)
&& *sid == var_sym
{
continue;
}
if !crate::base::walk::contains(arena, original_integrand, *sym_expr) {
continue;
}
let solutions = crate::transforms::solve::solve(arena, *denom, *sym_expr);
for sol in &solutions {
let degen_val = sol.value;
if contains_var(arena, degen_val, var_sym) || contains_non_finite(arena, degen_val)
{
continue;
}
if crate::base::walk::free_symbols(arena, degen_val)
.iter()
.any(|&s| !crate::base::walk::contains(arena, original_integrand, s))
{
continue;
}
if handled
.iter()
.any(|&(p, v)| p == *sym_expr && v == degen_val)
|| failed.contains(&(*sym_expr, degen_val))
{
continue;
}
if degenerate_value_covered(arena, &handled_by_denom, *denom, *sym_expr, degen_val)
{
continue;
}
if excluded_by_enclosing_case(arena, *sym_expr, degen_val) {
continue;
}
if nonreal_value_of_real_symbol(arena, *sym_expr, degen_val) {
continue;
}
let integrand_at_degen =
crate::transforms::subs::subs(arena, original_integrand, *sym_expr, degen_val);
let integrand_at_degen = crate::transforms::eval::eval(arena, integrand_at_degen);
if has_zero_denominator(arena, integrand_at_degen) {
continue;
}
tracing::debug!(
integrand = %arena.display(integrand_at_degen),
"try_piecewise_wrap: integrating a degenerate case"
);
let nonzero = NonzeroInCase::enter(arena, collect_denominators(arena, degen_val));
let _cap = (!real_at_generic_values(arena, degen_val))
.then(|| StepCap::enter(NONREAL_CASE_STEPS));
let degen_result = stage!(
arena,
"degenerate_case",
integrand_at_degen,
integrate(arena, integrand_at_degen, var)
);
drop(nonzero);
let degen_result = crate::transforms::eval::eval(arena, degen_result);
if crate::base::walk::has_unevaluated(arena, degen_result) {
failed.push((*sym_expr, degen_val));
continue;
}
let condition = arena.ne_(*sym_expr, degen_val);
let true_cond = arena.bool_true;
wrapped = arena.piecewise(&[(wrapped, condition), (degen_result, true_cond)]);
handled.push((*sym_expr, degen_val));
handled_by_denom.push((*denom, *sym_expr, degen_val));
}
}
}
wrapped
}
const NONREAL_CASE_STEPS: usize = 20;
fn real_at_generic_values(arena: &mut Arena, val: ExprId) -> bool {
let mut params: Vec<(String, ExprId)> = Vec::new();
for s in crate::base::walk::free_symbols(arena, val) {
let ExprNode::Symbol(sid) = *arena.node(s) else {
return true;
};
params.push((arena.symbol_name(sid).to_owned(), s));
}
if params.len() > 3 {
return true;
}
params.sort();
for signs in 0u32..(1 << params.len()) {
let mut bound = val;
for (k, &(_, s)) in params.iter().enumerate() {
let (p, q) = FTC_PARAMETER_VALUES[k];
let p = if signs & (1 << k) == 0 { p } else { -p };
let value = arena.rational(p, q);
bound = crate::transforms::subs::subs(arena, bound, s, value);
}
let bound = crate::transforms::eval::eval(arena, bound);
match crate::transforms::evalf::evalf_complex(arena, bound, 15) {
Ok(v) if !crate::transforms::evalf::is_real_to_digits(&v, 12) => {}
_ => return true,
}
}
false
}
thread_local! {
static NONZERO_IN_CASE: std::cell::RefCell<Vec<(usize, ExprId)>> =
const { std::cell::RefCell::new(Vec::new()) };
}
struct NonzeroInCase(usize);
impl NonzeroInCase {
fn enter(arena: &Arena, denominators: Vec<ExprId>) -> Self {
let key = std::ptr::from_ref::<Arena>(arena) as usize;
NONZERO_IN_CASE.with(|s| {
let mut s = s.borrow_mut();
let len = s.len();
s.extend(denominators.into_iter().map(|d| (key, d)));
NonzeroInCase(len)
})
}
}
impl Drop for NonzeroInCase {
fn drop(&mut self) {
let _ = NONZERO_IN_CASE.try_with(|s| {
if let Ok(mut s) = s.try_borrow_mut() {
s.truncate(self.0);
}
});
}
}
fn excluded_by_enclosing_case(arena: &mut Arena, sym: ExprId, val: ExprId) -> bool {
let key = std::ptr::from_ref::<Arena>(arena) as usize;
let nonzero: Vec<ExprId> = NONZERO_IN_CASE.with(|s| {
s.borrow()
.iter()
.filter(|&&(k, _)| k == key)
.map(|&(_, d)| d)
.collect()
});
nonzero.into_iter().any(|d| {
crate::base::walk::contains(arena, d, sym) && {
let at = crate::transforms::subs::subs(arena, d, sym, val);
let at = crate::transforms::expand::expand(arena, at);
let at = crate::transforms::eval::eval(arena, at);
arena.is_zero_structural(at)
}
})
}
fn nonreal_value_of_real_symbol(arena: &mut Arena, sym: ExprId, val: ExprId) -> bool {
let ExprNode::Symbol(sid) = *arena.node(sym) else {
return false;
};
crate::transforms::realness::symbol_declared_real(arena, sid)
&& crate::base::walk::free_symbols(arena, val).is_empty()
&& crate::transforms::evalf::evalf_complex64(arena, val)
.is_ok_and(|z| z.im.is_finite() && z.im.abs() > 1e-9 * z.norm().max(1.0))
}
fn degenerate_value_covered(
arena: &mut Arena,
handled: &[(ExprId, ExprId, ExprId)],
denom: ExprId,
sym: ExprId,
val: ExprId,
) -> bool {
let vanishes_at = |arena: &mut Arena, e: ExprId| {
let at = crate::transforms::subs::subs(arena, e, sym, val);
let at = crate::transforms::expand::expand(arena, at);
let at = crate::transforms::eval::eval(arena, at);
arena.is_zero_structural(at)
};
let mut same_denom: Vec<ExprId> = Vec::new();
for &(d, p, v) in handled {
let diff = arena.sub(p, v);
if vanishes_at(arena, diff) {
return true;
}
if d == denom {
same_denom.push(diff);
}
}
if same_denom.len() >= 2 {
let product = arena.mul(&same_denom);
return vanishes_at(arena, product);
}
false
}
fn collect_denominators(arena: &Arena, expr: ExprId) -> Vec<ExprId> {
let mut denoms = Vec::new();
let mut stack: Vec<ExprId> = vec![expr];
let mut visited: FxHashSet<ExprId> = FxHashSet::default();
while let Some(id) = stack.pop() {
if !visited.insert(id) {
continue;
}
if let ExprNode::Pow(base, exp) = arena.node(id).clone()
&& let Some(r) = arena.as_num(exp)
&& r.is_negative()
{
denoms.push(base);
}
arena.node(id).for_each_child(|c| stack.push(c));
}
denoms
}
fn has_zero_denominator(arena: &Arena, expr: ExprId) -> bool {
let mut stack: Vec<ExprId> = vec![expr];
let mut visited: FxHashSet<ExprId> = FxHashSet::default();
while let Some(id) = stack.pop() {
if !visited.insert(id) {
continue;
}
if matches!(arena.node(id), ExprNode::NaN | ExprNode::ComplexInfinity) {
return true;
}
if let ExprNode::Pow(base, exp) = arena.node(id).clone()
&& let Some(r) = arena.as_num(exp)
&& r.is_negative()
&& arena.is_zero_structural(base)
{
return true;
}
arena.node(id).for_each_child(|c| stack.push(c));
}
false
}
fn try_special_function_integral(
arena: &mut Arena,
dependent: &[ExprId],
constants: &[ExprId],
var: ExprId,
_var_sym: SymbolId,
) -> Option<ExprId> {
if dependent.len() != 2 {
return None;
}
let (func_factor, _inv_factor) = if is_inv_of_var(arena, dependent[0], var) {
(dependent[1], dependent[0])
} else if is_inv_of_var(arena, dependent[1], var) {
(dependent[0], dependent[1])
} else {
return None;
};
let sf_result = match arena.node(func_factor).clone() {
ExprNode::Sin(inner) if inner == var => Some(arena.si(var)),
ExprNode::Cos(inner) if inner == var => Some(arena.ci(var)),
ExprNode::Sinh(inner) if inner == var => Some(special_apply(arena, FN_SHI, var)),
ExprNode::Cosh(inner) if inner == var => Some(special_apply(arena, FN_CHI, var)),
ExprNode::Exp(inner) if inner == var => Some(arena.ei(var)),
ExprNode::Exp(inner) => {
if let ExprNode::Neg(neg_inner) = arena.node(inner).clone()
&& neg_inner == var
{
let neg_var = arena.neg(var);
let ei = arena.ei(neg_var);
let neg_ei = arena.neg(ei);
return Some(wrap_with_constants(arena, neg_ei, constants));
}
None
}
_ => None,
};
sf_result.map(|r| wrap_with_constants(arena, r, constants))
}
fn special_apply(arena: &mut Arena, name: &str, arg: ExprId) -> ExprId {
let sid = arena.symbols.intern(name);
arena.intern(ExprNode::Apply(sid, std::iter::once(arg).collect()))
}
fn is_inv_of_var(arena: &Arena, expr: ExprId, var: ExprId) -> bool {
if let ExprNode::Pow(base, exp) = arena.node(expr)
&& *base == var
&& let Some(r) = arena.as_num(*exp)
{
return *r == num_rational::Ratio::<num_bigint::BigInt>::from_integer((-1).into());
}
false
}
fn partition_factors(
arena: &mut Arena,
children: &[ExprId],
var_sym: SymbolId,
) -> (SmallVec<[ExprId; 4]>, SmallVec<[ExprId; 4]>) {
let mut constants: SmallVec<[ExprId; 4]> = SmallVec::new();
let mut dependent: SmallVec<[ExprId; 4]> = SmallVec::new();
for &child in children {
if contains_var(arena, child, var_sym) {
let flat = flatten_nested_pow(arena, child, var_sym)
.or_else(|| exp_power_as_exp(arena, child, var_sym))
.unwrap_or(child);
dependent.push(flat);
} else {
constants.push(child);
}
}
(constants, dependent)
}
fn exp_power_as_exp(arena: &mut Arena, expr: ExprId, var_sym: SymbolId) -> Option<ExprId> {
let ExprNode::Pow(base, c) = arena.node(expr).clone() else {
return None;
};
let ExprNode::Exp(u) = arena.node(base).clone() else {
return None;
};
if contains_var(arena, c, var_sym) {
return None;
}
let integer = arena.as_num(c).is_some_and(|q| q.is_integer());
if !integer && !real_for_real_variable(arena, u, var_sym) {
return None;
}
let cu = arena.mul(&[c, u]);
Some(arena.exp(cu))
}
fn real_for_real_variable(arena: &mut Arena, u: ExprId, var_sym: SymbolId) -> bool {
let mut reals = crate::base::assumptions::AssumptionCache::new();
let mut stack = vec![u];
let mut seen: FxHashSet<ExprId> = FxHashSet::default();
while let Some(id) = stack.pop() {
if !seen.insert(id) {
continue;
}
if !contains_var(arena, id, var_sym) {
if !constant_is_real(arena, id, &mut reals) {
return false;
}
continue;
}
let whole_line = match arena.node(id) {
ExprNode::Symbol(s) => *s == var_sym,
ExprNode::Add(_)
| ExprNode::Mul(_)
| ExprNode::Sin(_)
| ExprNode::Cos(_)
| ExprNode::Exp(_)
| ExprNode::Atan(_)
| ExprNode::Sinh(_)
| ExprNode::Cosh(_)
| ExprNode::Tanh(_)
| ExprNode::Abs(_) => true,
ExprNode::Pow(_, exp) => arena
.as_num(*exp)
.is_some_and(|r| r.is_integer() && !r.is_negative()),
_ => false,
};
if !whole_line {
return false;
}
arena.node(id).for_each_child(|c| stack.push(c));
}
true
}
const MAX_EXPANDED_DEGREE: u32 = 64;
fn expand_polynomial_product(
arena: &mut Arena,
dependent: &[ExprId],
var_sym: SymbolId,
) -> Option<ExprId> {
if dependent.len() < 2 {
return None;
}
let mut has_power_of_sum = false;
let mut degree: u32 = 0;
for &d in dependent {
let mut stack: Vec<(ExprId, u32)> = vec![(d, 1)];
let mut factor_degree: u32 = 0;
while let Some((id, mult)) = stack.pop() {
if !contains_var(arena, id, var_sym) {
continue;
}
match arena.node(id) {
ExprNode::Symbol(_) => factor_degree = factor_degree.max(mult),
ExprNode::Add(ch) => stack.extend(ch.iter().map(|&c| (c, mult))),
ExprNode::Mul(ch) => {
let k = u32::try_from(ch.len()).ok()?;
stack.extend(ch.iter().map(|&c| (c, mult.saturating_mul(k))));
}
ExprNode::Pow(b, n) => {
let k = arena
.as_num(*n)
.filter(|q| q.is_integer() && !q.is_negative())?;
let k = u32::try_from(k.to_integer()).ok()?;
if matches!(arena.node(*b), ExprNode::Add(_)) && k >= 2 {
has_power_of_sum = true;
}
stack.push((*b, mult.saturating_mul(k)));
}
_ => return None,
}
}
degree = degree.saturating_add(factor_degree);
}
if !has_power_of_sum || degree > MAX_EXPANDED_DEGREE {
return None;
}
let product = arena.mul(dependent);
let expanded = crate::transforms::expand::expand(arena, product);
(expanded != product).then_some(expanded)
}
fn merge_exponential_factors(
arena: &mut Arena,
dependent: &[ExprId],
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let mut exponents: SmallVec<[ExprId; 4]> = SmallVec::new();
let mut rest: SmallVec<[ExprId; 4]> = SmallVec::new();
for &d in dependent {
match arena.node(d) {
ExprNode::Exp(u) if is_polynomial_in(arena, *u, var, var_sym) => exponents.push(*u),
_ => rest.push(d),
}
}
if exponents.len() < 2 {
return None;
}
let sum = arena.add(&exponents);
let merged = arena.exp(sum);
rest.push(merged);
Some(arena.mul(&rest))
}
fn split_constant_factors(
arena: &mut Arena,
expr: ExprId,
var_sym: SymbolId,
) -> Option<(SmallVec<[ExprId; 4]>, ExprId)> {
let ExprNode::Mul(children) = arena.node(expr).clone() else {
return None;
};
let (constants, dependent) = partition_factors(arena, &children, var_sym);
if constants.is_empty() || dependent.is_empty() {
return None;
}
let f = if dependent.len() == 1 {
dependent[0]
} else {
arena.mul(&dependent)
};
Some((constants, f))
}
fn wrap_with_constants(arena: &mut Arena, result: ExprId, constants: &[ExprId]) -> ExprId {
if constants.is_empty() {
result
} else {
let mut all: SmallVec<[ExprId; 4]> = constants.iter().copied().collect();
all.push(result);
arena.mul(&all)
}
}
fn flatten_nested_pow(arena: &mut Arena, expr: ExprId, var_sym: SymbolId) -> Option<ExprId> {
let ExprNode::Pow(base, exp) = arena.node(expr).clone() else {
return None;
};
let ExprNode::Pow(inner_base, inner_exp) = arena.node(base).clone() else {
return None;
};
let (m, n) = (arena.as_num(inner_exp)?.clone(), arena.as_num(exp)?.clone());
let combined = &m * &n;
let two = num_bigint::BigInt::from(2);
let even_integer = |q: &crate::base::numeric::Q| q.is_integer() && (q.numer() % &two).is_zero();
let principal =
n.is_integer() || (m > crate::base::numeric::qi(-1) && m <= crate::base::numeric::qi(1));
let abs_form = even_integer(&m) && real_on_reals(arena, inner_base, var_sym);
let combined_id = arena.num_ratio(combined.clone());
if principal || (abs_form && even_integer(&combined)) {
Some(arena.pow(inner_base, combined_id))
} else if abs_form {
let abs_base = arena.abs(inner_base);
Some(arena.pow(abs_base, combined_id))
} else {
None
}
}
const FTC_SAMPLE_POINTS: [(i64, i64); 6] = [(1, 3), (-5, 7), (13, 11), (-17, 5), (7, 5), (-2, 1)];
const FTC_RETRY_POINTS: [(i64, i64); 6] =
[(1, 19), (-3, 23), (37, 4), (-29, 3), (211, 7), (-307, 9)];
const FTC_CHECK_DIGITS: u32 = 30;
const FTC_CHECK_REL_TOL: f64 = 1e-20;
const FTC_DECISIVE_FACTOR: f64 = 16.0;
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
enum FtcVerdict {
Verified,
Wrong,
Undecided,
Untestable,
}
fn antiderivative_rejected(
arena: &mut Arena,
f: ExprId,
big_f: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> bool {
let key = (std::ptr::from_ref::<Arena>(arena) as usize, f, big_f, var);
let verdict = match VERDICT_MEMO.with(|m| m.borrow().get(&key).copied()) {
Some(verdict) => verdict,
None => {
let verdict = stage!(
arena,
"verify",
big_f,
check_antiderivative(arena, f, big_f, var, var_sym)
);
VERDICT_MEMO.with(|m| m.borrow_mut().insert(key, verdict));
verdict
}
};
match verdict {
FtcVerdict::Verified | FtcVerdict::Untestable => false,
FtcVerdict::Wrong | FtcVerdict::Undecided => {
tracing::debug!(?verdict, "integrate: candidate antiderivative rejected");
true
}
}
}
thread_local! {
static VERDICT_MEMO: std::cell::RefCell<FxHashMap<(usize, ExprId, ExprId, ExprId), FtcVerdict>> =
std::cell::RefCell::new(FxHashMap::default());
}
const FTC_PARAMETER_VALUES: [(i64, i64); 8] = [
(17, 13),
(7, 11),
(23, 17),
(5, 19),
(31, 23),
(11, 29),
(41, 31),
(13, 37),
];
const FTC_PARAMETER_IMAG: [(i64, i64); 8] = [
(13, 41),
(19, 43),
(17, 47),
(23, 53),
(19, 59),
(29, 61),
(23, 67),
(31, 71),
];
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
enum ParameterValues {
Real,
Complex,
}
struct BoundCandidate {
f: ExprId,
big_f: ExprId,
may_be_complex: bool,
}
fn bind_parameters(
arena: &mut Arena,
f: ExprId,
big_f: ExprId,
var: ExprId,
mode: ParameterValues,
) -> Option<BoundCandidate> {
let mut params: Vec<(String, ExprId)> = Vec::new();
for root in [f, big_f] {
for s in crate::base::walk::free_symbols(arena, root) {
if s == var || params.iter().any(|&(_, p)| p == s) {
continue;
}
let ExprNode::Symbol(sid) = *arena.node(s) else {
return None;
};
params.push((arena.symbol_name(sid).to_owned(), s));
}
}
if params.len() > FTC_PARAMETER_VALUES.len() {
return None;
}
params.sort();
let (mut f, mut big_f) = (f, big_f);
let mut may_be_complex = false;
let mut facts = crate::base::assumptions::AssumptionCache::new();
for (k, &(_, s)) in params.iter().enumerate() {
let ExprNode::Symbol(sid) = *arena.node(s) else {
return None;
};
let mut declared = arena.symbol_assumptions(sid);
declared.normalize_declared();
let complex = declared.query(crate::base::assumptions::Props::REAL) != Some(true);
may_be_complex |= complex;
let (p, q) = FTC_PARAMETER_VALUES[k];
let re = arena.rational(p, q);
let mut candidates: SmallVec<[ExprId; 3]> = SmallVec::new();
if mode == ParameterValues::Complex && complex {
let (r, t) = FTC_PARAMETER_IMAG[k];
let im = arena.rational(r, t);
let i = arena.i_unit();
let im = arena.mul(&[im, i]);
candidates.push(arena.add(&[re, im]));
}
candidates.push(re);
candidates.push(arena.rational(-p, q));
let value = candidates
.into_iter()
.find(|&v| value_satisfies(arena, v, declared, &mut facts))?;
f = crate::transforms::subs::subs(arena, f, s, value);
big_f = crate::transforms::subs::subs(arena, big_f, s, value);
}
Some(BoundCandidate {
f,
big_f,
may_be_complex,
})
}
fn value_satisfies(
arena: &mut Arena,
value: ExprId,
declared: crate::base::assumptions::Assumptions,
facts: &mut crate::base::assumptions::AssumptionCache,
) -> bool {
declared
.known_true
.iter()
.all(|flag| facts.query(arena, value, flag) == Some(true))
&& declared
.known_false
.iter()
.all(|flag| facts.query(arena, value, flag) == Some(false))
}
fn check_antiderivative(
arena: &mut Arena,
f: ExprId,
big_f: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> FtcVerdict {
if crate::base::walk::has_unevaluated(arena, big_f) {
return FtcVerdict::Untestable;
}
let Some(bound) = bind_parameters(arena, f, big_f, var, ParameterValues::Real) else {
return FtcVerdict::Untestable;
};
let mut real_domain: Vec<ExprId> = Vec::new();
let verdict =
check_bound_antiderivative(arena, bound.f, bound.big_f, var, var_sym, &mut real_domain);
if verdict != FtcVerdict::Verified || !bound.may_be_complex {
return verdict;
}
if analytic_in_parameters(arena, f, var_sym) && analytic_in_parameters(arena, big_f, var_sym) {
return verdict;
}
let Some(complex) = bind_parameters(arena, f, big_f, var, ParameterValues::Complex) else {
return verdict;
};
if real_domain.is_empty() {
real_domain = real_domain_points(arena, bound.f, var);
}
if complex_round_differs(arena, complex.f, complex.big_f, var, var_sym, &real_domain) {
tracing::debug!(
"integrate: candidate antiderivative fails the FTC check at complex parameter values"
);
return FtcVerdict::Wrong;
}
verdict
}
fn analytic_in_parameters(arena: &Arena, e: ExprId, var_sym: SymbolId) -> bool {
let mut has_param: FxHashMap<ExprId, bool> = FxHashMap::default();
for id in crate::base::walk::post_order_ids(arena, e) {
let node = arena.node(id);
let mut below = matches!(node, ExprNode::Symbol(s) if *s != var_sym);
node.for_each_child(|c| below |= has_param.get(&c).copied().unwrap_or(false));
let analytic = matches!(
node,
ExprNode::Num(_)
| ExprNode::Symbol(_)
| ExprNode::Pi
| ExprNode::E
| ExprNode::ImaginaryUnit
| ExprNode::EulerGamma
| ExprNode::Catalan
| ExprNode::GoldenRatio
| ExprNode::Add(_)
| ExprNode::Mul(_)
| ExprNode::Pow(..)
| ExprNode::Neg(_)
| ExprNode::Sin(_)
| ExprNode::Cos(_)
| ExprNode::Tan(_)
| ExprNode::Exp(_)
| ExprNode::Ln(_)
| ExprNode::Asin(_)
| ExprNode::Acos(_)
| ExprNode::Atan(_)
| ExprNode::Sinh(_)
| ExprNode::Cosh(_)
| ExprNode::Tanh(_)
| ExprNode::Asinh(_)
| ExprNode::Acosh(_)
| ExprNode::Atanh(_)
| ExprNode::Gamma(_)
| ExprNode::Erf(_)
| ExprNode::Erfc(_)
| ExprNode::Si(_)
| ExprNode::Ci(_)
| ExprNode::Ei(_)
| ExprNode::Li(_)
| ExprNode::Zeta(_)
| ExprNode::RootSum(..)
);
if below && !analytic {
return false;
}
has_param.insert(id, below);
}
true
}
fn real_domain_points(arena: &mut Arena, f: ExprId, var: ExprId) -> Vec<ExprId> {
let mut out = Vec::new();
for &(p, q) in FTC_SAMPLE_POINTS.iter().chain(FTC_RETRY_POINTS.iter()) {
let point = arena.rational(p, q);
let f_at = crate::transforms::subs::subs(arena, f, var, point);
if numeric_value(arena, f_at)
.is_some_and(|v| crate::transforms::evalf::is_real_to_digits(&v, FTC_CHECK_DIGITS))
{
out.push(point);
}
}
out
}
const FTC_COMPLEX_AGREEMENTS: usize = 3;
fn complex_round_differs(
arena: &mut Arena,
f: ExprId,
big_f: ExprId,
var: ExprId,
var_sym: SymbolId,
real_domain: &[ExprId],
) -> bool {
let d_big_f = crate::transforms::diff::diff(arena, big_f, var);
if crate::base::walk::has_unevaluated(arena, d_big_f) {
return false;
}
let residual = arena.sub(d_big_f, f);
if !contains_var(arena, residual, var_sym) {
return residual_abs(arena, residual).is_some_and(|r| r > FTC_CHECK_REL_TOL);
}
let mut agreements = 0;
for &point in real_domain {
match ftc_point(arena, f, d_big_f, residual, var, point) {
FtcPoint::Differs { .. } => {
tracing::debug!(
point = %arena.display(point),
"integrate: F′ ≠ f at complex parameter values"
);
return true;
}
FtcPoint::Agrees { .. } => {
agreements += 1;
if agreements >= FTC_COMPLEX_AGREEMENTS {
return false;
}
}
FtcPoint::IntegrandUnevaluable | FtcPoint::DerivativeUnevaluable => {}
}
}
false
}
fn check_bound_antiderivative(
arena: &mut Arena,
f: ExprId,
big_f: ExprId,
var: ExprId,
var_sym: SymbolId,
real_domain: &mut Vec<ExprId>,
) -> FtcVerdict {
let d_big_f = crate::transforms::diff::diff(arena, big_f, var);
if crate::base::walk::has_unevaluated(arena, d_big_f) {
return FtcVerdict::Untestable;
}
let residual = arena.sub(d_big_f, f);
if !contains_var(arena, residual, var_sym)
&& let Some(r) = residual_abs(arena, residual)
{
if r > FTC_CHECK_REL_TOL {
return FtcVerdict::Wrong;
}
return if evaluates_at_a_sample_point(arena, big_f, var) {
FtcVerdict::Verified
} else if evaluates_at_a_sample_point(arena, f, var) {
FtcVerdict::Undecided
} else {
FtcVerdict::Untestable
};
}
let mut integrand_evaluates = false;
for points in [&FTC_SAMPLE_POINTS, &FTC_RETRY_POINTS] {
let mut agrees = false;
for &(p, q) in points {
let point = arena.rational(p, q);
match ftc_point(arena, f, d_big_f, residual, var, point) {
FtcPoint::IntegrandUnevaluable => {}
FtcPoint::Agrees { f_real } => {
integrand_evaluates = true;
if f_real {
real_domain.push(point);
}
if !agrees {
let big_f_at = crate::transforms::subs::subs(arena, big_f, var, point);
agrees = numeric_value(arena, big_f_at).is_some();
}
}
FtcPoint::DerivativeUnevaluable | FtcPoint::Differs { f_real: false } => {
integrand_evaluates = true;
}
FtcPoint::Differs { f_real: true } => {
tracing::debug!(
point = %arena.display(point),
"integrate: candidate antiderivative fails the numeric FTC check"
);
return FtcVerdict::Wrong;
}
}
}
if agrees {
return FtcVerdict::Verified;
}
}
if integrand_evaluates {
FtcVerdict::Undecided
} else {
FtcVerdict::Untestable
}
}
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
enum FtcPoint {
IntegrandUnevaluable,
DerivativeUnevaluable,
Agrees { f_real: bool },
Differs { f_real: bool },
}
fn ftc_point(
arena: &mut Arena,
f: ExprId,
d_big_f: ExprId,
residual: ExprId,
var: ExprId,
point: ExprId,
) -> FtcPoint {
let f_at = crate::transforms::subs::subs(arena, f, var, point);
let Some(f_val) = numeric_value(arena, f_at) else {
return FtcPoint::IntegrandUnevaluable;
};
let Some(f_abs) = crate::transforms::evalf::abs_to_f64(&f_val).filter(|v| v.is_finite()) else {
return FtcPoint::IntegrandUnevaluable;
};
let f_real = crate::transforms::evalf::is_real_to_digits(&f_val, FTC_CHECK_DIGITS);
let tol = FTC_CHECK_REL_TOL * f_abs.max(1.0);
let d_at = crate::transforms::subs::subs(arena, d_big_f, var, point);
let d_val = numeric_value(arena, d_at);
if let Some(d) = &d_val
&& let Some(gap) = crate::transforms::evalf::distance_to_f64(d, &f_val)
&& !gap.is_nan()
{
if gap * FTC_DECISIVE_FACTOR <= tol {
return FtcPoint::Agrees { f_real };
}
if gap >= tol * FTC_DECISIVE_FACTOR {
return FtcPoint::Differs { f_real };
}
}
let r_at = crate::transforms::subs::subs(arena, residual, var, point);
let gap = match residual_abs(arena, r_at) {
Some(gap) => gap,
None => {
let Some(d_val) = d_val else {
return FtcPoint::DerivativeUnevaluable;
};
match crate::transforms::evalf::distance_to_f64(&d_val, &f_val) {
Some(gap) if !gap.is_nan() => gap,
_ => return FtcPoint::DerivativeUnevaluable,
}
}
};
if gap <= tol {
FtcPoint::Agrees { f_real }
} else {
FtcPoint::Differs { f_real }
}
}
fn evaluates_at_a_sample_point(arena: &mut Arena, e: ExprId, var: ExprId) -> bool {
FTC_SAMPLE_POINTS
.iter()
.chain(FTC_RETRY_POINTS.iter())
.any(|&(p, q)| {
let point = arena.rational(p, q);
let e_at = crate::transforms::subs::subs(arena, e, var, point);
numeric_value(arena, e_at).is_some()
})
}
fn exact_numeric_root_sums(arena: &mut Arena, e: ExprId) -> ExprId {
let order = crate::base::walk::post_order_ids(arena, e);
let mut out = e;
for id in order {
if let ExprNode::RootSum(poly, body, sumvar) = *arena.node(id)
&& let Some(value) = crate::calculus::risch::log_to_real::numeric_rootsum_rational_body(
arena, poly, body, sumvar,
)
{
out = arena.subs_structural(out, id, value);
}
}
out
}
fn numeric_value(arena: &mut Arena, e: ExprId) -> Option<crate::base::bigcomplex::Complex> {
let e = exact_numeric_root_sums(arena, e);
let e = crate::transforms::eval::eval(arena, e);
if crate::base::walk::has_unevaluated(arena, e) {
return None;
}
crate::transforms::evalf::evalf_complex(arena, e, FTC_CHECK_DIGITS).ok()
}
fn residual_abs(arena: &mut Arena, residual: ExprId) -> Option<f64> {
let r = exact_numeric_root_sums(arena, residual);
let r = crate::transforms::eval::eval(arena, r);
if crate::base::walk::has_unevaluated(arena, r) {
return None;
}
let v = if let Some(q) = arena.as_num(r) {
crate::base::numeric::ratio_to_f64(q)?.abs()
} else {
let z = crate::transforms::evalf::evalf_complex(arena, r, FTC_CHECK_DIGITS).ok()?;
crate::transforms::evalf::abs_to_f64(&z)?
};
v.is_finite().then_some(v)
}
fn try_biquadratic_rational(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
var_sym: SymbolId,
) -> Option<ExprId> {
let (numer_id, denom_id) = crate::poly::polybridge::as_numer_denom(arena, expr);
if denom_id == arena.one {
return None;
}
let numer_id = crate::transforms::expand::expand(arena, numer_id);
let numer_id = crate::transforms::eval::eval(arena, numer_id);
let denom_id = crate::transforms::expand::expand(arena, denom_id);
let denom_id = crate::transforms::eval::eval(arena, denom_id);
let numer = crate::poly::polybridge::expr_to_poly(arena, numer_id, var)?;
let denom = crate::poly::polybridge::expr_to_poly(arena, denom_id, var)?;
if denom.degree() != Some(4) || numer.is_zero() {
return None;
}
if !denom.coeff(3).is_zero() || !denom.coeff(1).is_zero() {
return None;
}
let lead = denom.coeff(4);
if lead.is_zero() {
return None;
}
let (quotient, numer) = if numer.degree().unwrap_or(0) > 3 {
let (quot, rem) = numer.div_rem(&denom);
(Some(quot), rem)
} else {
(None, numer)
};
let poly_part = quotient.filter(|q| !q.is_zero()).map(|q| {
let deg = q.degree().unwrap_or(0);
let mut coeffs: Vec<Q> = vec![Q::zero(); deg + 2];
for (k, c) in coeffs.iter_mut().enumerate().skip(1) {
*c = q.coeff(k - 1) / Q::from_integer(k.into());
}
let antiderivative = Poly::from_coeffs(coeffs);
crate::poly::polybridge::poly_to_expr(arena, &antiderivative, var)
});
if numer.is_zero() {
return poly_part;
}
let p = denom.coeff(2) / &lead;
let q = denom.coeff(0) / &lead;
let n_ids: [ExprId; 4] = [
arena.num_ratio(numer.coeff(0) / &lead),
arena.num_ratio(numer.coeff(1) / &lead),
arena.num_ratio(numer.coeff(2) / &lead),
arena.num_ratio(numer.coeff(3) / &lead),
];
let disc = &p * &p - Q::from_integer(4.into()) * &q;
let quadratics: SmallVec<[(ExprId, ExprId, ExprId, ExprId); 2]> = if disc.is_negative() {
tracing::debug!(p = %p, q = %q, "integrate: biquadratic denominator, x² ± a·x + c");
biquadratic_conjugate_quadratics(arena, &p, &q, &n_ids)
} else if disc.is_positive() && p.is_positive() && q.is_positive() {
tracing::debug!(p = %p, q = %q, "integrate: biquadratic denominator, (x² + r₁)(x² + r₂)");
biquadratic_even_quadratics(arena, &p, &disc, &n_ids)
} else {
return None;
};
let two = arena.int(2);
let four = arena.int(4);
let half = arena.rational(1, 2);
let x_sq = arena.pow(var, two);
let two_x = arena.mul(&[two, var]);
let mut terms: SmallVec<[ExprId; 4]> = SmallVec::new();
for &(a, c, lin, cst) in quadratics.iter() {
let a_x = arena.mul(&[a, var]);
let quad = arena.add(&[x_sq, a_x, c]);
let ln_quad = arena.ln(quad);
let half_lin = arena.mul(&[half, lin]);
terms.push(arena.mul(&[half_lin, ln_quad]));
let four_c = arena.mul(&[four, c]);
let a_sq = arena.mul(&[a, a]);
let k_sq = arena.sub(four_c, a_sq);
let k = arena.sqrt(k_sq);
let two_cst = arena.mul(&[two, cst]);
let a_lin = arena.mul(&[a, lin]);
let atan_coeff_num = arena.sub(two_cst, a_lin);
let atan_coeff = arena.div(atan_coeff_num, k);
let atan_arg_num = arena.add(&[two_x, a]);
let atan_arg = arena.div(atan_arg_num, k);
let atan_val = arena.atan(atan_arg);
terms.push(arena.mul(&[atan_coeff, atan_val]));
}
if let Some(poly_part) = poly_part {
terms.push(poly_part);
}
let result = arena.add(&terms);
let result = crate::transforms::eval::eval(arena, result);
if candidate_rejected(arena, expr, result, var, var_sym) {
tracing::debug!("integrate: biquadratic closed form failed verification");
return None;
}
Some(result)
}
fn biquadratic_conjugate_quadratics(
arena: &mut Arena,
p: &Q,
q: &Q,
n: &[ExprId; 4],
) -> SmallVec<[(ExprId, ExprId, ExprId, ExprId); 2]> {
let p_id = arena.num_ratio(p.clone());
let q_id = arena.num_ratio(q.clone());
let two = arena.int(2);
let half = arena.rational(1, 2);
let c = arena.sqrt(q_id);
let two_c = arena.mul(&[two, c]);
let a_sq = arena.sub(two_c, p_id);
let a = arena.sqrt(a_sq);
let neg_a = arena.neg(a);
let sum_bd = arena.div(n[0], c); let sum_ac = n[3]; let diff_ac = {
let t = arena.sub(n[2], sum_bd);
arena.div(t, a)
}; let diff_bd = {
let c_n3 = arena.mul(&[c, n[3]]);
let t = arena.sub(n[1], c_n3);
arena.div(t, a)
}; let coeff_a = {
let s = arena.add(&[sum_ac, diff_ac]);
arena.mul(&[half, s])
};
let coeff_c = {
let s = arena.sub(sum_ac, diff_ac);
arena.mul(&[half, s])
};
let coeff_b = {
let s = arena.add(&[sum_bd, diff_bd]);
arena.mul(&[half, s])
};
let coeff_d = {
let s = arena.sub(sum_bd, diff_bd);
arena.mul(&[half, s])
};
let mut out = SmallVec::new();
out.push((neg_a, c, coeff_a, coeff_b));
out.push((a, c, coeff_c, coeff_d));
out
}
fn biquadratic_even_quadratics(
arena: &mut Arena,
p: &Q,
disc: &Q,
n: &[ExprId; 4],
) -> SmallVec<[(ExprId, ExprId, ExprId, ExprId); 2]> {
let p_id = arena.num_ratio(p.clone());
let disc_id = arena.num_ratio(disc.clone());
let half = arena.rational(1, 2);
let sqrt_disc = arena.sqrt(disc_id);
let r1 = {
let t = arena.sub(p_id, sqrt_disc);
arena.mul(&[half, t])
};
let r2 = {
let t = arena.add(&[p_id, sqrt_disc]);
arena.mul(&[half, t])
};
let coeff_a = {
let n3_r1 = arena.mul(&[n[3], r1]);
let t = arena.sub(n[1], n3_r1);
arena.div(t, sqrt_disc)
};
let coeff_c = arena.sub(n[3], coeff_a);
let coeff_b = {
let n2_r1 = arena.mul(&[n[2], r1]);
let t = arena.sub(n[0], n2_r1);
arena.div(t, sqrt_disc)
};
let coeff_d = arena.sub(n[2], coeff_b);
let zero = arena.zero;
let mut out = SmallVec::new();
out.push((zero, r1, coeff_a, coeff_b));
out.push((zero, r2, coeff_c, coeff_d));
out
}
#[cfg(test)]
mod tests {
use super::*;
use crate::base::arena::Arena;
fn sym(a: &mut Arena, name: &str) -> ExprId {
a.symbol(name)
}
fn display(a: &Arena, id: ExprId) -> String {
a.display(id).to_string()
}
#[test]
fn integrate_constant() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let five = a.int(5);
let result = integrate(&mut a, five, x);
assert_eq!(display(&a, result), "5*x");
}
#[test]
fn integrate_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let result = integrate(&mut a, x, x);
assert_eq!(display(&a, result), "1/2*x^2");
}
#[test]
fn integrate_x_squared() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let x2 = a.pow(x, two);
let result = integrate(&mut a, x2, x);
assert_eq!(display(&a, result), "1/3*x^3");
}
#[test]
fn integrate_x_inv() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let neg_one = a.int(-1);
let x_inv = a.pow(x, neg_one);
let result = integrate(&mut a, x_inv, x);
assert_eq!(display(&a, result), "ln(abs(x))");
}
#[test]
fn integrate_sin_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.sin(x);
let result = integrate(&mut a, expr, x);
assert_eq!(display(&a, result), "-cos(x)");
}
#[test]
fn integrate_cos_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.cos(x);
let result = integrate(&mut a, expr, x);
assert_eq!(display(&a, result), "sin(x)");
}
#[test]
fn integrate_exp_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.exp(x);
let result = integrate(&mut a, expr, x);
assert_eq!(display(&a, result), "exp(x)");
}
#[test]
fn integrate_sum() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let one = a.one;
let sum = a.add(&[x, one]);
let result = integrate(&mut a, sum, x);
let s = display(&a, result);
assert!(s.contains("x^2"), "should contain x^2: {s}");
assert!(s.contains("x"), "should contain x: {s}");
}
#[test]
fn integrate_constant_times_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let three = a.int(3);
let expr = a.mul(&[three, x]);
let result = integrate(&mut a, expr, x);
assert_eq!(display(&a, result), "3/2*x^2");
}
#[test]
fn integrate_other_symbol_is_constant() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let y = sym(&mut a, "y");
let result = integrate(&mut a, y, x);
assert_eq!(display(&a, result), "x*y");
}
#[test]
fn integrate_unevaluated_for_unknown() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.tan(x);
let result = integrate(&mut a, expr, x);
assert_eq!(display(&a, result), "-ln(abs(cos(x)))");
}
#[test]
fn integrate_x_sin_x_by_parts() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let sin_x = a.sin(x);
let expr = a.mul(&[x, sin_x]);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("sin(x)"), "should contain sin(x): {s}");
assert!(s.contains("cos(x)"), "should contain cos(x): {s}");
}
#[test]
fn integrate_x_exp_x_by_parts() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let exp_x = a.exp(x);
let expr = a.mul(&[x, exp_x]);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("exp(x)"), "should contain exp(x): {s}");
}
#[test]
fn integrate_x_cos_x_by_parts() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let cos_x = a.cos(x);
let expr = a.mul(&[x, cos_x]);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("sin(x)"), "should contain sin(x): {s}");
assert!(s.contains("cos(x)"), "should contain cos(x): {s}");
}
#[test]
fn integrate_sin_2x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let two_x = a.mul(&[two, x]);
let expr = a.sin(two_x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("cos"), "should contain cos: {s}");
assert!(
s.contains("1/2") || s.contains("2"),
"should have factor of 1/2: {s}"
);
}
#[test]
fn integrate_cos_3x_plus_1() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let three = a.int(3);
let one = a.one;
let three_x = a.mul(&[three, x]);
let inner = a.add(&[three_x, one]);
let expr = a.cos(inner);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("sin"), "should contain sin: {s}");
}
#[test]
fn integrate_exp_2x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let two_x = a.mul(&[two, x]);
let expr = a.exp(two_x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("exp"), "should contain exp: {s}");
}
#[test]
fn integrate_tan_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.tan(x);
let result = integrate(&mut a, expr, x);
assert_eq!(display(&a, result), "-ln(abs(cos(x)))");
}
#[test]
fn integrate_ln_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.ln(x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("ln(x)"), "should contain ln(x): {s}");
assert!(s.contains("x"), "should contain x: {s}");
assert!(
s.contains("x*ln(x)") || s.contains("ln(x)*x"),
"should contain x*ln(x): {s}"
);
}
#[test]
fn integrate_atan_form() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let x2 = a.pow(x, two);
let one = a.one;
let base = a.add(&[x2, one]);
let neg_one = a.int(-1);
let expr = a.pow(base, neg_one);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("atan(x)"), "should be atan(x), got: {s}");
}
#[test]
fn integrate_asin_form() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let x2 = a.pow(x, two);
let one = a.one;
let base = a.sub(one, x2); let neg_half = a.rational(-1, 2);
let expr = a.pow(base, neg_half);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("asin(x)"), "should be asin(x), got: {s}");
}
#[test]
fn integrate_tan_2x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let two_x = a.mul(&[two, x]);
let expr = a.tan(two_x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("ln") && s.contains("cos"),
"should contain ln and cos: {s}"
);
assert!(
s.contains("1/2") || s.contains("2"),
"should have factor involving 2: {s}"
);
}
#[test]
fn integrate_ln_3x_plus_1() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let three = a.int(3);
let one = a.one;
let three_x = a.mul(&[three, x]);
let inner = a.add(&[three_x, one]);
let expr = a.ln(inner);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("ln"), "should contain ln: {s}");
assert!(!s.contains("Integral"), "should not be unevaluated: {s}");
}
#[test]
fn integrate_tanh_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.tanh(x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("ln") && s.contains("cosh"),
"∫ tanh(x) dx should be ln(cosh(x)), got: {s}"
);
}
#[test]
fn integrate_tanh_2x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let two_x = a.mul(&[two, x]);
let expr = a.tanh(two_x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("ln") && s.contains("cosh"),
"∫ tanh(2x) dx should involve ln(cosh(2x)), got: {s}"
);
}
#[test]
fn integrate_sinh_2x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let two_x = a.mul(&[two, x]);
let expr = a.sinh(two_x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("cosh"),
"∫ sinh(2x) dx should involve cosh, got: {s}"
);
}
#[test]
fn integrate_cosh_3x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let three = a.int(3);
let three_x = a.mul(&[three, x]);
let expr = a.cosh(three_x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("sinh"),
"∫ cosh(3x) dx should involve sinh, got: {s}"
);
}
#[test]
fn integrate_2x_plus_1_cubed() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let one = a.one;
let two_x = a.mul(&[two, x]);
let inner = a.add(&[two_x, one]);
let three = a.int(3);
let expr = a.pow(inner, three);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(!s.contains("Integral"), "should not be unevaluated: {s}");
}
#[test]
fn integrate_x_plus_1_squared() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let one = a.one;
let inner = a.add(&[x, one]);
let two = a.int(2);
let expr = a.pow(inner, two);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(!s.contains("Integral"), "should not be unevaluated: {s}");
}
#[test]
fn integrate_asin_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.asin(x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("asin"), "should contain asin: {s}");
}
#[test]
fn integrate_atan_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.atan(x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("atan"), "should contain atan: {s}");
}
#[test]
fn integrate_acos_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let expr = a.acos(x);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(s.contains("acos"), "should contain acos: {s}");
}
#[test]
fn integrate_2x_exp_x_squared_u_sub() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let two_x = a.mul(&[two, x]);
let x2 = a.pow(x, two);
let exp_x2 = a.exp(x2);
let expr = a.mul(&[two_x, exp_x2]);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("exp"),
"∫ 2x·exp(x²) dx should contain exp, got: {s}"
);
assert!(
!s.contains("Integral"),
"∫ 2x·exp(x²) dx should not be unevaluated, got: {s}"
);
}
#[test]
fn integrate_cos_x_exp_sin_x_u_sub() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let sin_x = a.sin(x);
let cos_x = a.cos(x);
let exp_sin_x = a.exp(sin_x);
let expr = a.mul(&[cos_x, exp_sin_x]);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("exp"),
"∫ cos(x)·exp(sin(x)) dx should contain exp, got: {s}"
);
assert!(
!s.contains("Integral"),
"∫ cos(x)·exp(sin(x)) dx should not be unevaluated, got: {s}"
);
}
#[test]
fn integrate_x_over_x2_plus_1_u_sub() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let x2 = a.pow(x, two);
let one = a.one;
let x2_plus_1 = a.add(&[x2, one]);
let expr = a.div(x, x2_plus_1);
let result = integrate(&mut a, expr, x);
let s = display(&a, result);
assert!(
s.contains("ln"),
"∫ x/(x²+1) dx should contain ln, got: {s}"
);
assert!(
!s.contains("Integral"),
"∫ x/(x²+1) dx should not be unevaluated, got: {s}"
);
}
#[test]
fn integrate_complete_square() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let five = a.int(5);
let x2 = a.pow(x, two);
let two_x = a.mul(&[two, x]);
let quadratic = a.add(&[x2, two_x, five]);
let neg_one = a.int(-1);
let integrand = a.pow(quadratic, neg_one);
let result = integrate(&mut a, integrand, x);
let s = display(&a, result);
assert!(s.contains("atan"), "should use atan: {s}");
}
#[test]
fn integrate_complete_square_simple() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let two = a.int(2);
let x2 = a.pow(x, two);
let one = a.one;
let quadratic = a.add(&[x2, x, one]);
let neg_one = a.int(-1);
let integrand = a.pow(quadratic, neg_one);
let result = integrate(&mut a, integrand, x);
let s = display(&a, result);
assert!(s.contains("atan"), "should use atan: {s}");
}
#[test]
fn symbolic_linear_coeff_of_mul_a_x() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let param_a = sym(&mut a, "a");
let ax = a.mul(&[param_a, x]);
let var_sym = match a.node(x) {
ExprNode::Symbol(sid) => *sid,
_ => panic!("x should be a symbol"),
};
let result = super::symbolic_linear_coeff_of(&mut a, ax, x, var_sym);
assert!(
result.is_some(),
"a*x should be recognized as linear in x, node: {:?}",
a.node(ax)
);
let (coeff, constant) = result.unwrap();
assert_eq!(
coeff,
param_a,
"coefficient should be a, got {}",
display(&a, coeff)
);
assert_eq!(
constant,
a.zero,
"constant should be 0, got {}",
display(&a, constant)
);
}
#[test]
fn symbolic_linear_coeff_of_add_ax_b() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let param_a = sym(&mut a, "a");
let param_b = sym(&mut a, "b");
let ax = a.mul(&[param_a, x]);
let ax_plus_b = a.add(&[ax, param_b]);
let var_sym = match a.node(x) {
ExprNode::Symbol(sid) => *sid,
_ => panic!("x should be a symbol"),
};
let result = super::symbolic_linear_coeff_of(&mut a, ax_plus_b, x, var_sym);
assert!(
result.is_some(),
"a*x+b should be recognized as linear in x, expr: {}",
display(&a, ax_plus_b)
);
let (coeff, constant) = result.unwrap();
assert_eq!(
coeff,
param_a,
"coefficient should be a, got {}",
display(&a, coeff)
);
assert_eq!(
constant,
param_b,
"constant should be b, got {}",
display(&a, constant)
);
}
#[test]
fn symbolic_linear_coeff_of_bare_var() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let var_sym = match a.node(x) {
ExprNode::Symbol(sid) => *sid,
_ => panic!("x should be a symbol"),
};
let result = super::symbolic_linear_coeff_of(&mut a, x, x, var_sym);
assert!(result.is_some(), "x should be recognized as linear in x");
let (coeff, constant) = result.unwrap();
assert_eq!(coeff, a.one, "coefficient should be 1");
assert_eq!(constant, a.zero, "constant should be 0");
}
#[test]
fn integrate_sin_symbolic_coeff() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let param_a = sym(&mut a, "a");
let ax = a.mul(&[param_a, x]);
let sin_ax = a.sin(ax);
let result = integrate(&mut a, sin_ax, x);
let s = display(&a, result);
assert!(
!s.contains("Integral"),
"∫sin(a*x)dx should not be unevaluated: {s}"
);
assert!(s.contains("cos"), "should contain cos: {s}");
}
#[test]
fn integrate_exp_symbolic_coeff() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let param_a = sym(&mut a, "a");
let ax = a.mul(&[param_a, x]);
let exp_ax = a.exp(ax);
let result = integrate(&mut a, exp_ax, x);
let s = display(&a, result);
assert!(
!s.contains("Integral"),
"∫exp(a*x)dx should not be unevaluated: {s}"
);
assert!(s.contains("exp"), "should contain exp: {s}");
}
#[test]
fn integrate_cosh_symbolic_coeff() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let param_a = sym(&mut a, "a");
let ax = a.mul(&[param_a, x]);
let cosh_ax = a.cosh(ax);
let result = integrate(&mut a, cosh_ax, x);
let s = display(&a, result);
assert!(
!s.contains("Integral"),
"∫cosh(a*x)dx should not be unevaluated: {s}"
);
}
fn eval_at(
arena: &mut Arena,
expr: ExprId,
var: ExprId,
numer: i64,
denom: i64,
) -> Option<f64> {
let val = arena.rational(numer, denom);
let substituted = crate::transforms::subs::subs(arena, expr, var, val);
let evaled = crate::transforms::eval::eval(arena, substituted);
let s = crate::transforms::evalf::evalf(arena, evaled, 15).ok()?;
s.parse::<f64>().ok()
}
fn assert_ftc(
arena: &mut Arena,
integrand: ExprId,
antideriv: ExprId,
var: ExprId,
test_points: &[(i64, i64)],
tol: f64,
name: &str,
) {
let deriv = crate::transforms::diff::diff(arena, antideriv, var);
let deriv_simplified = crate::simplify::simplify_engine::smart_simplify(arena, deriv);
for &(n, d) in test_points {
let f_val = eval_at(arena, integrand, var, n, d);
let fp_val = eval_at(arena, deriv_simplified, var, n, d);
match (f_val, fp_val) {
(Some(f), Some(fp)) => {
assert!(
(f - fp).abs() < tol,
"FTC failed for {name} at x={n}/{d}: f(x)={f}, F'(x)={fp}, diff={}",
(f - fp).abs()
);
}
_ => {
}
}
}
}
#[test]
fn integrate_asinh_direct() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let asinh_x = a.asinh(x);
let result = integrate(&mut a, asinh_x, x);
assert!(
!matches!(a.node(result), ExprNode::Integral(_, _)),
"asinh integration should return a closed form, not Integral"
);
let s = display(&a, result);
assert!(s.contains("asinh"), "result should contain asinh: {s}");
assert!(s.contains("sqrt"), "result should contain sqrt: {s}");
assert_ftc(
&mut a,
asinh_x,
result,
x,
&[(1, 2), (3, 2), (5, 1)],
1e-8,
"∫asinh(x)dx",
);
}
#[test]
fn integrate_acosh_direct() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let acosh_x = a.acosh(x);
let result = integrate(&mut a, acosh_x, x);
assert!(
!matches!(a.node(result), ExprNode::Integral(_, _)),
"acosh integration should return a closed form, not Integral"
);
let s = display(&a, result);
assert!(s.contains("acosh"), "result should contain acosh: {s}");
assert_ftc(
&mut a,
acosh_x,
result,
x,
&[(3, 2), (2, 1), (5, 1)],
1e-8,
"∫acosh(x)dx",
);
}
#[test]
fn integrate_atanh_direct() {
let mut a = Arena::new();
let x = sym(&mut a, "x");
let atanh_x = a.atanh(x);
let result = integrate(&mut a, atanh_x, x);
assert!(
!matches!(a.node(result), ExprNode::Integral(_, _)),
"atanh integration should return a closed form, not Integral"
);
let s = display(&a, result);
assert!(s.contains("atanh"), "result should contain atanh: {s}");
assert!(s.contains("ln"), "result should contain ln: {s}");
assert_ftc(
&mut a,
atanh_x,
result,
x,
&[(1, 4), (1, 2), (3, 4)],
1e-8,
"∫atanh(x)dx",
);
}
fn ftc_verdict(f: &str, big_f: &str) -> FtcVerdict {
let k = "(sin(1)^2 + cos(1)^2 - 1)";
let ctx = crate::api::context::Context::new();
let f = ctx.parse(&f.replace('K', k)).unwrap().id();
let big_f = ctx.parse(&big_f.replace('K', k)).unwrap().id();
ctx.with_arena_mut(|a| {
let x = a.symbol("x");
let ExprNode::Symbol(x_sym) = *a.node(x) else {
unreachable!()
};
check_antiderivative(a, f, big_f, x, x_sym)
})
}
#[test]
fn ftc_check_vanishing_denominator_is_undecided_not_passed() {
assert_eq!(
ftc_verdict("x/(K*x + 1)", "(x - ln(abs(K*x + 1))/K)/K"),
FtcVerdict::Undecided
);
assert_eq!(
ftc_verdict("x*exp(K*x)", "x*exp(K*x)/K - exp(K*x)/K^2"),
FtcVerdict::Undecided
);
assert_eq!(ftc_verdict("x/K", "x^2/(2*K)"), FtcVerdict::Untestable);
}
#[test]
fn ftc_check_samples_parameters() {
assert_eq!(ftc_verdict("1/(x + a)", "ln(x + a)"), FtcVerdict::Verified);
assert_eq!(
ftc_verdict("1/(x + a)", "ln(abs(x + a))"),
FtcVerdict::Wrong
);
assert_eq!(
ftc_verdict("1/(x + a)", "ln(abs(x + a))/a"),
FtcVerdict::Wrong
);
}
#[test]
fn ftc_check_complex_integrand() {
let f = "1/(x - sqrt(1/2 - sqrt(5)/2))";
assert_eq!(
ftc_verdict(f, "ln(x - sqrt(1/2 - sqrt(5)/2))"),
FtcVerdict::Verified
);
assert_eq!(
ftc_verdict(f, "ln(abs(x - sqrt(1/2 - sqrt(5)/2)))"),
FtcVerdict::Undecided
);
assert_eq!(ftc_verdict("1/x", "ln(abs(x))"), FtcVerdict::Verified);
assert_eq!(ftc_verdict("1/x", "ln(abs(x)) + x"), FtcVerdict::Wrong);
assert_eq!(
ftc_verdict("sqrt(x - 5)", "2/3*(x - 5)^(3/2)"),
FtcVerdict::Verified
);
}
#[test]
fn ftc_check_samples_through_binders() {
let f = "1/(x - RootOf(x^5 - x + 1, 1))";
assert_eq!(
ftc_verdict(f, "ln(x - RootOf(x^5 - x + 1, 1))"),
FtcVerdict::Verified
);
assert_eq!(
ftc_verdict(f, "2*ln(x - RootOf(x^5 - x + 1, 1))"),
FtcVerdict::Undecided
);
}
#[test]
fn contains_var_respects_binders() {
let ctx = crate::api::context::Context::new();
let cases = [
("RootOf(x^5 - x + 1, 0)", false),
("x - RootOf(x^5 - x + 1, 0)", true),
("RootOf(x^3 + a*x + 1, x, 0)", false),
("RootOf(x^3 + a*x + 1, a, 0)", true),
("Sum(x^k, k, 0, 3)", true),
("Sum(k*x, x, 0, 3)", false),
("Sum(k*x, x, 0, x)", true),
("Integral(x, x)", true),
];
for (src, expected) in cases {
let e = ctx.parse(src).unwrap().id();
let got = ctx.with_arena_mut(|a| {
let x = a.symbol("x");
let ExprNode::Symbol(x_sym) = *a.node(x) else {
unreachable!()
};
contains_var(a, e, x_sym)
});
assert_eq!(got, expected, "contains_var({src}, x)");
}
}
#[test]
fn constant_realness_is_exact_or_certified() {
let ctx = crate::api::context::Context::new();
ctx.symbol_with("r", &[crate::base::assumptions::Assumption::Real])
.unwrap();
ctx.symbol_with("s", &[crate::base::assumptions::Assumption::Positive])
.unwrap();
let cases = [
("sqrt(2)", true),
("a*b + 3", false),
("abs(a) + 3", true),
("r*s + 3", true),
("ln(3)", true),
("sqrt(3 - sqrt(5))", true),
("sqrt(1/2 - sqrt(5)/2)", false),
("exp(pi*I/8)", false),
("RootOf(x^5 - x + 1, 1)", false),
("sqrt(r)", false),
("sqrt(-r)", false),
("sqrt(s)", true),
("sqrt(a)", false),
];
for (src, expected) in cases {
let e = ctx.parse(src).unwrap().id();
let got = ctx.with_arena_mut(|a| {
let mut reals = crate::base::assumptions::AssumptionCache::new();
constant_is_real(a, e, &mut reals)
});
assert_eq!(got, expected, "constant_is_real({src})");
}
}
#[test]
fn analytic_in_parameters_finds_real_parameter_forms() {
let ctx = crate::api::context::Context::new();
for (src, analytic) in [
("ln(a*x + 1)/a", true),
("ln(abs(x)) + a*x", true),
("atan(x/sqrt(a))/sqrt(a) + exp(b*x)*sin(x)", true),
("ln(abs(a*x + 1))/a", false),
("sign(a)*x", false),
("abs(a)*x^2/2", false),
("re(a)*x", false),
] {
let e = ctx.parse(src).unwrap().id();
let got = ctx.with_arena_mut(|a| {
let x = a.symbol("x");
let ExprNode::Symbol(x_sym) = *a.node(x) else {
unreachable!()
};
analytic_in_parameters(a, e, x_sym)
});
assert_eq!(got, analytic, "{src}");
}
}
#[test]
fn ftc_check_binds_unassumed_parameters_to_complex_values() {
use crate::base::assumptions::Assumption;
let ctx = crate::api::context::Context::new();
ctx.symbol_with("r", &[Assumption::Real]).unwrap();
ctx.symbol_with("n", &[Assumption::Integer]).unwrap();
ctx.symbol_with("z", &[Assumption::Complex]).unwrap();
ctx.symbol_with("m", &[Assumption::Negative]).unwrap();
let cases = [
("1/(z*x + 1)", "ln(abs(z*x + 1))/z", FtcVerdict::Wrong),
("1/(z*x + 1)", "ln(z*x + 1)/z", FtcVerdict::Verified),
("1/(x - m)", "ln(abs(x - m))", FtcVerdict::Verified),
("sqrt(m^2)", "m*x", FtcVerdict::Wrong),
("1/(a*x + 1)", "ln(abs(a*x + 1))/a", FtcVerdict::Wrong),
("1/(a*x + 1)", "ln(a*x + 1)/a", FtcVerdict::Verified),
("1/(r*x + 1)", "ln(abs(r*x + 1))/r", FtcVerdict::Verified),
("1/(x - 2)", "ln(abs(x - 2))", FtcVerdict::Verified),
("a*tan(x)", "-a*ln(abs(cos(x)))", FtcVerdict::Verified),
(
"a/sqrt(x^2 - 1)",
"a*ln(abs(x + sqrt(x^2 - 1)))",
FtcVerdict::Verified,
),
("a", "abs(a)*x", FtcVerdict::Wrong),
("exp(a*x)", "exp(a*x)/a", FtcVerdict::Verified),
("x^n", "x^(n + 1)/(n + 1)", FtcVerdict::Untestable),
];
for (f_src, big_f_src, expected) in cases {
let f = ctx.parse(f_src).unwrap().id();
let big_f = ctx.parse(big_f_src).unwrap().id();
let got = ctx.with_arena_mut(|a| {
let x = a.symbol("x");
let ExprNode::Symbol(x_sym) = *a.node(x) else {
unreachable!()
};
check_antiderivative(a, f, big_f, x, x_sym)
});
assert_eq!(got, expected, "∫ {f_src} = {big_f_src}");
}
}
}