Recurrences and Sequence Formula Guessing Functions in Sage

466 views
Skip to first unread message

Maxie Schmidt

unread,
Nov 19, 2016, 10:56:42 AM11/19/16
to sage-devel
Hello,

I have been working back and forth between Sage and Mathematica for a while now trying to learn how to use Sage to replace Mathematica's core functionality in my day to day use of it. One of the (sets of) functions haven't yet found a suitable open source alternative for is related to guessing formulas and generating functions for an input sequence (as in Mathematica's FindSequenceFunction and FindGeneratingFunction). Related functions are Mathematica's RSolve (https://reference.wolfram.com/language/ref/RSolve.html?q=RSolve) used to solve recurrence relations in closed-form formulas, and the closed-source RISC Sigma package (https://www.risc.jku.at/research/combinat/risc/software/Sigma/index.php) which is able to generate recurrences for many sums and simplify sums involving harmonic numbers.

The closest open source alternative I have come across so far is Rubey's software written in FriCAS described at http://axiom-wiki.newsynthesis.org/GuessingFormulasForSequences. However, looking at the source code for the package leaves quite a bit of work to do to reproduce the output of the Mathematica functions mentioned above. Is there a better open source alternative to the Mathematica functions FindSequenceFunction and RSolve that I'm missing in Sage?

If there isn't a good complete replacement for using Mathematica's routines for sequence formula guessing and solving recurrences, I'm interested in trying to write functions / packages that will replace Mathematica for these tasks. I have already done some work in Sage related to guessing formulas for special polynomial sequences in my Master's thesis at UIUC posted at https://arxiv.org/abs/1609.07301. Any suggestions or links to related open source software are appreciated.

Sincerely,

Maxie

Martin R

unread,
Nov 19, 2016, 12:30:09 PM11/19/16
to sage-devel
I'd be interested in what output you'd like to have.

The hard part in the FriCAS package was to get decent speed, changing output should be relatively straightforward.

(I guess that you are aware of the possibility of using the package from within sage)

Martin

Maxie Schmidt

unread,
Nov 19, 2016, 12:51:01 PM11/19/16
to sage-devel
I'm aware that the FriCAS package can be called from within Sage / Python, and I don't see any reason to rewrite what's already there, just extend it's current recognition capabilities and put a nice friendly wrapper function around all of these disparate routines. Something like find_sequence_formula([1,1,2,3,5,8,13]) or find_polynomial_generating_function([...], x) would be nice. I'm willing to work on writing the wrapper functions and help with extending the sequence recognition capabilities that are already there. It's been a little while since I've seriously thought about writing this sort of function in Sage, but my intuition is that the FriCAS package can already pick up a lot of what FindSequenceFunction will recognize out of the box.

I can / should put together a long list of sequence test cases that I'd like to see recognized by these routines in Sage. Anything that's in the Wolfram / Mathematica documentation is on the list. Also, sequences involving factors, or polynomial multiples of, special sequences like the Stirling, Bernoulli, and r-order harmonic numbers are something I'd like to see.

How much does the existing FriCAS rely on being able to find a homogeneous recurrence for the input sequences? There are a number of examples of sums involving the Stirling and Bernoulli numbers that cannot be expressed this way. Also, Mathematica tends to ignore closed-form sequence expressions involving harmonic numbers (say for known formulas satisfied by the Stirling numbers of the first kind), whereas if we have the FriCAS package's analog to the RISC Guess package along with an OSS version of Sigma and RSolve, these expressions are easy to obtain. These are just some suggestions at a starting point for replacing my having to switch over to Mathematica when I need sequence recognition functions.

Maxie

Dima Pasechnik

unread,
Nov 19, 2016, 1:23:13 PM11/19/16
to sage-devel
Needless to say perhaps, but Maxima has some capabilities in this respect, not sure how much of it is exposed in Sage.

E.g. ggf should be able to do rational generating functions:

And it does have things like Gosper and Zeilberger algorithms for finding generating functions
implemented too.

Fredrik Johansson

unread,
Nov 19, 2016, 2:00:27 PM11/19/16
to sage-devel

Open source alternatives to Mathematica would be most welcome.

The ore_algebra package allows guessing holonomic sequences (http://www.risc.jku.at/research/combinat/software/ore_algebra/doc/guessing.html) but there is nothing in Sage for closed forms, Stirling-type sequences or multi-sums as far as I know.

Fredrik

David Roe

unread,
Nov 19, 2016, 2:32:04 PM11/19/16
to sage-devel
Another thing that might be nice to tie in is the Online Encyclopedia of Integer Sequences.  Or at least include a link in the documentation, though most people looking to guess a sequence will already be aware of OEIS.
David

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscribe@googlegroups.com.
To post to this group, send email to sage-...@googlegroups.com.
Visit this group at https://groups.google.com/group/sage-devel.
For more options, visit https://groups.google.com/d/optout.

Thierry

unread,
Nov 19, 2016, 2:40:49 PM11/19/16
to sage-...@googlegroups.com
On Sat, Nov 19, 2016 at 02:31:42PM -0500, David Roe wrote:
> Another thing that might be nice to tie in is the Online Encyclopedia of
> Integer Sequences. Or at least include a link in the documentation, though
> most people looking to guess a sequence will already be aware of OEIS.
> David


sage: oeis?
> > email to sage-devel+...@googlegroups.com.
> > To post to this group, send email to sage-...@googlegroups.com.
> > Visit this group at https://groups.google.com/group/sage-devel.
> > For more options, visit https://groups.google.com/d/optout.
> >
>
> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

Martin R

unread,
Nov 19, 2016, 2:42:33 PM11/19/16
to sage-devel
out off the box, the FriCAS package can do the following (described in some detail in https://arxiv.org/abs/math/0702086, the journal published a shortened version):

1) Generating functions:

* rational (guessPade)
* algebraic (guessAlg)
* linear diffeq (guessHolo)
* polynomial diffeq (guessADE)

* Mahler type functional equations (guessFE)

* algebraic dependencies between a list of generating functions (guessAlgDep)

2) Recurrences:

* linear, polynomial  coefficients (guessPRec)
* polynomial, polynomial coefficients (guessRec)

3) Nested sums and products of the form sum_{k=1}^n f(k) and \prod_{k=1}^n f(k), the most useful case being where f is a rational function in k.  Note that f must not depend on n, unfortunately.

4) Moreover, some special forms (guessBinRat, guessExpRat), which were not really successful.

All of this works for lists of integers, algebraic numbers, polynomials, and over finite fields.
Also, there are q-analogues of all these (most interesting are probably q-differential equations).

There are some options, that can be used to restrict the search space, for example, restricting to homogeneous differential equations, restricting the degree of the polynomial coefficients, etc.  If there is sufficient interest, Jay Pantone found a (heuristically) better search space for polynomial diffeqs, which would be straightforward to implement.

Martin

Martin R

unread,
Nov 19, 2016, 2:49:39 PM11/19/16
to sage-devel
I forgot to add: once a diffeq for a generating function is found, the package can also give you the 1783-th coefficient of the generating function.

And, most importantly: the main thing to do to make it really useful for sage, is to implement a better conversion from sage lists to fricas lists, bypassing pexpect.  This is possible (Waldek Hebisch provided the necessary bits on the FriCAS side), but some work.

Martin

David Roe

unread,
Nov 19, 2016, 2:50:32 PM11/19/16
to sage-devel
On Sat, Nov 19, 2016 at 2:40 PM, Thierry <sage-goo...@lma.metelu.net> wrote:
On Sat, Nov 19, 2016 at 02:31:42PM -0500, David Roe wrote:
> Another thing that might be nice to tie in is the Online Encyclopedia of
> Integer Sequences.  Or at least include a link in the documentation, though
> most people looking to guess a sequence will already be aware of OEIS.
> David



        sage: oeis?

I just meant, in whatever functions Maxie is planning on writing, either use oeis as (part of) one of the guessing regimes, or at least include a link to Sage's oeis function in the docstring.
David


> > To post to this group, send email to sage-...@googlegroups.com.
> > Visit this group at https://groups.google.com/group/sage-devel.
> > For more options, visit https://groups.google.com/d/optout.
> >
>
> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscribe@googlegroups.com.

> To post to this group, send email to sage-...@googlegroups.com.
> Visit this group at https://groups.google.com/group/sage-devel.
> For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscribe@googlegroups.com.

Maxie Schmidt

unread,
Nov 19, 2016, 3:23:10 PM11/19/16
to sage-devel, roed...@gmail.com
Thanks for all of the suggestions. I think I will start by getting all of Martin's package routines working under a single wrapper function and try to improve from there. I will definitely include options to search the OEIS database.

The only thing I'm a little concerned about if this would eventually get bundled in Sage is getting the FriCAS source for the first package installed. My Debian box doesn't have FriCAS or Axiom bundled with my from source Sage install, and I remember a while back having multiple issues getting FriCAS installed on a Gentoo machine.

Maxie

> > To post to this group, send email to sage-...@googlegroups.com.
> > Visit this group at https://groups.google.com/group/sage-devel.
> > For more options, visit https://groups.google.com/d/optout.
> >
>
> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

> To post to this group, send email to sage-...@googlegroups.com.
> Visit this group at https://groups.google.com/group/sage-devel.
> For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

Martin R

unread,
Nov 19, 2016, 4:18:12 PM11/19/16
to sage-devel, roed...@gmail.com
Am Samstag, 19. November 2016 21:23:10 UTC+1 schrieb Maxie Schmidt:
Thanks for all of the suggestions. I think I will start by getting all of Martin's package routines working under a single wrapper function and try to improve from there.
 
I think that one question is how you want to represent the result of a "guess".  The fricas package always returns a list of elements of the fricas domain "Expression Integer", which is roughly equivalent to SymbolicRing in sage.  For example:

sage: list_g = fricas("guessADE([n^n/factorial n for n in 1..40])"); list_g

                                                           
2      3
     n          
,            3        2                   9x    32x       4
 
[[[x ]f(x): - f (x) + x f(x)  + 2f(x) = 0,f(x)= 1 + 2x + --- + ---- + O(x )]]
                                                           
2      3
sage
: g = list_g[0]
sage
: g.eval("[n=50]")

2635987571800400710022355920186541715689493871043967230302851194475838001
-------------------------------------------------------------------------
           
3353585829598799392450336149309868736512000000000000

The FriCAS "type" of g above is Expression Integer, but I don't think there is an exact equivalent in the current SymbolicRing, because g contains the differential equation, the initial values and (hidden in the display) a (slow, but generic) function that computes an arbitrary term.

I assume a convenient possibility would be to create a sage class that wraps this information...


The only thing I'm a little concerned about if this would eventually get bundled in Sage is getting the FriCAS source for the first package installed. My Debian box doesn't have FriCAS or Axiom bundled with my from source Sage install, and I remember a while back having multiple issues getting FriCAS installed on a Gentoo machine. 
 

Are you saying that "sage -i fricas" does/did not work?

Martin

Dima Pasechnik

unread,
Nov 19, 2016, 4:49:30 PM11/19/16
to sage-devel


On Saturday, November 19, 2016 at 8:23:10 PM UTC, Maxie Schmidt wrote:
Thanks for all of the suggestions. I think I will start by getting all of Martin's package routines working under a single wrapper function and try to improve from there. I will definitely include options to search the OEIS database.

The only thing I'm a little concerned about if this would eventually get bundled in Sage is getting the FriCAS source for the first package installed. My Debian box doesn't have FriCAS or Axiom bundled with my from source Sage install

Fricas is an experimental Sage package.
Just do

./sage -i fricas

Assuming that your ./sage is in good shape (say, compiled from source), you'll be all set.

Peter Luschny

unread,
Nov 19, 2016, 6:44:56 PM11/19/16
to sage-devel
> One of the (sets of) functions haven't yet found a suitable open source alternative
> for is related to guessing formulas and generating functions for an input sequence
> (as in Mathematica's FindSequenceFunction and FindGeneratingFunction).

I am amazed that no one has yet mentioned Gfun, the classic Maple package for 
the manipulation of generating functions by Bruno Salvy and others which is 
released under the GNU Lesser General Public License (LGPL v2.1). 


Peter
 

Ralf Stephan

unread,
Nov 20, 2016, 1:36:12 AM11/20/16
to sage-devel
1. Sage can guess rational generating functions and provide formulae for experimental math:

sage: C.<x> = CFiniteSequences(QQ)
sage
: C.guess([1,3,5,7,9,11,13])
C
-finite sequence, generated by (x + 1)/(x^2 - 2*x + 1)
sage
: _.closed_form()
2*n + 1
sage
: C.guess([0,1,1,2,3,5,8,13,21]).closed_form()
1/5*sqrt(5)*(1/2*sqrt(5) + 1/2)^n - 1/5*sqrt(5)*(-1/2*sqrt(5) + 1/2)^n

2.As to guessing holonomic recurrences there is the Ore algebra package. It would be nice to have a tutorial.

3. Subham Tibra has added a holonomic function package to SymPy:
which is however in a newer version than the SymPy in Sage

Regards,

Martin R

unread,
Nov 20, 2016, 2:13:25 AM11/20/16
to sage-devel
Did anybody ever compare the possibilities speed-wise?

The scope of these packages is quite different: ore_algebra is for the holonomic and q-holonomic universe, the guessing facility in fricas is for guessing only - but covers also polynomial recurrences and differential equations (and their q-analogues), functional equations, etc...

Two remarks:

* the ore_algebra package is currently broken, see #21344.
* gfun is free, but maple is not.

Martin

Martin R

unread,
Nov 20, 2016, 2:28:28 AM11/20/16
to sage-devel
2.As to guessing holonomic recurrences there is the Ore algebra package. It would be nice to have a tutorial.
 

Martin

Ralf Stephan

unread,
Nov 20, 2016, 3:38:42 AM11/20/16
to sage-devel
On Sunday, November 20, 2016 at 8:28:28 AM UTC+1, Martin R wrote:
2.As to guessing holonomic recurrences there is the Ore algebra package. It would be nice to have a tutorial.
 

Let me rephrase that. It would be nice to have ore-algebra as standard package and a Sage-internal tutorial.
 

Maxie Schmidt

unread,
Nov 25, 2016, 8:08:32 PM11/25/16
to sage-devel
I have a working draft of the sage code for sequence formula guessing posted at https://github.com/maxieds/sage-guess. I'm still very much working on adding the documentation / doc strings. Besides this, does anyone else have any suggestions for what else I could add to this package to improve it? It's obviously missing support for polynomial sequences right now. I'm planning on adding a tracking ticket and a tutorial once I have the documentation done and more tests written.

Maxie

Dima Pasechnik

unread,
Dec 7, 2018, 11:13:07 AM12/7/18
to sage-devel


On Saturday, November 26, 2016 at 1:08:32 AM UTC, Maxie Schmidt wrote:
I have a working draft of the sage code for sequence formula guessing posted at https://github.com/maxieds/sage-guess. I'm still very much working on adding the documentation / doc strings. Besides this, does anyone else have any suggestions for what else I could add to this package to improve it? It's obviously missing support for polynomial sequences right now. I'm planning on adding a tracking ticket and a tutorial once I have the documentation done and more tests written.

Maxie

I just did a pull request there to make that code run with current Sage. 
Reply all
Reply to author
Forward
0 new messages