ANN: Symbolic computation with Julia

[lapeyre....@gmail.com]

Hi,

Here is some code I have been working on. It is for symbolic computation based on pattern matching.
If you find it interesting, any feedback is welcome!

https://github.com/jlapeyre/SJulia

--John

[James Crist]

One of the SymPy devs here. Without looking too into your code, I'd like to comment on the pattern-dispatch system, with regards to associative-commutative expressions. This has been on my mind a lot lately, as I've been working on a rewrite of our pattern matcher to be more efficient/versatile.

With dispatch on patterns for symbolic computation, there seem to be two schools of thought. If you read papers from the theorem proving world, most use a sufficiently-smart-matcher(tm) that can handle patterns with certain properties (associative, commutative, identity, etc...). Most of these result in a single, heavily optimized data structure that then performs all pattern dispatches. In this case, you'd load all your rules, and then run the rewrite system.

In contrast, most of the computer algebra work use a more heuristic-y approach, with preprocessing of the patterns, and formation of many smaller matching structures. These are then applied manually in human-written logic code. By this I mean that you may have a structure for trig rewriting rules, and one for logarithm rules, etc... And then in your main routine (let's say a simplification routine), you would apply these manually (match(expr, trig_rules), match(expr, log_rules), etc...). The idea behind this is that in most cases you have some knowledge about the structure of the expression, which can reduce the number of rules you'd need to check.

For larger collections of patterns, a good way is to use a discrimination net (a trie like structure), which can aid in faster matching. This paper is a good reference: Associative-Commutative Discrimination Nets. This approach is great if you have a large group of patterns.

For smaller collections of patterns, Richard Jenk's paper "A Pattern Compiler" is a decent starting point. This system is much more similar to that of Mathematica's than the one dscribed above (describes the pattern compiler for SCRATCHPAD, which later became AXIOM). Unfortunately, it doesn't deeply discuss the implementation, but the basic idea behind it is given. The gist is that a collection of patterns is preprocessed, and all possible combinations of rules (associative/commutative/identity) are expanded. A matching program is then created and evaled, creating a specialized rewrite function for that set. For small patterns/collections of patterns (which for computer algebra is the norm) this is a great system, and one that I think could be implemented in Julia.

Sorry for the long response, this is something I've just been thinking about a lot lately.

- Jim

[lapeyre....@gmail.com]

Hello Jim,
Thanks for the long response. I looked at the papers; I think they are worth studying. I found that I can progress most quickly modeling on Mathematica. The design is fairly simple and consistent. Of course the documentation is vague in some places, but so far, it is not a problem. Of course, even if you know logically what the rules should do, there is still the matter of building data structures and applying them.  I am more concerned now in getting correct logical behavior than efficient application of rules. I want to see a cohesive system. A lot of risk is reduced, because we know Mathematica works. But, it is good to develop ideas about efficient implementation at the same time, and there is no avoiding it when doing commutative matching.

I installed sympy, but have not played with it yet, but I definitely will. I've looked just a bit at the code. If you have any experience optimizing pattern matching, or other resources, I'd like to hear about it.

Another resource is Richard Fateman's  mma4max, which I see has updates since I last built it. It aims more or less to implement mma in common lisp. In fact it has a Julia-compatible license and a complete parser written in common lisp. IIRC the comments show that the pattern matcher is not perfect, but I think it is complete.

--John

[Francesco Bonazzi]

On Monday, January 26, 2015 at 6:11:33 PM UTC+1, lapeyre....@gmail.com wrote:

I found that I can progress most quickly modeling on Mathematica. The design is fairly simple and consistent. Of course the documentation is vague in some places, but so far, it is not a problem. Of course, even if you know logically what the rules should do, there is still the matter of building data structures and applying them.  I am more concerned now in getting correct logical behavior than efficient application of rules. I want to see a cohesive system. A lot of risk is reduced, because we know Mathematica works. But, it is good to develop ideas about efficient implementation at the same time, and there is no avoiding it when doing commutative matching.

Have you had a look at Mathics?
https://github.com/mathics/Mathics

It is an open source Wolfram Mathematica interpreter written in Python, and uses SymPy as a callback for math algorithms.

You may try it from your browser:
http://www.mathics.net/

I installed sympy, but have not played with it yet, but I definitely will. I've looked just a bit at the code. If you have any experience optimizing pattern matching, or other resources, I'd like to hear about it.

If you wish to try sympy out, I may suggest to run " isympy -I " from the console.

Another resource is Richard Fateman's  mma4max, which I see has updates since I last built it. It aims more or less to implement mma in common lisp. In fact it has a Julia-compatible license and a complete parser written in common lisp. IIRC the comments show that the pattern matcher is not perfect, but I think it is complete.

That's better than Mathics concerning the license, as Mathics is GPL.

Anyways, I saw you use a lot of conventions from Mathematica, e.g. x_ with a trailing underscore to define parameters and x_Integer to restrict the types. But unlike Mathematica, _ is not an operator in Julia. Maybe there could be some confusion, what about using :: instead, or something like that?

Concerning the expression tree, Mathematica has an attribute "Orderless" to determine whether a type of node is commutative (e.g. Times, Plus). Julia's operators are dispatched according to the types of their arguments, so * may be commutative for some types, while not commutative for other types (most notably, for matrices it's not). In Mathematica you have to use the dot "." to multiply matrices, which creates a different kind of node (Dot). On the other hand in SymPy multiplicative commutativity is handled as an assumption on the single arguments of the "Times" node (Mul in SymPy).

In [2]: x = symbols('x')

In [3]: y, z = symbols('y z', commutative=False)

In [4]: x*y == y*x
Out[4]: True

In [5]: y*z == z*y
Out[5]: False

Here x is multiplicatively commutative, y and z are not.

I prefer this approach to that of Mathematica, unfortunately it is harder to implement.

[Francesco Bonazzi]

On Monday, January 26, 2015 at 6:11:33 PM UTC+1, lapeyre....@gmail.com wrote:

Another resource is Richard Fateman's  mma4max, which I see has updates since I last built it. It aims more or less to implement mma in common lisp. In fact it has a Julia-compatible license and a complete parser written in common lisp. IIRC the comments show that the pattern matcher is not perfect, but I think it is complete.

Maybe the parser could be useful if rewritten in Julia:
http://www.cs.berkeley.edu/~fateman/lisp/mma4max/parser.lisp

[James Crist]

Maybe the parser could be useful if rewritten in Julia:

Side tracking a bit, but one of the reasons I use SymPy instead of something like Maxima is that it meshes seamlessly into the language infrastructure. By writing a parser and creating another system means that it can't use julia's ecosystem without modification. I'd rather see a CAS work within Julia's parser framework. But that's just me.

I am more concerned now in getting correct logical behavior than efficient application of rules.

Making it just *work* is certainly a priority, but design decisions made now can make things tricky later on. Especially with the design of your internal representation. How you represent expression trees and assumptions can seriously affect the performance of matching/traversals/replacements.

[Francesco Bonazzi]

On Monday, January 26, 2015 at 6:57:48 PM UTC+1, James Crist wrote:

Side tracking a bit, but one of the reasons I use SymPy instead of something like Maxima is that it meshes seamlessly into the language infrastructure. By writing a parser and creating another system means that it can't use julia's ecosystem without modification. I'd rather see a CAS work within Julia's parser framework. But that's just me.

A Mathematica parser could be used as an optional addition, such as what Mathics is to SymPy.

Anyways, there are various options for Julia:

• use Expr and operate on them, this is similar to how Maxima works.
  
• declare symbolic types and use operator overloading: this is the approach by SymPy and Axiom.
• create a DSL invoked by a macro, this appears to be the case for SJulia.
  

[lapeyre....@gmail.com]

Hi Francesco,

 

It is an open source Wolfram Mathematica interpreter written in Python, and uses SymPy as a callback for math algorithms.

You may try it from your browser:
http://www.mathics.net/

 I looked at it, but didn't get a chance to use it.

If you wish to try sympy out, I may suggest to run " isympy -I " from the console.

    That's it! thanks. I didn't know where to start.

 
All of these approaches are interesting. I am most interested now in symbolic math based on rewriting expressions, rather than data types within a language representing mathematical objects. This would be good for solving integrals, but bad for manipulating polynomials in one variable with 1000's of terms.  Anyway,  I had trouble making the expression approach work via a Julia package, so it made sense to essentially make another language.  The particular syntax is not decided yet.

--John

[lapeyre....@gmail.com]

Maybe the parser could be useful if rewritten in Julia:
http://www.cs.berkeley.edu/~fateman/lisp/mma4max/parser.lisp

Another option is to rewrite it in Scheme, which is what Julia's parser is written in. I have hacked the mma4max parser some for a previous project, and I looked at the Julia parser. It looks like it might be a big undertaking.


btw. He has an mma5max now. Just a version number apparently.

[lapeyre....@gmail.com]

(My posts are not appearing where I want. Hope this one does!)

On Monday, January 26, 2015 at 6:57:48 PM UTC+1, James Crist wrote:

Maybe the parser could be useful if rewritten in Julia:

Side tracking a bit, but one of the reasons I use SymPy instead of something like Maxima is that it meshes seamlessly into the language infrastructure. By writing a parser and creating another system means that it can't use julia's ecosystem without modification. I'd rather see a CAS work within Julia's parser framework. But that's just me.

Yes, I played with both ideas. I think what I am looking for is a fundamentally different language from Julia. I can't see a way to make it an extension of Julia in the way that a normal package is. Although it may be possible to use it from within Julia in a convenient way. Writing classes to handle mathematical structures is more efficient in memory and runtime, but is less flexible than expressions.

 

I am more concerned now in getting correct logical behavior than efficient application of rules.

Making it just *work* is certainly a priority, but design decisions made now can make things tricky later on. Especially with the design of your internal representation. How you represent expression trees and assumptions can seriously affect the performance of matching/traversals/replacements.

I agree. But, the right decision is not always clear to me. So I am interested in the representation, not performance, because there is a good chance that parts will be thrown away.  I am fairly certain that I want to represent expressions in a way similar to a Julia, Mathematica, Maxima, or Lisp expression (arrays of pointers vs. linked list is another matter). Mathematica rules and patterns are also represented this way. Of course, I am sure they are used at some point to build data structures composed of types in some language. But, logically, and at the language level, everything is an expression.

[lapeyre....@gmail.com]

On Monday, January 26, 2015 at 7:20:54 PM UTC+1, Francesco Bonazzi wrote:

A Mathematica parser could be used as an optional addition, such as what Mathics is to SymPy.

Anyways, there are various options for Julia:

• use Expr and operate on them, this is similar to how Maxima works.
  

One experiment in the code I posted does this. And the main code does essentially the same thing, but it is not an Expr, but a Mxpr that
carries some additional fields.

> declare symbolic types and use operator overloading: this is the approach by SymPy and Axiom.

This is interesting and has advantages, but is not what I am looking for, at least not now. I  suppose hybridizing etc. is possible

 > create a DSL invoked by a macro, this appears to be the case for SJulia

Yes, the macro is just a way to allow more syntax than nested expressions. And this is really the same as your first item, as well.

Using the macro is much easier than writing or modifying a parser.

I expect that at some point,  all of these options will be tried with Julia.

[John Lapeyre]

To be more clear: First I think the ideas in this thread will be used
in Julia projects in interesting ways we aren't even thinking of---
it's well suited for all of them. At the moment, I'm only interested in
implementing a pattern and rule based system something like
Mathematica. An example of why, is the Rubi integration package

http://www.apmaths.uwo.ca/~arich/

The author claims that it "dramatically out-performs Maple and Mathematica". I have not tried
to verify this, but it is a serious project. It is written entirely with Mma
rules. It does not use Mma's integration system at all. Fateman's code
can employ at least a subset of these rules.

http://www.cs.berkeley.edu/~fateman/lisp/mma5max/rubi/

As far as I can see in mma5max, all expressions, including patterns and rules,
are maintained as singly linked lists during the whole procedure. It's
probably quite a bit slower than Mma running Rubi.  But, if you build
a CAS in Julia based on DataTypes representing mathematical objects,
it will have to be developed for a very, very, long time before it can
integrate like Rubi. Of course, like I already said, there are things
that Mathematica is inherently bad at. And other things, like tight
integration with the host language, might be more important for some
applications that being able to to 50K integrals.

Even if I or someone else, got the all the basic pattern matching
implemented, Rubi (for instance again) still uses Simplify and other
features. I imagine it would be very big project to develop to the
point of handling all the Rubi rules. Still, I am really surprised by
how quickly I was able to get as far as I have and how non-dismal the
performance is, even though I made no attempt to optimize the code at
any level.

--John

[Francesco Bonazzi]

On Tuesday, January 27, 2015 at 12:24:35 AM UTC+1, John Lapeyre wrote:

To be more clear: First I think the ideas in this thread will be used
in Julia projects in interesting ways we aren't even thinking of---
it's well suited for all of them. At the moment, I'm only interested in
implementing a pattern and rule based system something like
Mathematica. An example of why, is the Rubi integration package

http://www.apmaths.uwo.ca/~arich/

I read that the next version of Rubi will feature a decision tree, no longer pattern matching.
 

Even if I or someone else, got the all the basic pattern matching
implemented, Rubi (for instance again) still uses Simplify and other
features.

The easiest way would be to call SymPy's simplify from Julia via PyCall. Anyways, SymPy is BSD 3-clause licensed, you may copy from its code.

[John Lapeyre]

>     http://www.apmaths.uwo.ca/~arich/
>
> I read that the next version of Rubi will feature a decision tree, no longer pattern matching.

Interesting. I don't see it, do you have a link? I wonder if the switch is being made
for organizational, or efficiency reasons. I decision tree might be easier to implement in
more languages, which would be great. But, there still has to be an element of pattern
matching, via Mma or otherwise.

>     Even if I or someone else, got the all the basic pattern matching
>     implemented, Rubi (for instance again) still uses Simplify and other
>     features.

> The easiest way would be to call SymPy's simplify from Julia via PyCall. Anyways, SymPy is BSD 3-clause licensed, you may copy from its code.

Yes, except that they will behave differently, so it won't work directly. Of course you can see the
SymPy code, whereas the the Mma function is a black box.

--John

[Francesco Bonazzi]

On Tuesday, January 27, 2015 at 12:34:43 PM UTC+1, John Lapeyre wrote:

> I read that the next version of Rubi will feature a decision tree, no longer pattern matching.

Interesting. I don't see it, do you have a link?

https://github.com/sympy/sympy/issues/7749#issuecomment-54830230

By the way, I was thinking about the syntax, what about instead of

@ex f(x_, y_Integer, z) := ...

to use a more Julia-friendly one like

@ex f(x, y::Integer, :z) := ...

This does not need any string processing for underscore signs, colon prepended to the symbol (e.g. :z) to mark that z is not a wildcard pattern.

Anyways, it looks like on your patterns you are implicitly introducing an Mxpr <-> type correspondence, i.e. if in x_A the expression A were considered both as the head of an Mxpr and as a Julia type.

Any thoughts about interpolation of Julia variables in the @ex macro? I tried to call a=3; @ex f($a) but it doesn't work.

[lapeyre....@gmail.com]

On Wednesday, January 28, 2015 at 10:20:08 AM UTC+1, Francesco Bonazzi wrote:

On Tuesday, January 27, 2015 at 12:34:43 PM UTC+1, John Lapeyre wrote:

> I read that the next version of Rubi will feature a decision tree, no longer pattern matching.

Interesting. I don't see it, do you have a link?

https://github.com/sympy/sympy/issues/7749#issuecomment-54830230

 
   This looks great (at least based on the linked post) It no longer relies on the sophisticated part of the core Mma language, so it greatly broadens number of languages able to run Rubi. (Its only one now: Mma itself)  I guess it's not too much trouble to emit sympy code. I'm pretty sure this could be easily adapted to Maxima. In fact I wrote an Mma to Maxima translator a while ago (Hmm I should put it on github) so it might work immediately.  It uses RJF's  CL  Mma parser (uhh. too many acronyms). The translator can't handle Mma patterns (maybe just a liitle, I don't remember)  My translator could translate and pass tests for a fair size quantum information package. It also translated and correctly ran some symbolic quantum mechanics code used in research that was written by a colleague. But, it was unacceptably slow, because Maxima uses linked lists for expressions, which require fundamentally different algorithms, not just translation... I'm digressing.

 I guess Albert Rich still maintains the data in some other form ? Or is he editing these big nested conditional statements ...

 

By the way, I was thinking about the syntax, what about instead of

@ex f(x_, y_Integer, z) := ...

Well, I kind of tried different things, including :: at one point. I gradually moved away from a Julia extension to really another language: I just wanted the patterns, but they need support from expressions, but that requires a consistent evaluation sequence, etc. But, because I am relying on Julia's parser, the options are limited: the parser has to believe it's valid Julia code (although it won't eval). So, if I disallow underscores in identifiers and interpret them instead as signifying patterns, I don't consume any of the limited available Julia syntax. In Mma  you can have one, two, or three underscores before and/or after an identifier, each combination of which means something different. All of that has to be encoded in Julia syntax. By disallowing underscores, the whole problem is solved. I started to implement the Mma "Condition", which is yet another essential feature of the pattern matching. I tried :: for *that*, but it has high precedence, so you need to use parens, which I don't like. For now, I decided on just Condition(a,b). Furthermore, I stopped implementing Condition because it is relatively easy to do. AC matching is a b$&%h.  For me, there's no sense in continuing the experiment unless I'm convinced the harder stuff can be done.  Of course if you, or someone else want a particular feature to play with, I'd gladly consider implementing it. At present, I consider using the stock Julia parser a stopgap. I did spend some time with the RJF parser and even found a minor bug or two. I also looked at the Julia parser.  For me, it looks like a big job. But, if someone else wants to do it, great!

> Anyways, it looks like on your patterns you are implicitly introducing an Mxpr <-> type correspondence, i.e. if in x_A the expression A were considered both as the head of an Mxpr and as a Julia type.

Yeah, thats kinda it. Mma tries to maintain a fiction that non-expression types are really expressions with a particular head and zero length. So in your example, A might be a type or the head of an expression. I make the symbol 'Integer' Protected (except, I don't think I have yet put it on that list in the code), then under certain circumstances, it is evaled to get a DataType.

>  Any thoughts about interpolation of Julia variables in the @ex macro? I tried to call a=3; @ex f($a) but it doesn't work.

Yes, this works:

a = 3;  @ex f(JVar(a))

Not pretty, but I don't want to introduce new syntax/language features at this point. JVar is done very easily in the standard Mma way.

 apprules(mx::Mxpr{:JVar}) = eval(symname(mx.args[1]))

That's it. Note that this evaluates a, then evals the result, etc. until it stops at an SJulia symbol that evaluates to itself. This is the Mma behvavior, which
I am sticking with for now. If the attribute 'HoldAll' is set for 'JVar', then JVar(a) always evaluates the Julia symbol :a.  Both would be useful, I just didn't implement that yet.  I'll document this . Please ask any other questions you may have.
Note that I pollute the Julia side, by setting undefined Julia vars to true, so that the repl completes them in SJulia. You can
turn this off by commenting out a line in mxpr_type.jl .

--John

[lapeyre....@gmail.com]

Francesco, I just realized that Mathics is actually aiming to implement the Mma language. For some reason I thought it was only using the same syntax. So, this is indeed similar to what I am trying to do. Except that I consider the code that I have now to be an experiment.

On Wednesday, January 28, 2015 at 10:20:08 AM UTC+1, Francesco Bonazzi wrote:

[Francesco Bonazzi]

On Wednesday, January 28, 2015 at 5:53:00 PM UTC+1, lapeyre....@gmail.com wrote:

AC matching is a b$&%h.  For me, there's no sense in continuing the experiment unless I'm convinced the harder stuff can be done.

https://github.com/mathics/Mathics/blob/1a5cd204b76595730bdf4a91a6696a3292159d1a/mathics/core/util.py#L59

These are the util functions of Mathics for its own AC matching (GPL-licensed). They are lazy permutation generators, maybe you can find them useful: subranges for associative matching, subsets for associative-commutative matching. You can also use them for non-exhaustive matching, i.e. if you want just part of the node arguments to be matched. They are lazy python generators of subranges and subsets, they observe exactly the same order as Mathematica's matching priority.

Hope this can be useful.

Somewhere else there is also a key generator to sort DownValues with Python's sort( ) function, again observing the same order as Wolfram Mathematica.
 

Of course if you, or someone else want a particular feature to play with, I'd gladly consider implementing it. At present, I consider using the stock Julia parser a stopgap. I did spend some time with the RJF parser and even found a minor bug or two. I also looked at the Julia parser.  For me, it looks like a big job. But, if someone else wants to do it, great!

To get advanced feature working, such as integration, limits, equation and differential equation solving, I suggest to create an expression tree converter from and to SymPy using PyCall. It's a trick, but it's the only way you can achieve advanced features with little work.

[John Lapeyre]

> https://github.com/mathics/Mathics/blob/1a5cd204b76595730bdf4a91a6696a3292159d1a/mathics/core/util.py#L59
>
> These are the util functions of Mathics for its own AC matching (GPL-licensed). They are lazy permutation generators, maybe you can find them useful: subranges for associative matching, subsets for associative-commutative matching. You can also use them for non-exhaustive matching, i.e. if you want just part of the node arguments to be matched. They are lazy python generators of subranges and subsets, they observe exactly the same order as Mathematica's matching priority.
>

Thanks, this is indeed useful. The mmaXmax pattern matching code implements
the same two functions. The Julia Iterators package has subsets, but
I don't find the equivalent of subranges.

> Somewhere else there is also a key generator to sort DownValues with Python's sort( ) function, again observing the same order as Wolfram Mathematica.

For now at least, my code is licensed under the MIT license. Maybe I can find
something similar in mmaXmax, which has a less restrictive license. I wonder how Mathics got the information, reverse engineering maybe. The Wolfram docs that I saw were a bit vague.

>     Of course if you, or someone else want a particular feature to play with, I'd gladly consider implementing it. At present, I consider using the stock Julia parser a stopgap. I did spend some time with the RJF parser and even found a minor bug or two. I also looked at the Julia parser.  For me, it looks like a big job. But, if someone else wants to do it, great!
>
>
> To get advanced feature working, such as integration, limits, equation and differential equation solving, I suggest to create an expression tree converter from and to SymPy using PyCall. It's a trick, but it's the only way you can achieve advanced features with little work.

I was actually thinking of features closer to the core language at the moment. But,
I also think that down the road, calling python code would be the
easiest of available choices (Maxima, etc.), or maybe some specialized
C libraries.

--John

[Ondřej Čertík]

Hi John,

I just discovered this thread, some comments:

On Wed, Jan 28, 2015 at 9:53 AM, <lapeyre....@gmail.com> wrote:
>
>
> On Wednesday, January 28, 2015 at 10:20:08 AM UTC+1, Francesco Bonazzi
> wrote:
>>
>>
>>
>> On Tuesday, January 27, 2015 at 12:34:43 PM UTC+1, John Lapeyre wrote:
>>>
>>> > I read that the next version of Rubi will feature a decision tree, no
>>> > longer pattern matching.
>>>
>>> Interesting. I don't see it, do you have a link?
>>
>>
>> https://github.com/sympy/sympy/issues/7749#issuecomment-54830230
>>
>
> This looks great (at least based on the linked post) It no longer relies

I wrote that post. I also created a PR with actual SymPy code here:

https://github.com/sympy/sympy/pull/8036

It seems to work great so far. I am waiting for more code from Albert
(I CCed him) and also I need to polish the PR some more.

I also tested Rubi using Mathematica, and it seems faster and more
robust (it can also do more integrals) than Mathematica's Integrate[]
function. Albert said that on his benchmarks, the if/then/else form is
massively faster than the pattern matching (maybe up to 100x on some
examples).

> on the sophisticated part of the core Mma language, so it greatly broadens
> number of languages able to run Rubi. (Its only one now: Mma itself) I
> guess it's not too much trouble to emit sympy code. I'm pretty sure this
> could be easily adapted to Maxima. In fact I wrote an Mma to Maxima
> translator a while ago (Hmm I should put it on github) so it might work
> immediately. It uses RJF's CL Mma parser (uhh. too many acronyms). The
> translator can't handle Mma patterns (maybe just a liitle, I don't remember)
> My translator could translate and pass tests for a fair size quantum
> information package. It also translated and correctly ran some symbolic
> quantum mechanics code used in research that was written by a colleague.
> But, it was unacceptably slow, because Maxima uses linked lists for
> expressions, which require fundamentally different algorithms, not just
> translation... I'm digressing.
>
> I guess Albert Rich still maintains the data in some other form ? Or is he
> editing these big nested conditional statements ...

My understanding is that unfortunately he edits the nested conditional
statements. Apparently it is not easy to automatically translate the
pattern matching rules into a polished if/then/else form, though
ultimately it would be nice to have that. The if/then/else is nice,
that you can step through it easily via a debugger, and the only way
this might not finish is if the functions call each other recursively.

If you are interested in helping out with Rubi, please write to
Albert. He might really use some help with this, and lots of projects
would benefit. He maintains it in Mathematica, and then have automatic
translators to other codes, including Maxima and SymPy.

Finally I would mention that in SymPy, we are also developing a very
fast C++ core (https://github.com/sympy/csympy), which is actually
pure C++, the Python wrappers and seamless SymPy integration (using
the Python wrappers) is optional. We also have a C interface
(https://github.com/sympy/csympy/blob/master/src/cwrapper.h), though
it needs to be extended to cover all CSymPy's functionality. The goal
of the C interface is that it can now be trivially called from Julia,
here is some preliminary code from Isuru (also CCed):

https://gist.github.com/isuruf/cb4f16d05b7266de227b

If any of you would like to play with it and help us write a nice
Julia interface, that would be awesome. Then we can easily benchmark
it against any other computer algebra implementation in Julia. The
goal of CSymPy is match performance of any other computer algebra
system (or be faster). We are still working on it.

Ondrej

[lapeyre....@gmail.com]

(I am resending because, this my message to Ondřej did not go to the list)

Hi Ondřej,

The various sympy projects are very interesting, your work porting Rubi especially.

SJulia is using SymPy for a lot of features (Integrate, Solve, ...) . Here is some test code

 https://github.com/jlapeyre/SJulia/blob/master/test/sympy_test.jl

Francesco Bonazzi prompted me and started the code. From the little I've seen,
Sympy has impressive capabilities and the performance is good.

--John

[lapeyre....@gmail.com]

BTW.

Do you know much about the SymPy interface ? (or maybe the authors are reading this ?)
Francesco has already supplied a lot of information, but I am new to sympy.

Here, I ask for an integer to be converted with _to_mpmath

https://github.com/jlapeyre/SJulia/blob/master/src/sympy.jl#L89

And find that if the integer is actually a machine integer, then I get an Int64. So big
integers go to big integers and machine integers to machine integers.

Is there something similar I can do to get big float and machine floating point to map correctly
to Julia ?   I'm not sure how SymPy encodes  complex and rationals, but the code I have seems
to handle all of these correctly.

On Thursday, February 26, 2015 at 4:30:04 AM UTC+1, Ondřej Čertík wrote:

[Ondřej Čertík]

On Thu, Feb 26, 2015 at 3:44 AM, <lapeyre....@gmail.com> wrote:
> (I am resending because, this my message to Ondřej did not go to the list)

I can see my message right here:
https://groups.google.com/d/msg/julia-users/jP1VlOe80HQ/e4KzNARuI9AJ

> Do you know much about the SymPy interface ? (or maybe the authors are
> reading this ?)

I don't know who wrote SymPy interface to Julia.

> Francesco has already supplied a lot of information, but I am new to sympy.
>
> Here, I ask for an integer to be converted with _to_mpmath
>
> https://github.com/jlapeyre/SJulia/blob/master/src/sympy.jl#L89
>
> And find that if the integer is actually a machine integer, then I get an
> Int64. So big
> integers go to big integers and machine integers to machine integers.
>
> Is there something similar I can do to get big float and machine floating
> point to map correctly
> to Julia ? I'm not sure how SymPy encodes complex and rationals, but the
> code I have seems
> to handle all of these correctly.

SymPy uses mpmath to handle arbitrary precision floating point
numbers. SymPy can also use gmpy if it is installed, and that uses GMP
integers.
As to integers, we just use Python integers, so you just need to
figure out how to convert those. Rationals are defined here:

https://github.com/sympy/sympy/blob/8df828ea194f7a86dd900a80303e51c06ef3ec60/sympy/core/numbers.py#L1046

It has the "p" and "q" members for nominator and denominator. You can
find there the Integer and Float classes as well.

Complex numbers are represented symbolically using the "I" (=imaginary
unit) constant. In CSymPy, on the other hand, we represent complex
numbers using a Complex class with real and imaginary part (as
rationals):

https://github.com/sympy/csympy/blob/3fee2eb0ebd18099343b2f42670377aadf33e3cd/src/complex.h#L16

Ondrej

[j verzani]

I think PyCall already has what you need for conversion:

PyObject(big(pi)) will create a an mpf instance of a big float, like big(pi)

convert(BigFloat, PyObject(big(pi)))  will return a BigFloat from an mpf instance.

[lapeyre....@gmail.com]

Thanks for the info Ondrej and John.