Polynomials as Data
An Ex is a tree. When you know an expression is a polynomial in some symbols, you usually want a different view of it: a finite list of (monomial, coefficient) pairs that you can index, iterate, multiply, evaluate exactly, and lay out as a matrix. In 0.3 that view is Poly (symplex::poly_ex::Poly), reached from any expression with as_poly(&[&x, &y]).
Two things distinguish Poly from the rational-coefficient machinery in symplex::multipoly:
- Coefficients are
Ex. They may be exact rationals or symbolic parameters — anything free of the generators.a·x² + (a + b)·x + 3is a perfectly good polynomial inx. - Nothing is approximated. Every operation is exact, and
to_ex()rebuilds an expression equal to the (expanded) input.
Terms are always reported in descending lexicographic order of the exponent vectors — the order of SymPy’s Poly.terms().
Viewing an expression as a polynomial
as_poly (equivalently Poly::new) expands the expression and collects it by monomial. It returns None if a generator appears in a non-polynomial position — inside a function, under a negative or fractional power, or in an exponent.
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, y);
let e = (&x + &y * 2).powi(2) * &x - &y.powi(3);
let p = e.as_poly(&[&x, &y]).unwrap();
println!("{p}"); // Poly(x^3 + 4*x^2*y + 4*x*y^2 - y^3, x, y)
for (mono, coeff) in p.terms() {
println!("x^{} y^{} · {coeff}", mono[0], mono[1]);
}
// x^3 y^0 · 1
// x^2 y^1 · 4
// x^1 y^2 · 4
// x^0 y^3 · -1
println!("{:?}", p.monoms()); // [[3, 0], [2, 1], [1, 2], [0, 3]]
println!("{:?}", p.coeffs()); // [Ex(1), Ex(4), Ex(4), Ex(-1)]
println!("{}", p.coeff_monomial(&[1, 2]).unwrap()); // 4
println!("{}", p.coeff_monomial(&[5, 0]).unwrap()); // 0 (absent monomials are zero)
println!("{:?} {:?} {:?}", p.total_degree(), p.degree_in(&x), p.degree_list());
// Some(3) Some(3) [3, 3]
println!("{} {:?}", p.leading_coeff(), p.leading_monomial()); // 1 Some([3, 0])
println!("{} {} {}", p.num_terms(), p.is_homogeneous(), p.has_rational_coeffs());
// 4 true true
println!("{}", p.to_ex()); // x^3 + 4*x*y^2 - y^3 + 4*y*x^2
assert!(x.sin().as_poly(&[&x]).is_none()); // generator inside a function
assert!((ctx.int(1) / &x).as_poly(&[&x]).is_none()); // negative power
assert!(x.pow(&y).as_poly(&[&x]).is_none()); // generator in an exponent
}
Other structural queries: is_zero, is_ground (constant), is_univariate, is_linear, gens(), num_gens(), leading_term(), and equals(&other) (same generators, identical normalised coefficients). Poly::from_terms(&ctx, &gens, vec![(exps, coeff), …]), Poly::zero, Poly::one and Poly::constant build polynomials directly.
Symbolic coefficients
Any symbol that is not a generator becomes part of the coefficients. The same expression can be viewed with different generator lists, and the plain Ex methods degree, coeffs, coeff, leading_coeff and is_polynomial now accept parameter coefficients too (in 0.2 they required rational coefficients).
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, a, b);
let e = &a * &x.powi(2) + (&a + &b) * &x + &x.powi(2) + 3;
println!("{e}"); // a*x^2 + x^2 + x*(a + b) + 3
// Ex methods — ascending order, like 0.2:
println!("{:?}", e.degree(&x)); // Some(2)
let cs: Vec<String> = e.coeffs(&x).unwrap().iter().map(|c| c.to_string()).collect();
println!("{cs:?}"); // ["3", "a + b", "a + 1"]
println!("{}", e.coeff(&x, 1).unwrap()); // a + b
println!("{}", e.leading_coeff(&x).unwrap()); // a + 1
println!("{}", e.is_polynomial(&x)); // true
// Poly view in x alone — all_coeffs is dense and highest-degree first (SymPy order):
let p = e.as_poly(&[&x]).unwrap();
let dense: Vec<String> = p.all_coeffs().unwrap().iter().map(|c| c.to_string()).collect();
println!("{dense:?}"); // ["a + 1", "a + b", "3"]
println!("{}", p.has_rational_coeffs()); // false
// Promote a to a generator: now b is the only parameter.
let q = e.as_poly(&[&x, &a]).unwrap();
for (mono, coeff) in q.terms() {
println!("x^{} a^{} · {coeff}", mono[0], mono[1]);
}
// x^2 a^1 · 1
// x^2 a^0 · 1
// x^1 a^1 · 1
// x^1 a^0 · b
// x^0 a^0 · 3
println!("{:?}", x.pow(&a).degree(&x)); // None (x^a is not polynomial in x)
}
Note the two orderings: Ex::coeffs is ascending ([a₀, a₁, …], unchanged from 0.2), while Poly::all_coeffs is descending with zeros filled in, matching SymPy’s all_coeffs().
Exact evaluation and arithmetic
eval substitutes a value for every generator and evaluates; values may be rationals, radicals or expressions. eval_gen substitutes one generator (by a constant or a polynomial in the remaining generators) and returns a Poly with one generator fewer. Arithmetic (add, sub, mul, neg, scale, pow, derivative) requires identical generator lists and returns Result.
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, y);
let p = (&x.powi(2) * &y - &y / 2 + ctx.rational(1, 3)).as_poly(&[&x, &y]).unwrap();
println!("{}", p.eval(&[&ctx.rational(3, 2), &ctx.rational(-4, 5)]).unwrap()); // -16/15
println!("{}", p.eval(&[&ctx.int(2).sqrt(), &ctx.int(1)]).unwrap()); // 11/6
println!("{}", p.eval_gen(&x, &ctx.int(2)).unwrap()); // Poly(7/2*y + 1/3, y)
println!("{}", p.eval_gen(&x, &(&y + 1)).unwrap()); // Poly(y^3 + 2*y^2 + 1/2*y + 1/3, y)
println!("{}", p.derivative(&y).unwrap().to_ex()); // x^2 - 1/2
let s = (&x + &y).as_poly(&[&x, &y]).unwrap();
let d = (&x - &y).as_poly(&[&x, &y]).unwrap();
println!("{}", s.mul(&d).unwrap().to_ex()); // x^2 - y^2
println!("{}", s.pow(3).unwrap().to_ex()); // x^3 + 3*x*y^2 + y^3 + 3*y*x^2
println!("{}", s.add(&d).unwrap().to_ex()); // 2*x
println!("{}", s.sub(&d).unwrap().to_ex()); // 2*y
println!("{}", s.neg().to_ex()); // -x - y
println!("{}", s.scale(&ctx.rational(1, 2)).unwrap().to_ex()); // 1/2*x + 1/2*y
assert!(s.mul(&d).unwrap().equals(&(&x.powi(2) - &y.powi(2)).as_poly(&[&x, &y]).unwrap()));
assert!(s.add(&(&x + 1).as_poly(&[&x]).unwrap()).is_err()); // different generators
// Rational-coefficient helpers
let q = (&x.powi(2) * -4 + &x * 6).as_poly(&[&x]).unwrap();
let (c, prim) = q.content_and_primitive().unwrap();
println!("{c} · ({})", prim.to_ex()); // -2 · (2*x^2 - 3*x)
println!("{}", q.monic().unwrap().to_ex()); // x^2 - 3/2*x
}
Since 0.18 a Poly whose coefficients are all rational literals is held as an exact MultiPoly<Lex> rather than as one Ex per monomial, and add/sub/mul/pow/scale/derivative/eval/eval_gen between such polynomials run on rationals without touching the expression arena (products and powers accumulate integer numerators over a common denominator and reduce once). The representation is invisible: terms(), Display, to_ex() and to_multipoly() report exactly what they did before, a symbolic coefficient switches the polynomial to the Ex form, and mixed operations convert the exact side. Expect roughly 5–15× on dense products and 10–30× on exact evaluation compared with 0.17.
Coefficient matrices and exact linear systems
The reason Poly exists is to make questions like “is goal a linear combination of h₁, …, hₖ?” mechanical. Poly::monomial_basis collects every monomial that occurs in a family, and Poly::coefficient_matrix lays the family out with one row per monomial and one column per polynomial. The unknown multipliers λ then satisfy M·λ = coefficients of goal, which linsolve_matrix solves exactly — including the under- and over-determined cases.
use symplex::prelude::*;
use symplex::poly_ex::Poly;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x);
let h1 = (&x + 1).as_poly(&[&x]).unwrap();
let h2 = (&x.powi(2) - 1).as_poly(&[&x]).unwrap();
let goal = (&x + 1).powi(2).as_poly(&[&x]).unwrap();
let basis = Poly::monomial_basis(&[&h1, &h2, &goal]).unwrap();
println!("{basis:?}"); // [[2], [1], [0]]
let m = Poly::coefficient_matrix(&[&h1, &h2], &basis).unwrap();
println!("{m}");
let rhs: Vec<Ex> = basis.iter().map(|mono| goal.coeff_monomial(mono).unwrap()).collect();
match linsolve_matrix(&m, &Matrix::col_vector(rhs)).unwrap() {
LinearSolution::Unique(pairs) => {
for (var, val) in pairs {
println!("{var} = {val}"); // x1 = 2, x2 = 1
}
}
other => println!("{other:?}"),
}
// So (x + 1)² = 2·(x + 1) + 1·(x² − 1).
// x² + x + 1 is not in the span:
let goal2 = (&x.powi(2) + &x + 1).as_poly(&[&x]).unwrap();
let rhs2: Vec<Ex> = basis.iter().map(|mono| goal2.coeff_monomial(mono).unwrap()).collect();
println!("{:?}", linsolve_matrix(&m, &Matrix::col_vector(rhs2)).unwrap()); // Inconsistent
}
Output of the matrix:
[
[0, 1],
[1, 0],
[1, -1]
]
When the multipliers must be non-negative — the situation in every positivity certificate — feed the same matrix to the exact LP solver instead: Matrix::to_rational_rows() gives the rows in the form linprog::feasible_nonneg wants. The Polynomial Inequality Certificates cookbook entry does this end to end, and Exact Linear Programming describes the solver.
Rational normal form: ratsimp
ratsimp puts a rational expression into a canonical P/Q: one fraction, common factors cancelled by a multivariate GCD, integer-primitive numerator and denominator, and a positive leading coefficient in Q. Maximal non-rational subexpressions (sin x, π, √x) are treated as opaque indeterminates, exactly as SymPy’s cancel does. simplify_rational is now the same normal form.
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, y, r, j);
println!("{}", ((&x.powi(2) - &y.powi(2)) / (&x - &y)).ratsimp()); // x + y
println!("{}", (ctx.int(1) / &x + ctx.int(1) / &y).ratsimp()); // (x + y)/(x*y)
let nested = ctx.int(1) / (&x + ctx.int(1) / &y) + ctx.int(1) / (&y + ctx.int(1) / &x);
println!("{nested} → {}", nested.ratsimp());
// 1/(x + 1/y) + 1/(1/x + y) → (x + y)/(x*y + 1)
println!("{}", (x.sin().powi(2) / x.sin()).ratsimp()); // sin(x)
println!("{}", ((x.sin().powi(2) - x.cos().powi(2)) / (x.sin() - x.cos())).ratsimp());
// sin(x) + cos(x) — a difference of squares, no trig identity involved
let (n, d) = (ctx.int(1) / &x + ctx.int(1) / (&x + 1)).ratsimp().as_numer_denom();
println!("{n} / {d}"); // 2*x + 1 / x^2 + x
let mixed = (&x.powi(2) * 2 + &x * 4) / (&x * 6 + 12) + ctx.rational(1, 3);
println!("{mixed} → {}", mixed.ratsimp());
// (2*x^2 + 4*x)/(6*x + 12) + 1/3 → 1/3*x + 1/3
println!("{}", (&x + 1).ratsimp()); // x + 1 (already normal)
println!("{}", x.exp().ratsimp()); // exp(x) (unchanged)
// `solve` with parameter coefficients returns ratsimp'd solutions.
let eqn = (&r * 3 - 1) / (&j + 1) - (&r + 1) / (&j * 2);
println!("{}", eqn.solve(&r).unwrap()[0]); // (3*j + 1)/(5*j - 1)
println!("{}", eqn.simplify_rational()); // (5*j*r - 3*j - r - 1)/(2*j^2 + 2*j)
}
ratsimp is the right tool for checking an identity: (lhs − rhs).ratsimp() is structurally 0 exactly when the two sides agree as rational functions. cancel(&x) (single variable) and together() (no cancellation) remain available for lighter-weight jobs.
Sign of a polynomial on an interval
poly_is_nonnegative_on(&x, &lo, &hi) and poly_is_positive_on decide, exactly, whether a univariate polynomial with rational coefficients is ≥ 0 (resp. > 0) on the closed interval [lo, hi]. The method is a square-free decomposition (to find the roots where the sign can change), a Sturm count to check that none lies strictly inside the interval, and one sample point. Endpoints must be rationals or ±∞. The answer is None for non-polynomial input or symbolic coefficients.
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, a);
let (ninf, inf) = (ctx.neg_infinity(), ctx.infinity());
let sq = &x.powi(2) - &x * 2 + 1; // (x − 1)²
println!("{:?} {:?}",
sq.poly_is_nonnegative_on(&x, &ninf, &inf),
sq.poly_is_positive_on(&x, &ninf, &inf)); // Some(true) Some(false)
let cubic = &x.powi(3) - &x;
println!("{:?} {:?} {:?}",
cubic.poly_is_nonnegative_on(&x, &ctx.int(2), &inf),
cubic.poly_is_nonnegative_on(&x, &ctx.int(-2), &inf),
cubic.poly_is_nonnegative_on(&x, &ctx.int(-1), &ctx.int(0))); // Some(true) Some(false) Some(true)
let wobble = &x.powi(4) - &x.powi(2) * 5 + 4; // (x² − 1)(x² − 4)
println!("{:?} {:?} {:?}",
wobble.poly_is_nonnegative_on(&x, &ctx.rational(-1, 1), &ctx.rational(1, 1)),
wobble.poly_is_positive_on(&x, &ctx.rational(-1, 1), &ctx.rational(1, 1)),
wobble.poly_is_nonnegative_on(&x, &ctx.rational(3, 2), &ctx.rational(7, 4)));
// Some(true) Some(false) Some(false)
println!("{:?} {:?}",
(&x.powi(2) + 1).poly_is_positive_on(&x, &ninf, &inf),
(&a * &x + 1).poly_is_nonnegative_on(&x, &ninf, &inf)); // Some(true) None
println!("{:?}", x.sin().poly_is_nonnegative_on(&x, &ninf, &inf)); // None
println!("{:?}", wobble.count_real_roots(&x)); // Some(4)
}
For the roots themselves use count_real_roots, real_roots_isolate and nroots (see Algebra).
Numeric roots
Poly::nroots(digits) is Ex::nroots for the univariate, rational-coefficient case; anything else is an InvalidArgument error rather than a wrong answer. Roots come back as Vec<Complex64> (num_complex, re-exported in the prelude); real roots have im == 0.0 exactly.
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, a);
let p = (&x.powi(5) - &x - 1).as_poly(&[&x]).unwrap();
for Complex64 { re, im } in p.nroots(12).unwrap() {
println!("{re:.10} {im:+.10}i");
}
// -0.7648844336 -0.3524715460i
// -0.7648844336 +0.3524715460i
// 0.1812324445 -1.0839541013i
// 0.1812324445 +1.0839541013i
// 1.1673039783 +0.0000000000i
println!("{}", (&a * &x + 1).as_poly(&[&x]).unwrap().nroots(10).unwrap_err());
// Poly::nroots: invalid argument: polynomial must have rational coefficients
println!("{}", (&x * &a).as_poly(&[&x, &a]).unwrap().nroots(10).unwrap_err());
// Poly::nroots: invalid argument: polynomial must be univariate
}
Bridge to MultiPoly and Gröbner bases
symplex::groebner and symplex::polysys work on MultiPoly<GrevLex>, a rational-coefficient sparse polynomial indexed by variable position. to_multipoly() converts a Poly whose coefficients are all rational (None otherwise), and Poly::from_multipoly(&ctx, &gens, &mp) converts back, so you can move between the two worlds without touching the low-level representation.
use symplex::prelude::*;
use symplex::groebner::groebner_basis;
use symplex::poly_ex::Poly;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, y);
let f1 = (&x.powi(2) + &y.powi(2) - 1).as_poly(&[&x, &y]).unwrap();
let f2 = (&x - &y).as_poly(&[&x, &y]).unwrap();
let gb = groebner_basis(&[f1.to_multipoly().unwrap(), f2.to_multipoly().unwrap()]);
for g in &gb {
println!("{}", Poly::from_multipoly(&ctx, &[&x, &y], g).unwrap().to_ex());
}
// y^2 - 1/2
// x - y
// Normal form of x³ modulo the ideal:
let target = x.powi(3).as_poly(&[&x, &y]).unwrap().to_multipoly().unwrap();
let refs: Vec<_> = gb.iter().collect();
let rem = target.reduce(&refs);
println!("{}", Poly::from_multipoly(&ctx, &[&x, &y], &rem).unwrap().to_ex()); // 1/2*y
let a = ctx.symbol("a");
println!("{:?}", (&a * &x).as_poly(&[&x]).unwrap().to_multipoly().is_none()); // true
}
symplex::polysys::solve_system_ex (see Solving Equations) does the whole pipeline — Gröbner basis, triangularisation, algebraic back-substitution — when what you want is the solution set rather than the basis.
See cargo run --example polynomials for the complete program these snippets are drawn from.
Algebraic numbers and polynomial algebra (0.9)
0.9 adds a layer of Ex methods over the exact engines above, so the common polynomial-algebra questions no longer need a detour through Poly/MultiPoly. Every method takes and returns Ex; the ones that answer a query return Option (None when the input is not of the required shape), the ones that validate caller-supplied structure return Result.
Minimal polynomials
minimal_polynomial(&var) (SymPy minimal_polynomial) returns the minimal polynomial over ℚ of an algebraic constant built from rationals, radicals, i, φ, sums, products, integer powers and reciprocals. The result is integer-primitive with a positive leading coefficient, exactly as SymPy prints it; None means the number was not recognised as algebraic (π, e, free symbols, transcendental functions).
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
let x = ctx.symbol("x");
let a = ctx.int(2).sqrt() + ctx.int(3).sqrt();
println!("{}", a.minimal_polynomial(&x).unwrap()); // x^4 - 10*x^2 + 1
let cbrt2 = ctx.int(2).pow(&ctx.rational(1, 3));
println!("{}", cbrt2.minimal_polynomial(&x).unwrap()); // x^3 - 2
println!("{}", ctx.rational(3, 4).minimal_polynomial(&x).unwrap()); // 4*x - 3
let b = ctx.int(1) / (ctx.int(1) + ctx.int(2).sqrt());
println!("{}", b.minimal_polynomial(&x).unwrap()); // x^2 + 2*x - 1
assert!(ctx.pi().minimal_polynomial(&x).is_none());
}
Multivariate gcd and lcm without naming variables
gcd_all / lcm_all (SymPy gcd(f, g) / lcm(f, g)) treat both inputs as polynomials over ℚ in all of their free symbols. The result follows MultiPoly::gcd‘s normalisation: integer coefficients, positive leading coefficient (grevlex), integer content equal to the gcd of the inputs’ contents — for polynomials over ℤ that is the ordinary gcd over ℤ, gcd(2x, 4x) = 2x. Non-polynomial input (sin x, 1/x, π·x) gives None.
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, y, z);
let f = &x.powi(2) - &y.powi(2);
println!("{}", f.gcd_all(&(&x - &y)).unwrap()); // x - y
println!("{}", f.lcm_all(&(&x - &y)).unwrap()); // x^2 - y^2
println!("{}", (&x * &y * &z + &x * &y).gcd_all(&(&x * &z + &x)).unwrap()); // x*z + x
println!("{}", (&x * 2).gcd_all(&(&x * 4)).unwrap()); // 2*x
assert!(x.sin().gcd_all(&x).is_none());
}
Gröbner bases and normal forms from Ex
Ex::groebner(&polys, &vars, order) computes the reduced (monic) Gröbner basis in the given variables under MonomialOrder::Lex or MonomialOrder::GrevLex (symplex::multipoly::MonomialOrder); reduce_modulo(&basis, &vars, order) is the remainder of multivariate division — the unique normal form when basis is a Gröbner basis for that order, so it is zero exactly for members of the ideal. Variables must be distinct symbols and every polynomial must have rational coefficients; anything else is an InvalidArgument error naming the offending input.
use symplex::multipoly::MonomialOrder;
use symplex::prelude::*;
fn main() -> Result<(), SymplexError> {
let ctx = Context::new();
symplex::syms!(ctx; x, y);
let vars = [x.clone(), y.clone()];
let gens = [&x.powi(2) + &y.powi(2) - 1, &x - &y];
let lex = Ex::groebner(&gens, &vars, MonomialOrder::Lex)?;
for g in &lex {
println!("{g}");
}
// x - y
// y^2 - 1/2
let grevlex = Ex::groebner(&gens, &vars, MonomialOrder::GrevLex)?; // [y^2 - 1/2, x - y]
println!("{}", x.powi(3).reduce_modulo(&lex, &vars, MonomialOrder::Lex)?); // 1/2*y
let member = (&x.powi(2) - &y.powi(2)).reduce_modulo(&grevlex, &vars, MonomialOrder::GrevLex)?;
assert!(member.is_zero_structural());
Ok(())
}
Exact real roots
real_roots(&var) lists the distinct real roots of a rational-coefficient polynomial in increasing order (decided exactly with Sturm sequences): rational roots as numbers, every other root as the RootOf(g, k) node that solve already uses for degree ≥ 5, where g is the irreducible factor over ℤ and k its index among g’s complex roots. RootOf evaluates numerically (eval_f64, eval_decimal) and prints as such; root_of(&var, k) is the k-th real root (0-based). Unlike SymPy’s real_roots, a repeated root is listed once (as in count_real_roots).
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
let x = ctx.symbol("x");
let roots = (&x.powi(3) - &x * 2).real_roots(&x).unwrap();
for r in &roots {
println!("{r} ≈ {}", r.eval_f64().unwrap());
}
// RootOf(x^2 - 2, 0) ≈ -1.41421356…
// 0 ≈ 0
// RootOf(x^2 - 2, 1) ≈ 1.41421356…
let largest = (&x.powi(3) - &x * 2).root_of(&x, 2).unwrap(); // RootOf(x^2 - 2, 1)
assert_eq!(largest, roots[2]);
assert_eq!((&x.powi(2) + 1).real_roots(&x), Some(vec![])); // no real roots
}
Factoring modulo a prime
factor_mod(&var, p) (SymPy factor_list(f, modulus=p)) factors a rational-coefficient polynomial over GF(p) into (lc, [(monic irreducible factor, multiplicity)]) with coefficients in [0, p). A composite p is an InvalidArgument error, as is a coefficient whose denominator is divisible by p; the finite-field arithmetic supports odd primes below 2³¹.
use symplex::prelude::*;
fn main() -> Result<(), SymplexError> {
let ctx = Context::new();
let x = ctx.symbol("x");
let (lc, factors) = (&x.powi(2) + 1).factor_mod(&x, 5)?;
println!("{lc}: {:?}", factors.iter().map(|(f, m)| format!("({f})^{m}")).collect::<Vec<_>>());
// 1: ["(x + 2)^1", "(x + 3)^1"]
let (_, factors) = (&x.powi(2) + 1).factor_mod(&x, 3)?; // irreducible mod 3
assert_eq!(factors.len(), 1);
assert!((&x.powi(2) + 1).factor_mod(&x, 6).is_err());
Ok(())
}
Resultants and discriminants with symbolic coefficients
resultant and discriminant (0.2) require every coefficient to be rational. resultant_symbolic and discriminant_symbolic accept parameter coefficients: they build the Sylvester matrix over Ex entries, take its determinant with Matrix::det, and expand, so the answer is a polynomial in the parameters. discriminant_symbolic is division-free (the leading coefficient is eliminated with one row operation before the determinant).
use symplex::prelude::*;
fn main() {
let ctx = Context::new();
symplex::syms!(ctx; x, a, b, c, p, q);
let quad = &a * &x.powi(2) + &b * &x + &c;
println!("{}", quad.discriminant_symbolic(&x).unwrap()); // -4*a*c + b^2
let cubic = &x.powi(3) + &p * &x + &q;
println!("{}", cubic.discriminant_symbolic(&x).unwrap()); // -4*p^3 - 27*q^2
println!("{}", (&x - &a).resultant_symbolic(&(&x - &b), &x).unwrap()); // a - b
}