Solver API

alkahest.solve(equations, vars, *, numeric=False, method="groebner")

Solve a zero-dimensional polynomial — or rational — system.

An equation may be a ratio of polynomials in vars, which is the form nodal analysis produces: solve([(Vo - Vin)/R1 + Vo*s*C], [Vo]) → Vin/(1 + s*R1*C). Each equation is put over a common denominator and the numerator system is solved.

Clearing a denominator is not an equivalence — N/D = 0 means N = 0 and D != 0 — so every root is re-tested against the denominators that were multiplied through, and one that makes any of them vanish is dropped: solve([x/(x-1) - 1/(x-1)], [x]) returns [], not x = 1. Where the denominator still mentions a free parameter the question is not decidable, and the root is returned under a D != 0 hypothesis listed by solve_side_conditions() rather than under a silent assumption. Surviving roots are also substituted into the equations as you wrote them, before anything was cleared.

An equation whose denominator is identically zero (1/(x - x)) denotes no function at all and raises SolverError with .code == "E-SOLVE-005".

Only the unknowns must appear polynomially: a sub-expression free of them (sin(q1), cos(q1 - q2), sqrt(a)) is a coefficient, so solve(eom, [qdd1, qdd2]) solves equations of motion for the accelerations. Such sub-expressions are solved as opaque parameters and substituted back; every answer is then checked against the equations as written, and one that only exists because an identity between them was invisible (sin**2 + cos**2 == 1) is refused with E-SOLVE-001.

One equation in one unknown is solved in closed form, over the reals, when it reduces by taking logs (ka*exp(-ka*t) - ke*exp(-ke*t) → log(ke/ka)/(ke - ka); 2**x - 5), has a Lambert-W shape (x*exp(a*x) - b, a*x + b*log(x) - c, x**x - c), or is sin/cos/tan of an affine argument equal to a constant. A trigonometric equation returns the whole periodic family in an integer symbol k (asin(c) + 2*pi*k and pi - asin(c) + 2*pi*k), listed as k ∈ Integer by solve_side_conditions(); with numeric=True the k = 0 member of each family is returned and the side conditions say so. Hypotheses such as ke/ka > 0 are listed there too. A Lambert equation with two real roots returns both — lambert_w (W₀) and lambert_wm1 (W₋₁): x*exp(x) + 1/4 gives W₀(-1/4) and W₋₁(-1/4); a symbolic argument known to be negative records c + exp(-1) > 0. Anything else keeps refusing with E-SOLVE-001 — with numeric=True the message points at interval=(a, b) and nsolve().

equationslist[Expr]

Each expression represents p(vars) = 0.

varslist[Expr]

Variables to solve for (symbols). Free symbols that appear in equations but are omitted from vars are treated as parameters; solutions may be expressions in those symbols (e.g. solve([x**2 - y], [x])` → `±sqrt(y)).

numericbool, default False

Used when method="groebner": symbolic Expr values vs numbers. The domain is ℂ, so a real value is a float and a non-real one a Python complex (solve([s**2 + 1], [s], numeric=True) → ±1j). When Lex back-substitution refuses a degree ≥ 3 univariate (E-SOLVE-002, see below), numeric=True falls back to complex homotopy continuation, which returns every finite root or raises HomotopyError (E-HOMOTOPY-004) if a path is lost. A solution that depends on a free parameter has no numeric value, so it is returned as its exact symbolic Expr (solve([x - y], [x], numeric=True) → [{x: y}]) and solve_side_conditions() says so; where no symbolic solution exists (the homotopy fallback) a free parameter still raises SolverError with .code == "E-SOLVE-006". A parameter-free value that cannot be evaluated raises E-SOLVE-007.

methodstr, default "groebner"

"groebner" — Lex basis + triangular back-substitution. "homotopy" — total-degree homotopy continuation in ℂⁿ followed by Newton projection to real tuples (always float dict values).

interval(float, float), optional

With numeric=True and one equation in one unknown: return every real root in [a, b], found numerically and certified (see nsolve() with all_in) — the route for equations with no closed form, e.g. solve([sin(x) - x/10], [x], numeric=True, interval=(-12, 12)) → seven roots. Each root's enclosure is listed in solve_side_conditions(). When the root list cannot be proved complete (a double root, a stretch outside the equation's domain) the call raises E-SOLVE-012 instead.

list[dict]

Each dict maps a variable Expr to Expr (symbolic Groebner), float or complex (Groebner with numeric=True), or float (method="homotopy", which reports real roots only). A univariate of degree ≥ 3 is factored over ℚ (or ℚ[params]) first; each irreducible factor is solved by formula, by an even / palindromic / pure-cubic reduction to radicals, or as exact RootOf(f, k) values (see root_of() for their order and evaluation). An irreducible factor of degree ≥ 3 that still involves a parameter, or has algebraic-number coefficients, raises E-SOLVE-002. Solutions are a set: a double root is one entry, not two. Every parameter-free tuple has been substituted back into the equations you passed and could not be shown to violate them. A tuple containing a free parameter is the generic solution: it has been substituted back and its residual is identically zero in the parameters (equations that factor are solved factor by factor, so equilibrium systems are complete). It holds under the non-vanishing hypotheses reported by solve_side_conditions() — solve([a*x - b], [x]) is b/a for a ≠ 0, and that condition is listed there. An unknown declared Positive, NonNegative, NonZero or Integer gets no value provably outside its domain.

GroebnerBasis

When method="groebner" and no finite solution list could be produced — usually a positive-dimensional ideal, but also when the Lex basis admits no complete triangular elimination in vars. It is "here is the ideal" rather than a claim that the solutions are infinite.

alkahest.nsolve(equation, x, *, bracket=None, x0=None, all_in=None)

Numeric real root(s) of equation = 0 in the symbol x, for equations with no closed form — implied volatility BS(sigma) - price, cos(x) - x, sin(x) - x/10. Every other symbol must already have a value (substitute it first). Pass exactly one of:

bracket : (float, float) (a, b) with f(a), f(b) of opposite sign. Brent's method; always converges, to a root in [a, b]. Returns a NumericRoot. x0 : float A starting point. Damped Newton; which root it reaches depends on x0. Returns a NumericRoot. all_in : (float, float) Every real root in [a, b], by subdivision with ball-arithmetic exclusion and an interval-derivative (monotonicity) test. Returns a RootIsolation; its roots are all "verified", and it is complete unless some sub-interval (a double root, a stretch outside the domain) could not be settled — those are listed in unresolved, never guessed at. Needs an equation built from smooth functions (no abs/floor/Piecewise of x).

Each root carries an enclosure [lo, hi] and a status: the search runs in f64, then the enclosure is proved in 128-bit outward-rounded ball arithmetic — "verified" means f changes sign across [lo, hi] and its derivative's enclosure excludes zero there, so exactly one root lies inside; "exists" drops uniqueness; "heuristic" is the bare iterate (a double root, a non-smooth equation). Loops honour an active Budget.

SolverError

E-SOLVE-009 — a free parameter, an empty/non-finite interval, a bracket without a sign change, all_in on a non-smooth equation; E-SOLVE-012 — Newton/Brent did not converge or the equation has no real value where the search needed one.

BudgetExceededError

The active budget stopped the search.

Example:

PythonCopy

>>> r = alkahest.nsolve(alkahest.cos(x) - x, x, bracket=(0, 1))
>>> r.value, r.status
(0.7390851332151607, 'verified')
>>> len(alkahest.nsolve(alkahest.sin(x) - x / 10, x, all_in=(-12, 12)))
7

A numeric root returned by nsolve(), with an enclosure and a statement of what the enclosure proves.

status is "verified" (exactly one root lies in [lo, hi]: a sign change proved in 128-bit ball arithmetic and an interval enclosure of the derivative that excludes zero), "exists" (a root lies in [lo, hi], uniqueness not proved) or "heuristic" (the f64 iterate only — lo == hi == value). float(root) is root.value.

value

lo

hi

enclosure

status

verified

Every real root of an equation in an interval, from alkahest.nsolve(eq, x, all_in=(a, b)).

roots lists NumericRoot values in increasing order, each "verified" with a disjoint enclosure. unresolved lists the sub-intervals (lo, hi) that could neither be excluded nor shown to hold exactly one root — a multiple root ((x - 1)**2), a tangency, a stretch outside the equation's domain. complete is True exactly when unresolved is empty, i.e. when roots is provably the whole list. Indexing, len and iteration go over roots.

roots

unresolved

complete

The index-th distinct root of the polynomial poly in var, exactly.

solve returns these for irreducible factors of degree ≥ 3 that have no radical form it uses (an irreducible quintic, a cubic with three real roots, …); this function builds the same objects directly.

Roots are numbered in one fixed order: the real roots in ascending order, then the non-real roots by ascending real part, ties broken by ascending imaginary part (a conjugate pair is listed lower half-plane first). Repeated roots count once.

polyExpr

A non-constant polynomial in var with rational coefficients.

varExpr

The polynomial variable.

indexint

0 <= index < number of distinct roots.

Expr

A rational number when the root is rational; otherwise RootOf(f, k), where f is the irreducible factor of poly the root belongs to and k its index among f's roots. A RootOf has no free symbols; it is a number to diff, subs and simplify, and it evaluates with eval_expr (real roots), evaluate(..., mode="complex"), and to any precision with interval_eval / evaluate(..., precision_bits=...) (real roots, rigorous enclosures from FLINT's root isolation).

SolverError

E-SOLVE-008 when poly is not a non-constant polynomial in var over ℚ, or index is out of range.

Example:

pool = ExprPool() x = pool.symbol("x") r = root_of(x**3 - x - 1, x, 0) # RootOf(x^3 - x - 1, 0) eval_expr(r, {}) # 1.3247179572447...

A computed Gröbner basis for a polynomial ideal.

The basis is a sequence: it supports len(), integer indexing and iteration, yielding its generators as GbPoly. Read the generators back as expressions with to_exprs(), or one at a time with GbPoly.to_expr() — that is how an elimination ideal (a "lex" basis, or the differential elimination performed by rosenfeld_groebner()) is turned into readable relations.

orderstr

Monomial order the generators are reduced under: "lex", "grlex" or "grevlex".

Example:

gb = alkahest.GroebnerBasis.compute([x2 + y2 - one, x - y], [x, y]) len(gb) # number of generators [g.to_expr() for g in gb] # generators as Expr

order

proof_status

compute(polys, vars, order=None, params=None)

compute_f5(polys, vars, order=None)

compute_raw(gb_polys, order=None)

eliminate(vars)

# Implicitize (t, t**2): eliminate t from {x - t, y - t**2}
gb = alkahest.GroebnerBasis.compute([x - t, y - t**2], [t, x, y])
gb.eliminate([t]).to_exprs()   # [((y * -1) + x^2)]  i.e. y = x**2

variables() → list[Expr]

polynomials() → list[GbPoly]

to_exprs(vars=None)

reduce(p)

contains(p)

A sparse multivariate polynomial over ℚ, as used by the Gröbner machinery.

A GbPoly stores exponent vectors, not variable names, so reading one back needs the variable list its exponent slots refer to. Polynomials handed out by Alkahest carry that list with them (see variables()), so to_expr() normally takes no arguments; pass vars explicitly only for a polynomial built against a different list.

is_zerobool

n_varsint

n_termsint

is_zero

n_vars

n_terms

variables() → list[Expr]

terms() → list[tuple[tuple[int, ...], int | Fraction]]

p = alkahest.expr_to_gbpoly(x**2 - 3*y, [x, y])
p.terms()   # [((0, 1), -3), ((2, 0), 1)]

to_expr(vars=None)

alkahest.expr_to_gbpoly(expr, vars) — convert a polynomial Expr into the GbPoly representation.

The inverse of GbPoly.to_expr(). Exponent slot i of the result refers to vars[i], and the polynomial remembers vars, so it can be fed straight to GroebnerBasis.reduce(), GroebnerBasis.contains() or GroebnerBasis.compute_raw().

Raises ValueError if expr is not polynomial in vars — a free symbol outside vars, a negative or symbolic exponent, or a transcendental call.

Example:

p = alkahest.expr_to_gbpoly(x**2 + y**2 - pool.integer(1), [x, y])
gb = alkahest.GroebnerBasis.compute_raw([p])

No docstring.

alkahest.solve_side_conditions() -> list[str]

The hypotheses the most recent solve() on this thread assumed in order to return the solutions it did — one string per condition, e.g. "a ≠ 0".

Two things land here. A leading coefficient the back-substitution divided by without being able to decide it, and — for a rational equation — a denominator that was cleared and whose vanishing at the returned root could not be decided either. solve([(Vo - Vin)/R1 + Vo*s*C], [Vo]) lists R1 ≠ 0: at R1 = 0 the equation you wrote has no value at all, so the answer is the answer for R1 ≠ 0. Declaring the symbol positive (pool.symbol("R1", "positive")) discharges that condition instead of reporting it.

solve([a*x - b], [x]) returns b/a, which is the solution for a ≠ 0: at a = 0 the equation reads -b = 0, so there is either no solution (b ≠ 0) or every x (b = 0), and neither of those is b/a. That generic-parameter reading is deliberate and useful: the tuple is checked to be a solution for generic parameter values, and the hypotheses listed here — those of the branches that were kept, factored — say where it stops being one. A value of an unknown declared Positive (etc.) whose membership depends on a parameter is listed too (a > 0).

This mirrors DerivedResult.verification["side_conditions"] and ZeilbergerCertificate.side_conditions; solve returns plain dict s, which cannot carry the attribute, so it is reported beside the result instead.

An empty list means the solver proved every coefficient it divided by to be non-zero — not that it did not look. Reset by each solve call, so read it before the next one; repeated reads of the same call agree.

One homotopy endpoint with optional Smale diagnostics.

coordinates

max_residual

smale_alpha

smale_certified

to_dict()

enclosures()

Lex-basis triangular decomposition (possibly split on factored univariates).

One triangular component of a decomposition, as returned by triangularize().

The component's polynomials are GbPoly; each carries the variable ordering triangularize() was called with, so p.to_expr() reads it back as an Expr.

n_vars

polys()

variables()

to_exprs()

No docstring.

No docstring.

primary()

associated_prime()

Radical √I of the ideal generated by polys.

If vars is omitted, free symbols across all polys are used (sorted by internal id). At least one free symbol must be present.

alkahest.diophantine(equation, vars) — integer solutions of a binary Diophantine equation.

Solves equation = 0 over the integers in exactly two unknowns. Supported patterns:

• Linear — returns a parametric family x(t), y(t) (kind == "parametric_linear").
• Sum of two squares — x² + y² = n (finite list of points).
• Pell / generalized Pell — x² − D·y² = 1 or = N.

Raises DiophantineError when the equation is not a polynomial in vars, has non-integer coefficients, is an unsupported pattern, or has no integer solution.

Example:

sol = diophantine(x + 2*y - 5, [x, y])
# sol.kind == "parametric_linear"; sol.parametric = [x(t), y(t)]

Result of diophantine: a parametric linear family, a finite list of integer points, or a Pell-type description.

kind

parameter

parametric

points

pell_d

fundamental

pell_n

pell_particular

pell_unit

alkahest.pell_classes(d, n) — every solution class of x² − d·y² = n.

diophantine reports one particular solution of a generalized Pell equation, but n can have several fundamental solutions (inequivalent under the unit group). This returns all of them.

d

positive non-square integer.

n

nonzero integer.

Returns a dict with d, n, unit (the fundamental solution (ux, uy) of x² − d·y² = 1), negative_unit (fundamental solution of x² − d·y² = −1 or None) and classes — a list of (x, y), one per class, sorted by (y, x); empty when the equation has no solution. Every integer solution is x + y√d = ±(x_c + y_c√d)·(ux + uy√d)^k for exactly one class c and one k ∈ ℤ.

Raises DiophantineError (E-DIOPH-003) for a square or non-positive d, n = 0, or an instance whose complete search bound is too large; BudgetExceededError when the active budget runs out.

Example:

pell_classes(10, 9)["classes"]   # [(3, 0), (7, 2), (13, 4)]

∃ var . body (same pool).

∀ var . body (same pool).

Decide a closed polynomial sentence over ℝ (one outer \forall/\exists; purely polynomial body in the bound symbol with integer coefficients).

Returns (truth, witness_or_none) where witness maps symbol names to rational decimal strings when an existential sentence is deduced satisfied.

Brown-style CAD projection polynomials after eliminating elim_var.

No docstring.

alkahest.prove_nonneg(expr, vars, *, constraints=(), basis_degree=None, level=2)

Prove p ≥ 0 on {x : g_i(x) ≥ 0} with a Handelman-style certificate p = Σ_α c_α·Π g_i^{α_i} (c_α ≥ 0). With no constraints this is sos_decompose().

alkahest.sos_decompose(expr, vars, *, basis_degree=None) -> PositivityCertificate

Exact rational sum-of-squares decomposition p = Σ_j σ_j·q_j².

Raises SosError rather than guessing: E-SOS-003 with a witness point when p is negative somewhere, E-SOS-002 when no certificate of the searched shape exists at this basis degree (which is not a proof that none exists — see the docs).

An exact, re-checkable proof that a polynomial is non-negative.

Returned by sos_decompose() and prove_nonneg(). The identity it carries has already been re-expanded and compared against the target exactly; verify() re-runs that check on demand, so a downstream consumer never has to trust the search that produced it.

kind

degree

num_squares

identity

claim

log

expression

multipliers

verify()

to_lean()