Maple and Python Comparison¶
Doron Zeilberger’s original Pisot.txt was implemented in the computer algebra system Maple. Maple has a very robust symbolic computation system, but its programming language is clunky. I have chosen to implement Pisot.txt in the programming language Python to enhance its readability and show that Maple (and its ilk) are not the only languages that can perform symbolic computation.
Python was designed to be a readable and general purpose language. Applications specific to mathematicians are not very general purpose, so symbolic computation is not a native feature of Python. Fortunately the library SymPy implements symbolic computation in Python and along with it most features found in Maple. By performing all of the “behind the scenes” calculation with SymPy, we obtain both the readability of Python and the symbolic power of Maple.
As a short example of the differences between Maple and Python, we will compare
two snippets from Zeilberger’s Maple procedure Pis
and the corresponding
Python version. The procedure Pis
computes the absolute value of the second
largest root of the characteristic equation for a C-finite recurrence. As a
part of this procedure, Zeilberger discards 1
if it is a root. If lu
contains the roots, then the following Maple snippet accomplishes this task:
if member(1,lu) then
lu:=convert({op(lu)} minus {1},list):
if lu=[] then
RETURN(FAIL):
fi:
fi:
Perhaps this is an obvious design pattern to an experienced Maple programmer, but to my eyes, there is a cognitive barrier between what we want to do (remove the root 1) and what is being done (convert a list to a set, perform set subtraction, and then converting the result back to a list).
If roots
is a dictionary of roots, then in Python the same task can be
accomplished as:
if 1 in roots:
del roots[1]
if not roots:
return None
The Python version communicates the intent of the code more clearly and without as many technical details.
Detailed Comparison¶
As an example of the differences between the two, let us compare the full
implementations of Zeilberger’s Maple procedure Pis
. His procedure is as
follows:
#Pis(C): Inputs a C-finite sequence and outputs the absolute value of the second-largest root
#It is a Pisot number if it is less than 1.
#Fis([[1,1],[1,1]]);
#Pis([[10,219,4796,105030],[22,-3,18,-11]]);
Pis:=proc(C) local x,lu,i,aluf,mu:
lu:=[solve(x^nops(C[2])-add(C[2][i]*x^(nops(C[2])-i),i=1..nops(C[2])))]:
if nops(lu)<>nops(C[1]) then
RETURN(FAIL):
fi:
if member(1,lu) then
lu:=convert({op(lu)} minus {1},list):
if lu=[] then
RETURN(FAIL):
fi:
fi:
aluf:=1:
for i from 2 to nops(lu) do
if abs(evalf(lu[i]))>abs(evalf(lu[aluf])) then
aluf:=i:
fi:
od:
mu:=evalf([op(1..aluf-1,lu),op(aluf+1..nops(lu),lu)]):
max(seq(abs(mu[i]),i=1..nops(mu))):
end:
Now, taking our implementation of cfinite.CFinite
for granted, the
Python
implementation is as follows:
def pisot_root(c_seq):
"""
Compute the absolute value of the second-largest root of the characteristic
equation for the C-finite sequence.
:c_seq: :class:`.CFinite` instance.
:returns: Floating point evaluation of the absolute value of the root, or
None.
"""
roots = c_seq.characteristic_roots()
n_roots = len(roots.keys())
if n_roots != c_seq.degree:
return None
if 1 in roots:
del roots[1]
if not roots:
return None
root_norms = [abs(root) for root in roots.keys()]
root_norms = [sympy.re(sympy.N(norm)) for norm in root_norms]
max_index = root_norms.index(max(root_norms))
del root_norms[max_index]
return max(root_norms)
The procedures are about the same length, accounting for blank lines and comments.
The first feature is the documentation. Both versions document their inputs and
outputs, but the Python version is written in a standard way that allows for
automatic documentation generation. In fact, the documentation at
pisot.pisot_root()
is automatically generated, as is every other piece of
documentation on this site.
The next part computes the roots of the characteristic polynomial. Maple:
lu:=[solve(x^nops(C[2])-add(C[2][i]*x^(nops(C[2])-i),i=1..nops(C[2])))]:
Python (taking advantage of cfinite.CFinite.characteristic_roots()
):
roots = c_seq.characteristic_roots()
Using Python’s classes, it is clear that the characteristic roots are a property of the C-finite sequence, and that these roots are what we are computing.
Next, we look to see if there are any repeated roots. This is true if the number of distinct roots is less than the degree of the sequence. Maple:
if nops(lu)<>nops(C[1]) then
RETURN(FAIL):
fi:
Python:
n_roots = len(roots.keys())
if n_roots != c_seq.degree:
return None
The Maple version counts the number of coefficients in the C-finite sequence
and calls this the degree. The structure of the Maple C-finite sequences is
[[coeffs], [initials]]
, so it counts the number of elements in the first
element of a nested list. The Python version, relying on the class
cfinite.CFinite
, simply asks for the sequence’s degree. Again,
Python’s classes allow us to embed information into an object, rather than
relying on its actual representation.