Opened 5 years ago
Closed 5 years ago
#15382 closed enhancement (fixed)
Implementing Macaulay Resultant (sage days 55)
Reported by:  haochen_uw  Owned by:  

Priority:  major  Milestone:  sage6.3 
Component:  algebraic geometry  Keywords:  Macaulay Resultant sagedays55 
Cc:  Merged in:  
Authors:  Soli Vishkautsan, Hao Chen  Reviewers:  Ben Hutz 
Report Upstream:  N/A  Work issues:  
Branch:  db145e3 (Commits)  Commit:  db145e34be5d524b10f1504e5892e256dd921d05 
Dependencies:  Stopgaps: 
Description (last modified by )
Trying to Implement Macaulay Resultant for universal polynomials as well as polynomials with coefficients in a field k. Newest version of patch is .3.patch
Attachments (3)
Change History (34)
Changed 5 years ago by
comment:1 Changed 5 years ago by
 Keywords sagedays55 added
comment:2 Changed 5 years ago by
 Description modified (diff)
comment:3 Changed 5 years ago by
 Status changed from new to needs_review
Changed 5 years ago by
comment:4 Changed 5 years ago by
 Description modified (diff)
comment:5 Changed 5 years ago by
http://trac.sagemath.org/attachment/ticket/15382/trac15382macaulayresultant.3.patch changed function names to lower case and put names into global name space.
comment:6 Changed 5 years ago by
 Status changed from needs_review to needs_work
I haven't had a chance to look at the algorithm yet, but I browsed through the code in the patch and there a couple things.
1) Instead of adding .n.patch for each new version, unless the new version is drastically different simply replace the previous patch with the new one.
2) add a comment "apply patch.name" to specify which patch should be applied
Actually in the code:
3)Replace the 'assert' which are raising AssertionError?, with a 'raise' statement and something more meaning full like TypeError? or ValueError?
4) The monomials() function takes as input the ring it is creating the monomials over. This seems like it should be a member function of the appropriate polynomial ring class instead of a standalone function (assuming it isn't already implemented).
5) I'm still debating whether I think macalauy_resultant should also be a member function of some polynomial ring. I'm inclined to say it should be (since resultant() is) but I'm open to arguments against.
6) It appears you can pass in a bunch of polynomials from all different polynomial rings and it won't detect the problem (you just get the parent of the first flist element).
That's all for now. I still need actually look at the algorithm...
comment:7 Changed 5 years ago by
 Branch set to u/wishcow/ticket/15382
 Created changed from 11/08/13 19:46:30 to 11/08/13 19:46:30
 Modified changed from 11/15/13 15:30:12 to 11/15/13 15:30:12
comment:8 Changed 5 years ago by
 Commit set to 32fa6c14411947b49bf32c3b6be176575bade538
Added a branch including all the code in trac15382macaulayresultant.3.patch. I am unable to delete extraneous patches, as we do not have permission to do so.
New commits:
32fa6c1  trac 15382 
comment:9 Changed 5 years ago by
ok. I looked at the algorithm today and it looks fine. I also ran some tests comparing against the primes_of_bad_reduction
function for projective morphisms and the results match (well, at least the primes match).
In additional to my previous comments, I have two suggestions to make this faster and one other comment
1) It seems like you could combine these two loops:
 for mon in mon_d:
 for f in newflist:
2) It also seems like you could always just work with the vector of exponents instead of converting back and forth to monomials with mon_d; i.e. don't use mon_d at all, just use degs.
3) You have a monomials() function, yet you don't use it in the main macaulay_resultant function. Although with (2) you don't need it anyway. You may however want to use some kind of iterator instead of making lists of monomials since you don't need to keep them around after you've used it. I'm sure WeightedIntegerVectors
has an iterator of some kind.
It would be nice if there were a few comments in the code for macaulay_resultant
so that it is easier to follow.
comment:10 Changed 5 years ago by
Another issue is that macaulay_resultant() returns an element of the fraction field of the polynomial ring instead of a polynomial.
comment:11 Changed 5 years ago by
 Commit changed from 32fa6c14411947b49bf32c3b6be176575bade538 to 94bcf0ec174e1ba5f697d4bf2e446d8b57874a32
comment:12 Changed 5 years ago by
I merged the for loops as much as possible. This produced about a 25% improvement in performance, not much considering the work. The main problem is that for now we can't avoid working with actual monomials yet (I got rid of the bidirectional translation between vectors and monomials as much as I possibly could). The ring of polynomials is not equipped to work with vectors instead of monomials, for example in retrieving the coefficients.
The doctests fail on one test in a very strange way. Please let me know if you have any idea how to fix it.
comment:13 Changed 5 years ago by
The doctest failure is simply a typo. You have "equal to number" instead of "equal number".
Actually, polynomial rings *are* equipped to work with exponent vectors in place of monomials. You can get a coefficient with
R.<x,y,z>=PolynomialRing(QQ) f=2*x^3 + 3*x*y*z + 4*y^3 + 5*y^2+z f[[0,3,0]]
You should also look back up at comment 6 again as issues (3)(6) still need to be addressed.
comment:14 Changed 5 years ago by
 Thanks! You don't know how much time I spent trying to figure out why the test failed :P
 thanks for the pointer about the exponent vectors!
 I am working on the other comments, don't worry :)
comment:15 Changed 5 years ago by
 Commit changed from 94bcf0ec174e1ba5f697d4bf2e446d8b57874a32 to a9da35924d555a5d79e0d5d36587adebcdeffff6
comment:16 Changed 5 years ago by
OK, I made the change to exponentvectors for extracting the coefficients of the polynomial. The code is much cleaner now, but sadly, it actually works a little slower now. I cannot really explain why this is happening, since doing f[[0,3,0]]
is faster than doing f.monomial_coefficient(y^3)
. In any case, if you have more performanceimproving suggestions, let me know.
comment:17 Changed 5 years ago by
 Milestone changed from sage6.1 to sage6.2
comment:18 Changed 5 years ago by
 Milestone changed from sage6.2 to sage6.3
comment:19 Changed 5 years ago by
 Commit changed from a9da35924d555a5d79e0d5d36587adebcdeffff6 to e6e898513b51b865e7574f5534c0c9bc0051addc
Branch pushed to git repo; I updated commit sha1. New commits:
a81f728  trac 15382

e08b29c  improve performance of macaulay_resultant

95639b7  some fixes to performance issues

d258bf8  minor changes

ca92141  finished implementing minor performance improvements. There is still a failed test to take care of.

3ece4d9  removed all usage of actual monomials, only using degreelists now.

c99a01a  minor changes to comments

b206fb2  moved macaulay_resultant to MPolynomialRing_generic class; deleted macaulay_resultant.py

e6e8985  Merged all.py to fix conflict

comment:20 Changed 5 years ago by
 Commit changed from e6e898513b51b865e7574f5534c0c9bc0051addc to c27d8026caf3bfcac9e7a442c841db6dbf432c40
Branch pushed to git repo; I updated commit sha1. New commits:
c27d802  deleted macaulay_resultant.py

comment:21 Changed 5 years ago by
 Commit changed from c27d8026caf3bfcac9e7a442c841db6dbf432c40 to 2a15c39fe2d10a6eace8d5250291dabae008822c
Branch pushed to git repo; I updated commit sha1. New commits:
2a15c39  fixed failed doctest

comment:22 Changed 5 years ago by
 Commit changed from 2a15c39fe2d10a6eace8d5250291dabae008822c to 8baa0984e422ce04713bda0f5be32d2bdf7d12d6
comment:23 Changed 5 years ago by
 Status changed from needs_work to needs_review
comment:24 Changed 5 years ago by
 Branch changed from u/wishcow/ticket/15382 to u/bhutz/ticket/15382
 Modified changed from 07/09/14 18:46:58 to 07/09/14 18:46:58
comment:25 Changed 5 years ago by
 Commit changed from 8baa0984e422ce04713bda0f5be32d2bdf7d12d6 to 97ee6f284baa43bc348c42a8b0a4ca254024e5f5
 Reviewers set to Ben Hutz
 Status changed from needs_review to needs_work
The docbuild issue was from making duplicate references. You only need to declare the references once. I fixed that and some minor formatting/doc issues.
As for testing: The functionality is working successfully, but there are two issues to address:
1) the assert
statements need to be switched to raiseerror
so that is raises the appropriate error on failure (not just asserterror)
2) I think the speed can be greatly improved on average by changing the way the matrix is built. Currently, all possible monomials are generated and then each entry in the matrix is filled in by checking the coefficient of that monomial in the polynomials. Since I would expect 'most' polynomials to be sparse, it should be more efficient to start with a matrix of all 0s and fill in only those entries corresponding to appropriate monomials in the polynomials.
New commits:
97ee6f2  15382: fixed some minor doc issues

comment:26 Changed 5 years ago by
 Branch changed from u/bhutz/ticket/15382 to u/wishcow/ticket/15382
comment:27 Changed 5 years ago by
 Commit changed from 97ee6f284baa43bc348c42a8b0a4ca254024e5f5 to f1ec984c85ef10089392f32385745f629a4fb29f
 Status changed from needs_work to needs_review
New commits:
f1ec984  replaced assert with raise. cleaned up for loops. added reverse lookup of monomials

comment:28 Changed 5 years ago by
 Branch changed from u/wishcow/ticket/15382 to u/bhutz/ticket/15382
 Modified changed from 07/11/14 09:33:59 to 07/11/14 09:33:59
comment:29 Changed 5 years ago by
 Commit changed from f1ec984c85ef10089392f32385745f629a4fb29f to db145e34be5d524b10f1504e5892e256dd921d05
 Status changed from needs_review to positive_review
This looks good to me. I just removed some trailing whitespaces and reformatted the author note.
New commits:
db145e3  removed trailing whitespaces

comment:30 Changed 5 years ago by
 Type changed from task to enhancement
comment:31 Changed 5 years ago by
 Branch changed from u/bhutz/ticket/15382 to db145e34be5d524b10f1504e5892e256dd921d05
 Resolution set to fixed
 Status changed from positive_review to closed
patch for macaulay resultant