Thread: [q-lang-users] Problem using function "isprime"
Brought to you by:
agraef
From: <we...@gm...> - 2004-10-07 18:20:54
|
Hello, I got a problem using function "isprime" from GMP: : : ==> isprime 999907; true ==> isprime 1000003; isprime 1000003 ==> : : It shouldn't be a problem to determine that 1000003 is a prime? Best regards, Willi |
From: <Dr....@t-...> - 2004-10-09 12:06:50
|
we...@gm... wrote: > ==> isprime 1000003; > isprime 1000003 > > It shouldn't be a problem to determine that 1000003 is a prime? Well, this uses GMP's probabilistic prime test (Miller-Rabin algorithm), which may fail (in which case the Q function also fails). As a remedy, you can try to increase the clib::ISPRIME_REP value, or add your own test which catches the few cases when the Miller-Rabin test fails. See the chapter on "Clib", Section "Additional Integer Functions" in the manual for details. Here's a suitable definition which you can add to your program, copied straight from the docs: isprime N:Int = true if N = 2; = false if N and 1 = 0; // even number = try_div 3 (intsqrt N) N otherwise; try_div L M N = true if L > M; = false if N mod L = 0; = try_div (L+2) M N otherwise; -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-09 12:53:10
|
Dr....@t-... (Albert Graef) writes: [snip] > Here's a suitable definition which you can add to your program, copied > straight from the docs: [snip] While we're here, consider this pair of mutually-recursive functions for the job, as well: isprime N = not (any ((=0).(N mod))) (primes (intsqrt N)) ; primes N = [2] if N=2; = filter isprime [2..N] otherwise; Not quite as pretty to look at, but I think it may run faster as it only uses the primes between 2..intsqrt N rather than all odds, if I read your algorithm right... HTH, ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-09 21:15:21
|
Tim Haynes wrote: > While we're here, consider this pair of mutually-recursive functions for > the job, as well: > > isprime N = not (any ((=0).(N mod))) (primes (intsqrt N)) ; > > primes N = [2] if N=2; > = filter isprime [2..N] otherwise; > > Not quite as pretty to look at, but I think it may run faster as it only > uses the primes between 2..intsqrt N rather than all odds, if I read your > algorithm right... Hi, Tim! From the algorithm I posted: > = try_div 3 (intsqrt N) N otherwise; ^^^^^^^ See? :) While your algorithm looks quite neat, I doubt that it will be any faster, since it uses a list, while the one from the docs is just plain tail-recursive number crunching (which has the added benefit that it only needs O(1) space). I'm also a bit worried about the recursive isprime call, which might slow down things quite a bit if enough members of the [2..N] list fail the Miller-Rabin test. Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-09 22:39:31
|
Dr....@t-... (Albert Graef) writes: [snip] >> Not quite as pretty to look at, but I think it may run faster as it only >> uses the primes between 2..intsqrt N rather than all odds, if I read your >> algorithm right... > > Hi, Tim! From the algorithm I posted: > > > =3D try_div 3 (intsqrt N) N otherwise; > ^^^^^^^ > > See? :) Yes... I had that too... > While your algorithm looks quite neat, I doubt that it will be any > faster, since it uses a list, while the one from the docs is just plain > tail-recursive number crunching (which has the added benefit that it only > needs O(1) space). I'm also a bit worried about the recursive isprime > call, which might slow down things quite a bit if enough members of the > [2..N] list fail the Miller-Rabin test. Very interesting concerns, but I've just been very rude and benchmarked the two... :) For smallish N<1E4, it's evenly matched (3% more occasions on which mine is faster than yours).=20 However, for N approaching 1E9 or so, the number of occasions on which my algorithm is faster pulls away - around 2.1E9 I've been seeing about 3x as many occasions on which mine's been faster compared to yours. It seems to be taking 2-3% less time, when it's better. The recursive isprime call is actually a win; effectively it's calling intsqrt(intsqrt(intsqrt(intsqrt.....))) a few times (per prime in range, not per-odd-number in range) which tends towards 1 quicker than your concern about implementation slows it down. Bear in mind that only 1/12 numbers in the range [1..1E6] is a prime, too... It can of course be considerably quicker: just by amending=20 primes N =3D [2] if N=3D2;=20=20=20=20=20=20=20=20=20 =3D [2,3] if N=3D3; //cheating! =3D filter isprime [2..N] otherwise; it's faster than yours on 10x as many ocassasion than where yours was faster than mine, around N<=3D2E9, so it's hitting the bottom few cases of the recursion quite a *lot*. ~~~~ What I'd like to know is whether there's any memoization going on - by the time I've evaluated `isprime 13' once for intsqrt(N), I don't need to do it all over again for intsqrt(intsqrt(N)), etc. Actually, how would I write a memoize function in Q, anyway? I don't see that, or the word `cache', in the docs anywhere :) ~Tim =2D-=20 <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-11 11:35:35
|
Tim Haynes wrote: > Very interesting concerns, but I've just been very rude and benchmarked the > two... :) > For smallish N<1E4, it's evenly matched (3% more occasions on which mine is > faster than yours). The trouble with theory is that it's theoretically true but practically useless. :) I stand corrected... > It can of course be considerably quicker: just by amending > > primes N = [2] if N=2; > = [2,3] if N=3; //cheating! I don't consider that cheating, IIRC the Milner Rabin algorithm itself (at least the GMP version) also uses a (fairly large) table of known primes. I'd even go farther: primes N = [2] if N=2; = [2,3] if N=3; = [2,3|filter isprime [5,7..N]] otherwise; Does that speed up your algorithm even further? I haven't tried. > it's faster than yours on 10x as many ocassasion than where yours was > faster than mine, around N<=2E9, so it's hitting the bottom few cases of > the recursion quite a *lot*. Ok, ok, I'm convinced. An update of the isprime section is on the TODO list. ;-) > What I'd like to know is whether there's any memoization going on - by the > time I've evaluated `isprime 13' once for intsqrt(N), I don't need to do it > all over again for intsqrt(intsqrt(N)), etc. Maybe the GMP function does some memoization internally? The Q interpreter doesn't do it, that I know for sure. > Actually, how would I write a memoize function in Q, anyway? I don't see > that, or the word `cache', in the docs anywhere :) At your request :) http://q-lang.sourceforge.net/qdoc/qdoc_12.html#SEC118 Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-11 22:17:16
|
Dr....@t-... (Albert Graef) writes: [snip] > I don't consider that cheating, IIRC the Milner Rabin algorithm itself (at > least the GMP version) also uses a (fairly large) table of known primes. > I'd even go farther: > > primes N = [2] if N=2; > = [2,3] if N=3; > = [2,3|filter isprime [5,7..N]] otherwise; > > Does that speed up your algorithm even further? I haven't tried. Yes, a little, but not that much, I think. [snip] >> Actually, how would I write a memoize function in Q, anyway? I don't see >> that, or the word `cache', in the docs anywhere :) > > At your request :) > > http://q-lang.sourceforge.net/qdoc/qdoc_12.html#SEC118 Oh, *excellent* <rubs hands in glee>. Thanks for that. Refs as well? I have much Q to be learning. :) I've implemented memoization on `prime' and `isprime' functions - works first time ;) It emphasizes the extreme contrast between composite and prime (or few-large-factors not-quite-prime) inputs; whereas it was always close (3% ish) previously, with both algorithms suffering longer computation-times for larger nearly-prime input N, I'm now seeing improvements of 50% or more in individual runs - my algorithm is almost constant/linear by comparison. This one wins :) Cheers, ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-12 02:09:28
|
Tim Haynes wrote: > Oh, *excellent* <rubs hands in glee>. Thanks for that. Refs as well? I have > much Q to be learning. :) Also soon to be unleashed upon the unexpecting world: SWIG support! This will turn the creation of library interfaces into a no-brainer. But don't tell anyone just yet, I'm still testing. ;-) C already works very well, I wish I had that when I wrote the clib module ... As soon as I'm satisfied with the result of auto-wrapping a huge C++ library like the entire Qt it's time for a release again. > I've implemented memoization on `prime' and `isprime' functions - works > first time ;) Care to post your final algorithm here? Cheers, Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-12 09:59:14
|
Dr....@t-... (Albert Graef) writes: > Tim Haynes wrote: >> Oh, *excellent* <rubs hands in glee>. Thanks for that. Refs as well? I have >> much Q to be learning. :) > > Also soon to be unleashed upon the unexpecting world: SWIG support! This > will turn the creation of library interfaces into a no-brainer. But don't > tell anyone just yet, I'm still testing. ;-) C already works very well, I > wish I had that when I wrote the clib module ... As soon as I'm satisfied > with the result of auto-wrapping a huge C++ library like the entire Qt > it's time for a release again. Oh wow. If it's really capable of coping with something that size, it'll be amazing. I'm already wondering about writing a tiny graphical app in Q (an RSS aggregator with a difference) as it is :) >> I've implemented memoization on `prime' and `isprime' functions - works >> first time ;) > > Care to post your final algorithm here? Sure. :) I think this is the best of all worlds so far: #!/usr/bin/env q #! -q #! -cmain ARGS || quit def HASH = ref emptyhdict; public hashed F X; hashed F X = get HASH!(F,X) if member (get HASH) (F,X); = put HASH (update (get HASH) (F,X) Y) || Y where Y = F X otherwise; hashed isprime = hisprime; hisprime N = not (any ((=0).(N mod))) (primes (intsqrt N)) ; hashed primes = hprimes; hprimes N = [2] if N=2; = [2,3] if N=3; = [2,3|filter isprime [5,7..N]] otherwise; def P=val (ARGS!1); dump N = strcat [str N, " is ", (str (isprime N)), "\n"]; main ARGS = dump P; (The reason for `dump' and main being separate is that I was running dump, and therefore isprime, across a list [(P-10)..(P+10)], for benchmarking purposes.. :) You're welcome to use it for docs / consider it GPL'd / whatever's easiest, as you wish. :) ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-21 00:48:25
|
Hi, Tim, > Oh wow. If it's really capable of coping with something that size, it'll be > amazing. I'm already wondering about writing a tiny graphical app in Q (an > RSS aggregator with a difference) as it is :) You could also use the tk module for that. I've already written some substantial applications with that module, so it's fairly well-tested. I also use vtcl as a GUI builder, the necessary boilerplate code to load the created Tcl scripts from Q is quite simple. While Tcl/Tk looks quite old and ugly when compared to Qt, it's easy to program it and you have the added benefit of cross-platform portability. >>Care to post your final algorithm here? > > Sure. :) I think this is the best of all worlds so far: Thanks a lot, I'll add this as an example to the next Q release. Which is imminent, as I'm almost finished with Q-SWIG. Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-21 08:31:49
|
Dr....@t-... (Albert Graef) writes: > Hi, Tim, > >> Oh wow. If it's really capable of coping with something that size, it'll >> be amazing. I'm already wondering about writing a tiny graphical app in >> Q (an RSS aggregator with a difference) as it is :) > > You could also use the tk module for that. I've already written some > substantial applications with that module, so it's fairly well-tested. I > also use vtcl as a GUI builder, the necessary boilerplate code to load > the created Tcl scripts from Q is quite simple. While Tcl/Tk looks quite > old and ugly when compared to Qt, it's easy to program it and you have > the added benefit of cross-platform portability. Yeah, I'll bear it in mind. Actually, I was rather hoping to run it against GTK.. :) Meanwhile, I'm working on the core algorithms for the RSS-reader... the integration with libxml2 is coming into its own :) ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-22 16:28:54
|
Tim Haynes wrote: > Yeah, I'll bear it in mind. Actually, I was rather hoping to run it against > GTK.. :) Well, while SWIG-Q seems to work pretty well now, and I've written my first Qt hello world program in Q, my first experiments with wrapping the entire Qt and GTK libraries have been rather sobering. :( The problem is not with SWIG-Q, but with SWIG itself. Neither the GTK nor the Qt headers are SWIG-friendly at all, so I'd have to manually extract all relevant declarations from that gazillion of header files. While this is already a lot better than writing wrapper modules by hand, it's still a big amount of work. And I haven't been able to find ready-made SWIG interfaces for GTK and Qt on the web either. Maybe I should just release the SWIG-enabled Q version as it is now (it's already in cvs, I just need to build/test it on Windows and I plan to do that tonight), and sit back and see what you and others can do with it. ;-) > Meanwhile, I'm working on the core algorithms for the RSS-reader... the > integration with libxml2 is coming into its own :) Feedback on that module is appreciated, of course. Also, I'm awaiting some nifty new Q-CGI examples to be included in a future Q release ... :) Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: <Dr....@t-...> - 2004-10-23 09:20:58
Attachments:
primetest.q
|
Tim Haynes wrote: > def HASH = ref emptyhdict; > > public hashed F X; > hashed F X = get HASH!(F,X) if member (get HASH) (F,X); > = put HASH (update (get HASH) (F,X) Y) || Y > where Y = F X otherwise; > > hashed isprime = hisprime; > > hisprime N = not (any ((=0).(N mod))) (primes (intsqrt N)) ; > > hashed primes = hprimes; > > hprimes N = [2] if N=2; > = [2,3] if N=3; > = [2,3|filter isprime [5,7..N]] otherwise; Hmm, when I run this code over here, the hashed function is apparently never invoked. In fact, with filter isprime [5000000..5000100] I even get a failed condition: ! Error in conditional 0> stdlib.q, line 134: filter isprime [5000011,5000012,5000013,...] ==> [5000011|filter isprime [...]] if isprime 5000011 I think the applications of hashed are the wrong way round (it's easy to fall into this trap, it has happened to me, too). I've attached my (hopefully corrected) version of the script to this post. This one *does* invoke hashed now (verified with tbreak). I'm ready to add this version to the examples now (if you don't complain -- better do that quickly, as I'm currently doing the final touches for Q 6.0 ;-). Cheers, Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-26 16:23:29
|
Dr....@t-... (Albert Graef) writes: > Tim Haynes wrote: >> def HASH = ref emptyhdict; >> public hashed F X; >> hashed F X = get HASH!(F,X) if member (get HASH) (F,X); >> = put HASH (update (get HASH) (F,X) Y) || Y >> where Y = F X otherwise; >> hashed isprime = hisprime; >> hisprime N = not (any ((=0).(N mod))) (primes (intsqrt N)) ; >> hashed primes = hprimes; >> hprimes N = [2] if N=2; = [2,3] if N=3; >> = [2,3|filter isprime [5,7..N]] otherwise; > > Hmm, when I run this code over here, the hashed function is apparently > never invoked. Very strange! THe performance was completely different from the unhashed version - I observed it staying all-but linear when the all-the-odds algorithm, or even the unhashed one, were showing increases, here.. And yes it seemed to be producing the right answers, for a few random spot-checks, too :) > In fact, with filter isprime [5000000..5000100] I even get a failed > condition: > > ! Error in conditional > 0> stdlib.q, line 134: filter isprime [5000011,5000012,5000013,...] ==> > [5000011|filter isprime [...]] if isprime 5000011 Very weird. I did actually see a glitch with one large number, so I'll have to go back and investigate it more... > I think the applications of hashed are the wrong way round (it's easy to > fall into this trap, it has happened to me, too). I've attached my > (hopefully corrected) version of the script to this post. This one *does* > invoke hashed now (verified with tbreak). > > I'm ready to add this version to the examples now (if you don't complain > -- better do that quickly, as I'm currently doing the final touches for Q > 6.0 ;-). If it's not too late, I'd like to hold back on that - I'm travelling atm and will be back home tomorrowish. Can run yours and mine through the mill again then. I'm particularly worried about those error-in-conditional problems... :/ ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-30 19:31:21
|
Tim Haynes wrote: > If it's not too late, I'd like to hold back on that - I'm travelling atm > and will be back home tomorrowish. Can run yours and mine through the mill > again then. I'm particularly worried about those error-in-conditional > problems... :/ Hi, Tim, unfortunately Q 6.0 was already released last weekend (didn't have the time to announce this on the ML this time, I'm so busy now that the semester has begun). My corrected version of your script is included in the examples. That versions seems to work ok, though. I'd suggest that if you have a new version I'll update it in cvs immediately and then release it with the next Q version. Ok? Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Tim H. <q...@st...> - 2004-10-30 21:29:44
|
Dr....@t-... (Albert Graef) writes: > Tim Haynes wrote: >> If it's not too late, I'd like to hold back on that - I'm travelling atm >> and will be back home tomorrowish. Can run yours and mine through the mill >> again then. I'm particularly worried about those error-in-conditional >> problems... :/ > > Hi, Tim, > > unfortunately Q 6.0 was already released last weekend (didn't have the > time to announce this on the ML this time, I'm so busy now that the > semester has begun). That is not unfortunate ... :) > My corrected version of your script is included in the examples. That > versions seems to work ok, though. I'd suggest that if you have a new > version I'll update it in cvs immediately and then release it with the > next Q version. Ok? I've just got some initial benchmark results - it's still good, it's still obviously quicker than the naiive baseline (faster about 3.5x as often as naiive is faster than mine). And when it's faster, it's by a factor of 2-3, while I've never seen it being over 12% slower at worst. Graphs of times taken show a couple of curves much lower than the naiive algorithm, too. Neat :) ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-10-31 09:19:55
|
Tim Haynes wrote: > I've just got some initial benchmark results - it's still good, it's > still obviously quicker than the naiive baseline (faster about 3.5x as > often as naiive is faster than mine). And when it's faster, it's by a > factor of 2-3, while I've never seen it being over 12% slower at worst. Great! If you have some benchmark results or other comments ready to be added to the script, just send them to me. > Graphs of times taken show a couple of curves much lower than the naiive > algorithm, too. Just remember that all of this is achieved at the expense of O(N log N) memory requirements (whereas the naive algorithm uses only O(1) data). ;-) > Neat :) I agree. :) Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |
From: Larry G. <lg...@ac...> - 2004-10-31 22:03:42
|
People: I pretty much copied this from the help screen: def F = fopen c:\test\mytest.txt "w" ; fwrites F "Hello, world." ; fclose F; Why does it give me this error in Qpad? % C:\Program Files\QPAD\TEST\Filew01.q Error Filew01.q, line 1: parse error at or near symbol `:' Error Filew01.q, line 2: parse error at or near symbol `;' Error Filew01.q, line 3: parse error at or near symbol `;' The error message is not too helpful. Should the filename be in double quotes? Are there some simple, complete, examples of file I/O? Larry Gregg |
From: Tim H. <q...@st...> - 2004-10-31 22:20:49
|
Larry Gregg <lg...@ac...> writes: > People: > > I pretty much copied this from the help screen: > > def F = fopen c:\test\mytest.txt "w" ; > fwrites F "Hello, world." ; > fclose F; > > Why does it give me this error in Qpad? > > % C:\Program Files\QPAD\TEST\Filew01.q > Error Filew01.q, line 1: parse error at or near symbol `:' > Error Filew01.q, line 2: parse error at or near symbol `;' > Error Filew01.q, line 3: parse error at or near symbol `;' > > The error message is not too helpful. Should the filename be in double > quotes? Yes, it's a string field. Also I think the \ should be doubled-up, as well. > Are there some simple, complete, examples of file I/O? <http://q-lang.sourceforge.net/qdoc/qdoc_12.html#SEC91> might be a reasonable place to start. HTH, ~Tim -- <http://spodzone.org.uk/> |
From: <Dr....@t-...> - 2004-11-01 08:26:42
|
Larry Gregg wrote: > def F = fopen c:\test\mytest.txt "w" ; > fwrites F "Hello, world." ; > fclose F; > > Why does it give me this error in Qpad? > > % C:\Program Files\QPAD\TEST\Filew01.q > Error Filew01.q, line 1: parse error at or near symbol `:' > Error Filew01.q, line 2: parse error at or near symbol `;' > Error Filew01.q, line 3: parse error at or near symbol `;' As Tim already pointed out, you have to put the filename in double quotes and escape the backslashes (or replace them with slashes). Probably you'll also want a newline at the end of the string. You can enter this directly in the interpreter by running an empty script and then entering these commands at the interpreter prompt: ==> def F = fopen "c:/test/mytest.txt" "w" ==> fwrites F "Hello, world.\n" ==> fclose F You also have to distinguish between the *programming language* in which you write a Q script (the "upper pane" in Qpad), and the *command language* of the interpreter which allows you (among other things) to type expressions to be evaluated (the "lower pane" in Qpad). In the manual, something to be entered in the interpreter is usually shown with the ==> prompt in front of it, just like I did above. If you want to put the above into a script, you'll have to wrap it up in a function definition (i.e., a collection of equations). Scripts can only contain variable definitions and equations, not just an expression by itself. So you can write something like: test = fwrites F "Hello, world.\n" where F = fopen "mytest.txt" "w"; Once you enter this script, save it to a file and run it with the interpreter, you can call the test function in the interpreter: ==> test Some important differences between the programming language and the interpreter's command language are: - The programming language is free-format, while the command language is line-oriented, to facilitate interactive usage. - The programming language allows you to define functions, while the command language allows you to evaluate expressions. - The programming language has a lot of additional constructs to declare functions, variables and types, while in the command language you have some special commands for interactive usage (such as debugger and profiling commands). I know that this can be a little confusing for beginners since a subset of the Q language (expression and variable definition syntax) is used both in the programming and in the command language. A closer description of the command language of the interpreter can be found in Appendix B.2 of the manual, see http://q-lang.sourceforge.net/qdoc/qdoc_14.html#SEC143. HTH Albert -- Dr. Albert Gr"af Dept. of Music-Informatics, University of Mainz, Germany Email: Dr....@t-..., ag...@mu... WWW: http://www.musikwissenschaft.uni-mainz.de/~ag |