Opened 6 years ago

Closed 5 years ago

#18824 closed enhancement (fixed)

Add MultiPolynomialIdeal.groebner_basis("libgiac")

Reported by: malb Owned by:
Priority: major Milestone: sage-6.8
Component: interfaces: optional Keywords:
Cc: frederichan Merged in:
Authors: Martin Albrecht, Frederic Han Reviewers: Martin Albrecht, Frederic Han
Report Upstream: N/A Work issues:
Branch: 02cf665 (Commits) Commit: 02cf665c79729f856ae8d43408caa1f6d1506973
Dependencies: #18841 Stopgaps:

Description (last modified by frederichan)

Giac is an optional package and supports fast computation of Gröbner bases. Let's make that easily accessible.

This trac depends on #18841 in the following points: 1) a sig_count error detection in the doctest. 2) the output of the Katsura3 doctest in multi_polynomial_ideal.py was generated with giac-1.2.0-19 as provided in #18841. (Earlier versions provides a reversed output)

Change History (50)

comment:1 Changed 6 years ago by frederichan

  • Cc frederichan added

Thank you. NB: Currently giacpy is just an optional spkg not a sage library so there is no sage.lib.giac

I don't know if there is a suitable place on sage git to expand the so called upstream spkg, I would be happy to put it in public so that sage developper could also modify it.

comment:2 Changed 6 years ago by frederichan

  • Branch set to public/giacpyGB

comment:3 Changed 6 years ago by malb

giacpy being optional is no problem. I'd suggest to add sage.libs.libgiac which imports the right stuff from giacpy and raises an error if that fails, e.g.

def groebner_basis_libgiac(gens, epsilon=0, prot=False, <other args>):
   from giacpy import libgiac
   F = libgiac(gens)
   giacsettings.proba_epsilon=epsilon
   # get ring from gens, check term order etc.
   B = F.gbasis([P.gens()],'revlex') 
   # convert result to PolynomialSequence and return

comment:4 Changed 6 years ago by frederichan

  • Branch public/giacpyGB deleted

Sorry I am struggling with git (I am failing to create a branch), could you please start a branch with your function so that I can pull?

NB: you need to import giacsettings also.

comment:5 Changed 6 years ago by frederichan

  • Branch set to public/giacpyGB
  • Commit set to 6c54aacce07d0885c49b4942a46a1bf47d2c86e9

ouf! this time I have succeed to push a trivial thing. I can start to work


New commits:

6c54aacstart groebner basis branch with libgiac

comment:6 Changed 6 years ago by git

  • Commit changed from 6c54aacce07d0885c49b4942a46a1bf47d2c86e9 to 0b695a051582653dfd9d67fa3f47a43db2786a17

Branch pushed to git repo; I updated commit sha1. New commits:

9e02904Merge branch 'public/giacpyGB' of trac.sagemath.org:sage into public/giacpy-groebner
0988f8drename sage.libs.libgiac to sage.libs.giac
0b695a0poc of groebner_basis_libgiac

comment:7 Changed 6 years ago by malb

I added a proof of concept implementation

comment:8 Changed 6 years ago by git

  • Commit changed from 0b695a051582653dfd9d67fa3f47a43db2786a17 to fab096e4658fb1300ec56943f9601b7cf9c915ed

Branch pushed to git repo; I updated commit sha1. New commits:

f57602fMerge branch 'develop' of trac.sagemath.org:sage into t/public/giacpyGB
2fa05baMerge branch 'public/giacpyGB' of trac.sagemath.org:sage into t/public/giacpyGB
fab096etest import, gens, rings. change conversion to sage

comment:9 Changed 6 years ago by frederichan

I have noticed that with giac1.2.0-13 (spkg version) it may happen that the giac answer contains a 0 term, but it is fixed in giac1.2.0-15. May be updating the spkg is the best solution?

comment:10 Changed 6 years ago by git

  • Commit changed from fab096e4658fb1300ec56943f9601b7cf9c915ed to 86bd0a445cee318d3cca0fef87babf4118c910d8

Branch pushed to git repo; I updated commit sha1. New commits:

55fff07drop unused code and imports
5f95d87don't check for list type explicitly
b8725bcobey PEP8 style conventions
86bd0a4use decorators for preserving global settings

comment:11 Changed 6 years ago by malb

I cleaned it up a bit more. It still needs tests, documentation and a way to shut Giac up. Also, does Giac have an option to print a protocol of the GB computation (most systems do).

Finally, we get this doctest failure:

File "src/sage/libs/giac.py", line 19, in sage.libs.giac
Failed example:
    sig_on_count()
Expected:
    0
Got:
    1

which suggests there's a bug somewhere.

comment:12 Changed 6 years ago by malb

Updating Giac should be a separate ticket.

comment:13 Changed 6 years ago by frederichan

Thank you. I am still working on this. I found a way for the info level.

comment:14 follow-up: Changed 6 years ago by frederichan

The sig_count pb is from giacpy 0.4.8 Pygen __iter__, it could/should? be fixed in a separate trac? we can avoid it with: return PolynomialSequence(P, tuple(gb_giac), immutable=True) instead of just gb_giac. I don't remark time loss with tuple, should we adopt this?

is there a way to substitute efficiently variables in PolynomialRing? (ex: 'i' will give use troubles in giac) if not I will just add some protection with an error message

Last edited 6 years ago by frederichan (previous) (diff)

comment:15 Changed 6 years ago by frederichan

sorry I wanted to say: I don't remark time loss with tuple. (so I have edited my comment:14)

comment:16 in reply to: ↑ 14 Changed 6 years ago by malb

Replying to frederichan:

The sig_count pb is from giacpy 0.4.8 Pygen __iter__, it could/should? be fixed in a separate trac?

yes, this should be a separate ticket but this ticket IMHO should depend on that ticket. That is, this ticket should wait until the fix is in.

we can avoid it with: return PolynomialSequence(P, tuple(gb_giac), immutable=True) instead of just gb_giac. I remark time loss with tuple, should we adopt this?

You mean instead of PolynomialSequence(gb_giac, P, immutable=True)? I don't think we should work around a bug which should be fixed instead.

is there a way to substitute efficiently variables in PolynomialRing? (ex: 'i' will give use troubles in giac) if not I will just add some protection with an error message

You can try but it's not very fast:

sage: P.<a,b,c,i> = QQ[]
sage: R.<x,y,z,w> = QQ[]
sage: f = P.random_element()
sage: R(f)
-2*y^2 - 1/4*y*z + 6*z*w - 2*z - 59/5
sage: %timeit R(f)
10000 loops, best of 3: 77.3 µs per loop

Maybe check P.gens_dict().keys() for problematic names and only switch over if there is a problem?

comment:17 Changed 6 years ago by frederichan

  • Dependencies set to 18841

comment:18 Changed 6 years ago by frederichan

  • Dependencies changed from 18841 to #18841

comment:19 Changed 6 years ago by git

  • Commit changed from 86bd0a445cee318d3cca0fef87babf4118c910d8 to 2f653024a1db32c326b87f8dca48503b87007acb

Branch pushed to git repo; I updated commit sha1. New commits:

6d810f6add prot and restore global proba
f725910Merge branch 'public/giacpyGB' of trac.sagemath.org:sage into t/public/giacpyGB
cc81303add prot feature + fixes in doctests
2f65302add and change doctests. Add a threads option

comment:20 Changed 6 years ago by frederichan

giac gbasis have also a rur option to obtain a Rational Univariate Representation of a 0 dimensional ideal. Are you also interested by a more convient interface to this or should we leave it for later?

Last edited 6 years ago by frederichan (previous) (diff)

comment:21 Changed 6 years ago by git

  • Commit changed from 2f653024a1db32c326b87f8dca48503b87007acb to 425b8ddc68ca59c0f24e39995aebcdff0d9a263c

Branch pushed to git repo; I updated commit sha1. New commits:

425b8ddadd doctests

comment:22 Changed 6 years ago by malb

I'm a one-feature-per-ticket kind of guy, so I'd put that rur enhancement in a separate ticket.

comment:23 Changed 6 years ago by git

  • Commit changed from 425b8ddc68ca59c0f24e39995aebcdff0d9a263c to 28e2db98e53e67ce9b45544d59283a133f109c41

Branch pushed to git repo; I updated commit sha1. New commits:

756d223wrap long lines as recommended by style guide
7076dd0follow PEP8 style guide
28e2db9make error message more concise

comment:24 Changed 6 years ago by git

  • Commit changed from 28e2db98e53e67ce9b45544d59283a133f109c41 to 061f2779431eceb921bedaecf6810c3a04808f04

Branch pushed to git repo; I updated commit sha1. New commits:

061f277Merge branch 'develop' of trac.sagemath.org:sage into giacpyGB

comment:25 Changed 6 years ago by git

  • Commit changed from 061f2779431eceb921bedaecf6810c3a04808f04 to 8bf1320413e63fe2f71d4164dea93ae3d27d73bd

Branch pushed to git repo; I updated commit sha1. New commits:

81b68faadd groebner_basis("giac") option
8bf1320we have namespaces, use them

comment:26 Changed 6 years ago by frederichan

Dear martin, thank you for the modifs. About the doctests I don't trust even the first 2 characters: // of the message // Giac ... I have not yet tested this giac spkg 1.2.0-19 on os X to see but last time I got pbs with the begining of this message on OSX versus linux. so I think that in a doctest it is safer to do it in 2 steps: sage: gb = J.groebner_basis('giac') # optional - giacpy, random sage: gb # optional - giacpy

Last edited 6 years ago by frederichan (previous) (diff)

comment:27 Changed 6 years ago by git

  • Commit changed from 8bf1320413e63fe2f71d4164dea93ae3d27d73bd to 5d23482104cf6a10406cc23f81874b93973650b0

Branch pushed to git repo; I updated commit sha1. New commits:

5d23482fix-up doctest

comment:28 Changed 6 years ago by malb

I guess this ticket is done then (?)

comment:29 Changed 6 years ago by git

  • Commit changed from 5d23482104cf6a10406cc23f81874b93973650b0 to b590a50900bccd1e9d5fe6d62c5636d946175332

Branch pushed to git repo; I updated commit sha1. New commits:

b590a50blacklist 'e' for variable names

comment:30 Changed 6 years ago by frederichan

may be I should also add some doc in multi_polynomial_ideal.py or pointer to the sage.libs.giac.groebner_basis doc

comment:31 Changed 6 years ago by git

  • Commit changed from b590a50900bccd1e9d5fe6d62c5636d946175332 to 5603124d85096ba10d0f8dd525105b3b697961af

Branch pushed to git repo; I updated commit sha1. New commits:

5603124add doc in sage.rings.polynomial.multi_polynomial_ideal.groebner_basis

comment:32 Changed 6 years ago by frederichan

  • Description modified (diff)

comment:33 Changed 6 years ago by malb

Perhaps

sage: I9.groebner_basis("giac",proba_epsilon=1e-7) # optional - giacpy

should be tagged long?

comment:34 Changed 6 years ago by frederichan

On a 2012 (2core i5-2435M CPU @ 2.40GHz) notebook with 4Go of RAM:

sage: A9=PolynomialRing(QQ,9,'x') # optional - giacpy
sage: I9=sage.rings.ideal.Katsura(A9) # optional - giacpy
sage: time I9.groebner_basis("giac",proba_epsilon=1e-7) # optional - giacpy
// Giac share root-directory:/home/fred/dev/sage/develop/sage.run/local/share/giac/
// Unable to find keyword file /home/fred/dev/sage/develop/sage.run/local/share/giac/doc/fr/keywords
// Giac share root-directory:/home/fred/dev/sage/develop/sage.run/local/share/giac/
Help file /home/fred/dev/sage/develop/sage.run/local/share/giac/doc/fr/aide_cas not found
Added 0 synonyms
Running a probabilistic check for the reconstructed Groebner basis. If successfull, error probability is less than 1e-07 and is estimated to be less than 10^-133. Use proba_epsilon:=0 to certify (this takes more time).
CPU times: user 3.29 s, sys: 19.9 ms, total: 3.31 s
Wall time: 3.23 s
Polynomial Sequence with 143 Polynomials in 9 Variables

comment:35 Changed 6 years ago by malb

Three seconds is long by sage doctest standards. We're running thousands of doctests.

comment:36 Changed 6 years ago by frederichan

OK so some others in giac.py need it also. I see that sometime there is an indication of the time. would # optional giacpy, long time (3s on a 2012 notebook)

be OK or is it too long.

comment:37 Changed 6 years ago by malb

It's fine but not mandatory to add this (I'd say)

comment:38 Changed 6 years ago by git

  • Commit changed from 5603124d85096ba10d0f8dd525105b3b697961af to 36620a9e988a7453bf8655d950ae842924e9f8bf

Branch pushed to git repo; I updated commit sha1. New commits:

36620a9add long time tags to some doctests

comment:39 Changed 6 years ago by frederichan

to change I have tried it on OSX(10.7) with #18841 and got no problems.

comment:40 Changed 6 years ago by frederichan

  • Authors set to Martin Albrecht, Frederic Han

comment:41 Changed 6 years ago by malb

  • Status changed from new to needs_review

comment:42 Changed 5 years ago by git

  • Commit changed from 36620a9e988a7453bf8655d950ae842924e9f8bf to f3a235b3265710b7f35b9e86d96f2ac84d5ff6b6

Branch pushed to git repo; I updated commit sha1. New commits:

24a7254Merge branch 'develop' of trac.sagemath.org:sage into t/public/giacpyGB
f3a235bupdate docstring in cachefunc.pyx (groebner_basis call)

comment:43 Changed 5 years ago by frederichan

  • Status changed from needs_review to needs_work

The zero ideal gives a segfault

comment:44 Changed 5 years ago by git

  • Commit changed from f3a235b3265710b7f35b9e86d96f2ac84d5ff6b6 to 02cf665c79729f856ae8d43408caa1f6d1506973

Branch pushed to git repo; I updated commit sha1. New commits:

02cf665add an is_zero ideal test before calling giac gbasis

comment:45 Changed 5 years ago by frederichan

  • Status changed from needs_work to needs_review

comment:46 Changed 5 years ago by malb

I think this ticket is done, i.e. I'm positively reviewing your contributions. If you review mine and agree with them, we can set this to positive review.

comment:47 Changed 5 years ago by frederichan

OK, I thought we had to wait for a third person. I will look if I can reproduce this doctest error: http://patchbot.sagemath.org/log/18824/Ubuntu/14.04/x86_64/3.13.0-61-generic/librae/2015-08-17%2014:47:28%20+0000?short do you think it could be related to this patch?

comment:48 Changed 5 years ago by frederichan

I am positively reviewing your changes.

I was not able to reproduce the previous patchbot doctest error with sage 6.9.beta2.

comment:49 Changed 5 years ago by malb

  • Reviewers set to Martin Albrecht, Frederic Han
  • Status changed from needs_review to positive_review

comment:50 Changed 5 years ago by vbraun

  • Branch changed from public/giacpyGB to 02cf665c79729f856ae8d43408caa1f6d1506973
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.