analytic combinatorics: new code for computing asymptotics for multivariate generating functions

[3b09ef35-8a33-427f-b580-f9c14c161f29]

This code is a collection of functions designed to compute asymptotics of Maclaurin coefficients of certain classes of multivariate generating functions.

The main function asymptotics() returns the first N terms of the asymptotic expansion of the Maclaurin coefficients F_{n\alpha} of the multivariate meromorphic function F=G/H as n\to\infty. It assumes that F is holomorphic in a neighborhood of the origin, that H is a polynomial, and that asymptotics in the direction of \alpha (a tuple of positive integers) are controlled by smooth or multiple points.

CC: @sagetrac-cyril-banderier @zafeirakopoulos @sagetrac-tmonteil @dkrenn @behackl

Component: combinatorics

Keywords: analytic combinatorics, multivariate generating functions

Author: Daniel Krenn, Alex Raichev

Branch/Commit: a951f08

Reviewer: Daniel Krenn, David Loeffler, Travis Scrimshaw

Issue created by migration from https://trac.sagemath.org/ticket/10519

[3b09ef35-8a33-427f-b580-f9c14c161f29]

Attachment: mgf.sage.gz

[1d7ec08f-60ae-4512-91a6-8324c06eab9f]

comment:2

Hi Alex,

I can't address the mathematics on this one, but can give some hints about getting this reviewed.

1. You'll want to submit this as a patch - then folks can run the doctests, check the documentation formatting. There is a lot of checking that can't happen easily with a *.sage file. (See below.) Read the "Walk Through" section of the developer guide for help.

2. Your examples need slightly different formatting, like

    EXAMPLES::

        sage: input 1
        output-1

    Something nice to say. ::

        sage: input-2
        output-2

Double-colons start verbatim blocks and sometimes print as colons.

Sage code samples needs indentation.

3. Sage is mostly object-oriented. "Internal" method names should begin with an underscore, then they won't show in tab-completion etc (Python convention).

All your functions won't be available in teh global namespace (where they could conflict with other). You should find an object (like maybe some kind of polynomial, especially if it is already int eh combinatorial library?) and make your functions be methods on that object.

Thanks for your contribution - I hope the above is helpful for you.

Rob

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:3

Question for you experienced Sage developers: how can i better package my code?

I think it makes sense to move the major functions such as 'asymptotics' to methods in the class /Applications/sage/devel/sage/sage/rings/polynomial/multi_polynomial.pyx, but what about the subsidiary functions that don't concern polynomials in particular such as 'diff_prod' which calculates derivatives of products and solves equations?

[seblabbe]

comment:4

Question for you experienced Sage developers: how can i better package my code?

Hi, I am currently at Sage Days 28. Right now, there is a discussion about Analytic combinatorics in Sage and your code was mentionned in the discussion. I am not an expert of the domain, but I am coding oriented object Python since some time now. So below are just some of my thoughts to, I hope, help you turn your set of functions into an oriented object structure.

How to structure a bunch of functions into classes? How to find which objects (python classes) you need? Here is the trick I personaly use. Consider each of your functions as a question you ask. Then, ask yourself to who are you asking each of your questions? Answers often gives you a good hint about the objects you need to implement. EXAMPLE. Suppose I code the function determinant. Question : To who do I ask the determinant?. Answer: To a matrix. Hence, matrix might be a good object (a python class) to implement.

You are the best person to answer to these questions. You might have 30 functions in you file, but only two or three different answers to the above question. Regroup the similar functions together: they will become the methods of a same class.

The sage file you uploaded starts with :

> This code relates to analytic combinatorics.
> More specifically, it is a collection of functions designed
> to compute asymptotics of Maclaurin coefficients of certain classes of
> multivariate generating functions.

> The main function asymptotics() returns the first `N` terms of
> the asymptotic expansion of the Maclaurin coefficients `F_{n\alpha}`
> of the multivariate meromorphic function `F=G/H` as `n\to\infty`.
> It assumes that `F` is holomorphic in a neighborhood of the origin,
> that `H` is a polynomial, and that asymptotics in the direction of
> `\alpha` (a tuple of positive integers) are controlled by smooth
> or multiple points.

Reading only these lines, I imagine the following structure:

class HolomorphicMultivariateMeromorphicFunction(object):

    # Constructor of the object
    def __init__(self, F, G):
        #stores important information on the object as attributes of self
        self._F = F
        self._G = G

    def maclaurin_coefficients(self, n, alpha):
        r"""
        Return the maclaurin coefficients of self.

        INPUT:

        - ``alpha`` - tuple of positive integers

        OUTPUT:

        a python list of the first terms 

        OR 

        maybe an object of a class you implement if there exists pertinent
        questions to ask to it.
        """
        #Do some computations based (I guess) on self._F and self._G
        intermediate_result1 = self.some_intermediate_computations_1()
        #Do more computations
        return something

    def asymptotics(self, N, alpha):
        r"""
        Returns the asymptotics of Maclaurin coefficients.
        """
        #Do some computations based (I guess) on self._F and self._G
        intermediate_result2 = self.some_intermediate_computations_2()
        intermediate_result3 = self.some_intermediate_computations_3()
        return something

    #put here all the others functions needed to compute the asymptotics
    def some_intermediate_computations_1(self):
        pass
    def some_intermediate_computations_2(self):
        pass
    def some_intermediate_computations_3(self):
        pass

    ...

It also looks like you need some robustness somehow. But I need to know more information about what means

that asymptotics in the direction of \alpha (a tuple of positive integers) are controlled by smooth or multiple points.

to decide whether this is checked at the creation of the object or before returning the asymptotics. But these hypothesis should be checked somewhere.

Hope this helps.

Cheers,

Sébastien Labbé, Montréal, (but currently at Sage Days 28, Orsay, France)

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:6

Hi Sebastien:

Thanks very much for your helpful comments. I've been working on rewriting my code in an object-oriented way and will incorporate your suggestions. Some functions in mgf.sage, such as deciding whether a list of polynomials is algebraically dependent, are general purpose, and so i'm trying to put them into existing Sage classes. The others are more specific to asymptotics and deserve their own class similar to the one you outlined above. Robustness is also something i need to improve upon.

I'll post my modifications here after i write and test them more.

Thanks again.

Alex

[3b09ef35-8a33-427f-b580-f9c14c161f29]

Attachment: amgf_alpha_0_5.sage.gz

New object oriented version

[3b09ef35-8a33-427f-b580-f9c14c161f29]

Attachment: amgf.sage.gz

Fixed a bug in cohomologous_integrand()

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:7

Hi Sebastien et al:

I rewrote my code in object-oriented style. Please have a look if you have the time. Do i need to make the file into a patch, or is it useable as is?

Thanks for your help. Alex

[seblabbe]

comment:9

Cool. I now see that a lot of the functions are now inside a class, which looks better. So, as I said earlier, this is not really my field of study, so I don't know how far I can go on that review. But, I can suggest at least two things:

1. It would be better if you post a Mercurial patch instead of a Sage file so that people can apply it and test it. Are you able to do this? Take a look at the Sage Developper's Guide, where it should be explain how to create a patch.

To add your file to the documentation, you might want to look here.

This winter, I gave a presentation on How to Contribute to Sage. This may helps you to start with Mercurial for instance.

Also, right now, I am working on #11379 where I am adding a new file to Sage. You might find it helpfull to see how I proceed.

2. Some the Python convention are not always followed. See also pep-0008. Especially for spaces :

      Yes:

          i = i + 1

      No:

          i= i+1

[seblabbe]

comment:10

To add a file to a patch, you must do:

sage -hg add your_file.py

I was not able to find this information anywhere which is strange...

Tell me if you have questions.

Sébastien

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:11

Thanks for your help, Se'bastien. I'm trying to put my Sage file into a patch by following your presentation slides, but am running into difficulties. I cloned Sage and copied my Sage file as 'amgf.sage' into

/Applications/sage/devel/sage-araichev/sage/combinat/

I rebuilt Sage etc. and made it to step 10 "Test the changes" of your slides. But when i run Sage and test a few commands, Sage doesn't seem to see the QuasiRationalExpression class i created in /Applications/sage/devel/sage-araichev/sage/combinat/amgf.sage.

I thought this might be because amgf.sage is a Sage file and not a Python file. So i renamed it amgf.py instead and got syntax errors. What's up with that? Do i have to rewrite all my code in Python instead of Sage?

Thanks for your attention. Alex

[hivert]

comment:12

Hi Alex,

I thought this might be because amgf.sage is a Sage file and not a Python file. So i renamed it amgf.py instead and got syntax errors. What's up with that? Do i have to rewrite all my code in Python instead of Sage?

The basic answer is yes, but Rewriting is a big word for what is really needed. There is little work to do since Sage mostly follows Python syntax. The two main difference are handling of integer (see http://www.sagemath.org/doc/tutorial/afterword.html), and the necessity to import what you need.

Handling of integer

• Notation for exponentiation: In Python ** means exponentiation and ^ means “xor”.

• If you need to return an integer for the user, write it return Integer(1) instead of return 1. In Python 1 is a machine integer int (32 or 64 bits depending on your machine) and Integer(1) is a Sage/Gmp arbitrary precision integer. Also Integer are much more powerful than int, for example they know about prime and factorization.

Importing stuff

The second big change is the necessity to import all what you need. More precisely, each time you use some Sage function, you need to import it at the beginning of the file. for example if you want to you PolynomialRing, you need to write

from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing

You can ask Sage where to find PolynomialRing using:

sage: PolynomialRing.__module__
'sage.rings.polynomial.polynomial_ring_constructor'

This also correspond to the path starting after site-packages given when you are asking Sage for PolynomialRing help:

sage: PolynomialRing?
Type:           function
[...]
File:           /home/florent/src/Sage/sage/local/lib/python2.6/site-packages/sage/rings/polynomial/polynomial_ring_constructor.py
[...]

I hope this helps,

Florent

[seblabbe]

comment:13

Just answering as a complement to Florent's answer.

But when i run Sage and test a few commands, Sage doesn't seem to see the QuasiRationalExpression class i created in /Applications/sage/devel/sage-araichev/sage/combinat/amgf.sage.

To test your code in a Sage session, you will first need to import the good class (or function). For instance :

sage: from sage.combinat.amgf import QuasiRationalExpression

If you would like QuasiRationalExpression to be already imported once Sage is started, you must add the line from sage.combinat.amgf import QuasiRationalExpression to the file sage/combinat/all.py. Note that it is not "à la mode" to add new stuff to all.py files because Sage is very slow to start because of all the things that are imported. If you think a lot of people will use the code, then it might be defendable.

Sébastien

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:14

Thanks for the tips, Florent and Se'bastien. I'll get to work on converting amgf.sage to amgf.py.

Alex

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:15

Hang on, can i somehow use the Sage pre-parser to automatically convert amgf.sage to amgf.py?

Alex

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:16

Aha! The command 'sage amgf.sage' pre-parses amgf.sage and saves it as amgf.py. Thank you ask.sagemath.org!

Alex

[3b09ef35-8a33-427f-b580-f9c14c161f29]

Pythonified amgf.sage

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:17

Attachment: amgf.py.gz

[seblabbe]

comment:18

The patch is empty. I hope you did not lose anything?

Sébastien

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:19

Attachment: trac_10519.patch.gz

Hmm, i couldn't figure out how to add the code to the patch. I tried 'sage -hg add amgf.py' as you suggested, Se'bastien, but that didn't work. So i finally copied and pasted it in. I hope that works. Let me know the proper way, yo.

Sorry for the delay.

Alex

[11d1fc49-71a1-44e1-869f-76be013245a0]

Attachment: trac_10519-fixed.patch.gz

Patch against 5.0

[11d1fc49-71a1-44e1-869f-76be013245a0]

Reviewer: David Loeffler

[11d1fc49-71a1-44e1-869f-76be013245a0]

comment:22

Apply trac_10519-fixed.patch

(for the patchbot).

Alex's last "patch" wasn't a patch at all, just a Mercurial header with some Python code arbitrarily pasted into it. I've just uploaded a non-broken version of the patch. I also added the necessary import statements -- a file that's supposed to be part of the Sage library can't import the Sage library (with from all_cmdline import *), since that would lead to an infinite loop, so the necessary imports have to be done one-by-one.

The patch now passes doctests on my machine, but there are some worrying errors that came up when I ran a syntax checker on the new file: {{{amgf.py:908: local variable 'R' is assigned to but never used amgf.py:909: local variable 'd' is assigned to but never used amgf.py:1070: local variable 'Ht' is assigned to but never used amgf.py:1251: local variable 'Ht' is assigned to but never used amgf.py:1591: undefined name 'n' amgf.py:2199: local variable 'H' is assigned to but never used amgf.py:2242: undefined name 'verdict' }}} The "undefined name" errors tend to suggest that some cases have bugs in them, and those cases are not checked in any doctest, which is a double whammy. Hence "needs work".

[11d1fc49-71a1-44e1-869f-76be013245a0]

comment:23

Sorry, that list got mangled. It should be

amgf.py:908: local variable 'R' is assigned to but never used
amgf.py:909: local variable 'd' is assigned to but never used
amgf.py:1070: local variable 'Ht' is assigned to but never used
amgf.py:1251: local variable 'Ht' is assigned to but never used
amgf.py:1591: undefined name 'n'
amgf.py:2199: local variable 'H' is assigned to but never used
amgf.py:2242: undefined name 'verdict'

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:24

Thanks for your review, David. This patch was my first on Sage Trac and i've since learned a lot, thanks to the constructive feedback of fellow Sagers such as yourself.

Since amgf.py was somewhat of a code bomb, i decided to split it up into simpler pieces, submit those pieces, and, once they get accepted, recombine them. If you or anyone else would like to help me in this endeavor, check out the latest amgf piece at https://github.com/sagemath/sage-prod/issues/12417.

Thanks!

I'll keep you all updated on the progress of amgf with these comments.

[3b09ef35-8a33-427f-b580-f9c14c161f29]

Attachment: trac_10519.2.patch.gz

Major rewrite.

[fchapoton]

comment:26

Apply only trac_10519.2.patch

[fchapoton]

Work Issues: coverage, trailing whitespace

[fchapoton]

comment:27

• Please ensure that every function is doctested

• Please remove all trailing whitespaces

[fchapoton]

comment:28

Attachment: trac_10519-asymptotic_multiple-v3.patch.gz

I have removed all whitespaces in the new patch. All test pass.

BUT there remains 6 procedures with no doctest. Remember that it is neccesary to have 100% doctesting.

Also the name amgf is not good. It should be changed to something more descriptive.

[fchapoton]

Description changed:

--- 
+++ 
@@ -9,3 +9,5 @@
 that `H` is a polynomial, and that asymptotics in the direction of
 `\alpha` (a tuple of positive integers) are controlled by smooth
 or multiple points.
+
+Apply [attachment: trac_10519-asymptotic_multiple-v3.patch](https://github.com/sagemath/sage-prod/files/10651709/trac_10519-asymptotic_multiple-v3.patch.gz)

[fchapoton]

Changed work issues from coverage, trailing whitespace to coverage

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:30

Thanks, folks. I'll document those outstanding 6 procedures within a week.

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:31

Hmm, i might take longer, because now that i've upgraded to Sage 5.3, i can't create a Sage clone :-( See http://ask.sagemath.org/question/1821/cant-create-sage-clone-again.

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:32

Still can't create a Sage clone. In the meantime, here's my patch as a Python file (which i tested on Sage 5.3 with sage -t --long) for anyone who wants to review it: https://www.dropbox.com/s/fcxsbuw119qclxt/trac_10519v4.py.

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:33

Attachment: trac_10519v4.patch.gz

[fchapoton]

comment:34

apply only trac_10519v4.patch

[fchapoton]

Description changed:

--- 
+++ 
@@ -10,4 +10,4 @@
 `\alpha` (a tuple of positive integers) are controlled by smooth
 or multiple points.

-Apply [attachment: trac_10519-asymptotic_multiple-v3.patch](https://github.com/sagemath/sage-prod/files/10651709/trac_10519-asymptotic_multiple-v3.patch.gz)
+Apply [attachment:trac_10519v4.patch ]

[fchapoton]

Changed work issues from coverage to failing tests, trailing whitespaces

[fchapoton]

comment:35

please remove all trailing whitespaces

please correct the failing doctests

[fchapoton]

comment:36

apply trac_10519-v5.patch

[fchapoton]

Description changed:

--- 
+++ 
@@ -10,4 +10,4 @@
 `\alpha` (a tuple of positive integers) are controlled by smooth
 or multiple points.

-Apply [attachment:trac_10519v4.patch ]
+Apply [attachment:trac_10519-v5.patch ]

[fchapoton]

comment:37

Attachment: trac_10519-v5.patch.gz

apply trac_10519-v5.patch

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:38

Attachment: trac_10519-v6.patch.gz

Hmm, i don't understand why i didn't get any test errors when i submitted v4. Anyway, fixed the errors.

[3b09ef35-8a33-427f-b580-f9c14c161f29]

Changed work issues from failing tests, trailing whitespaces to none

[fchapoton]

Description changed:

--- 
+++ 
@@ -10,4 +10,4 @@
 `\alpha` (a tuple of positive integers) are controlled by smooth
 or multiple points.

-Apply [attachment:trac_10519-v5.patch ]
+Apply [attachment:trac_10519-v6.patch ]

[fchapoton]

comment:40

Apply trac_10519-v6.patch

[fchapoton]

comment:41

well, now version V5 works and version V6 does not pass the tests.

I had already made all the necessary corrections in version V5.

Can we just forget version V6 ?

[3b09ef35-8a33-427f-b580-f9c14c161f29]

comment:42

OK, let's forget V6. I bet the problem is one or two docstring lines that output a dictionary. The random order of the keys is causing a test failure when it doesn't match up with the docstring key order. We can fix that by sorting the keys before printing them.

[fchapoton]

comment:43

apply trac_10519-v5.patch