fredrikj.net / blog /
Things I would like to see in a computer algebra system
April 20, 2022
If I were to design a computer algebra system (CAS) from scratch today, I would try to achieve the following goals. Several ideas have already been implemented in one system or another, but certainly not all at once. Some ideas could be adopted by existing CASes; others would require fundamentally different designs or new programming languages. Some ideas are possibly unrealistic or naive.
Contents
Free and open source
This immediately rules out using or building on certain non-free systems. No matter how many other boxes they check.
Inert expressions
I would like a CAS to distinguish clearly between syntactical symbolic expressions and values. Most CASes perform some kind of automatic rewriting or "canonicalization": for example, if I input the expression $2(b-a)$, the CAS might return $-2a + 2b$. This is fine if I want an expanded normal form of a value in $\mathbb{Q}[a,b]$, but sometimes I really want to preserve the exact expression or work with some canonical form other than the one the CAS author thought of. Working around builtin canonicalization is a nuisance! A far better solution is to have syntactical symbolic expressions that are completely inert by default and require explicit commands for rewriting them. (This is what I did for the Fungrim backend.) This is simpler to implement too, so it's a double win.
Here, I really consider symbolic expressions separate from data structures for explicit computation in particular algebraic domains like $\mathbb{Q}[a,b]$, where canonicalization indeed is crucial for performance.
Based on math (even in analysis!)
I would like to have a CAS that is based on math. Real numbers should actually be (isomorphic to) Cauchy sequences over $\mathbb{Q}$, rings should actually satisfy the ring axioms, sets should actually be able to hold infinitely many elements (but not themselves), holomorphic functions should not occasionally have branch cuts on their domain, "solving an equation" should give me exactly the set of solutions (neither more nor less) over some well-defined domain, $x = 0$ should be true if and only if $x$ is zero...
Pretty much the only systems that more or less get math right are the mathematical libraries for proof assistants (Lean, Coq, etc.). The mathematical structures and interfaces in Lean's MathLib are actually pretty close to what I'd like to see in a CAS.
The algebraically oriented CASes like Sage, Magma and GAP are pretty good as long as you stick to pure algebra involving finite objects. They generally start to cheat when it comes to mathematical analysis. A CAS worth its name ought to be able to represent structures like transcendental number fields, rings of holomorphic functions, differential fields, transseries, and so on. Unfortunately, to the extent current-generation CASes do mathematical analysis at all, it is either numerically or using ad-hoc "symbolics" without solid algebraic foundations. (More on this in a future blog post, if I have time.)
Based on typed math
The programming language of a CAS should have a type system suitable for doing mathematics. This presumably means having fully dependent types in some form. You can implement mathematical structures on top of a poor type system (Sage does this, for example), but proper language support would be infinitely better. I have no idea what the ideal CAS type system might look like, and I believe this is a hard research problem; the following are just some scattered thoughts on the matter.
I will note that "type system" does not necessarily mean "static type system", nor "decidable", although statically decidable properties are useful (especially for reusable library code). A core ability of a CAS is to handle untyped expressions and map them to different algebraic structures, often requiring type information to be computed at runtime, and sometimes depending on properties that are not guaranteed to be computable. (Make everything decidable and you get a formal theorem prover, but a practical CAS also requires the ability to work with unproven math!)
There are tradeoffs involved in how much information you want to encode in types. On one extreme, you have a system where the type signature of a function completely specifies the function. On the other extreme (most mainstream programming languages), the type system does not even guarantee that functions are total functions. I do think that the type system in a CAS should be expressive enough to encode arbitrary mathematical relations as invariants/preconditions/postconditions (perhaps as "type restrictions"), even if you don't end up using this all the time in practice.
When it comes to actually modelling mathematical structures with types, it is important to understand that mathematics requires the ability to consider the same abstract value as an element of different algebraic structures. The number 6 can be an element of the ring of integers $\mathbb{Z}$, but it can also be an element of the semiring $\mathbb{N}$, or the set $\{n \in \mathbb{Z} : n \text{ composite}\}$, the rational field $\mathbb{Q}$, the complex field $\mathbb{C}$, the polynomial ring $\mathbb{Z}[[x]]$, and so on. You might want to distinguish between "$\mathbb{Z}$" as a type representing the integers, as a set containing all the values of this type, or as a subset or subtype of something else (like the set of rational numbers or the type of rational numbers).
The most useful type or domain in a given situation typically depends on what kind of closure properties you need rather than on the values you start with. This means, e.g. rings tend to be more useful than semirings. For example, 6 has an additive inverse in $\mathbb{Z}$ but not in $\mathbb{N}$, so if you ever want to do a subtraction, chances are that you want to represent 6 as an element of $\mathbb{Z}$ rather than as an element of $\mathbb{N}$ (this is the reason why int often makes more sense than unsigned int in C, even when you have an unsigned value). In this situation, I'd still like to be able to attach restrictions like "$n \ge 1$" in some form.
Classically "analysis-oriented" CASes like Mathematica just implicitly work over a single mathematical universe containing the field of the complex numbers and also different objects like formal indeterminates and infinities. This maximizes closure properties (every function can be applied to everything!) and avoids a lot of hairy type conversion issues, but gives very weak invariants (you cannot assume properties like associativity... though the system often does so anyway, as a practical necessity, resulting in bugs).
Modern "algebra-oriented" CASes like Sage give elements "parents" corresponding explicit algebraic structures. This has pros and cons, but on balance I think it is the better approach. The main danger is that the user gets bogged down in types, but this is perhaps a language/interface issue that can be solved. You certainly need a robust system for automatic coercions; Sage has the right idea here with its coercion system based on category theory, even if the implementation isn't perfect.
I will end this section with a more concrete example. Consider a function that divides two integers: a / b -> c. What is its type signature? Here are some possibilities:
- (Integer, Integer) -> Integer. Inexact division floors the result, and division by zero returns a dummy value like 0. (This is the preferred behavior in proof assistants.)
- (Integer, Integer) -> Union(Integer, DivisionByZeroError). As above, but division by zero signals an error.
- (Integer, Integer) -> ProjectiveExtendedInteger or (Integer, Integer) -> Union(Integer, Infinity). As above, but division by zero returns an infinity in a structure corresponding to $\mathbb{Z} \cup \{ \infty \}$. (What about $0 / 0$ though?)
- (Integer, Integer) -> Union(Integer, InexactError, DivisionByZeroError). Division by zero as well as inexact division signal an error.
- (Integer, Integer) -> Union(Integer, Rational, DivisionByZeroError). Returns an Integer when the division is exact, Rational when it is a fraction, and signals an error on division by zero.
- (Integer, Integer) -> Union(Rational, DivisionByZeroError). Always returns a fraction except when dividing by zero.
- (Integer, Integer) -> Rational. Always returns a fraction, with a dummy value like 0 on division by zero.
- (Integer, Integer) -> ProjectiveExtendedRational. Returns a fraction, or an infinity on division by zero.
- (Integer, NonZeroInteger) -> Integer. Floor division; division by zero excluded by type-level contract.
- (Integer, NonZeroInteger) -> Rational. Division producing fraction; division by zero excluded by type-level contract.
- (Integer x, NonZeroInteger y [y divides x]) -> Integer. Inexact division also excluded by type-level contract.
- (Integer, Integer) -> Union(Integer, DivisionByZeroError). Division presumed to be exact; inexact division returns a garbage value. (This is GMP's divexact.)
- (Integer, Integer) -> IEEEFloat64. This one makes no sense whatsoever but it is the way some mainstream programming languages define division.
- (Integer, Integer) -> Union(IEEEFloat64, DivisionByZeroError). Ditto.
I can think of situations where any of these would be useful. One possible takeway is that it is not sufficient for a CAS to infer the target type of an abstract mathematical function from the types of the arguments, and certainly not from their values. A CAS should, accordingly:
- Make it easy to choose both the domain and the codomain of a mathematical function that can have multiple "concretizations".
- Base any automatic coercions on mathematically sensible structure (e.g. preserving ring homomorphisms).
- Provide robust error handling (mismatched "concretizations" of abstract operations are a leading source of bugs).
Type-integrated symbolics and enclosures
A CAS should be able to deal with objects that are known precisely ("the integer 3") as well as objects that are known partially ("some integer $n \ge 3$") or approximately ("some real number near 3.14159"). My earlier post computing with metavalues discusses this among other issues.
There are several poor ways to deal with this. The Mathematica approach is to make everything a symbolic expression, where symbolic expressions can return unevaluated or partially unevaluated. For example, if you attempt to evaluate $x = 1 + 1$ where nothing is known about $x$, you simply get the expression $x = 2$ back. This is an elegant solution, but you don't get the benefits of a proper type system as discussed above. Another issue is that symbolic expressions easily blow up in size when fed through algorithms.
The Sage approach is to have a nice type system designed for precise objects but then throw inexact objects into "inexact rings" (which don't really behave like rings) and all symbolics into a Mathematica-like "symbolic ring" (which is certainly not a mathematical ring). (You then add a buggy "assumptions system" as a band-aid to make the symbolic ring vaguely useful.)
I believe a better approach is to have a type system with a clean separation between "mathematical types" and "implementation types/data" with uniform support for alternative representations (canonical, exact, symbolic, partial, inexact) on the implementation level. For example:
- An Integer can hold a binary blob representing the integer 3 in canonical form, but it can alternatively hold a non-canonical, symbolic or partial description of an integer like "$10^{10^{1000}}$" (unexpanded!), "$n \ge 3$" or "$n \equiv 15 \bmod 60$" or "$n$ is-prime", or "$n$ not equal to this other symbolic integer $m$". Or some combination of such attributes.
- A Matrix instance can hold a canonical array of coefficients, but it might also just hold information about the number of rows and columns (and these numbers may themselves be symbolic integers), the values on the main diagonal, the symbolic description $A = B C$ where $B$ and $C$ are other (potentially symbolic) matrices, etc. A "matrix over $\mathbb{Z}$" can be $[[1, 1], [0, 1]]$, but it can also be $[[1, n], [0, 1]]$ where $n$ is a symbolic integer, an unknown $2 \times 2$ integer matrix known to have determinant 1, or a completely unspecified integer matrix.
- An "element of a a ring" can be an element of the specific ring $\mathbb{Z}$, but it can also be an element of a symbolic ring $R$. This symbolic ring $R$ itself can have additional information associated with it, e.g. the property "commutative", or "characteristic $n$" (where $n$ can be symbolic). I should be able to represent the specific element $1 \in R$ of a symbolic ring $R$.
Enclosures in the sense of interval arithmetic ($x \in [3.14, 3.15]$) also fit into this framework. A key point here is that a "symbolic integer", "symbolic function" "symbolic matrix", "symbolic element of a ring" or "real number represented inexactly by an interval" should not be a second-class citizen; I should be able to create a symbolic instance or value enclosure of any mathematical type, and these instances should work with the usual mathematical type system (unlike Sage's symbolics). I should also be able (indeed, encouraged) to write code that is agnostic about the underlying representation or level of specification of mathematical objects, at least until I specifically need access to the internals.
This design does lead to two difficulties, though I think they can be solved. The first difficulty is that even basic programming elements (logic, collections and control flow) need to fit into the same framework:
- Predicates must be able to return booleans that can be unknown or symbolic.
- A loop or recursion with respect to a symbolic integer $n$ must be evaluated abstractly or symbolically.
- A sorted list or associative array needs to deal correctly with elements whose representations are not comparable or hashable.
The second issue is the potential implementation complexity and unpredictable performance characteristics. It would be difficult to hardcode data structures and algorithms for all combinations of internal representations. You'd rather want to define a mathematical type in declarative form as a set of possible properties (e.g. for a Matrix: domain, shape, actual elements, elements on the diagonal, determinant, trace, eigenvalues, has-$x$-as-eigenvalue, is-invertible, is-symmetric, is-real, ...) with a schema for types and invariants (e.g. shape is a pair of natural numbers; is-symmetric implies is-square, ...), have a compiler that automatically generates compact binary data representations, and then provide user-configurable rules for evaluation (e.g. "reduce all booleans to True/False/Unknown", or "only track the shapes of matrices" or "enclose real numbers within 53 bits of precision").
Some formal theorem provers have the right ideas here already; it is mostly a question of making things efficient and versatile for calculation.
Good for inequalities
This is a rather specific feature request. Any CAS will let me compute some terms of Stirling's series for $\Gamma(z)$. How come I can't get a truncation of Stirling's series with an explicit error term? Analysis-oriented CASes are generally good at manipulating equalities and limits, but strangely poor at manipulating inequalities. This is bizarre since inequalities are at least as important as equalities in analysis!
Large but lean
The most advanced CASes are multi-gigabyte behemoths, and even more specialized systems start to weigh a bit heavy when you count all dependencies. They take time to download, they are complex to compile, and they are difficult to develop. Even if we get around the build/installation hurdles (let's just put everything into clouds! or containers! or containers in the clouds!), we have startup time and interactive latency to deal with. There are times where I want to run a multi-hour job in a CAS, but most of the time I just want to do a very quick calculation. I will not launch sage in the terminal (≈ 3 seconds startup) if I can get away with python or gp (≈ 0.01 seconds startup). I will certainly not reach for a cloud notebook (perhaps more than a minute to navigate a website, log in, wait for some JS framework to load, and spin up a remote VM) if I have a better option. Add to this the latency for loading any extra libraries/modules and compiling code at runtime (Julia currently leaves a lot to be desired here).
Now, size isn't something you can escape from, at least entirely. A useful CAS needs to know a great deal of mathematics, how to translate that mathematics to machine instructions, and how to interact with the user and other software; this knowledge will not fit on a floppy disk. However, I think the user-unfriendliness of large CASes is mostly of question of design. Consider Pari/GP for a counterpoint: it is self-contained (except for an optional dependency on GMP), the entire binary is something like 10 MB, starts instantly (it even loads pretty quickly in a web browser; try it here), and it does a ton of advanced number theory very efficiently.
Pie-in-the-sky: I would like to see a CAS deliberately designed around having a small language/interface kernel together with a large database of mathematics (types, functions, constants, theorems, properties, algorithms), with automatic, efficient lazy loading.
I say "database" rather than "library" because I think the usual library/module way to organize software leaves something to be desired. Computational mathematics is inherently highly integrated; a numerical calculation might invoke some number-theoretic subroutine that transforms the problem to a calculation over finite fields aided by a lookup in a database of Conway polynomials. Even if a computer algebra system is neatly organized into libraries and modules, a single calculation may end up using a sparse subset of functionality of several components. This shouldn't require downloading, possibly compiling and then loading into memory all the contents of several huge libraries. (And as a user, I should certainly not have to do any of these things manually!)
I propose that the mathematical library of a CAS should be some kind of semantic database with code-as-data, fine-grained automatic dependency tracking, efficient loading/installation/caching, and tiered interpretation/JIT compilation for dynamic code. With this design, you could have a kernel that fits in a small WASM or native executable, making it possible to launch a local instance in a terminal or a web browser where you can do factor(12345) or plot(sin(x)) instantly while also having convenient access to efficient implementations of modular forms on quaternion algebras (or whatever) that will be loaded quickly on demand.
I could see computer algebra moving server-side entirely, but it will be a tall order to to solve the usability and performance problems. Even the free-as-in-beer web services offered by the multimillion dollar company Wolfram (Alpha and Cloud) have notable issues. Surprisingly, SymPy Live and Gamma are much more usable; that is, at least until you try to calculate something nontrivial and the backend chokes.
Integer performance
The core datatype in numerical computing is the machine-precision floating-point number. High-performance computing is virtually synonymous with turning tasks into cache-friendly and vectorized operations on homogeneous arrays of floats. Hardware, compilers, programming languages and numerical libraries are strongly designed around this idea.
The core datatype in computer algebra is the arbitrary-size integer. Hardware, compilers, programming languages are not at all designed for performant computation with this datatype. Indeed, applications that support arbitrary-size integers typically just wrap GMP's mpz_t (sometimes behind extra layers of indirection) and implement mathematical operations naively on top. This results in massive amounts of overhead for small integers and often in algorithms with poor asymptotic complexity, giving users the misleading impression that exact integer arithmetic is inherently inefficient (often 10 to 1000 times slower than necessary), which in turn leads to a vicious cycle of neglect.
Performance-oriented computer algebra libraries like FLINT take pains packing integers and arrays of integers into arrays of bits or words, using modular arithmetic where appropriate to avoid intermediate coefficient swell and facilitate vectorization. This requires a lot of manual coding: you have one standard representation for the coefficients of polynomials, several other internal representations for the coefficients when multiplying or computing GCDs, others for exponents of multivariate polynomials... with lots of fiddling for conversions, resizing, overflow detection, modular arithmetic, precomputed inverses, delayed reduction...
I imagine that a CAS created by an advanced extraterrestrial civilization would be written in a programming language with a native notion of arbitrary-size integer with the ability to generate context-optimized low-level packed representations, bit-level magic and parallelization/vectorization from high-level descriptions of data structures, morphisms and algorithms. Barring alien technology, a good start would be some kind of compiler/language infrastructure where doing all these things manually is easier than what it is right now!
Good math display
Superficial but important: a CAS should be able to produce publication-quality formulas, graphics and tables without special effort by the user. In practice, the results are often less than impressive. These should not even be particularly hard problems to solve nowadays; they are just not prioritized, or don't get worked on by persons with a sufficiently developed eye for aesthetics or detail (sorry).
Text-friendly (and human-friendly)
Tools that work well with plain human-readable ASCII input and output are wonderful. I would like every object in a CAS to print to a plain text expression in the syntax of the language which evaluates back to precisely the same object, independent of context. Serialization will, as a consequence, be trivial. Other forms (pretty-printed, condensed, binary) can by all means be supported as well, but a verbose and faithful text representation should be trivial to access (if not the default form of output). This also calls for having a simple, regular syntax so that it is easy to extract and compose subexpressions, and so that the syntax is comprehensible for new users.
(I'm on the fence about Unicode math. It seems to occupy an uncanny valley between easily-typed ASCII math and proper typesetting.)
Concerns about readability and serializability may also influence semantic design decisions. It seems desirable to have strongly functional semantics that avoid implicit context and state. I also think a CAS should have good support for decimal floating-point numbers, without forced conversions to binary. (Lossy binary-decimal conversions are a frequent interface pain point with Arb.)
Well-named
Since mathematics is highly integrated, a CAS should provide a unified namespace with a consistent and convenient naming convention. You can by all means have sub-namespaces for specialized functionality, and the library implementation can be hierarchical, but for the love of things don't make users guess which modules they need to import. Granted, naming things is the hardest problem in computer science, so this is perhaps wishful thinking!
Actually, I lied: the hardest problem in computer science is backward compatibility. A simple, feature-frozen kernel with a versioned standard library could be a way forward (but easier said than done).
fredrikj.net | Blog index | RSS feed | Follow me on Mastodon | Become a sponsor