How Python could improve to better support something like Sympy?

[Henri Tuhola]

Hello,

I would directly want to reach people who have been developing sympy, so I posted on your mailing list.

Development of sympy on top of python has not been frictionless, and I think you know that. Unfortunately only few people know what's even wrong there.

I am working on a new programming language and found out that the choice of how to allow extension of arithmetic is a really tough problem to figure out. Simply implementing arbitrary operator overloading and calling it a day doesn't seem to be so good choice.

I wonder that if you had a choice of redesigning arithmetic in Python to better support sympy, what would you attempt to solve? If you know what you would do, I would gladly read that.

Thank you ahead of the time.

[Richard Fateman]

It might be useful to list the points of friction that you are aware of in the design of python with

respect to use in sympy,

to see if people have work-arounds that you might have missed, or to encourage

people to point out solutions that exist in alternative languages that already exist.

As long as you are considering other programming languages, maybe you should look at  Lisp,

which is used in a number of well-known symbolic math systems, and is in some respects

python-friendly.

RJF

[Aaron Meurer]

You might look at multiple dispatch as an improved method of operator
overloading, such as what is implemented in Julia. There are also
potentially more advanced things which can be useful, such as pattern
matching (Mathematica, Haskell).

You have to consider the tradeoffs, however, with expresibility. The
most general possible system is a full macro system, which lets you
effectively define whatever syntax you want. But the downside to that
is that you no longer have a fully consistent language. Someone
reading the code for a library must first learn the syntax. To
contrast something like Python, which does not have a macro system and
has (relatively) limited operator overloading, someone who already
knows the language from other uses can read the SymPy code and have a
good idea of what things mean (it also helps that SymPy and Python use
a fairly strict adherence to standard mathematical notation).

To give an example of this, we sometimes wish that Python allowed
changing the operator precedence. For instance, we can't use ^ for
exponentiation, because even though Python allows overriding that
operator, it is fixed at a precedence lower than +, so a + b ^ 2
without parentheses will give (a + b)^2 (^ is XOR). The downside is
that new users to Python must learn to use ** instead of ^ for
exponents. The upside is that anyone already familiar with Python can
see any SymPy expression and know how it will be interpreted, because
operator precedence is uniform across all Python code, even if it
overloads the operators to mean different things.

So I think there's a sweet spot. We have run into limitations of
Python's operator overloading mechanism. I personally think that
multiple dispatch and algebraic pattern matching are both interesting,
but having not used a language that implements either extensively, I
can't say for sure if either are that sweet spot.

Aaron Meurer

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
> To post to this group, send email to sy...@googlegroups.com.
> Visit this group at https://groups.google.com/group/sympy.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/2ab5c94f-97f2-46d9-aa92-39839a68566a%40googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

[Henri Tuhola]

The idea about changing the precedence and meaning of ^ from

XOR to exponentiation is a bit provoking. I've been

programming so long that I often just forget that ^

means exponentiation for most people.

I did a quick search and it reveals that the current Python

symbol for exponentiation (**) comes from Fortran. It

used that symbol because the ^ did not exist yet. The ^ as

bitwise XOR was introduced by C language. Even the search

results pointed out that this is confusing for newcomers.

If I did what you propose, everyone who has been programming

C would do the mix up from time to time. But these guys are

actually more capable of dealing with the problem, compared

to people who might enter programming through symbolic

algebra. Replacing the XOR ^ with POW ^ could make a lot of

sense even if I did not really consider it at first.

Another, slightly wilder idea that I've had is that what if

you actually made ordinary numbers written by users bignums,

and numbers such as 1.2 would be fractionals by default.

This could slow down some programs considerably but the way

I see it the dynamic language is usually an excellent

starting point and when you go further from that it might be

better to instead have tools to optimize the dynamically

typed programs into something else.

I did follow Richard's advice and trying to list and/or

identify the points of friction.

To help this off I made myself a model of what a typical

operator overloading looks like in Python at this moment:

https://gist.github.com/cheery/228b6651fb6a460b91f26195fe58e397

To resolve a+b, Python first calls a.__add__(b)  for the

left-side object. If it returns NotImplemented, then it

tries b.__radd__(a). If it still returns NotImplemented then

it fails.

Operator overloading in Python looks better than I remember.

It's probably because many guides and memos fail to state

the important part of how Python specifically resolves the

call and they do not even include the 'isinstance' -part to

identify the other side on the expression.

In my language I have implemented the '+' as a multimethod

and added another multimethod to coerce. I were happy to

that for a while but then I realised that it may make some

things worse because the multimethod resolution ends up

being more complex rather than simpler, and not that many

problems appear to come from the dispatch mechanism after

all.

Another thing that seemed important but is maybe not that

crucial flaw is that Python's approach ends up dominating

with the left side of __add__. If the operation is found in

the __add__ then it won't be searched from the __radd__ and

this can cause conflicts between different libraries that

extend the arithmetic.

In my opinion the real problem is, and will likely always

be in getting multiple libraries interoperate. I've been

looking that statically typed languages seem to handle this

better than dynamically typed ones but in the end it turns

out to be the same ways by which you successfully cope with

the challenges.

If you have M different kinds of things, then in worst case

you need M*M different implementations for arithmetic. In an

one project it's likely that you cope well, but in presence

of multiple systems it becomes harder.

I think these ideas boil down to this: When you create new

behavior for arithmetic you're creating a new numerical

system that extends from the base types you have. This new

"number system" describes how the new values behave with the

existing values.

Languages like Haskell seem acknowledge that people create

new numerical systems when they overload operators. The

numbers are implicitly wrapped with (fromInteger N), and

you are supposed to define those conversion functions from

numeric values when you define the + operation. Different

numerical systems defined this way do not interact with each

other in Haskell.

I've been studying subtyping in order to optimize

dynamically typed programs directly from their source code

and get them to translate and run in C or GLSL -like

environments alongside the dynamic portions of the programs.

I think that this interacts with the way how you overload

the arithmetic.

The another, kind of a potential starting point I've been

thinking about is related to how numerical libraries deal

efficiently with computing.

The kind of an obvious thing that a programming language

implementation could provide would be tight memory layouts,

or even in-memory relational tables where you can fill up

the numerical data. The tight memory layouts would appear to

matter the most for efficient computing.

Numerical data itself seems to always be bundles of numbers.

You either got quaternions or matrix data with single or

double precision floating point.

I've been thinking that actually when you have quaternions

or matrices, it looks like neither one of those actually

need to take precedence and the number bundle end up staying

as an array. Perhaps those bundles could be treated as sums

of terms*constants, where the terms -part is a parametric

type? 

Lisp is often mentioned being a family of various languages

such as CLisp, scheme, racket. If you people provided specific

examples it would help but I will try to find out those

myself as well.

[Aaron Meurer]

On Mon, Jan 15, 2018 at 3:38 PM, Henri Tuhola <henri....@gmail.com> wrote:
> The idea about changing the precedence and meaning of ^ from
> XOR to exponentiation is a bit provoking. I've been
> programming so long that I often just forget that ^
> means exponentiation for most people.
>
> I did a quick search and it reveals that the current Python
> symbol for exponentiation (**) comes from Fortran. It
> used that symbol because the ^ did not exist yet. The ^ as
> bitwise XOR was introduced by C language. Even the search
> results pointed out that this is confusing for newcomers.
>
> If I did what you propose, everyone who has been programming
> C would do the mix up from time to time. But these guys are
> actually more capable of dealing with the problem, compared
> to people who might enter programming through symbolic
> algebra. Replacing the XOR ^ with POW ^ could make a lot of
> sense even if I did not really consider it at first.
>
>
> Another, slightly wilder idea that I've had is that what if
> you actually made ordinary numbers written by users bignums,
> and numbers such as 1.2 would be fractionals by default.
> This could slow down some programs considerably but the way
> I see it the dynamic language is usually an excellent
> starting point and when you go further from that it might be
> better to instead have tools to optimize the dynamically
> typed programs into something else.

Making floats rational numbers has a ton of pitfalls. The reason
floating point numbers work the way they do is that floats always
"truncate" the result of any calculation to fit in the floating point
format (8 bytes for a double precision float). With rationals, the
numerator and denominator invariable end up increasing, sometimes
massively. For instance, if you add n rational numbers together, the
denominator can be at least as big as the product of the denominators
of the original numbers. Guido has also talked about this here
https://python-history.blogspot.com/search?q=rational

Floats also have an advantage of being very fast, since they are
supported directly in the hardware.

One thing that would be interesting, though, would be the ability to
recover the original floating point number as a user typed it. In
Python there's no way to distinguish 0.1 and 0.10000000000000001, for
instance. In SymPy, if you want to create a Float or Rational from a
float, the best practice is to input it as a string, so that you don't
lose any precision from Python's interpretation of the float literal.

In fact, as an aside, Python's eager evaluation is both a boon and a
bane for SymPy. On the one hand, it makes the language very easy to
reason about. On the other, we run into issues like the gotcha that
expressions like x**(1/2) create x**0.5, when the version with the
rational number x**Rational(1, 2) is almost always much preferred.
There's no way around this at the language level. We just have to tell
users to be careful when doing integer/integer (things are even worse
in Python 2 where 1/2 produces 0, but fortunately, this is fixed in
Python 3 and more and more people are using it).

>
>
> I did follow Richard's advice and trying to list and/or
> identify the points of friction.
>
> To help this off I made myself a model of what a typical
> operator overloading looks like in Python at this moment:
> https://gist.github.com/cheery/228b6651fb6a460b91f26195fe58e397
>
> To resolve a+b, Python first calls a.__add__(b) for the
> left-side object. If it returns NotImplemented, then it
> tries b.__radd__(a). If it still returns NotImplemented then
> it fails.

That's the basic idea, though it does get more complicated than that
when you consider multiple inheritance.

>
> Operator overloading in Python looks better than I remember.
> It's probably because many guides and memos fail to state
> the important part of how Python specifically resolves the
> call and they do not even include the 'isinstance' -part to
> identify the other side on the expression.
>
> In my language I have implemented the '+' as a multimethod
> and added another multimethod to coerce. I were happy to
> that for a while but then I realised that it may make some
> things worse because the multimethod resolution ends up
> being more complex rather than simpler, and not that many
> problems appear to come from the dispatch mechanism after
> all.
>
> Another thing that seemed important but is maybe not that
> crucial flaw is that Python's approach ends up dominating
> with the left side of __add__. If the operation is found in
> the __add__ then it won't be searched from the __radd__ and
> this can cause conflicts between different libraries that
> extend the arithmetic.

That is a problem, especially for things like SymPy that want to
gobble things up into their own objects (like sympy.Add). In SymPy we
haven't yet properly figured out how to allow objects to define their
behavior in an Add object (things are complicated by the fact that Add
is an n-ary function).

Aaron Meurer

> https://groups.google.com/d/msgid/sympy/05adcf16-ef85-4101-87da-aa421d9dd4f5%40googlegroups.com.

[Henri Tuhola]

The nice thing about real numbers seem to also be the nasty thing about them. I suppose there should be a way to select what integers and floating point values mean in a module.

I implemented the representational floats, they are float subtypes that retain their string representation. They are produced by parse_float, so anything that takes string floating point items as an input and does not change them will return these. You are able to recognize them with isinstance(x, float_repr). The input string can be retrieved with x.to_string().

The idea of annotating input into the results is nice enough that I'll also consider it with string values.

It didn't require as many changes I thought. The changeset introducing these objects: https://github.com/cheery/lever/commit/c7af1fce2bcf0035fc6752eae7dd126742455420

I also removed the ^ from the syntax and put the users do xor(a, b) when they want bitwise xor. I did it this way because I've written quite lot of code with my language despite that it's not finished and some of that code is using xor. The intent is that I will recognize the old code that uses it and later add the symbol back in as an exponentiation. This ought alleviate the need to change the precedence rules.

I feel I'm still not done here, but I've been wondering about these things for months so there's no need to rush.

Henri Tuhola

[Aaron Meurer]

I personally believe that Python gets more right than wrong, so if you
want more inspiration I would look at the things Python has done. For
SymPy, it's nice that Python allows full mathematical syntax in more
or less the exact way you'd expect. Little touches like complex number
literals (more useful for numeric libraries than SymPy since they are
floating point only), uniary +, and x**-y precedence are missing in
many other languages. For numerics, you should look at the history of
NumPy and what has allowed it to exist, especially the buffer
protocol.

Aaron Meurer

> https://groups.google.com/d/msgid/sympy/91325214-00da-4bae-9008-47d09539a2fa%40googlegroups.com.

[Henri Tuhola]

The things that Python gets wrong started getting on my way and I certainly wouldn't like to get things wrong that Python got right, but I guess I can't avoid it more than getting things wrong that Python got wrong.

I did look into NumPy a bit, but like Richard Fateman did before you, you're giving me a really big target to look into.

At first I thought that NumPy being solely an optimization wouldn't be very interesting for this discussion. There's some point in pointing to NumPy though.

NumPy appears to be a "parametric array", exactly like I mentioned earlier on. I slightly disregarded it because it's not parametric in every way, doesn't handle quaternions and I haven't seen it used on computer graphics, and the big thing about it is that it represents large arrays of values and the fact that it's parametric has less importance.

My project sits on top of PyPy and they've probably done some NumPy stuff that I could pick up as a starting point. I think I'll look there.

Henri Tuhola

[Richard Fateman]

There's actually a long history of people designing languages for computer algebra,

going back to at least 1963. One burning question has been whether the implementation

language and the user-visible language should be the same.  Or mostly the same.

So Maple is mostly implemented in Maple, plus some kernel in a C-like language.

Maxima is mostly implemented in Lisp, but some is in Maxima's top-level language

which can be compiled to Lisp/assembler.  Mathematica is a mix of some

Objective-C (or something similar) and the wolfram language.

Axiom is (I think) mostly in its own language.

Reduce has Lisp plus an infix-version of something that

is pretty much Lisp, I think.

Various other front-ends exists (menu/graph/etc) as alternatives to

languages.

  I think the emphasis placed on syntax in most peoples' minds

concerning programming languages is a barrier to understanding

what matters in terms of semantics.  One reason to learn Lisp

is that it frees you from that viewpoint, and then you can go back

to whatever other language you were using and write better

programs.  Or stick with Lisp, and be enjoy its other advantages :)

There are systems entirely written in Lisp where the user

language is also Lisp. (conveniently, these days, using

the Common Lisp Object System.  Which if you have not

studied, you might consider doing, since it is sort of the

mother ship of object systems.  )

RJF

[Henri Tuhola]

I did a study, which means that I searched google with "designing languages for computer algebra" and went along for 5 pages of results.

Designing a new language for computer algebra used to be how it was done. Axiom CAS dates to 1971 and was apparently called scratchpad back then. The name is confusingly similar to "sketchpad" which was an early CAD in 1963 when human computer interaction through GUI was novel

Maple's paper appeared to give large value for list-like records that contain expressions. Maxima resembles SymPy a lot in the sense that it's got the scheme-like evaluation and then simplification before display.

Mathematica has nice documentation that mentions all the algorithms they use, but the overall design they use appears to be even more minimal than lisp usually is. It's just hash-tagged expression trees with numbers, strings and symbols as leaves in it. The system evaluates that tree through iteration of several algorithms and it's got a large database of subroutines to support it.

I also stumbled upon Redberry and Magma and read some papers about them. The first one seems to be focused on tensor algebra, the second one is obsessed with representation of mathematical domains. Overall I think they aren't bringing that much to the table that'd be important now.

There is a common theme. Computer algebra system designers do not seem to worry much about the representations of their numbers. They settle to integers, reals and fractionals of the system they build upon and consistently use the same numerical system everywhere.

After seeing that I decided to start over and look into numerical analysis. Now there I found them caring about number representations. Looking into MATLAB really gave some much needed new ideas.

Now I think I may have the pieces of the puzzle to continue on my own. Those came up pretty quickly when I knew what to search for!

NumPy is related to MATLAB by quite many pieces. I think that the MATLAB one is a better study piece as they have allowed the language details to change in order to make it work for numerical computing. To make the language support numerical analysis it appears that the packed arrays will be the important subject. I guess that in order to learn more about this I have to play with these things a bit. I'll maybe read through few courses and play out with the ideas a bit and see how they fit my language.

Object systems are themselves also relevant and perhaps important subject of this all. I've been looking at my clone-python-object-system for several months now. I think I see points where it is outright failing to do what I need. Also the illusion it gives of a free customization is outright treacherous and that has driven me to think that the object oriented programming is better for selling books and ruining expensive enterprise projects than for anything else.

Inheritance with single parent only seems to be quite useful with some pattern matching and dispatch tools, but it looks like I am going to need several different ways to extend the language, each with their own qualities. One of them (python object style extensible records) is not special and should not have any distinctions of the other ways to extend.

I already knew about the CLOS. Knowing about it has not helped much because I used to think of it as just yet another class-object system.

Henri Tuhola

[oldk1331]

Hi Henri Tuhola:

> I am working on a new programming language and found
> out that the choice of how to allow extension of arithmetic
> is a really tough problem to figure out. Simply implementing
> arbitrary operator overloading and calling it a day doesn't
> seem to be so good choice.

I agree this is the key question: polymorphism.  Overloading
is ad-hoc polymorphism and template (C++) is parametric
polymorphism.

> Designing a new language for computer algebra used to be
> how it was done. Axiom CAS dates to 1971 and was
> apparently called scratchpad back then.

You mentioned the history of ScratchPad/AXIOM, but didn't
comment on it.  It's a very interesting system. It's the only CAS
with a staticly typed compiled language, like Julia/Rust/Haskell.
And it is also built on Lisp, but uses indentation like python.

Just like Julia's abstract type, Rust's trait, Haskell's typeclass,
(C++'s concept), AXIOM has the notation 'Category'.
That is used to solve the problem of polymorphism.
Category in AXIOM is very mathematical, like
Ring/Field/IntegralDomain.

FriCAS is an active fork of AXIOM, you can learn more at
http://fricas.sourceforge.net/
and
https://fricas.github.io/

[Henri Tuhola]

Without your mention I was going to skip AXIOM because I didn't get to grasp the ideas it presented immediately and the claims felt quite pompous.

I read about those categories you mentioned and it is still baffling by which logic that system works. Reading further on it actually seems interesting. It proposes it provides compositions of objects that restrict by their mathematical properties.

I'm going to read deeper and play out with that stuff.

Henri Tuhola