|
From: William H. <ha...@ya...> - 2008-08-20 17:08:10
|
Hi Peter,
This is very interesting. There is a fairly dramatic
jump at 1122004669633. This is the point where 5
pseudoprime tests are used instead of 4, and not 10^16
where Miller Rabin is used.
The reason for the jump is that SPRP64 has to be used
instead of SPRP since the numbers can exceed the size
of an unsigned long.
The Lucas test seems to be very good, and given that a
counterexample to the conjecture that it plus a 2-SPRP
test is a primality test is worth money, we may as
well use it in FLINT for our probably prime testing up
to 2^64.
We'll use Pocklington-Lehmer for z_isprime after
10^16, I think. I'll be very interested to see what
the time comparison is like for that, not that it
matters, because we currently don't have any test in
FLINT which is guaranteed to separate composites from
primes after 10^16.
Bill.
--- Peter Shrimpton <ps...@gm...> wrote:
> Hello
>
> Attached are spreadsheets (one excel the other
> openoffice) showing the
> relative speeds of z_isprime and a Baillie-PSW test.
> In contradiction to my
> earlier conversation with Bill (I was getting the
> colours confused) the
> Baillie-PSW test does perform a lot better then then
> z_isprime as soon as
> z_isprime has to do 5 pseudprime tests. However
> whist testing this I found
> that there is a bug in my lucas test code which
> makes it return a handfull
> of composites as prime. It only seems to affect
> about 10 or so numbers which
> are all 'quite small' (about 8 digits long).
>
> Also I have managed to get mongomery powering
> working for unsigned longs
> (see attached). I realised that as you only need the
> first few bits for some
> of the multiplications and so the whole opperation
> can fit into 2 unsigned
> longs. Currently it is slightly slower than z_powmod
> in most places although
> I can occasionally get improvements. I think It
> depends on the relative
> sizes of the parameters and/or how many 1's or 0's
> they have but as there
> are 3 of them I have no idea how to test it
> completely to see which cases it
> is faster.
>
> Peter
> > #include <gmp.h>
> #include "flint.h"
> #include "long_extras.h"
> #include <stdio.h>
> #include "longlong.h"
>
> unsigned long mprod_precomp(unsigned long x,
> unsigned long y, unsigned long n, unsigned long
> nprime, unsigned long r1, unsigned long bits)
> {
> unsigned long z,high,low,xyh,xyl;
>
> umul_ppmm(xyh,xyl,x,y);
> z=(xyl)&r1;
> z=(z*(nprime))&r1;
> umul_ppmm(high,low,z,n);
> add_ssaaaa(high,low,xyh,xyl,high,low);
> z=(low>>bits)+(high<<(FLINT_BITS-bits));
>
>
> if (z>n){
> return z-n;
> }
> return z;
> }
>
>
> unsigned long mpow(unsigned long const x, unsigned
> long const y, unsigned long const n)
> {
> int bits=FLINT_BIT_COUNT(n);
> unsigned long r,xhat,phat,power,nprime;
> pre_inv2_t ninv;
> r=(1L<<bits);
> power=1L<<FLINT_BIT_COUNT(y);
> ninv=z_precompute_inverse(n);
> nprime=r-z_invert(n,r);
> phat=z_mod2_precomp(r,n,ninv);
> xhat=z_mulmod2_precomp(x,phat,n,ninv);
> r=r-1;
>
> while (power)
> {
> phat=mprod_precomp(phat,phat,n,nprime,r,bits);
> if ((y&power)>0){
> phat=mprod_precomp(phat,xhat,n,nprime,r,bits);
> }
> power=power>>1;
> }
> return mprod_precomp(phat,1L,n,nprime,r,bits);
>
> }
>
> int main()
> {
> unsigned long x,y,n;
> int i,lengthx,lengthy,lengthn,minlen,maxlen;
> minlen=60;
> maxlen=63;
> //x=2L;
>
> for (lengthx=minlen; lengthx<maxlen; lengthx++)
> {
> for (lengthy=minlen; lengthy<maxlen; lengthy++)
> {
> for (lengthn=minlen; lengthn<maxlen; lengthn++)
> {
> //printf("%d^%d mod
> %d\n",lengthx,lengthy,lengthn);
> for (i=0; i<100000; i++)
> {
> x=z_randbits(lengthx)+1;
> //y=z_randbits(lengthy)+1;
> n=2*z_randbits(lengthn-1)+1;
> y=n-1;
> //mpow(x,y,n);
> z_powmod2(x,y,n);
> }
> }
> }
> }
> }
>
> >
-------------------------------------------------------------------------
> This SF.Net email is sponsored by the Moblin Your
> Move Developer's challenge
> Build the coolest Linux based applications with
> Moblin SDK & win great prizes
> Grand prize is a trip for two to an Open Source
> event anywhere in the world
>
http://moblin-contest.org/redirect.php?banner_id=100&url=/>
_______________________________________________
> flint-devel mailing list
> fli...@li...
>
https://lists.sourceforge.net/lists/listinfo/flint-devel
>
|