FastECPP software and 50000 digit primality proof (reposted from NMBRTHRY) - Page 13

paulunderwood · 2022-05-29, 10:20

Quote:

Originally Posted by Luminescence

Upon the direct implementation each step was very roughly about 10% faster compared to calling the external executable. Sometimes more, sometimes less.

I got the entire thing to compile by adding the necessary headers to nt.c and putting the gwnum headers into the cm lib folder. Somehow gwthreads.h also ended up in there, not sure if required. I also put them into /usr/local/include, though I am not sure if that's actually necessary.

Then I grabbed the gwnum.a file, renamed it to libgwnum.a and put it in my /usr/local/lib folder.

In the cm lib folder, open Makefile.am, add the gwnum headers + gwthread.h to include_HEADERS and add -lgwnum -ldl -lpthread -lstdc++ to libcm_la_LDFLAGS

Open a terminal in the top level directory, run make distclean, automake and then ./configure --enable_mpi --enable-shared=no. It should (hopefully) compile.

In nt.c I have:

Code:

#include <stdio.h>
#include "gwthread.h"
#include "giants.h"
#include "gwnum.h"
#include "gwcommon.h"

plus the "transplant" of my code.

In Makefile.am I have:

Code:

ACLOCAL_AMFLAGS = -I m4

SUBDIRS = data doc lib src tests

include_HEADERS = gwcommon.h gwthread.h gwcommon.h gwnum.h

libcm_la_LDFLAGS = -lgwnum -ldl -lpthread -lstdc++

Code:

ls lib/g*.h      
lib/giants.h  lib/gwcommon.h  lib/gwnum.h  lib/gwthread.h

Code:

ls /usr/local/lib/libgwnum.a
/usr/local/lib/libgwnum.a

Code:

echo $LD_LIBRARY_PATH                                                        
/usr/local/lib

I ran:

Code:

make distclean
automake
./configure --enable-mpi --enable-shared=no
make clean
make

but get these errors:

Code:

libtool: link: gcc -std=c99 -D CM_MODPOLDIR=\"/usr/local/share/cm\" -I../lib -g -O2 -o cm params.o cm.o  ../lib/.libs/libcm.a -lpari -lmpfrcx -lmpc -lmpfr -lgmp -lz -lm
/usr/bin/ld: ../lib/.libs/libcm.a(nt.o): in function `cm_nt_is_prime':
/home/paul/Downloads/cm-0.4.1dev/lib/nt.c:177: undefined reference to `allocgiant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:178: undefined reference to `allocgiant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:178: undefined reference to `allocgiant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:179: undefined reference to `allocgiant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:181: undefined reference to `ctog'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:184: undefined reference to `gwinit2'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:186: undefined reference to `gwsetup_general_mod_giant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:187: undefined reference to `gwalloc'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:187: undefined reference to `gwalloc'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:191: undefined reference to `bitlen'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:192: undefined reference to `binarytogw'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:193: undefined reference to `binarytogw'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:198: undefined reference to `gwfftfftmul'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:199: undefined reference to `gwsmalladd'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:200: undefined reference to `gwfftfftmul'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:201: undefined reference to `gwsmalladd'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:195: undefined reference to `gwfft'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:196: undefined reference to `gwfft'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:203: undefined reference to `gwfftfftmul'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:204: undefined reference to `gwsmalladd'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:205: undefined reference to `gwfftfftmul'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:206: undefined reference to `gwsmalladd'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:209: undefined reference to `gwtogiant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:210: undefined reference to `gwtogiant'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:211: undefined reference to `modg'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:212: undefined reference to `itog'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:213: undefined reference to `subg'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:214: undefined reference to `modg'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:215: undefined reference to `itog'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:216: undefined reference to `subg'
/usr/bin/ld: /home/paul/Downloads/cm-0.4.1dev/lib/nt.c:217: undefined reference to `gwdone'
collect2: error: ld returned 1 exit status
make[2]: *** [Makefile:419: cm] Error 1
make[2]: Leaving directory '/home/paul/Downloads/cm-0.4.1dev/src'
make[1]: *** [Makefile:460: all-recursive] Error 1
make[1]: Leaving directory '/home/paul/Downloads/cm-0.4.1dev'
make: *** [Makefile:371: all] Error 2

What am I doing wrong?

Luminescence · 2022-05-29, 11:56

Quote:

Originally Posted by paulunderwood

In Makefile.am I have:

Code:

ACLOCAL_AMFLAGS = -I m4

SUBDIRS = data doc lib src tests

include_HEADERS = gwcommon.h gwthread.h gwcommon.h gwnum.h

libcm_la_LDFLAGS = -lgwnum -ldl -lpthread -lstdc++

I ran:

Code:

make distclean
automake
./configure --enable-mpi --enable-shared=no
make clean
make

What am I doing wrong?

Seems like you pasted the headers and LDFLAGS things in the top level directory Makefile.am. Go into the libs/Makefile.am and append the headers and flags to the existing lines.

I forgot that I ran autoreconf -vfi as the very first step in the TLD.

R. Gerbicz · 2022-05-29, 12:01

Quote:

Originally Posted by frmky

This is great! It solves the issue I was having getting it to run on the 32GB Arm system. Thanks!

Great, don't know what is the memory usage in the 2nd phase, if that would be too large then you can complete that phase on another computer with more memory. That phase is much faster than the 1st phase, so at least in time that is not critical. Not touched that part of the code.

There was also an ask to me the version (number) I've modified, used the tarball from https://www.multiprecision.org/cm/
with name cm-0.4.1dev-41c4bce.tar.gz . With that it should compile, anyway there was no real improvement to the ecpp code in the recent weeks, so the concrete version does not matter.

frmky · 2022-05-29, 18:44

Please make use of gwnum optional with an option to configure. The current code supports Arm and we shouldn't break that.

paulunderwood · 2022-05-30, 06:14

Make gwnum.a from the p95 source (in the gwnum directory) and copy gwnum.a to /usr/local/lib/libgwnum.a

Copy gwnum.h gwcommon.h gwthread.h giants.h to CM's lib directory.

Add this line (red) to CM's configure.ac

Code:

LIBS="$LIBS -lgwnum -lpthread -ldl -lstdc++"

AC_CONFIG_FILES([Makefile doc/Makefile lib/Makefile data/Makefile data/af/Makefile data/df/Ma>
AC_OUTPUT

Add this line to lib/Makefile.am

Code:

include_HEADERS = cm.h cm-arith.h

include_HEADERS += gwnum.h gwcommon.h gwthread.h giants.h

Add to lib/nt.c:

Code:

#include <stdio.h>
#include "giants.h"
#include "gwnum.h"
#include "gwcommon.h"
#include "gwthread.h"

and alter the function:

Code:

int cm_nt_is_prime (mpz_t a) {

        char string[100000];
        float fA;
        int i, k, LEN;
        long D;
        unsigned init_va[1] = {2}, init_vb[1], A;
        giant  n, gva, gvb, z;
        gwnum  wva, wvb;
        gwhandle *gwdata;

        if ( mpz_sizeinbase ( a, 2 ) < 1000 ) return ( mpz_probab_prime_p ( a, 0 ) > 0 );

        if ( mpz_even_p ( a ) ) return ( 0 );

        mpz_get_str( string, 10, a );
        LEN = strlen ( string ); LEN = ( LEN >> 2 ) + 8;
        n = newgiant ( LEN );
        gva = newgiant ( LEN );  gvb = newgiant ( LEN );
        z = newgiant ( LEN );

        ctog ( string, n );

        gwdata = (gwhandle*) malloc (sizeof ( gwhandle ) );
        gwinit ( gwdata );

        gwsetup_general_mod_giant( gwdata, n );
        wva = gwalloc ( gwdata ); wvb = gwalloc ( gwdata );

        init_vb[0] = 2;
        binarytogw ( gwdata, init_vb, 1, wvb );
        for ( i = mpz_sizeinbase ( a, 2 ) - 2; i >= 0; i-- ) {
                gwsquare ( gwdata, wvb );
                if ( mpz_tstbit ( a, i ) ) {
                        gwsmallmul ( gwdata, 2.0, wvb );
                }
        }
        gwtogiant ( gwdata, wvb, gvb );
        modg ( n, gvb );
        itog ( 2, z );
        subg ( z, gvb );
        if ( !(isZero ( gvb )) ) return ( 0 );

        if ( mpz_perfect_square_p ( a ) ) return ( 0 );

        A = 2;
        do {
                A++;
                if ( A == 0 ) return ( 0 );
                D = A * A - 4;
                k = mpz_ui_kronecker ( D, a );
                if ( k == 0 ) return ( 0 );
        } while ( k == 1 );

        init_vb[0] = A;
        fA = (float) A;
        binarytogw ( gwdata, init_va, 1, wva );
        binarytogw ( gwdata, init_vb, 1, wvb );
        for ( i = mpz_sizeinbase ( a, 2 ) - 1; i >= 0; i-- ) {
                gwfft ( gwdata, wva, wva );
                gwfft ( gwdata, wvb, wvb );
                if ( mpz_tstbit ( a, i ) ) {
                        gwfftfftmul ( gwdata, wvb , wva, wva );
                        gwsmalladd ( gwdata, -fA, wva );
                        gwfftfftmul ( gwdata, wvb , wvb, wvb );
                        gwsmalladd ( gwdata, -2.0, wvb );
                } else {
                        gwfftfftmul ( gwdata, wvb , wva, wvb );
                        gwsmalladd ( gwdata, -fA, wvb );
                        gwfftfftmul ( gwdata, wva , wva, wva );
                        gwsmalladd ( gwdata, -2.0, wva );
                }
        }
        gwtogiant ( gwdata, wva, gva );
        gwtogiant ( gwdata, wvb, gvb );
        modg ( n, gva );
        itog ( A, z );
        subg ( z, gva );
        modg ( n, gvb );
        itog ( 2, z );
        subg ( z, gvb );
        gwdone ( gwdata );
        if ( isZero ( gva ) && isZero ( gvb ) )
                return ( 1 );
        else
                return ( 0 );

}

Then run these commands

Code:

autoreconf -vfi  
make distclean
./configure --enable-mpi --enable-shared=no
make clean
make
sudo make install

Sorry about there being no conditional GWNUM compile yet.

The 26000 threshold might be able to be reduced.

Edit: I have changed the threshold in the above function to 1000

Luminescence · 2022-05-31, 04:07

Quote:

Originally Posted by paulunderwood

Edit: I have changed the threshold in the above function to 1000

That's what I have noticed as well. At any bitlevel the gwnum test was atleast as fast or faster than the mpz_probab_prime_p function. It was much faster on longer steps. (On a dual Xeon Platinum 8375C box for reference)

I will dig out the number and specific step when I am back from work, but I had FastECPP take different steps with gwnum on low bitlevels.
The math done here is a mystery to me. What PRP test is actually stronger? The one with gwnum or the GMP function?

paulunderwood · 2022-05-31, 05:27

Quote:

Originally Posted by Luminescence

That's what I have noticed as well. At any bitlevel the gwnum test was atleast as fast or faster than the mpz_probab_prime_p function. It was much faster on longer steps. (On a dual Xeon Platinum 8375C box for reference)

I will dig out the number and specific step when I am back from work, but I had FastECPP take different steps with gwnum on low bitlevels.
The math done here is a mystery to me. What PRP test is actually stronger? The one with gwnum or the GMP function?

Both my gwnum function and GMP's mpz_probab_prime_p are some variation of Fermat+Lucas tests (a la BPSW).The Fermat PRP test could be made "strong" as could the Lucas component. Think of a Lucas PRP test as a two dimensional Fermat PRP test. The root of Lucas numbers goes back to Fibonacci and before.

mathwiz · 2022-05-31, 07:00

Quote:

Originally Posted by paulunderwood

Both my gwnum function and GMP's mpz_probab_prime_p are some variation of Fermat+Lucas tests (a la BPSW).The Fermat PRP test could be made "strong" as could the Lucas component. Think of a Lucas PRP test as a two dimensional Fermat PRP test. The root of Lucas numbers goes back to Fibonacci and before.

One question: mprime does a bunch of internal error checking (Jacobi / gw_square_carefully / etc) to be resilient to hardware errors. Is that done here?

paulunderwood · 2022-05-31, 07:39

Quote:

Originally Posted by mathwiz

One question: mprime does a bunch of internal error checking (Jacobi / gw_square_carefully / etc) to be resilient to hardware errors. Is that done here?

A good question. Unless some beats me to it I will endeavour to:

use gw_square_carefully on intro and outro bits -- 50 bits each should be enough.
make the Lucas test more robust by not using "illegal" functions like gwfftfftmul ( gwdata, wva, wva, wva );
turn the Fermat test into a strong PRP test.
Anything else?

paulunderwood · 2022-05-31, 12:32

Quote:

Originally Posted by paulunderwood

Edit: I have changed the threshold in the above function to 1000

That figure was ambitious! 4096 is realistic, maybe not quite optimal.

Here is the latest cm_nt_is_prime function (also attached), Follow previous instructions and use this code.

Code:

int cm_nt_is_prime (mpz_t a) {

        char string[100000];
        float fA;
        int i, j, k, LEN, r, A;
        long D;
        unsigned init_va[1] = {2}, init_vb[1] = {2};
        giant  n, gva, gvb;
        gwnum  wva, wvb, wtmp;
        gwhandle *gwdata;
        mpz_t gmp_var;

	     if ( mpz_sizeinbase ( a, 2 ) < 4096 ) return ( mpz_probab_prime_p ( a, 0 ) > 0 );
        
        if ( mpz_even_p ( a ) ) return ( 0 );

        mpz_get_str( string, 10, a );
        LEN = strlen ( string ); LEN = ( LEN >> 2 ) + 8;
        n = newgiant ( LEN ); gva = newgiant ( LEN );  gvb = newgiant ( LEN );

        ctog ( string, n );

        gwdata = (gwhandle*) malloc (sizeof ( gwhandle ) );
        gwinit ( gwdata );

        gwsetup_general_mod_giant( gwdata, n );
        wva = gwalloc ( gwdata ); wvb = gwalloc ( gwdata ); wtmp = gwalloc ( gwdata );
        // Base 2 strong Fermat PRP test
        mpz_init ( gmp_var );
        mpz_sub_ui ( gmp_var, a, 1);
        r = 0;
        do {
           r++;
           mpz_tdiv_q_2exp ( gmp_var, gmp_var, 1);
        } while mpz_even_p ( gmp_var );
        mpz_sub_ui ( gmp_var, a, 1); 
        binarytogw ( gwdata, init_vb, 1, wvb );
        j = mpz_sizeinbase ( gmp_var, 2 ) - 2;
        k = j - 50;
        for ( i = j; i > k; i-- ) {
                gwsquare2_carefully ( gwdata, wvb, wvb );
                if ( mpz_tstbit ( gmp_var, i ) ) {
                        gwsmallmul ( gwdata, 2.0, wvb );
                }
        }
        for ( i = k; i > 50; i-- ) {
                gwsquare2 ( gwdata, wvb, wvb );
                if ( mpz_tstbit ( gmp_var, i ) ) {
                        gwsmallmul ( gwdata, 2.0, wvb );
                }                
        }
        for ( i = 50; i >= r; i-- ) {
                gwsquare2_carefully ( gwdata, wvb, wvb );
                if ( mpz_tstbit ( gmp_var, i ) ) {
                        gwsmallmul ( gwdata, 2.0, wvb );
                }                
        }
        mpz_clear ( gmp_var );
        gwtogiant ( gwdata, wvb, gvb );
        modg ( n, gvb );
        if ( !( isone ( gvb ) ) ) {
                itog ( 1, gva );
                for ( i = r; i > 0; i-- ) {
                        gwtogiant ( gwdata, wvb, gvb );
                        addg ( gva, gvb );
                        modg ( n, gvb );
                        if ( isZero ( gvb ) ) break;
                        gwsquare2_carefully ( gwdata, wvb, wvb );
                }
                if ( !isZero ( gvb ) ) return ( 0 );
        }

        if ( mpz_perfect_square_p ( a ) ) return ( 0 );

        // Find strong discriminant
        A = 2;
        do {
                A++;
                if ( A == 0 ) return ( 0 );
                D = A * A - 4;
                k = mpz_ui_kronecker ( D, a );
                if ( k == 0 ) return ( 0 );
        } while ( k == 1 );
        
        // Lucas PRP test over x^2-A*x+1
        init_vb[0] = A;
        fA = (float) A;
        binarytogw ( gwdata, init_va, 1, wva );
        binarytogw ( gwdata, init_vb, 1, wvb );
        j++;
        k = j - 50;
        for ( i = j; i > k; i-- ) {
                if ( mpz_tstbit ( a, i ) ) {
                        gwcopy ( gwdata, wvb, wtmp );
                        gwmul_carefully ( gwdata, wvb, wva );
                        gwsmalladd ( gwdata, -fA, wva );
                        gwsquare2_carefully ( gwdata, wtmp, wvb );
                        gwsmalladd ( gwdata, -2.0, wvb );
                } else {
                        gwcopy ( gwdata, wva, wtmp );
                        gwmul_carefully ( gwdata, wva , wvb );
                        gwsmalladd ( gwdata, -fA, wvb );
                        gwsquare2_carefully ( gwdata, wtmp, wva );
                        gwsmalladd ( gwdata, -2.0, wva );
                }
        }
        for ( i = k; i > 50; i-- ) {
                if ( mpz_tstbit ( a, i ) ) {
                        gwcopy ( gwdata, wvb, wtmp );
                        gwmul ( gwdata, wvb, wva );
                        gwsmalladd ( gwdata, -fA, wva );
                        gwsquare2 ( gwdata, wtmp, wvb );
                        gwsmalladd ( gwdata, -2.0, wvb );
                } else {
                        gwcopy ( gwdata, wva, wtmp );
                        gwmul ( gwdata, wva, wvb );
                        gwsmalladd ( gwdata, -fA, wvb );
                        gwsquare2 ( gwdata, wtmp, wva );
                        gwsmalladd ( gwdata, -2.0, wva );
                }
        }
        for ( i = 50; i >= 0; i-- ) {
                if ( mpz_tstbit ( a, i ) ) {
                        gwcopy ( gwdata, wvb, wtmp );
                        gwmul_carefully ( gwdata, wvb, wva );
                        gwsmalladd ( gwdata, -fA, wva );
                        gwsquare2_carefully ( gwdata, wtmp, wvb );
                        gwsmalladd ( gwdata, -2.0, wvb );
                } else {
                        gwcopy ( gwdata, wva, wtmp );
                        gwmul_carefully ( gwdata, wva, wvb );
                        gwsmalladd ( gwdata, -fA, wvb );
                        gwsquare2_carefully ( gwdata, wtmp, wva );
                        gwsmalladd ( gwdata, -2.0, wva );
                }
        }
        gwsmalladd ( gwdata, -fA, wva );
        gwtogiant ( gwdata, wva, gva );
        gwtogiant ( gwdata, wvb, gvb );
        gwdone ( gwdata );
        modg ( n, gva );
        modg ( n, gvb );
        if ( isZero ( gva ) && istwo ( gvb ) )
                return ( 1 );
        else
                return ( 0 );
}

This need eyeballing and testing.

R. Gerbicz · 2022-05-31, 13:51

Quote:

Originally Posted by mathwiz

One question: mprime does a bunch of internal error checking (Jacobi / gw_square_carefully / etc) to be resilient to hardware errors. Is that done here?

There is no Jacobi check on general numbers, we have error check but that is somewhat slow on these general numbers and you really don't need any check on these tiny numbers.

Quote:

Originally Posted by paulunderwood

A good question. Unless some beats me to it I will endeavour to:

use gw_square_carefully on intro and outro bits -- 50 bits each should be enough.
make the Lucas test more robust by not using "illegal" functions like gwfftfftmul ( gwdata, wva, wva, wva );
turn the Fermat test into a strong PRP test.
Anything else?

Any fancy Lucas test is also pointless, you would never get a 200 digits composite prp number. And this slows down the project a lot/somewhat: if you find a prp number after the first Miller Rabin test, then in roughly at the same time all other cores are also done (with their own prp test), so you'd start a Lucas test on a single core and all other cores are doing nothing, and it is independent from the code: you can't move to the next N value, because you haven't done the test, and even if you'd do any other computation (on this step) then you would need to throw it, because you have already found a (prp) prime.
This is also problematic in my code, because from gmp-6.2.1(?) it is doing a mpz_stronglucas after a successful Miller Rabin hit. I will change this.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
For which types of primes is GPU primality test software available?	bur	GPU Computing	6	2020-08-28 06:20
Fastest software for Mersenne primality test?	JonathanM	Information & Answers	25	2020-06-16 02:47
APR-CL as primality proof	f1pokerspeed	FactorDB	14	2014-01-09 21:06
Proof of Primality Test for Fermat Numbers	princeps	Math	15	2012-04-02 21:49
PRIMALITY PROOF for Wagstaff numbers!	AntonVrba	Math	96	2009-02-25 10:37

2022-05-29, 18:44	#136
frmky Jul 2003 So Cal 5×479 Posts	Please make use of gwnum optional with an option to configure. The current code supports Arm and we shouldn't break that.