20221127, 19:51  #210  
Sep 2002
Database er0rr
2^{3}×3×11×17 Posts 
Quote:
#include "/home/paul/Downloads/gmp6.2.1/gmpimpl.h" I am looking forward to seeing a nice speed up in production runs. Last fiddled with by paulunderwood on 20221127 at 20:59 

20221129, 01:04  #211 
"Serge"
Mar 2008
Phi(4,2^7658614+1)/2
3^{2}×1,117 Posts 
Thank you, Robert.
The two hack together cut down the run time by 2/3 in two categories: Cornacchia, is_prime. The next contender for a 2/3cut would be probably qroot? Code:
... discriminants: 99.6 (99.8) 17174 qroot: 561090.6 (19718.9) < now the largest spender 58575671 Cornacchia: 238202.6 (8039.3) 380884 trial div: 314247.4 (9941.6) 98055 is_prime: 522238.2 (18328.6)  Size [92]: 68343 bits 
20221129, 01:57  #212  
"Robert Gerbicz"
Oct 2005
Hungary
3^{2}·179 Posts 
Quote:
on nt.c see line 268, 300, 308, you don't need more change, the others are trivial or expected O(1) number of mulmod, just use the gwnum library for powm. Paul could make the changes. 

20221129, 14:25  #213  
Sep 2002
Database er0rr
2^{3}·3·11·17 Posts 
Quote:
I am finding it difficult to debug. If someone else want to have a crack at it then please go ahead. Last fiddled with by paulunderwood on 20221129 at 14:26 

20221129, 14:38  #214  
"Robert Gerbicz"
Oct 2005
Hungary
1611_{10} Posts 
Quote:
And on line 268 the b value is quite small (fits in 32 bits), for the other two lines it is just random number from [0,n). Btw even in general it would not be a problem because mpz_mod(b,b,n) reduces b into the [0,n) interval, and then you could do the exponentiation. Last fiddled with by R. Gerbicz on 20221129 at 14:40 

20221129, 20:20  #215 
Jul 2003
So Cal
7^{2}·53 Posts 
Can this be relatively easily rewritten with mpz_size(), mpz_limbs_modify(), and mpz_limbs_finish() to remove the dependency on gmpimpl.h, or is there something I'm missing in my quick review of the code?

20221129, 21:16  #216  
"Robert Gerbicz"
Oct 2005
Hungary
3113_{8} Posts 
Quote:
OK, then we could copy those definitions, but it could be not the same in different gmp versions and could change in the future... 

20221129, 21:35  #217 
Jul 2003
So Cal
7^{2}×53 Posts 
Thanks. I agree that copying those definitions isn't a good idea. I was hoping that it could be done using the documented public GMP interface relatively easily without loss of speed. Too much to ask perhaps!

20221130, 06:28  #218  
Sep 2002
Database er0rr
2^{3}·3·11·17 Posts 
Quote:
Code:
/*****************************************************************************/ unsigned int cm_nt_mpz_tonelli_generator (mpz_ptr q, mpz_ptr z, mpz_srcptr p) /* For p an odd prime, compute and return e such that p1 = 2^e * q with q odd. If e>1 (that is, p != 3 mod 4), also q and a generator z of the 2Sylow subgroup of (Z/pZ)^* are returned. Otherwise q and z are not modified. */ { unsigned int e = mpz_scan1 (p, 1); if (e > 1) { mpz_tdiv_q_2exp (q, p, e); for (mpz_set_ui (z, 2ul); mpz_legendre (z, p) != 1; mpz_add_ui (z, z, 1ul)); if ( mpz_sizeinbase ( q, 2 ) < 120 ) mpz_powm (z, z, q, p); else { #define SAFE 50 char string [110000]; int i, length; unsigned fft_size = 0; giant gz, gb, gp; gwnum wr, wb, wc; gwhandle *gwdata; length = ( strlen ( string ) >> 2 ) + 8; gp = newgiant ( length ); mpz_get_str ( string, 10, p ); ctog ( string, gp ); gwdata = (gwhandle*) malloc ( sizeof ( gwhandle ) ); gwinit ( gwdata ); gwset_larger_fftlen_count ( gwdata , fft_size ); gwsetup_general_mod_giant ( gwdata, gp ); wr = gwalloc ( gwdata ); wb = gwalloc ( gwdata ); wc = gwalloc ( gwdata ); mpz_get_str ( string, 10, z ); gb = newgiant ( length ); ctog ( string, gb ); gianttogw ( gwdata, gb, wb ); int j = mpz_sizeinbase ( q, 2 )  2; int k = j  SAFE; gwcopy ( gwdata, wb, wr ); for ( i = j; i > k; i ) { gwsquare2_carefully ( gwdata, wr, wr ); if ( mpz_tstbit ( q, i ) ) gwmul ( gwdata, wb, wr ); } for ( i = k; i > SAFE; i ) { gwcopy ( gwdata, wr, wc ); gwmul ( gwdata, wr, wr ); if ( mpz_tstbit ( q, i ) ) gwmul ( gwdata, wb, wr ); /* if ( gw_get_maxerr ( gwdata ) < 0.44 ) { gwswap ( wr, wc ); } else { printf ( "*** Maximum error %lf excessive at bit %d. Increasing FFT size.\n", gw_get_maxerr( gwdata ), i ); gwtogiant ( gwdata, wr, gz ); gwdone ( gwdata ); ++fft_size; gwdata = (gwhandle*) malloc ( sizeof ( gwhandle ) ); gwinit ( gwdata ); gwset_larger_fftlen_count ( gwdata, fft_size ); gwsetup_general_mod_giant ( gwdata, gp ); wr = gwalloc ( gwdata ); wb = gwalloc ( gwdata ); wc = gwalloc ( gwdata ); gwerror_checking ( gwdata, 1 ); gianttogw ( gwdata, gz, wr ); gianttogw ( gwdata, gb, wb ); i++; } */ } for ( i = SAFE; i >= 0; i ) { gwsquare2_carefully ( gwdata, wr, wr ); if ( mpz_tstbit ( q, i ) ) gwmul ( gwdata, wb, wr ); } gz = newgiant ( length ); gwtogiant ( gwdata, wr, gz ); gtompz( gz, z ); //mpz_out_str (NULL,10,z);printf("z\n"); gwdone ( gwdata ); } } return e; } /*****************************************************************************/ Last fiddled with by paulunderwood on 20221130 at 07:25 

20221201, 23:40  #219 
Sep 2002
Database er0rr
1000110001000_{2} Posts 
tonelli beta
I found the bug in my code. Attached is beta code. I think I have "length" right at 2^17 characters. Sorry for the ugliness of the code  sheesh, two functions when one would suffice. The kickin threshold is set at 4096 bits in 3 places. I shall inflict this code on my own production run.
(To install cm (ecppmpi) on Debian based systems follow EdH's instructions.) Follow the instructions here and install Robert's code for some nice speedups. Then just overwrite the two sections in lib/nt.c with the attached. The beta3 file uses faster gwswap rather than gwcopy. It also presumes no round off error is thrown by gwsmallmul. (beta and beta2 have been pulled because error checking was not turned on and because gwmul3_carefully was not there too  thanks to forumite Luminescence for pointing out the errors.) Last fiddled with by paulunderwood on 20221203 at 03:54 
20221204, 03:17  #220 
Sep 2002
Database er0rr
4488_{10} Posts 
New working code
(To install cm (ecppmpi) on Debian based systems follow EdH's instructions.) After the instructions given here to get GWNUM working, install Robert's code for some nice speedups. To cut along story short, attached is a C file that can be dropped into CM's lib directory. Then make the following changes to lib/nt.c: Code:
#include <stdio.h> #include "gw_utility5.c" #include "cmimpl.h" Code:
int cm_nt_is_prime (mpz_t a) { return ( gw_prp ( a ) ); } Also in nt/curve.c add #include gw_utility5.c and change mpz_powm to gw_powm. Let me know if you have any trouble with the C code. (Please redownload gw_utility5.c because I have added a mpz_clear.) Last fiddled with by paulunderwood on 20221208 at 19:48 
Thread Tools  
Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
what is the best primality software?  bbb120  Proth Prime Search  75  20221008 11:45 
Fastest software for Mersenne primality test?  JonathanM  Information & Answers  25  20200616 02:47 
APRCL as primality proof  f1pokerspeed  FactorDB  14  20140109 21:06 
Proof of Primality Test for Fermat Numbers  princeps  Math  15  20120402 21:49 
PRIMALITY PROOF for Wagstaff numbers!  AntonVrba  Math  96  20090225 10:37 