[sage-devel] Re: computing the number of partitions of an integer (fixes in latest pari)
This problem seems to be fixed in pari-2.4.2 (not yet realead!, code in CVS) ?numbpart((5*10535+4)) time = 16 ms. %2 = 13277569485341661441070935908285030462857350709714871167223602317834599507871542657403046634712636713000828040255875556477097762416076415269139010200121379629789951459033537585723875697364877024644629549129563676618917774238138926865673007168971398575 The CHANGES file says: Fixed 12- numbpart inaccurate (and slow). E.g numbpart(52602) off by 1 [#557]. This version of pari also includes some new functions that could be useful to Sage. For example, it includes the function polhermite and pollegendre, for computing the Hermite and Lengendre orthogonal polynomials, that currently sage computes using maxima and I think that using libpari could be much faster. Pablo On 7/26/07, Bill Hart [EMAIL PROTECTED] wrote: Pari implements it, but incorrectly: sage: number_of_partitions(5*10535+4,algorithm=pari) 132775694853416614410709359082850304628573507097 148711672236023178345995078715426574030466347126 367130008280402558755564770977624160764152691390 102001213796297899514590335375857238756973648770 246446295491295636766189177742381389268656730071 68971398574 But: sage: number_of_partitions(5*10535+4) 132775694853416614410709359082850304628573507097 148711672236023178345995078715426574030466347126 367130008280402558755564770977624160764152691390 102001213796297899514590335375857238756973648770 246446295491295636766189177742381389268656730071 68971398575 --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~--~~~~--~~--~--~---
[sage-devel] Re: computing the number of partitions of an integer
I've tested it but seems to be buggy:: it works up to 156 ./part 156 p(156)= 73232243759 but for 157 gives a floating point exception error (and a gdb tracing says it is in the line 176 of the source code r2 = r0 % r1; in function g I've compiled it using gcc part.c -g -lmpfr -lm -DTEST_CODE -o part and my mpfr library is 2.2.1_p5 and gmp is 4.2.1 The code theres seems indeed to be GPL-ed, it has a copyright notice that says: (C) 2002 Ralf Stephan ([EMAIL PROTECTED]). * See part.pdf on the same path for a summary of the math. * Distributed under GPL (see gnu.org). Version 2002-12. The pdf is: http://www.ark.in-berlin.de/part.pdf Pablo On /26/07, Bill Hart [EMAIL PROTECTED] wrote: Alternatively you could just use the implementation here: http://www.ark.in-berlin.de/part.c which seems to only rely on mpfr and GMP. However, since the Pari version was based on it, I suspect it may also need correction. He does seem to adjust the precision as he goes, but I've no idea how far above or below the required precision this is. However there is no need to use a precision of 55 from a certain point on. That probably forces it to use two doubles instead of one in the floating point computations, which I think would be unnecesary. You can look on mathworld for how well the Rademacher series converges. It would probably be quicker to do a fresh implementation from scratch given the application in mind. The code there may not be GPL'd also. To check it works, one should at least look at the congruences modulo 5, 7, 11 and 13. The author of the mpfr version does seem to check them mod 11 but for very small n and only for a very few values of n. If the gcds prove to be a bottleneck, I have some code that is roughly twice as fast as GMP for computing these. Certainly the code in the mpfr version is suboptimal and needs some serious optimisation when the terms in the dedekind sum fit into a single limb, which on a 64 bit machine is probably always, when n is of a reasonable size. Another suggestion is that for sufficiently small values of h and or k, it may be quicker to compute the dedekind sum directly rather than use the gcd formula. A lookup table would also be useful for the tail end of the gcd/ dedekind sum computation. I'd be shocked if the computation for n = 10^9 couldn't be done in under 10s on sage.math. Bill. --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~--~~~~--~~--~--~---
[sage-devel] 2.7.1 and gcc-2.4.1 linbox problem
SAGE-2.7.1 sadly does NOT build if you use the latest stable release of the gcc compiler (gcc-4.2.1). It fails in linbox for probably the same reason as was given in the thread: 2.5.3 and gcc-2.4.0 linbox problem SAGE-2.7.1 does build and pass its test suite on my x86-Linux (pentium4) if built with gcc-4.1.2. -- Kate Minola University of Maryland, College Park --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~--~~~~--~~--~--~---
[sage-devel] Re: computing the number of partitions of an integer
On 7/27/07, Pablo De Napoli [EMAIL PROTECTED] wrote: I've tested it but seems to be buggy:: it works up to 156 ./part 156 p(156)= 73232243759 but for 157 gives a floating point exception error (and a gdb tracing says it is in the line 176 of the source code r2 = r0 % r1; in function g I've attached a slightly modified version of part.c that seems to work (it doesn't crash and gives correct answers in all the cases I tested). I just changed the code in g slighlty so if r1 is 0, then r2 is also set to 0. I compile it using gcc part.c -O3 -g -lmpfr -lm -o part -I/home/was/sage/local/include -L/home/was/sage/local/lib/ -lgmp where /home/was/sage is the path to my SAGE install. Interestingly, when I time part.c it is not much better than CVS pari (see below). Bill Hart wrote: This, with the other ideas I gave above, will definitely let us beat Mathematica convincingly, as a simple back of the envelope calculation shows. It will have the added benefit of allowing us to beat Magma at something. Magma is already terrible at computing NumberOfPartitions, though not as bad as GAP (which is even worse): [EMAIL PROTECTED]:~$ magma Magma V2.13-5 Fri Jul 27 2007 10:54:07 on ubuntu [Seed = 1720324423] Type ? for help. Type Ctrl-D to quit. time k := NumberOfPartitions(10^5); Time: 3.070 With the latest CVS pari: ? gettime; n=numbpart(10^5); gettime/1000.0 %1 = 0.01200 ? gettime; n=numbpart(10^6); gettime/1000.0 %2 = 0.2080 ? gettime; n=numbpart(10^7); gettime/1000.0 %3 = 6.973 ? gettime; n=numbpart(10^8); gettime/1000.0 With part.c: I use time ./part 10 out and record system + user time: 10^60.04 seconds 10^74.584 seconds 10^8 47.483 seconds So, to clarify or summarize, Bill is there one thing that we could change in part.c that would speed it up? I realized you sort of answered this before, but the sequence of previous emails yesterday was rather long and may have got confused. Thanks, william The code theres seems indeed to be GPL-ed, it has a copyright notice that says: (C) 2002 Ralf Stephan ([EMAIL PROTECTED]). * See part.pdf on the same path for a summary of the math. * Distributed under GPL (see gnu.org). Version 2002-12. The pdf is: http://www.ark.in-berlin.de/part.pdf Pablo On /26/07, Bill Hart [EMAIL PROTECTED] wrote: Alternatively you could just use the implementation here: http://www.ark.in-berlin.de/part.c which seems to only rely on mpfr and GMP. However, since the Pari version was based on it, I suspect it may also need correction. He does seem to adjust the precision as he goes, but I've no idea how far above or below the required precision this is. However there is no need to use a precision of 55 from a certain point on. That probably forces it to use two doubles instead of one in the floating point computations, which I think would be unnecesary. You can look on mathworld for how well the Rademacher series converges. It would probably be quicker to do a fresh implementation from scratch given the application in mind. The code there may not be GPL'd also. To check it works, one should at least look at the congruences modulo 5, 7, 11 and 13. The author of the mpfr version does seem to check them mod 11 but for very small n and only for a very few values of n. If the gcds prove to be a bottleneck, I have some code that is roughly twice as fast as GMP for computing these. Certainly the code in the mpfr version is suboptimal and needs some serious optimisation when the terms in the dedekind sum fit into a single limb, which on a 64 bit machine is probably always, when n is of a reasonable size. Another suggestion is that for sufficiently small values of h and or k, it may be quicker to compute the dedekind sum directly rather than use the gcd formula. A lookup table would also be useful for the tail end of the gcd/ dedekind sum computation. I'd be shocked if the computation for n = 10^9 couldn't be done in under 10s on sage.math. Bill. -- William Stein Associate Professor of Mathematics University of Washington http://www.williamstein.org --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~--~~~~--~~--~--~--- /* Implementation of the unrestricted partition function using mpfr (see * mpfr.org). (C) 2002 Ralf Stephan ([EMAIL PROTECTED]). * See part.pdf on the same path for a summary of the math. * Distributed under GPL (see gnu.org). Version 2002-12. * * The unrestricted partition function counts the nonnegative integer * solutions to a+2b+3c+...=n.
[sage-devel] Re: computing the number of partitions of an integer
Your version did work for me! Pablo On 7/27/07, William Stein [EMAIL PROTECTED] wrote: On 7/27/07, Pablo De Napoli [EMAIL PROTECTED] wrote: I've tested it but seems to be buggy:: it works up to 156 ./part 156 p(156)= 73232243759 but for 157 gives a floating point exception error (and a gdb tracing says it is in the line 176 of the source code r2 = r0 % r1; in function g I've attached a slightly modified version of part.c that seems to work (it doesn't crash and gives correct answers in all the cases I tested). I just changed the code in g slighlty so if r1 is 0, then r2 is also set to 0. I compile it using gcc part.c -O3 -g -lmpfr -lm -o part -I/home/was/sage/local/include -L/home/was/sage/local/lib/ -lgmp where /home/was/sage is the path to my SAGE install. Interestingly, when I time part.c it is not much better than CVS pari (see below). Bill Hart wrote: This, with the other ideas I gave above, will definitely let us beat Mathematica convincingly, as a simple back of the envelope calculation shows. It will have the added benefit of allowing us to beat Magma at something. Magma is already terrible at computing NumberOfPartitions, though not as bad as GAP (which is even worse): [EMAIL PROTECTED]:~$ magma Magma V2.13-5 Fri Jul 27 2007 10:54:07 on ubuntu [Seed = 1720324423] Type ? for help. Type Ctrl-D to quit. time k := NumberOfPartitions(10^5); Time: 3.070 With the latest CVS pari: ? gettime; n=numbpart(10^5); gettime/1000.0 %1 = 0.01200 ? gettime; n=numbpart(10^6); gettime/1000.0 %2 = 0.2080 ? gettime; n=numbpart(10^7); gettime/1000.0 %3 = 6.973 ? gettime; n=numbpart(10^8); gettime/1000.0 With part.c: I use time ./part 10 out and record system + user time: 10^60.04 seconds 10^74.584 seconds 10^8 47.483 seconds So, to clarify or summarize, Bill is there one thing that we could change in part.c that would speed it up? I realized you sort of answered this before, but the sequence of previous emails yesterday was rather long and may have got confused. Thanks, william The code theres seems indeed to be GPL-ed, it has a copyright notice that says: (C) 2002 Ralf Stephan ([EMAIL PROTECTED]). * See part.pdf on the same path for a summary of the math. * Distributed under GPL (see gnu.org). Version 2002-12. The pdf is: http://www.ark.in-berlin.de/part.pdf Pablo On /26/07, Bill Hart [EMAIL PROTECTED] wrote: Alternatively you could just use the implementation here: http://www.ark.in-berlin.de/part.c which seems to only rely on mpfr and GMP. However, since the Pari version was based on it, I suspect it may also need correction. He does seem to adjust the precision as he goes, but I've no idea how far above or below the required precision this is. However there is no need to use a precision of 55 from a certain point on. That probably forces it to use two doubles instead of one in the floating point computations, which I think would be unnecesary. You can look on mathworld for how well the Rademacher series converges. It would probably be quicker to do a fresh implementation from scratch given the application in mind. The code there may not be GPL'd also. To check it works, one should at least look at the congruences modulo 5, 7, 11 and 13. The author of the mpfr version does seem to check them mod 11 but for very small n and only for a very few values of n. If the gcds prove to be a bottleneck, I have some code that is roughly twice as fast as GMP for computing these. Certainly the code in the mpfr version is suboptimal and needs some serious optimisation when the terms in the dedekind sum fit into a single limb, which on a 64 bit machine is probably always, when n is of a reasonable size. Another suggestion is that for sufficiently small values of h and or k, it may be quicker to compute the dedekind sum directly rather than use the gcd formula. A lookup table would also be useful for the tail end of the gcd/ dedekind sum computation. I'd be shocked if the computation for n = 10^9 couldn't be done in under 10s on sage.math. Bill. -- William Stein Associate Professor of Mathematics University of Washington http://www.williamstein.org --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~--~~~~--~~--~--~---
[sage-devel] Re: computing the number of partitions of an integer
On 7/27/07, Pablo De Napoli [EMAIL PROTECTED] wrote: Your version did work for me! Pablo Great. Is there any chance you could work on any of the optimization or timing ideas Bill Hart suggested in the previous post? I won't be able to work on this myself for a while. -- William --~--~-~--~~~---~--~~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~--~~~~--~~--~--~---
