Implement arithmetic product of cycle index series

[sagetrac-agd]

In a 2008 paper (see below), Maia and Méndez define an operation on combinatorial species which they dub the "arithmetic product". This operation corresponds to a nice combinatorial operation: the arithmetic product F⊡G corresponds to the species of "F-assemblies of cloned G-structures", which are structures of the partitional composite species F∘G with the additional requirement that all the G-structures be isomorphic.

As shown in the paper, the cycle index of the arithmetic product F⊡G can be computed in terms of the cycle indices of the species F and G. The attached patch adds code to calculate the result of this operation. It includes a doctest which verifies a nontrivial computation related to the species of "regular octopuses".

• Maia, Manuel and Méndez, Miguel. On the arithmetic product of combinatorial species. 1 February 2008. http://arxiv.org/abs/math/0503436

Apply:

• attachment: trac_14542_cycle_index_arithmetic_product.patch
• attachment: trac_14542-review-dg.patch​

CC: @sagetrac-sage-combinat

Component: combinatorics

Keywords: species, cycle index

Author: Andrew Gainer-Dewar

Reviewer: Darij Grinberg

Merged: sage-5.12.beta3

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

[darijgr]

comment:2

Huh -- can you explain the very first edit in the patch? This looks like an inadvertent page-down or a case of trying to use search but forgetting to press Ctrl+F...

[sagetrac-agd]

comment:3

Yikes! Sorry about that. I've replaced the patch with a version without that line.

Thanks for looking it over!

[darijgr]

comment:4

line 602: "g" should be inside double backticks.

line 620: Where does the "1 +" come from? Is there a regular octopus on 0 vertices? (I don't think so.)

line 623: xrange will be deprecated in Python 3.0; just saying.

line 643: Please explain what the arithmetic product of partition is supposed to mean. The Maia-Mendez paper defines the arithmetic product of polynomials; as far as I understand, by the arithmetic product of two partitions \lambda and \mu you mean the partition obtained by writing down gcd(\lambda_i, \mu_j) times the integer lcm(\lambda_i, \mu_j) for every i and every j, and then sorting the resulting sequence in nonincreasing order. That's fine, but you should clarify the relation between

• this operation on partitions,
• the Maia-Mendez operation on polynomials, and
• the operation on symmetric functions that you have actually implemented.
  (You are interpreting the x_i of Maia-Mendez as the i-th power sum, as far as I understand, and the partition \lambda corresponds to the product p_\lambda = p_{\lambda_1} p_{\lambda_2} ....)

About your implementation of arith_prod_of_partitions: I'm not sure how fast it is. It computes each lcm and each gcd very often (maxval times), whereas the definition I gave would only have to compute it once. On the other hand, the definition I gave requires sorting (or, better, insort) and probably handles repeated entries suboptimally.

line 685: I disagree with this, as on line 620.

line 690: Why do you put a p.zero() after res? (I'm new to the CycleIndexSeries class, so maybe it's important...)

On another note rather unrelated to Sage: The species viewpoint on the arithmetic product shows that the arithmetic product on the ring of symmetric functions (the one given by

p_\lambda \boxdot p_\mu = \prod_{i,j} p_{\lcm(\lambda_i, \mu_j)}^{\gcd(\lambda_i, \mu_j)}

for all partitions \lambda and \mu) is defined over the integers, not just over the rationals (despite the p_\lambda not forming a Z-basis of Symm). Is there an algebraic proof of this? It looks like a nice application of species to proving a nontrivial algebraic result. Is there a species-theoretical proof of the integrality of \Delta_3 in http://mathoverflow.net/questions/120924/is-the-renormalized-third-comultiplication-on-mathbfsymm-integral as well?

EDIT: This patch made me aware that I have no idea what the sum_generator method in sage/combinat/species/series.py does. Is it possible to have a docstring for it that actually explains its behavior? The examples given only confuse me.

[sagetrac-agd]

comment:5

Replying to @darijgr:

Thanks so much for your thorough review! I'm currently in the middle of writing and grading final exams, but I'll dig into your comments and the code next week.

[sagetrac-agd]

Patch to implement arithmetic product of cycle indices

[sagetrac-agd]

comment:6

Attachment: trac_14542_cycle_index_arithmetic_product.patch.gz

Replying to @darijgr:

Thanks again for your review! I've had a chance now to go back and look over your comments more closely.

line 602: "g" should be inside double backticks.

Fixed.

line 620: Where does the "1 +" come from? Is there a regular octopus on 0 vertices? (I don't think so.)

I've gotten into the rather unfortunate habit of ignoring the situation at n=0, but you're right about this one. I've updated the code to reflect this.

line 623: xrange will be deprecated in Python 3.0; just saying.

Since the numbers involved here are small, just using range instead won't hurt anything, so I've switched it over.

line 643: Please explain what the arithmetic product of partition is supposed to mean. The Maia-Mendez paper defines the arithmetic product of polynomials; as far as I understand, by the arithmetic product of two partitions \lambda and \mu you mean the partition obtained by writing down gcd(\lambda_i, \mu_j) times the integer lcm(\lambda_i, \mu_j) for every i and every j, and then sorting the resulting sequence in nonincreasing order. That's fine, but you should clarify the relation between

• this operation on partitions,
• the Maia-Mendez operation on polynomials, and
• the operation on symmetric functions that you have actually implemented.
  (You are interpreting the x_i of Maia-Mendez as the i-th power sum, as far as I understand, and the partition \lambda corresponds to the product p_\lambda = p_{\lambda_1} p_{\lambda_2} ....)

About your implementation of arith_prod_of_partitions: I'm not sure how fast it is. It computes each lcm and each gcd very often (maxval times), whereas the definition I gave would only have to compute it once. On the other hand, the definition I gave requires sorting (or, better, insort) and probably handles repeated entries suboptimally.

On review, I think your way is much better; if nothing else, it results in much cleaner and (to my eye) more mathematical code. I've added comments explaining what happens in each method and how they connect to the paper.

line 685: I disagree with this, as on line 620.

Also fixed.

line 690: Why do you put a p.zero() after res? (I'm new to the CycleIndexSeries class, so maybe it's important...)

The sum_generator method takes finite lists, but they need to represent infinite objects, so it assumes that the last element of the list is meant to repeat. Thus, the zero.

On another note rather unrelated to Sage: The species viewpoint on the arithmetic product shows that the arithmetic product on the ring of symmetric functions (the one given by

p_\lambda \boxdot p_\mu = \prod_{i,j} p_{\lcm(\lambda_i, \lambda_j)}^{\gcd(\lambda_i, \lambda_j)}

for all partitions \lambda and \mu) is defined over the integers, not just over the rationals (despite the p_\lambda not forming a Z-basis of Symm). Is there an algebraic proof of this? It looks like a nice application of species to proving a nontrivial algebraic result. Is there a species-theoretical proof of the integrality of \Delta_3 in http://mathoverflow.net/questions/120924/is-the-renormalized-third-comultiplication-on-mathbfsymm-integral as well?

I'll need to give this some more thought, but I wanted to go ahead and get the revised code up. I'm much more a combinatorialist than an algebraist, so I may not be the best person for the job, but I will definitely give it a thinking-over.

EDIT: This patch made me aware that I have no idea what the sum_generator method in sage/combinat/species/series.py does. Is it possible to have a docstring for it that actually explains its behavior? The examples given only confuse me.

This would definitely be helpful. My understanding of it is very limited; I seem to have managed to bend it to my will, but I'm quite far from mastery. The basic idea seems to be that it formally constructs a LazyPowerSeries by adding over an infinite iterable of "generators", but the details are buried deep in Stream and LazyPowerSeries.

[darijgr]

comment:7

Thanks for the update! I'll have to take a closer look at it.

Meanwhile, are you sure you want those asserts in lines 606-607? I don't know what a typical use case of this method should be, but I can imagine applying it in a case where the base ring is not discrete (i. e., you cannot check for equality) and these asserts will throw a NotImplemented. Think of a lazy power series ring in 2 variables implemented as a lazy power series ring over a lazy power series ring; then the base ring itself is lazy and this stuff breaks. I have a similar issue with your #14543 patch. Maybe add a check_input=True keyword which can optionally be set to False to avoid these asserts?

In other news:

sage: L = LazyPowerSeriesRing(QQ)                             
sage: L.sum_generator([L([0]),L([1])]).coefficients(10)       
[0, 1, 1, 1, 1, 1, 1, 1, 1, 1]
sage: L.sum_generator([L([0]),L([1]),L([0])]).coefficients(10)
[0, 1, 1, 1, 1, 1, 1, 1, 1, 1]

I don't think the 1 repeats because it's the last element -- it is not, in the second case. I think it repeats because of the "sum" in "sum_generator". But like you, I don't have a clue what's going on with sum_generator really (I actually think I have less of a clue than you do). I'll try to ask around in Orsay...

[sagetrac-agd]

Branch: u/agd/cis_arith_prod

[darijgr]

review patch

[darijgr]

comment:9

Attachment: trac_14542-review-dg.patch.gz

I've attached a review patch. Sorry for it taking that long; I was waiting for a thorough explanation of sum_generator, which now has been added to the otherwise unrelated ticket #13433. The explanation shows that it is not necessary to pad a sequence with zeroes to use sum_generators on it; sum_generator([a,b]) does precisely the same as sum_generator([a,b,0]). But this isn't really relevant to this ticket here. What is correct that L([a,b]) is not the same as L([a,b,0]), and the latter syntax is the correct one in most realistic cases. This has nothing to do with sum_generator. Sorry for the confusion.

I have added some more info to the docstring and the comments, replaced the arXiv link by the proper syntax for arXiv identifiers (hopefully correctly), improved PEP 8 compliance (it's def f(n), not def f (n)), and replaced n/d for n//d (both do the same when d divides n, but the latter is faster). If you agree with these all, please mark this ticket as positive_review.

For the patchbot:

Apply trac_14542_cycle_index_arithmetic_product.patch​ trac_14542-review-dg.patch​

[sagetrac-agd]

comment:10

Replying to @darijgr:

Thanks for the review!

I've attached a review patch. Sorry for it taking that long; I was waiting for a thorough explanation of sum_generator, which now has been added to the otherwise unrelated ticket #13433. The explanation shows that it is not necessary to pad a sequence with zeroes to use sum_generators on it; sum_generator([a,b]) does precisely the same as sum_generator([a,b,0]). But this isn't really relevant to this ticket here. What is correct that L([a,b]) is not the same as L([a,b,0]), and the latter syntax is the correct one in most realistic cases. This has nothing to do with sum_generator. Sorry for the confusion.

Yes, the semantics of sum_generator and LazyPowerSeriesRing are quite mysterious. I did some investigating while working on #14906 and learned a bit about working with them. I now think the correct way to handle this kind of construction is with LazyPowerSeriesRing.term, but unfortunately I am away from my build machine for the rest of the week (I'm at MathFest), so I won't be able to update this until I get back.

Alternately, if you'd be willing to put together a new patch, the correction is simple. The lines

    # Build a list which has res in the nth slot and 0's before and after
    # to feed to sum_generator
    res_in_seq = [p.zero()]*n + [res, p.zero()]
            
    return self.parent(res_in_seq)

should be replaced by the single line

    return self.parent.term(res, n)    

I have added some more info to the docstring and the comments, replaced the arXiv link by the proper syntax for arXiv identifiers (hopefully correctly), improved PEP 8 compliance (it's def f(n), not def f (n)), and replaced n/d for n//d (both do the same when d divides n, but the latter is faster). If you agree with these all, please mark this ticket as positive_review.

These changes all look great. Thanks for the help! I'm not sure whether it's appropriate to mark the patch as positive_review with the issue above still outstanding, so for now I'll just say it:

Positive review!

[darijgr]

comment:11

You want self.parent().term, not self.parent.term. The fact that self.parent(res_in_seq) works (and is equivalent to self.parent()(res_in_seq)) is due to syntactic sugar.

You are correct that self.parent.term(res, n) works and makes more sense, but it's also 3 times slower than self.parent(res_in_seq) (not sure why), so I'd keep what I've written.

[sagetrac-agd]

comment:14

Replying to @darijgr:

You want self.parent().term, not self.parent.term. The fact that self.parent(res_in_seq) works (and is equivalent to self.parent()(res_in_seq)) is due to syntactic sugar.

Whoops! This is what I get from trying to write code from a browser after a day of workshops.

You are correct that self.parent.term(res, n) works and makes more sense, but it's also 3 times slower than self.parent(res_in_seq) (not sure why), so I'd keep what I've written.

How strange! I'll add this to my pile of things to investigate once I get back on my feet. Using res_in_seq for this is fine for now, though.

Once again, thanks for the feedback and the review!

[jdemeyer]

comment:15

Please fill in the real names of the Authors and Reviewers.

[darijgr]

Author: Andrew Gainer-Dewar

[darijgr]

Reviewer: Darij Grinberg

[jdemeyer]

comment:17

I am assuming the patches should be merged and the git branch is just there for reference.

[jdemeyer]

Merged: sage-5.12.beta3

[sagetrac-agd]

Changed branch from u/agd/cis_arith_prod to none