From c3216aec66f9756265edc370163f42f0525697c9 Mon Sep 17 00:00:00 2001 From: czurnieden Date: Thu, 4 Apr 2019 22:34:17 +0200 Subject: [PATCH] A prime sieve --- demo/test.c | 178 +++++++++++++++++++++++++- doc/bn.tex | 170 ++++++++++++++++++++++++ helper.pl | 2 +- libtommath_VS2008.vcproj | 60 +++++++++ makefile | 34 ++--- makefile.mingw | 35 ++--- makefile.msvc | 35 ++--- makefile.shared | 34 ++--- makefile.unix | 34 ++--- mp_error_to_string.c | 2 + mp_is_small_prime.c | 96 ++++++++++++++ mp_next_small_prime.c | 35 +++++ mp_prec_small_prime.c | 38 ++++++ mp_sieve_clear.c | 26 ++++ mp_sieve_init.c | 21 +++ s_mp_eratosthenes.c | 25 ++++ s_mp_eratosthenes_init.c | 28 ++++ s_mp_eratosthenes_segment.c | 37 ++++++ s_mp_eratosthenes_segment_init.c | 36 ++++++ s_mp_init_single_segment_with_start.c | 34 +++++ s_mp_isqrt.c | 44 +++++++ s_mp_sieve_clear.c | 17 +++ s_mp_sieve_get.c | 16 +++ s_mp_sieve_nextset.c | 18 +++ s_mp_sieve_setall.c | 33 +++++ sources.cmake | 15 +++ tommath.def | 5 + tommath.h | 67 +++++++++- tommath_class.h | 77 +++++++++++ tommath_private.h | 32 +++++ 30 files changed, 1200 insertions(+), 84 deletions(-) create mode 100644 mp_is_small_prime.c create mode 100644 mp_next_small_prime.c create mode 100644 mp_prec_small_prime.c create mode 100644 mp_sieve_clear.c create mode 100644 mp_sieve_init.c create mode 100644 s_mp_eratosthenes.c create mode 100644 s_mp_eratosthenes_init.c create mode 100644 s_mp_eratosthenes_segment.c create mode 100644 s_mp_eratosthenes_segment_init.c create mode 100644 s_mp_init_single_segment_with_start.c create mode 100644 s_mp_isqrt.c create mode 100644 s_mp_sieve_clear.c create mode 100644 s_mp_sieve_get.c create mode 100644 s_mp_sieve_nextset.c create mode 100644 s_mp_sieve_setall.c diff --git a/demo/test.c b/demo/test.c index 16fef5570..e1ca74a52 100644 --- a/demo/test.c +++ b/demo/test.c @@ -811,7 +811,7 @@ static int test_mp_prime_rand(void) /* test for size */ for (ix = 10; ix < 128; ix++) { - printf("Testing (not safe-prime): %9d bits \n", ix); + printf("\rTesting (not safe-prime): %9d bits ", ix); fflush(stdout); DO(mp_prime_rand(&a, 8, ix, (rand_int() & 1) ? 0 : MP_PRIME_2MSB_ON)); EXPECT(mp_count_bits(&a) == ix); @@ -1529,6 +1529,177 @@ static int test_mp_decr(void) return EXIT_FAILURE; } +static int test_mp_is_small_prime(void) +{ + mp_sieve sieve; + mp_err e; + int i, test_size; + + const mp_sieve_prime to_test[] = { + 52, 137, 153, 179, 6, 153, 53, 132, 150, 65, + 27414, 36339, 36155, 11067, 52060, 5741, + 29755, 2698, 52572, 13053, 9375, 47241,39626, + 207423, 128857, 37419, 141696, 189465, + 41503, 127370, 91673, 8473, 479142414, 465566339, + 961126169, 1057886067, 1222702060, 1017450741, + 1019879755, 72282698, 2048787577, 2058368053 + }; + const bool tested[] = { + false, true, false, true, false, false, true, false, false, + false, false, false, false, false, false, true, false, false, + false, false, false, false, false, false, true, false, false, + false, false, false, true, false, false, false, true, false, + false, false, false, false, true, false + }; + bool result; + + mp_sieve_init(&sieve); + + test_size = (int)(sizeof(to_test)/sizeof(mp_sieve_prime)); + + for (i = 0; i < test_size; i++) { + if ((e = mp_is_small_prime(to_test[i], &result, &sieve)) != MP_OKAY) { + fprintf(stderr,"mp_is_small_prime failed: \"%s\"\n",mp_error_to_string(e)); + goto LTM_ERR; + } + if (result != tested[i]) { + fprintf(stderr,"mp_is_small_prime failed for %lu. Said %lu but is %lu \n", + (unsigned long)to_test[i], (unsigned long)result, (unsigned long)tested[i]); + goto LTM_ERR; + } + } + mp_sieve_clear(&sieve); + return EXIT_SUCCESS; +LTM_ERR: + mp_sieve_clear(&sieve); + return EXIT_FAILURE; +} + +static int test_mp_next_small_prime(void) +{ + mp_sieve sieve; + mp_sieve_prime ret = 0lu, p; + mp_int primesum, t; + mp_err e; + int i, test_size; + + /* Jumping wildly to and fro */ + const mp_sieve_prime to_test[] = { + 52, 137, 153, 179, 6, 153, 53, 132, 150, 65, + 27414, 36339, 36155, 11067, 52060, 5741, + 29755, 2698, 52572, 13053, 9375, 47241, + 39626, 207423, 128857, 37419, 141696, 189465, + 41503, 127370, 91673, 8473, 479142414, 465566339, + 961126169, 1057886067, 1222702060, 1017450741, + 1019879755, 72282698, 2048787577, 2058368053 + }; + const mp_sieve_prime tested[] = { + 53, 137, 157, 179, 7, 157, 53, 137, 151, 67, + 27427, 36341, 36161, 11069, 52067, 5741, + 29759, 2699, 52579, 13063, 9377, 47251, + 39631, 207433, 128857, 37423, 141697, 189467, + 41507, 127373, 91673, 8501, 479142427, 465566393, + 961126169, 1057886083, 1222702081, 1017450823, + 1019879761, 72282701, 2048787577, 2058368113 + }; + const char *primesum_32 = "202259606268580"; + + mp_sieve_init(&sieve); + + test_size = (int)(sizeof(to_test)/sizeof(mp_sieve_prime)); + + for (i = 0; i < test_size; i++) { + if ((e = mp_next_small_prime(to_test[i], &ret, &sieve)) != MP_OKAY) { + fprintf(stderr,"mp_next_small_prime failed with \"%s\" at index %d\n", + mp_error_to_string(e), i); + goto LBL_ERR; + } + if (ret != tested[i]) { + fprintf(stderr,"mp_next_small_prime failed for %lu. Said %lu but is %lu \n", + (unsigned long)to_test[i], (unsigned long)ret, (unsigned long)tested[i]); + goto LBL_ERR; + } + } + + DOR(mp_init_multi(&primesum, &t, NULL)); + + for (p = 4293918720lu; ret < (mp_sieve_prime)MP_SIEVE_BIGGEST_PRIME;) { + DO(mp_next_small_prime(p, &ret, &sieve)); + p = ret + 1u; +#ifdef MP_64BIT + mp_set_u64(&t, ret); +#else + mp_set_u32(&t, ret); +#endif + DO(mp_add(&primesum, &t, &primesum)); + } + printf("Primesum computed: "); + DO(mp_fwrite(&primesum, 10, stdout)); + puts(""); + DO(mp_read_radix(&t, primesum_32, 10)); + EXPECT(mp_cmp(&primesum, &t) == MP_EQ); + + mp_sieve_clear(&sieve); + mp_clear_multi(&primesum, &t, NULL); + return EXIT_SUCCESS; +LBL_ERR: + mp_clear_multi(&primesum, &t, NULL); + mp_sieve_clear(&sieve); + return EXIT_FAILURE; +} + +static int test_mp_prec_small_prime(void) +{ + mp_sieve sieve; + mp_sieve_prime ret; + mp_err e; + int i, test_size; + + const mp_sieve_prime to_test[] = { + 52, 137, 153, 179, 6, 153, 53, 132, 150, 65, + 27414, 36339, 36155, 11067, 52060, 5741, + 29755, 2698, 52572, 13053, 9375, 47241, + 39626, 207423, 128857, 37419, 141696, 189465, + 41503, 127370, 91673, 8473, 479142414, 465566339, + 961126169, 1057886067, 1222702060, 1017450741, + 1019879755, 72282698, 2048787577, 2058368053 + }; + const mp_sieve_prime tested[] = { + 47, 137, 151, 179, 5, 151, 53, 131, 149, 61, + 27409, 36319, 36151, 11059, 52057, 5741, + 29753, 2693, 52571, 13049, 9371, 47237, + 39623, 207409, 128857, 37409, 141689, 189463, + 41491, 127363, 91673, 8467, 479142413, 465566323, + 961126169, 1057886029, 1222702051, 1017450739, + 1019879717, 72282697, 2048787577, 2058368051 + }; + + mp_sieve_init(&sieve); + + test_size = (int)(sizeof(to_test)/sizeof(mp_sieve_prime)); + + for (i = 0; i < test_size; i++) { + if ((e = mp_prec_small_prime(to_test[i], &ret, &sieve)) != MP_OKAY) { + fprintf(stderr,"mp_prec_small_prime failed with \"%s\" at index %d\n", + mp_error_to_string(e), i); + goto LTM_ERR; + } + if (ret != tested[i]) { + fprintf(stderr,"mp_prec_small_prime failed for %lu. Said %lu but is %lu \n", + (unsigned long)to_test[i], (unsigned long)ret, (unsigned long)tested[i]); + goto LTM_ERR; + } + } + + mp_sieve_clear(&sieve); + return EXIT_SUCCESS; +LTM_ERR: + mp_sieve_clear(&sieve); + return EXIT_FAILURE; +} + + + /* Cannot test mp_exp(_d) without mp_root_n and vice versa. So one of the two has to be tested from scratch. @@ -2240,6 +2411,11 @@ static int unit_tests(int argc, char **argv) T1(mp_prime_next_prime, MP_PRIME_NEXT_PRIME), T1(mp_prime_rand, MP_PRIME_RAND), T1(mp_rand, MP_RAND), + + T1(mp_is_small_prime, MP_IS_SMALL_PRIME), + T1(mp_next_small_prime, MP_NEXT_SMALL_PRIME), + T1(mp_prec_small_prime, MP_PREC_SMALL_PRIME), + T1(mp_read_radix, MP_READ_RADIX), T1(mp_read_write_ubin, MP_TO_UBIN), T1(mp_read_write_sbin, MP_TO_SBIN), diff --git a/doc/bn.tex b/doc/bn.tex index 566b3be32..6913cbf02 100644 --- a/doc/bn.tex +++ b/doc/bn.tex @@ -2347,6 +2347,173 @@ \section{Random Primes} \label{fig:primeopts} \end{figure} + +\section{Prime Sieve} +The prime sieve is implemented as a simple segmented Sieve of Eratosthenes. + +This library needs some small sequential amounts for divisibility tests starting at two and +some random small primes for the Miller--Rabin test. That means we need a small base sieve +for quick results with a cold start and small segments to get random small primes fast. + +The base sieve holds odd values only and even with a size of \texttt{4096} bytes it is +small enough to get build quickly, in about 50 $\mu$sec on the author's old machine. + +The choice of the size for the individual segments varies more in the results. See table + \ref{table:sievebenchmarks} for some data. Size must be of the form $-1 + 2^n$. +Default size is \texttt{0xFFF}. It might be of interest that the biggest prime gap +below $2^{32}$ is $335$. + +\begin{table}[h] + \begin{center} + \begin{tabular}{r c l} + \textbf{Segment Size (bits)} & \textbf{Random Access in $\mu$sec} & \textbf{Full Primesum} \\ + \texttt{ 0x1F} & (60) & - \\ + \texttt{ 0x3F} & (75) & - \\ + \texttt{ 0x7F} & 63 & - \\ + \texttt{ 0xFF} & 70 & 26 min \\ + \texttt{ 0x1FF} & 73 & 13 min \\ + \texttt{ 0x3FF} & 75 & 7 min 17 sec \\ + \texttt{ 0x7FF} & 85 & 4 min 10 sec \\ + \texttt{ 0xFFF} & 100 & 2 min 52 sec \\ + \texttt{ 0x1FFF} & 120 & 1 min 45 sec \\ + \texttt{ 0x3FFF} & 145 & 1 min 19 sec \\ + \texttt{ 0x7FFF} & 200 & 1 min 04 sec \\ + \texttt{ 0xFFFF} & 350 & 0 min 57 sec \\ + \texttt{ 0xFFFFFFFF} & 300 & 0 min 35 sec + \end{tabular} + \caption{Average access times (warm) with the default base sieve and \texttt{MP-64BIT}} \label{table:sievebenchmarks} + \end{center} +\end{table} + +The default sizes are in \texttt{tommath\_private.h}: \texttt{MP\_SIEVE\_PRIME\_MAX\_SQRT} is +the size of the base sieve and \texttt{MP\_SIEVE\_RANGE\_A\_B} is the size of the segment. + +\subsection{Initialization and Clearing} +Initializing. It cannot fail because it only sets some default values. Memory is allocated later according to needs. +\index{mp\_sieve\_init} +\begin{alltt} +void mp_sieve_init(mp_sieve *sieve); +\end{alltt} +The function \texttt{mp\_sieve\_init} is equivalent to +\begin{alltt} +mp_sieve sieve = {NULL, NULL, 0}; +\end{alltt} + +Free the memory used by the sieve with +\index{mp\_sieve\_clear} +\begin{alltt} +void mp_sieve_clear(mp_sieve *sieve); +\end{alltt} + +\subsection{Primality Test of Small Numbers} +Individual small numbers can be tested for primality with +\index{mp\_is\_small\_prime} +\begin{alltt} +int mp_is_small_prime(mp_sieve_prime n, bool *result, + mp_sieve *sieve); +\end{alltt} +It sets the variable \texttt{result} to \texttt{false} if the given number in \texttt{n} is a composite and +\texttt{true} if it is a prime. + +It will return \texttt{MP\_SIEVE\_MAX\_REACHED} to flag the content of \texttt{result} as the last valid one. + +The implementation of this function does all the heavy lifting--the building of the base sieve and the segment +if one is necessary. The base sieve stays, so this function can be used to ``warm up'' the sieve and, if +\texttt{n} is slightly larger than the upper limit of the base sieve, ``warm up'' the first segment, too. + +\subsection{Find Adjacent Primes} +To find the prime bigger than a number \texttt{n} use +\index{mp\_next\_small\_prime} +\begin{alltt} +int mp_next_small_prime(mp_sieve_prime n, mp_sieve_prime *result, + mp_sieve *sieve); +\end{alltt} +and to find the one smaller than \texttt{n} +\index{mp\_prec\_small\_prime} +\begin{alltt} +int mp_prec_small_prime(mp_sieve_prime n, mp_sieve_prime *result, + mp_sieve *sieve); +\end{alltt} +Both functions return \texttt{n} if \texttt{n} is prime. Use \texttt{n + 1} to get +the prime after \texttt{n} if \texttt{n} is prime and \texttt{n - 1} to get the +the prime preceding \texttt{n} if \texttt{n} is prime. + +\subsection{Useful Constants} +\begin{description} +\item[\texttt{MP\_SIEVE\_BIGGEST\_PRIME}] \texttt{read-only} The biggest prime the sieve can offer. +It is $4\,294\,967\,291$ for all \texttt{MP\_xBIT}. + +\item[\texttt{mp\_sieve\_prime}] \texttt{read-only} The basic type for the numbers in the sieve. It + is \texttt{uint32\_t} for \texttt{MP\_16BIT}, \texttt{MP\_32BIT}; and + \texttt{uint64\_t} for \texttt{MP\_64BIT}.. + +\item[\texttt{MP\_SIEVE\_PR\_UINT}] Choses the correct macro from \texttt{inttypes.h} to print a\\ + \texttt{mp\_sieve\_prime}. The header \texttt{inttypes.h} must be included before\\ + \texttt{tommath.h} for it to work. +\end{description} +\subsection{Examples}\label{sec:spnexamples} +\subsubsection{Initialization and Clearing} +Using a sieve follows the same procedure as using a bigint: +\begin{alltt} +/* Declaration */ +mp_sieve sieve; + +/* + Initialization. + Cannot fail, only sets a handful of default values. + Memory allocation is done in the library itself + according to needs. + */ +mp_sieve_init(&sieve); + +/* use the sieve */ + +/* Clean up */ +mp_sieve_clear(&sieve); +\end{alltt} +\subsubsection{Primality Test} +A small program to test the input of a small number for primality. +\begin{alltt} +#include +#include +#include +/*inttypes.h is needed for printing and must be included before tommath.h*/ +#include +#include "tommath.h" +int main(int argc, char **argv) +{ + mp_sieve_prime number; + mp_sieve sieve; + int e; + /* variable holding the result of mp_is_small_prime */ + mp_sieve_prime result; + + if (argc != 2) { + fprintf(stderr,"Usage %s number\textbackslash{}n", argv[0]); + exit(EXIT_FAILURE); + } + number = (mp_sieve_prime) strtoul(argv[1],NULL, 10); + if (errno == ERANGE) { + fprintf(stderr,"strtoul(l) failed: input out of range\textbackslash{}n"); + goto LTM_ERR; + } + mp_sieve_init(&sieve); + if ((e = mp_is_small_prime(number, &result, &sieve)) != MP_OKAY) { + fprintf(stderr,"mp_is_small_prime failed: \textbackslash{}"%s\textbackslash{}"\textbackslash{}n", + mp_error_to_string(e)); + goto LTM_ERR; + } + printf("The number %" LTM_SIEVE_PR_UINT " is %s prime\textbackslash{}n", + number,(result)?"":"not"); + + mp_sieve_clear(&sieve); + exit(EXIT_SUCCESS); +LTM_ERR: + mp_sieve_clear(&sieve); + exit(EXIT_FAILURE); +} +\end{alltt} + \chapter{Random Number Generation} \section{PRNG} \index{mp\_rand} @@ -2417,6 +2584,9 @@ \subsection{To ASCII} mp_err mp_fwrite(const mp_int *a, int radix, FILE *stream); \end{alltt} + + + \subsection{From ASCII} \index{mp\_read\_radix} \begin{alltt} diff --git a/helper.pl b/helper.pl index e5ad27ac5..c7c82c29b 100755 --- a/helper.pl +++ b/helper.pl @@ -58,7 +58,7 @@ sub check_source { push @{$troubles->{unwanted_free}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bfree\s*\(/; # and we probably want to also avoid the following push @{$troubles->{unwanted_memcpy}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bmemcpy\s*\(/ && $file !~ /s_mp_copy_digs.c/; - push @{$troubles->{unwanted_memset}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bmemset\s*\(/ && $file !~ /s_mp_zero_buf.c/ && $file !~ /s_mp_zero_digs.c/; + push @{$troubles->{unwanted_memset}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bmemset\s*\(/ && $file !~ /s_mp_zero_buf.c/ && $file !~ /s_mp_zero_digs.c/ && $file !~ /s_mp_sieve_setall.c/; push @{$troubles->{unwanted_memmove}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bmemmove\s*\(/; push @{$troubles->{unwanted_memcmp}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bmemcmp\s*\(/; push @{$troubles->{unwanted_strcmp}}, $lineno if $file =~ /^[^\/]+\.c$/ && $l =~ /\bstrcmp\s*\(/; diff --git a/libtommath_VS2008.vcproj b/libtommath_VS2008.vcproj index e27aa9860..072c54ceb 100644 --- a/libtommath_VS2008.vcproj +++ b/libtommath_VS2008.vcproj @@ -532,6 +532,10 @@ RelativePath="mp_invmod.c" > + + @@ -596,6 +600,10 @@ RelativePath="mp_neg.c" > + + @@ -608,6 +616,10 @@ RelativePath="mp_pack_count.c" > + + @@ -740,6 +752,14 @@ RelativePath="mp_shrink.c" > + + + + @@ -820,6 +840,22 @@ RelativePath="s_mp_div_small.c" > + + + + + + + + @@ -832,6 +868,10 @@ RelativePath="s_mp_get_bit.c" > + + @@ -840,6 +880,10 @@ RelativePath="s_mp_invmod_odd.c" > + + @@ -904,6 +948,22 @@ RelativePath="s_mp_rand_platform.c" > + + + + + + + + diff --git a/makefile b/makefile index 666274378..0dc77b2b7 100644 --- a/makefile +++ b/makefile @@ -33,23 +33,27 @@ mp_error_to_string.o mp_exch.o mp_expt_n.o mp_exptmod.o mp_exteuclid.o mp_fread. mp_from_ubin.o mp_fwrite.o mp_gcd.o mp_get_double.o mp_get_i32.o mp_get_i64.o mp_get_l.o mp_get_mag_u32.o \ mp_get_mag_u64.o mp_get_mag_ul.o mp_grow.o mp_hash.o mp_init.o mp_init_copy.o mp_init_i32.o mp_init_i64.o \ mp_init_l.o mp_init_multi.o mp_init_set.o mp_init_size.o mp_init_u32.o mp_init_u64.o mp_init_ul.o \ -mp_invmod.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o mp_mod_2d.o \ -mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o mp_mul_2.o \ -mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_or.o mp_pack.o mp_pack_count.o mp_prime_fermat.o \ -mp_prime_frobenius_underwood.o mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o \ -mp_prime_rabin_miller_trials.o mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o \ -mp_radix_size_overestimate.o mp_rand.o mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o \ -mp_reduce_2k_l.o mp_reduce_2k_setup.o mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o \ -mp_reduce_setup.o mp_root_n.o mp_rshd.o mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o \ -mp_set_l.o mp_set_u32.o mp_set_u64.o mp_set_ul.o mp_shrink.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ +mp_invmod.o mp_is_small_prime.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o \ +mp_mod_2d.o mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o \ +mp_mul_2.o mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_next_small_prime.o mp_or.o mp_pack.o \ +mp_pack_count.o mp_prec_small_prime.o mp_prime_fermat.o mp_prime_frobenius_underwood.o \ +mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o mp_prime_rabin_miller_trials.o \ +mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o mp_radix_size_overestimate.o mp_rand.o \ +mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o mp_reduce_2k_l.o mp_reduce_2k_setup.o \ +mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o mp_reduce_setup.o mp_root_n.o mp_rshd.o \ +mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o mp_set_l.o mp_set_u32.o mp_set_u64.o \ +mp_set_ul.o mp_shrink.o mp_sieve_clear.o mp_sieve_init.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ mp_sqrtmod_prime.o mp_sub.o mp_sub_d.o mp_submod.o mp_to_radix.o mp_to_sbin.o mp_to_ubin.o mp_ubin_size.o \ mp_unpack.o mp_xor.o mp_zero.o s_mp_add.o s_mp_copy_digs.o s_mp_div_3.o s_mp_div_recursive.o \ -s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_exptmod_fast.o s_mp_get_bit.o s_mp_invmod.o \ -s_mp_invmod_odd.o s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o \ -s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \ -s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o \ -s_mp_radix_size_overestimate.o s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o \ -s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o s_mp_zero_digs.o +s_mp_div_school.o s_mp_div_small.o s_mp_eratosthenes.o s_mp_eratosthenes_init.o \ +s_mp_eratosthenes_segment.o s_mp_eratosthenes_segment_init.o s_mp_exptmod.o s_mp_exptmod_fast.o \ +s_mp_get_bit.o s_mp_init_single_segment_with_start.o s_mp_invmod.o s_mp_invmod_odd.o s_mp_isqrt.o \ +s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o s_mp_mul_balance.o \ +s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o s_mp_mul_toom.o \ +s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_radix_size_overestimate.o \ +s_mp_rand_platform.o s_mp_sieve_clear.o s_mp_sieve_get.o s_mp_sieve_nextset.o s_mp_sieve_setall.o \ +s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \ +s_mp_zero_digs.o #END_INS diff --git a/makefile.mingw b/makefile.mingw index 7aeb564e7..e59d6fea9 100644 --- a/makefile.mingw +++ b/makefile.mingw @@ -35,27 +35,32 @@ mp_error_to_string.o mp_exch.o mp_expt_n.o mp_exptmod.o mp_exteuclid.o mp_fread. mp_from_ubin.o mp_fwrite.o mp_gcd.o mp_get_double.o mp_get_i32.o mp_get_i64.o mp_get_l.o mp_get_mag_u32.o \ mp_get_mag_u64.o mp_get_mag_ul.o mp_grow.o mp_hash.o mp_init.o mp_init_copy.o mp_init_i32.o mp_init_i64.o \ mp_init_l.o mp_init_multi.o mp_init_set.o mp_init_size.o mp_init_u32.o mp_init_u64.o mp_init_ul.o \ -mp_invmod.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o mp_mod_2d.o \ -mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o mp_mul_2.o \ -mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_or.o mp_pack.o mp_pack_count.o mp_prime_fermat.o \ -mp_prime_frobenius_underwood.o mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o \ -mp_prime_rabin_miller_trials.o mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o \ -mp_radix_size_overestimate.o mp_rand.o mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o \ -mp_reduce_2k_l.o mp_reduce_2k_setup.o mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o \ -mp_reduce_setup.o mp_root_n.o mp_rshd.o mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o \ -mp_set_l.o mp_set_u32.o mp_set_u64.o mp_set_ul.o mp_shrink.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ +mp_invmod.o mp_is_small_prime.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o \ +mp_mod_2d.o mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o \ +mp_mul_2.o mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_next_small_prime.o mp_or.o mp_pack.o \ +mp_pack_count.o mp_prec_small_prime.o mp_prime_fermat.o mp_prime_frobenius_underwood.o \ +mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o mp_prime_rabin_miller_trials.o \ +mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o mp_radix_size_overestimate.o mp_rand.o \ +mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o mp_reduce_2k_l.o mp_reduce_2k_setup.o \ +mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o mp_reduce_setup.o mp_root_n.o mp_rshd.o \ +mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o mp_set_l.o mp_set_u32.o mp_set_u64.o \ +mp_set_ul.o mp_shrink.o mp_sieve_clear.o mp_sieve_init.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ mp_sqrtmod_prime.o mp_sub.o mp_sub_d.o mp_submod.o mp_to_radix.o mp_to_sbin.o mp_to_ubin.o mp_ubin_size.o \ mp_unpack.o mp_xor.o mp_zero.o s_mp_add.o s_mp_copy_digs.o s_mp_div_3.o s_mp_div_recursive.o \ -s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_exptmod_fast.o s_mp_get_bit.o s_mp_invmod.o \ -s_mp_invmod_odd.o s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o \ -s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \ -s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o \ -s_mp_radix_size_overestimate.o s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o \ -s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o s_mp_zero_digs.o +s_mp_div_school.o s_mp_div_small.o s_mp_eratosthenes.o s_mp_eratosthenes_init.o \ +s_mp_eratosthenes_segment.o s_mp_eratosthenes_segment_init.o s_mp_exptmod.o s_mp_exptmod_fast.o \ +s_mp_get_bit.o s_mp_init_single_segment_with_start.o s_mp_invmod.o s_mp_invmod_odd.o s_mp_isqrt.o \ +s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o s_mp_mul_balance.o \ +s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o s_mp_mul_toom.o \ +s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_radix_size_overestimate.o \ +s_mp_rand_platform.o s_mp_sieve_clear.o s_mp_sieve_get.o s_mp_sieve_nextset.o s_mp_sieve_setall.o \ +s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \ +s_mp_zero_digs.o HEADERS_PUB=tommath.h HEADERS=tommath_private.h tommath_class.h tommath_superclass.h tommath_cutoffs.h $(HEADERS_PUB) + #The default rule for make builds the libtommath.a library (static) default: $(LIBMAIN_S) diff --git a/makefile.msvc b/makefile.msvc index 59a365e1c..052d4775f 100644 --- a/makefile.msvc +++ b/makefile.msvc @@ -31,27 +31,32 @@ mp_error_to_string.obj mp_exch.obj mp_expt_n.obj mp_exptmod.obj mp_exteuclid.obj mp_from_ubin.obj mp_fwrite.obj mp_gcd.obj mp_get_double.obj mp_get_i32.obj mp_get_i64.obj mp_get_l.obj mp_get_mag_u32.obj \ mp_get_mag_u64.obj mp_get_mag_ul.obj mp_grow.obj mp_hash.obj mp_init.obj mp_init_copy.obj mp_init_i32.obj mp_init_i64.obj \ mp_init_l.obj mp_init_multi.obj mp_init_set.obj mp_init_size.obj mp_init_u32.obj mp_init_u64.obj mp_init_ul.obj \ -mp_invmod.obj mp_is_square.obj mp_kronecker.obj mp_lcm.obj mp_log_n.obj mp_lshd.obj mp_mod.obj mp_mod_2d.obj \ -mp_montgomery_calc_normalization.obj mp_montgomery_reduce.obj mp_montgomery_setup.obj mp_mul.obj mp_mul_2.obj \ -mp_mul_2d.obj mp_mul_d.obj mp_mulmod.obj mp_neg.obj mp_or.obj mp_pack.obj mp_pack_count.obj mp_prime_fermat.obj \ -mp_prime_frobenius_underwood.obj mp_prime_is_prime.obj mp_prime_miller_rabin.obj mp_prime_next_prime.obj \ -mp_prime_rabin_miller_trials.obj mp_prime_rand.obj mp_prime_strong_lucas_selfridge.obj mp_radix_size.obj \ -mp_radix_size_overestimate.obj mp_rand.obj mp_rand_source.obj mp_read_radix.obj mp_reduce.obj mp_reduce_2k.obj \ -mp_reduce_2k_l.obj mp_reduce_2k_setup.obj mp_reduce_2k_setup_l.obj mp_reduce_is_2k.obj mp_reduce_is_2k_l.obj \ -mp_reduce_setup.obj mp_root_n.obj mp_rshd.obj mp_sbin_size.obj mp_set.obj mp_set_double.obj mp_set_i32.obj mp_set_i64.obj \ -mp_set_l.obj mp_set_u32.obj mp_set_u64.obj mp_set_ul.obj mp_shrink.obj mp_signed_rsh.obj mp_sqrmod.obj mp_sqrt.obj \ +mp_invmod.obj mp_is_small_prime.obj mp_is_square.obj mp_kronecker.obj mp_lcm.obj mp_log_n.obj mp_lshd.obj mp_mod.obj \ +mp_mod_2d.obj mp_montgomery_calc_normalization.obj mp_montgomery_reduce.obj mp_montgomery_setup.obj mp_mul.obj \ +mp_mul_2.obj mp_mul_2d.obj mp_mul_d.obj mp_mulmod.obj mp_neg.obj mp_next_small_prime.obj mp_or.obj mp_pack.obj \ +mp_pack_count.obj mp_prec_small_prime.obj mp_prime_fermat.obj mp_prime_frobenius_underwood.obj \ +mp_prime_is_prime.obj mp_prime_miller_rabin.obj mp_prime_next_prime.obj mp_prime_rabin_miller_trials.obj \ +mp_prime_rand.obj mp_prime_strong_lucas_selfridge.obj mp_radix_size.obj mp_radix_size_overestimate.obj mp_rand.obj \ +mp_rand_source.obj mp_read_radix.obj mp_reduce.obj mp_reduce_2k.obj mp_reduce_2k_l.obj mp_reduce_2k_setup.obj \ +mp_reduce_2k_setup_l.obj mp_reduce_is_2k.obj mp_reduce_is_2k_l.obj mp_reduce_setup.obj mp_root_n.obj mp_rshd.obj \ +mp_sbin_size.obj mp_set.obj mp_set_double.obj mp_set_i32.obj mp_set_i64.obj mp_set_l.obj mp_set_u32.obj mp_set_u64.obj \ +mp_set_ul.obj mp_shrink.obj mp_sieve_clear.obj mp_sieve_init.obj mp_signed_rsh.obj mp_sqrmod.obj mp_sqrt.obj \ mp_sqrtmod_prime.obj mp_sub.obj mp_sub_d.obj mp_submod.obj mp_to_radix.obj mp_to_sbin.obj mp_to_ubin.obj mp_ubin_size.obj \ mp_unpack.obj mp_xor.obj mp_zero.obj s_mp_add.obj s_mp_copy_digs.obj s_mp_div_3.obj s_mp_div_recursive.obj \ -s_mp_div_school.obj s_mp_div_small.obj s_mp_exptmod.obj s_mp_exptmod_fast.obj s_mp_get_bit.obj s_mp_invmod.obj \ -s_mp_invmod_odd.obj s_mp_log.obj s_mp_log_2expt.obj s_mp_log_d.obj s_mp_montgomery_reduce_comba.obj s_mp_mul.obj \ -s_mp_mul_balance.obj s_mp_mul_comba.obj s_mp_mul_high.obj s_mp_mul_high_comba.obj s_mp_mul_karatsuba.obj \ -s_mp_mul_toom.obj s_mp_prime_is_divisible.obj s_mp_prime_tab.obj s_mp_radix_map.obj \ -s_mp_radix_size_overestimate.obj s_mp_rand_platform.obj s_mp_sqr.obj s_mp_sqr_comba.obj s_mp_sqr_karatsuba.obj \ -s_mp_sqr_toom.obj s_mp_sub.obj s_mp_zero_buf.obj s_mp_zero_digs.obj +s_mp_div_school.obj s_mp_div_small.obj s_mp_eratosthenes.obj s_mp_eratosthenes_init.obj \ +s_mp_eratosthenes_segment.obj s_mp_eratosthenes_segment_init.obj s_mp_exptmod.obj s_mp_exptmod_fast.obj \ +s_mp_get_bit.obj s_mp_init_single_segment_with_start.obj s_mp_invmod.obj s_mp_invmod_odd.obj s_mp_isqrt.obj \ +s_mp_log.obj s_mp_log_2expt.obj s_mp_log_d.obj s_mp_montgomery_reduce_comba.obj s_mp_mul.obj s_mp_mul_balance.obj \ +s_mp_mul_comba.obj s_mp_mul_high.obj s_mp_mul_high_comba.obj s_mp_mul_karatsuba.obj s_mp_mul_toom.obj \ +s_mp_prime_is_divisible.obj s_mp_prime_tab.obj s_mp_radix_map.obj s_mp_radix_size_overestimate.obj \ +s_mp_rand_platform.obj s_mp_sieve_clear.obj s_mp_sieve_get.obj s_mp_sieve_nextset.obj s_mp_sieve_setall.obj \ +s_mp_sqr.obj s_mp_sqr_comba.obj s_mp_sqr_karatsuba.obj s_mp_sqr_toom.obj s_mp_sub.obj s_mp_zero_buf.obj \ +s_mp_zero_digs.obj HEADERS_PUB=tommath.h HEADERS=tommath_private.h tommath_class.h tommath_superclass.h tommath_cutoffs.h $(HEADERS_PUB) + #The default rule for make builds the tommath.lib library (static) default: $(LIBMAIN_S) diff --git a/makefile.shared b/makefile.shared index 4938da820..10714571a 100644 --- a/makefile.shared +++ b/makefile.shared @@ -30,23 +30,27 @@ mp_error_to_string.o mp_exch.o mp_expt_n.o mp_exptmod.o mp_exteuclid.o mp_fread. mp_from_ubin.o mp_fwrite.o mp_gcd.o mp_get_double.o mp_get_i32.o mp_get_i64.o mp_get_l.o mp_get_mag_u32.o \ mp_get_mag_u64.o mp_get_mag_ul.o mp_grow.o mp_hash.o mp_init.o mp_init_copy.o mp_init_i32.o mp_init_i64.o \ mp_init_l.o mp_init_multi.o mp_init_set.o mp_init_size.o mp_init_u32.o mp_init_u64.o mp_init_ul.o \ -mp_invmod.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o mp_mod_2d.o \ -mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o mp_mul_2.o \ -mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_or.o mp_pack.o mp_pack_count.o mp_prime_fermat.o \ -mp_prime_frobenius_underwood.o mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o \ -mp_prime_rabin_miller_trials.o mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o \ -mp_radix_size_overestimate.o mp_rand.o mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o \ -mp_reduce_2k_l.o mp_reduce_2k_setup.o mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o \ -mp_reduce_setup.o mp_root_n.o mp_rshd.o mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o \ -mp_set_l.o mp_set_u32.o mp_set_u64.o mp_set_ul.o mp_shrink.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ +mp_invmod.o mp_is_small_prime.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o \ +mp_mod_2d.o mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o \ +mp_mul_2.o mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_next_small_prime.o mp_or.o mp_pack.o \ +mp_pack_count.o mp_prec_small_prime.o mp_prime_fermat.o mp_prime_frobenius_underwood.o \ +mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o mp_prime_rabin_miller_trials.o \ +mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o mp_radix_size_overestimate.o mp_rand.o \ +mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o mp_reduce_2k_l.o mp_reduce_2k_setup.o \ +mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o mp_reduce_setup.o mp_root_n.o mp_rshd.o \ +mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o mp_set_l.o mp_set_u32.o mp_set_u64.o \ +mp_set_ul.o mp_shrink.o mp_sieve_clear.o mp_sieve_init.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ mp_sqrtmod_prime.o mp_sub.o mp_sub_d.o mp_submod.o mp_to_radix.o mp_to_sbin.o mp_to_ubin.o mp_ubin_size.o \ mp_unpack.o mp_xor.o mp_zero.o s_mp_add.o s_mp_copy_digs.o s_mp_div_3.o s_mp_div_recursive.o \ -s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_exptmod_fast.o s_mp_get_bit.o s_mp_invmod.o \ -s_mp_invmod_odd.o s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o \ -s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \ -s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o \ -s_mp_radix_size_overestimate.o s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o \ -s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o s_mp_zero_digs.o +s_mp_div_school.o s_mp_div_small.o s_mp_eratosthenes.o s_mp_eratosthenes_init.o \ +s_mp_eratosthenes_segment.o s_mp_eratosthenes_segment_init.o s_mp_exptmod.o s_mp_exptmod_fast.o \ +s_mp_get_bit.o s_mp_init_single_segment_with_start.o s_mp_invmod.o s_mp_invmod_odd.o s_mp_isqrt.o \ +s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o s_mp_mul_balance.o \ +s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o s_mp_mul_toom.o \ +s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_radix_size_overestimate.o \ +s_mp_rand_platform.o s_mp_sieve_clear.o s_mp_sieve_get.o s_mp_sieve_nextset.o s_mp_sieve_setall.o \ +s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \ +s_mp_zero_digs.o #END_INS diff --git a/makefile.unix b/makefile.unix index be0467480..3b4ea1008 100644 --- a/makefile.unix +++ b/makefile.unix @@ -36,23 +36,27 @@ mp_error_to_string.o mp_exch.o mp_expt_n.o mp_exptmod.o mp_exteuclid.o mp_fread. mp_from_ubin.o mp_fwrite.o mp_gcd.o mp_get_double.o mp_get_i32.o mp_get_i64.o mp_get_l.o mp_get_mag_u32.o \ mp_get_mag_u64.o mp_get_mag_ul.o mp_grow.o mp_hash.o mp_init.o mp_init_copy.o mp_init_i32.o mp_init_i64.o \ mp_init_l.o mp_init_multi.o mp_init_set.o mp_init_size.o mp_init_u32.o mp_init_u64.o mp_init_ul.o \ -mp_invmod.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o mp_mod_2d.o \ -mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o mp_mul_2.o \ -mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_or.o mp_pack.o mp_pack_count.o mp_prime_fermat.o \ -mp_prime_frobenius_underwood.o mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o \ -mp_prime_rabin_miller_trials.o mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o \ -mp_radix_size_overestimate.o mp_rand.o mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o \ -mp_reduce_2k_l.o mp_reduce_2k_setup.o mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o \ -mp_reduce_setup.o mp_root_n.o mp_rshd.o mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o \ -mp_set_l.o mp_set_u32.o mp_set_u64.o mp_set_ul.o mp_shrink.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ +mp_invmod.o mp_is_small_prime.o mp_is_square.o mp_kronecker.o mp_lcm.o mp_log_n.o mp_lshd.o mp_mod.o \ +mp_mod_2d.o mp_montgomery_calc_normalization.o mp_montgomery_reduce.o mp_montgomery_setup.o mp_mul.o \ +mp_mul_2.o mp_mul_2d.o mp_mul_d.o mp_mulmod.o mp_neg.o mp_next_small_prime.o mp_or.o mp_pack.o \ +mp_pack_count.o mp_prec_small_prime.o mp_prime_fermat.o mp_prime_frobenius_underwood.o \ +mp_prime_is_prime.o mp_prime_miller_rabin.o mp_prime_next_prime.o mp_prime_rabin_miller_trials.o \ +mp_prime_rand.o mp_prime_strong_lucas_selfridge.o mp_radix_size.o mp_radix_size_overestimate.o mp_rand.o \ +mp_rand_source.o mp_read_radix.o mp_reduce.o mp_reduce_2k.o mp_reduce_2k_l.o mp_reduce_2k_setup.o \ +mp_reduce_2k_setup_l.o mp_reduce_is_2k.o mp_reduce_is_2k_l.o mp_reduce_setup.o mp_root_n.o mp_rshd.o \ +mp_sbin_size.o mp_set.o mp_set_double.o mp_set_i32.o mp_set_i64.o mp_set_l.o mp_set_u32.o mp_set_u64.o \ +mp_set_ul.o mp_shrink.o mp_sieve_clear.o mp_sieve_init.o mp_signed_rsh.o mp_sqrmod.o mp_sqrt.o \ mp_sqrtmod_prime.o mp_sub.o mp_sub_d.o mp_submod.o mp_to_radix.o mp_to_sbin.o mp_to_ubin.o mp_ubin_size.o \ mp_unpack.o mp_xor.o mp_zero.o s_mp_add.o s_mp_copy_digs.o s_mp_div_3.o s_mp_div_recursive.o \ -s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_exptmod_fast.o s_mp_get_bit.o s_mp_invmod.o \ -s_mp_invmod_odd.o s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o \ -s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \ -s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o \ -s_mp_radix_size_overestimate.o s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o \ -s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o s_mp_zero_digs.o +s_mp_div_school.o s_mp_div_small.o s_mp_eratosthenes.o s_mp_eratosthenes_init.o \ +s_mp_eratosthenes_segment.o s_mp_eratosthenes_segment_init.o s_mp_exptmod.o s_mp_exptmod_fast.o \ +s_mp_get_bit.o s_mp_init_single_segment_with_start.o s_mp_invmod.o s_mp_invmod_odd.o s_mp_isqrt.o \ +s_mp_log.o s_mp_log_2expt.o s_mp_log_d.o s_mp_montgomery_reduce_comba.o s_mp_mul.o s_mp_mul_balance.o \ +s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o s_mp_mul_toom.o \ +s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_radix_size_overestimate.o \ +s_mp_rand_platform.o s_mp_sieve_clear.o s_mp_sieve_get.o s_mp_sieve_nextset.o s_mp_sieve_setall.o \ +s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \ +s_mp_zero_digs.o HEADERS_PUB=tommath.h diff --git a/mp_error_to_string.c b/mp_error_to_string.c index 39adcd124..8d587f201 100644 --- a/mp_error_to_string.c +++ b/mp_error_to_string.c @@ -21,6 +21,8 @@ const char *mp_error_to_string(mp_err code) return "Buffer overflow"; case MP_OVF: return "Integer overflow"; + case MP_SIEVE_MAX_REACHED: + return "Upper limit of prime sieve reached"; default: return "Invalid error code"; } diff --git a/mp_is_small_prime.c b/mp_is_small_prime.c new file mode 100644 index 000000000..78b8fcfb0 --- /dev/null +++ b/mp_is_small_prime.c @@ -0,0 +1,96 @@ +#include "tommath_private.h" +#ifdef MP_IS_SMALL_PRIME_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + +/* + * Sets "result" to true if n is prime or false respectively. + * Also sets "result" to zero in case of error. + * Worst case runtime is: building a base sieve and a segment and + * search the segment + */ +/* TODO: add ability to keep segments? */ +mp_err mp_is_small_prime(mp_sieve_prime n, bool *result, mp_sieve *sieve) +{ + mp_err err = MP_OKAY; + mp_sieve_prime a = 0uL, b = 0uL; + + if (n < 2uL) { + *result = false; + return err; + } + + if (n == 2uL) { + *result = true; + return err; + } + + if ((n & 1uL) == 0uL) { + *result = false; + return err; + } + + /* neither of 2^16-x, 2^32-x, or 2^64-x are prime for 0<=x<=4 */ + if (n >= (mp_sieve_prime)(MP_SIEVE_PRIME_MAX - 3u)) { + *result = false; + return err; + } + + if (sieve->base.content == NULL) { + if ((err = s_mp_eratosthenes_init((mp_sieve_prime)MP_SIEVE_BASE_SIEVE_SIZE, &(sieve->base))) != MP_OKAY) { + *result = false; + return err; + } + } + + /* No need to generate a segment if n is in the base sieve */ + if (n < MP_SIEVE_BASE_SIEVE_SIZE) { + /* might be a small sieve, so check size of sieve first */ + if (n < sieve->base.size) { + *result = (s_mp_sieve_get(&(sieve->base), (n - 1uL) / 2uL) == 1ull); + return err; + } else { + return MP_VAL; + } + } + + /* no further shortcuts to apply, build and search a segment */ + + /* we have a segment and may be able to use it */ + if (sieve->segment.content != NULL) { + a = sieve->single_segment_a; + /* last segment may not fit into range_a_b */ + if (a > (MP_SIEVE_PRIME_MAX - MP_SIEVE_RANGE_A_B)) { + b = (mp_sieve_prime)(MP_SIEVE_PRIME_MAX); + } else { + b = a + MP_SIEVE_RANGE_A_B; + } + /* check if n is inside the bounds of the segment */ + if (n >= a && n <= b) { + *result = (s_mp_sieve_get(&(sieve->segment), n - a) == 1ull); + return err; + } else { + a = 0ul; + b = 0ul; + } + } + + if (n > a) { + if (n > (MP_SIEVE_PRIME_MAX - MP_SIEVE_RANGE_A_B)) { + a = (mp_sieve_prime)(MP_SIEVE_PRIME_MAX - MP_SIEVE_RANGE_A_B); + b = (mp_sieve_prime) MP_SIEVE_PRIME_MAX; + } else { + a = n; + b = a + (mp_sieve_prime)MP_SIEVE_RANGE_A_B; + } + } + if ((err = s_mp_init_single_segment_with_start(a, + &(sieve->base), &(sieve->segment), &(sieve->single_segment_a))) != MP_OKAY) { + *result = false; + return err; + } + /* finally, check for primality */ + *result = (s_mp_sieve_get(&(sieve->segment), n - a) == 1ull); + return err; +} +#endif diff --git a/mp_next_small_prime.c b/mp_next_small_prime.c new file mode 100644 index 000000000..5106bbc40 --- /dev/null +++ b/mp_next_small_prime.c @@ -0,0 +1,35 @@ +#include "tommath_private.h" +#ifdef MP_NEXT_SMALL_PRIME_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + +/* + * Mimics behaviour of Pari/GP's nextprime(n) + * If n is prime set *result to n else set *result to first prime > n + * and 0 in case of error + */ +mp_err mp_next_small_prime(mp_sieve_prime n, mp_sieve_prime *result, mp_sieve *sieve) +{ + bool ret = false; + mp_err e = MP_OKAY; + + if (n < 2) { + *result = 2; + return e; + } + + for (; (ret == false) && (n != 0lu); n++) { + if (n > MP_SIEVE_BIGGEST_PRIME) { + return MP_SIEVE_MAX_REACHED; + } + /* just call mp_is_small_prime(), it does all of the heavy work */ + if ((e = mp_is_small_prime(n, &ret, sieve)) != MP_OKAY) { + *result = 0; + return e; + } + } + *result = n - 1; + return e; +} + +#endif diff --git a/mp_prec_small_prime.c b/mp_prec_small_prime.c new file mode 100644 index 000000000..d73411122 --- /dev/null +++ b/mp_prec_small_prime.c @@ -0,0 +1,38 @@ +#include "tommath_private.h" +#ifdef MP_PREC_SMALL_PRIME_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + +/* + * Mimics behaviour of Pari/GP's precprime(n) + * If n is prime set *result to n else set *result to first prime < n + * and 0 in case of error + */ +mp_err mp_prec_small_prime(mp_sieve_prime n, mp_sieve_prime *result, mp_sieve *sieve) +{ + bool ret = false; + mp_err e = MP_OKAY; + + if (n == 2lu) { + *result = 2lu; + return e; + } + + if (n < 2lu) { + *result = 0lu; + return e; + } + + for (; ret == false; n--) { + if ((e = mp_is_small_prime(n, &ret, sieve)) != MP_OKAY) { + *result = 0lu; + return e; + } + } + *result = n + 1lu; + + return e; +} + +#endif diff --git a/mp_sieve_clear.c b/mp_sieve_clear.c new file mode 100644 index 000000000..e96186935 --- /dev/null +++ b/mp_sieve_clear.c @@ -0,0 +1,26 @@ +#include "tommath_private.h" +#ifdef MP_SIEVE_CLEAR_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + +/* Free memory of one sieve */ +static void s_clear_one(mp_single_sieve *sieve) +{ + if (sieve->content != NULL) { + MP_FREE(sieve->content, (size_t)segment->alloc); + } +} + +void mp_sieve_clear(mp_sieve *sieve) +{ + s_clear_one(&(sieve->base)); + s_clear_one(&(sieve->segment)); +} + + +#endif diff --git a/mp_sieve_init.c b/mp_sieve_init.c new file mode 100644 index 000000000..87bc30da6 --- /dev/null +++ b/mp_sieve_init.c @@ -0,0 +1,21 @@ +#include "tommath_private.h" +#ifdef MP_SIEVE_INIT_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + +/* Set the default values for the sieve */ +void mp_sieve_init(mp_sieve *sieve) +{ + sieve->base.content = NULL; + sieve->segment.content = NULL; + sieve->single_segment_a = 0uL; +} + + + +#endif diff --git a/s_mp_eratosthenes.c b/s_mp_eratosthenes.c new file mode 100644 index 000000000..3f78c7557 --- /dev/null +++ b/s_mp_eratosthenes.c @@ -0,0 +1,25 @@ +#include "tommath_private.h" +#ifdef S_MP_ERATOSTHENES_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + +/* + * Simple Eratosthenes' sieve, starting at zero + * Keeping odd primes only as the single optimization + */ +void s_mp_eratosthenes(mp_single_sieve *bst) +{ + mp_sieve_prime n, k, r, j; + n = bst->size; + r = s_mp_isqrt(n); + s_mp_sieve_setall(bst); + for (k = 1uL; k < ((r - 1uL) / 2uL); k += 1uL) { + if (s_mp_sieve_get(bst, k)) { + for (j = k * 2uL * (k + 1uL); j < (n - 1uL) / 2uL; j += 2uL * k + 1uL) { + s_mp_sieve_clear(bst, j); + } + } + } +} + +#endif diff --git a/s_mp_eratosthenes_init.c b/s_mp_eratosthenes_init.c new file mode 100644 index 000000000..0686e1fe5 --- /dev/null +++ b/s_mp_eratosthenes_init.c @@ -0,0 +1,28 @@ +#include "tommath_private.h" +#ifdef S_MP_ERATOSTHENES_INIT_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* + * Initiate a sieve that stores the odd numbers only: + * allocate memory, set actual size and allocated size and fill it + */ +mp_err s_mp_eratosthenes_init(mp_sieve_prime n, mp_single_sieve *bst) +{ + size_t bits_in_sp; + + bits_in_sp = ((((size_t)n + MP_SIEVE_PRIME_NUM_BITS - 1u)/MP_SIEVE_PRIME_NUM_BITS)/2u); + bst->content = MP_MALLOC(bits_in_sp * sizeof(mp_sieve_prime) + sizeof(mp_sieve_prime)); + if (bst->content == NULL) { + return MP_MEM; + } + bst->alloc = (mp_sieve_prime)((bits_in_sp) * sizeof(mp_sieve_prime) + sizeof(mp_sieve_prime)); + bst->size = n; + s_mp_eratosthenes(bst); + return MP_OKAY; +} + + +#endif diff --git a/s_mp_eratosthenes_segment.c b/s_mp_eratosthenes_segment.c new file mode 100644 index 000000000..796512fcf --- /dev/null +++ b/s_mp_eratosthenes_segment.c @@ -0,0 +1,37 @@ +#include "tommath_private.h" +#ifdef S_MP_ERATOSTHENES_SEGMENT_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + +/* Fill sieve "segment" of the range [a,b] from the basic sieve "base" */ +void s_mp_eratosthenes_segment(mp_sieve_prime a, mp_sieve_prime b, + mp_single_sieve *base, mp_single_sieve *segment) +{ + mp_sieve_prime r, j, i; + mp_sieve_prime p; + r = s_mp_isqrt(b); + s_mp_sieve_setall(segment); + p = 0u; + for (i = 0u; p < r; i++) { + if (p <= 1ul) { + p = 2ul; + } else { + p = s_mp_sieve_nextset(base, ((p - 1uL) / 2uL) + 1uL); + p = (2uL * p) + 1uL; + } + j = (a / p) * p; + if (j < a) { + j += p; + } + for (; j < b; j += p) { + /* j+=p can overflow */ + if (j >= a) { + s_mp_sieve_clear(segment, j - a); + } else { + break; + } + } + } +} + +#endif diff --git a/s_mp_eratosthenes_segment_init.c b/s_mp_eratosthenes_segment_init.c new file mode 100644 index 000000000..84934084f --- /dev/null +++ b/s_mp_eratosthenes_segment_init.c @@ -0,0 +1,36 @@ +#include "tommath_private.h" +#ifdef S_MP_ERATOSTHENES_SEGMENT_INIT_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* + * Init sieve "segment" of the range [a,b]: + * allocate memory, set actual size and allocated size and fill it from sieve "base" + * TODO: merge with s_mp_eratosthenes_init() + */ +mp_err s_mp_eratosthenes_segment_init(mp_sieve_prime a, mp_sieve_prime b, + mp_single_sieve *segment, mp_single_sieve *base) +{ + mp_sieve_prime n; + size_t bits_in_sp; + + n = b - a; + bits_in_sp = ((size_t)n + MP_SIEVE_PRIME_NUM_BITS - 1u)/MP_SIEVE_PRIME_NUM_BITS; + bits_in_sp = bits_in_sp * sizeof(mp_sieve_prime) + sizeof(mp_sieve_prime); + /* All segments besides the last one have the same size--we can keep the memory */ + if (segment->content == NULL) { + segment->content = MP_MALLOC(bits_in_sp); + if (segment->content == NULL) { + return MP_MEM; + } + } + segment->alloc = (mp_sieve_prime)bits_in_sp; + segment->size = n; + s_mp_eratosthenes_segment(a, b, base, segment); + return MP_OKAY; +} + + +#endif diff --git a/s_mp_init_single_segment_with_start.c b/s_mp_init_single_segment_with_start.c new file mode 100644 index 000000000..fe5793c84 --- /dev/null +++ b/s_mp_init_single_segment_with_start.c @@ -0,0 +1,34 @@ +#include "tommath_private.h" +#ifdef S_MP_INIT_SINGLE_SEGMENT_WITH_START_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* + * Build a segment sieve with the largest reasonable size. "a" is the start of + * the sieve Size is MIN(range_a_b,MP_SIEVE_PRIME_MAX-a) + */ +mp_err s_mp_init_single_segment_with_start(mp_sieve_prime a, + mp_single_sieve *base_sieve, + mp_single_sieve *single_segment, + mp_sieve_prime *single_segment_a) +{ + mp_sieve_prime b; + mp_err err = MP_OKAY; + + /* last segment might not fit, depending on size of range_a_b */ + if (a > (MP_SIEVE_PRIME_MAX - MP_SIEVE_RANGE_A_B)) { + b = (mp_sieve_prime)MP_SIEVE_PRIME_MAX - 1uL; + } else { + b = a + (mp_sieve_prime)MP_SIEVE_RANGE_A_B; + } + if ((err = s_mp_eratosthenes_segment_init(a, b, single_segment, base_sieve)) != MP_OKAY) { + return err; + } + *single_segment_a = a; + return err; +} + + +#endif diff --git a/s_mp_isqrt.c b/s_mp_isqrt.c new file mode 100644 index 000000000..7dccd84f4 --- /dev/null +++ b/s_mp_isqrt.c @@ -0,0 +1,44 @@ +#include "tommath_private.h" +#ifdef S_MP_ISQRT_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* + * Integer square root, hardware style + * Wikipedia calls it "Digit-by-digit calculation" + * https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Digit-by-digit_calculation + * This is the base 2 method described at + * https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Binary_numeral_system_(base_2) + */ +mp_sieve_prime s_mp_isqrt(mp_sieve_prime n) +{ + mp_sieve_prime s, rem, root; + if (n < 1uL) { + return 0uL; + } + + /* highest power of four <= n */ + s = (mp_sieve_prime)(1uLL << (MP_SIEVE_PRIME_NUM_BITS - 2)); + rem = n; + root = 0uL; + while (s > n) { + s >>= 2uL; + } + while (s != 0uL) { + if (rem >= (s | root)) { + rem -= (s | root); + root >>= 1uL; + root |= s; + } else { + root >>= 1uL; + } + s >>= 2uL; + } + + return root; +} + + +#endif diff --git a/s_mp_sieve_clear.c b/s_mp_sieve_clear.c new file mode 100644 index 000000000..299186346 --- /dev/null +++ b/s_mp_sieve_clear.c @@ -0,0 +1,17 @@ +#include "tommath_private.h" +#ifdef S_MP_SIEVE_CLEAR_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + + +/* Clear bit */ +void s_mp_sieve_clear(mp_single_sieve *bst, mp_sieve_prime n) +{ + ((*((bst)->content+(n/MP_SIEVE_PRIME_NUM_BITS)) + &= ~(1llu<<(n % MP_SIEVE_PRIME_NUM_BITS)))); +} + + +#endif diff --git a/s_mp_sieve_get.c b/s_mp_sieve_get.c new file mode 100644 index 000000000..e3c70a937 --- /dev/null +++ b/s_mp_sieve_get.c @@ -0,0 +1,16 @@ +#include "tommath_private.h" +#ifdef S_MP_SIEVE_GET_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* Read bit */ +mp_sieve_prime s_mp_sieve_get(mp_single_sieve *bst, mp_sieve_prime n) +{ + return (((*((bst)->content+(n/MP_SIEVE_PRIME_NUM_BITS)) + & (1llu<<(n % MP_SIEVE_PRIME_NUM_BITS))) != 0)); +} + + +#endif diff --git a/s_mp_sieve_nextset.c b/s_mp_sieve_nextset.c new file mode 100644 index 000000000..97128b3a9 --- /dev/null +++ b/s_mp_sieve_nextset.c @@ -0,0 +1,18 @@ +#include "tommath_private.h" +#ifdef S_MP_SIEVE_NEXTSET_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + +/* Searches next set bit independent of the kind of sieve */ +mp_sieve_prime s_mp_sieve_nextset(mp_single_sieve *bst, mp_sieve_prime n) +{ + while ((n < ((bst)->size)) && (!s_mp_sieve_get(bst, n))) { + n++; + } + return n; +} + + +#endif diff --git a/s_mp_sieve_setall.c b/s_mp_sieve_setall.c new file mode 100644 index 000000000..96744d9b9 --- /dev/null +++ b/s_mp_sieve_setall.c @@ -0,0 +1,33 @@ +#include "tommath_private.h" +#ifdef S_MP_SIEVE_SETALL_C +/* LibTomMath, multiple-precision integer library -- Tom St Denis */ +/* SPDX-License-Identifier: Unlicense */ + + + + +#ifdef MP_USE_MEMOPS +# include +#endif +void s_mp_sieve_setall(mp_single_sieve *bst) +{ + /* + TODO: To set all even bits to zero + u32: 5 * 17*257 * 65537 + u64: u32 * 4294967297 + + Useful or problems with endianess? + */ +#ifdef MP_USE_MEMOPS + memset((bst)->content, 0xff, bst->alloc); +#else + mp_sieve_prime i, bs_size; + bs_size = bst->alloc / sizeof(mp_sieve_prime); + for (i = 0; i < bs_size; i++) { + (bst)->content[i] = MP_SIEVE_FILL; + } +#endif +} + + +#endif diff --git a/sources.cmake b/sources.cmake index 797d461e2..04e1c4c1b 100644 --- a/sources.cmake +++ b/sources.cmake @@ -57,6 +57,7 @@ mp_init_u32.c mp_init_u64.c mp_init_ul.c mp_invmod.c +mp_is_small_prime.c mp_is_square.c mp_kronecker.c mp_lcm.c @@ -73,9 +74,11 @@ mp_mul_2d.c mp_mul_d.c mp_mulmod.c mp_neg.c +mp_next_small_prime.c mp_or.c mp_pack.c mp_pack_count.c +mp_prec_small_prime.c mp_prime_fermat.c mp_prime_frobenius_underwood.c mp_prime_is_prime.c @@ -109,6 +112,8 @@ mp_set_u32.c mp_set_u64.c mp_set_ul.c mp_shrink.c +mp_sieve_clear.c +mp_sieve_init.c mp_signed_rsh.c mp_sqrmod.c mp_sqrt.c @@ -129,11 +134,17 @@ s_mp_div_3.c s_mp_div_recursive.c s_mp_div_school.c s_mp_div_small.c +s_mp_eratosthenes.c +s_mp_eratosthenes_init.c +s_mp_eratosthenes_segment.c +s_mp_eratosthenes_segment_init.c s_mp_exptmod.c s_mp_exptmod_fast.c s_mp_get_bit.c +s_mp_init_single_segment_with_start.c s_mp_invmod.c s_mp_invmod_odd.c +s_mp_isqrt.c s_mp_log.c s_mp_log_2expt.c s_mp_log_d.c @@ -150,6 +161,10 @@ s_mp_prime_tab.c s_mp_radix_map.c s_mp_radix_size_overestimate.c s_mp_rand_platform.c +s_mp_sieve_clear.c +s_mp_sieve_get.c +s_mp_sieve_nextset.c +s_mp_sieve_setall.c s_mp_sqr.c s_mp_sqr_comba.c s_mp_sqr_karatsuba.c diff --git a/tommath.def b/tommath.def index 1af4e9e79..2d0bffbbe 100644 --- a/tommath.def +++ b/tommath.def @@ -60,6 +60,7 @@ EXPORTS mp_init_u64 mp_init_ul mp_invmod + mp_is_small_prime mp_is_square mp_kronecker mp_lcm @@ -76,9 +77,11 @@ EXPORTS mp_mul_d mp_mulmod mp_neg + mp_next_small_prime mp_or mp_pack mp_pack_count + mp_prec_small_prime mp_prime_fermat mp_prime_frobenius_underwood mp_prime_is_prime @@ -112,6 +115,8 @@ EXPORTS mp_set_u64 mp_set_ul mp_shrink + mp_sieve_clear + mp_sieve_init mp_signed_rsh mp_sqrmod mp_sqrt diff --git a/tommath.h b/tommath.h index 44c92b22a..9faee3807 100644 --- a/tommath.h +++ b/tommath.h @@ -95,13 +95,14 @@ typedef enum { } mp_ord; typedef enum { - MP_OKAY = 0, /* no error */ - MP_ERR = -1, /* unknown error */ - MP_MEM = -2, /* out of mem */ - MP_VAL = -3, /* invalid input */ - MP_ITER = -4, /* maximum iterations reached */ - MP_BUF = -5, /* buffer overflow, supplied buffer too small */ - MP_OVF = -6 /* mp_int overflow, too many digits */ + MP_OKAY = 0, /* no error */ + MP_ERR = -1, /* unknown error */ + MP_MEM = -2, /* out of mem */ + MP_VAL = -3, /* invalid input */ + MP_ITER = -4, /* maximum iterations reached */ + MP_BUF = -5, /* buffer overflow, supplied buffer too small */ + MP_OVF = -6, /* mp_int overflow, too many digits */ + MP_SIEVE_MAX_REACHED = -7 /* Upper limit of prime sieve reached */ } mp_err; typedef enum { @@ -540,6 +541,7 @@ mp_err mp_prime_frobenius_underwood(const mp_int *N, bool *result) MP_WUR; */ mp_err mp_prime_is_prime(const mp_int *a, int t, bool *result) MP_WUR; + /* finds the next prime after the number "a" using "t" trials * of Miller-Rabin. * @@ -562,6 +564,57 @@ mp_err mp_prime_next_prime(mp_int *a, int t, bool bbs_style) MP_WUR; */ mp_err mp_prime_rand(mp_int *a, int t, int size, int flags) MP_WUR; +/* ---> full prime sieve <--- */ + +/* Segmented version of an Eratosthenes sieve, optimized for random access */ + +#ifdef MP_64BIT +typedef uint64_t mp_sieve_prime; +# define MP_SIEVE_PR_UINT PRIu64 +#else +typedef uint32_t mp_sieve_prime; +# define MP_SIEVE_PR_UINT PRIu32 +#endif + +#define MP_SIEVE_BIGGEST_PRIME 4294967291lu + +typedef struct mp_single_sieve_t { + mp_sieve_prime *content; /* bitset holding the sieve */ + mp_sieve_prime size; /* number of entries (which is a slightly misleading description) */ + mp_sieve_prime alloc; /* size in bytes */ +} mp_single_sieve; + +typedef struct mp_sieve_t { + mp_single_sieve base; /* base sieve (0..MP_SIEVE_PRIME_MAX_SQRT) */ + mp_single_sieve segment; /* segment (range_a_b) */ + mp_sieve_prime single_segment_a; /* startpoint of segment */ +} mp_sieve; + + +/* Init a sieve. Sets the necessary defaults */ +void mp_sieve_init(mp_sieve *sieve); + +/* Clear a sieve. Frees the memory for the contents of base and segment */ +void mp_sieve_clear(mp_sieve *sieve); + +/* + Deterministicaly checks if a small prime (< MP_SIEVE_PRIME_MAX) is prime. + Optimized for random access. + */ +mp_err mp_is_small_prime(mp_sieve_prime n, bool *result, mp_sieve *sieve) MP_WUR; +/* + Puts the next prime >= n in "result". May return MP_SIEVE_MAX_REACHED to flag the content + of "result" as the last valid one. +*/ +mp_err mp_next_small_prime(mp_sieve_prime n, mp_sieve_prime *result, mp_sieve *sieve) MP_WUR; +/* + Puts the prime preceeding "n" in "result" +*/ +mp_err mp_prec_small_prime(mp_sieve_prime n, mp_sieve_prime *result, mp_sieve *sieve) MP_WUR; + + + + /* ---> radix conversion <--- */ int mp_count_bits(const mp_int *a) MP_WUR; diff --git a/tommath_class.h b/tommath_class.h index becbcb193..a70322f7c 100644 --- a/tommath_class.h +++ b/tommath_class.h @@ -66,6 +66,7 @@ # define MP_INIT_U64_C # define MP_INIT_UL_C # define MP_INVMOD_C +# define MP_IS_SMALL_PRIME_C # define MP_IS_SQUARE_C # define MP_KRONECKER_C # define MP_LCM_C @@ -82,9 +83,11 @@ # define MP_MUL_D_C # define MP_MULMOD_C # define MP_NEG_C +# define MP_NEXT_SMALL_PRIME_C # define MP_OR_C # define MP_PACK_C # define MP_PACK_COUNT_C +# define MP_PREC_SMALL_PRIME_C # define MP_PRIME_FERMAT_C # define MP_PRIME_FROBENIUS_UNDERWOOD_C # define MP_PRIME_IS_PRIME_C @@ -118,6 +121,8 @@ # define MP_SET_U64_C # define MP_SET_UL_C # define MP_SHRINK_C +# define MP_SIEVE_CLEAR_C +# define MP_SIEVE_INIT_C # define MP_SIGNED_RSH_C # define MP_SQRMOD_C # define MP_SQRT_C @@ -138,11 +143,17 @@ # define S_MP_DIV_RECURSIVE_C # define S_MP_DIV_SCHOOL_C # define S_MP_DIV_SMALL_C +# define S_MP_ERATOSTHENES_C +# define S_MP_ERATOSTHENES_INIT_C +# define S_MP_ERATOSTHENES_SEGMENT_C +# define S_MP_ERATOSTHENES_SEGMENT_INIT_C # define S_MP_EXPTMOD_C # define S_MP_EXPTMOD_FAST_C # define S_MP_GET_BIT_C +# define S_MP_INIT_SINGLE_SEGMENT_WITH_START_C # define S_MP_INVMOD_C # define S_MP_INVMOD_ODD_C +# define S_MP_ISQRT_C # define S_MP_LOG_C # define S_MP_LOG_2EXPT_C # define S_MP_LOG_D_C @@ -159,6 +170,10 @@ # define S_MP_RADIX_MAP_C # define S_MP_RADIX_SIZE_OVERESTIMATE_C # define S_MP_RAND_PLATFORM_C +# define S_MP_SIEVE_CLEAR_C +# define S_MP_SIEVE_GET_C +# define S_MP_SIEVE_NEXTSET_C +# define S_MP_SIEVE_SETALL_C # define S_MP_SQR_C # define S_MP_SQR_COMBA_C # define S_MP_SQR_KARATSUBA_C @@ -449,6 +464,12 @@ # define S_MP_INVMOD_ODD_C #endif +#if defined(MP_IS_SMALL_PRIME_C) +# define S_MP_ERATOSTHENES_INIT_C +# define S_MP_INIT_SINGLE_SEGMENT_WITH_START_C +# define S_MP_SIEVE_GET_C +#endif + #if defined(MP_IS_SQUARE_C) # define MP_CLEAR_C # define MP_CMP_MAG_C @@ -566,6 +587,10 @@ # define MP_COPY_C #endif +#if defined(MP_NEXT_SMALL_PRIME_C) +# define MP_IS_SMALL_PRIME_C +#endif + #if defined(MP_OR_C) # define MP_CLAMP_C # define MP_GROW_C @@ -582,6 +607,10 @@ # define MP_COUNT_BITS_C #endif +#if defined(MP_PREC_SMALL_PRIME_C) +# define MP_IS_SMALL_PRIME_C +#endif + #if defined(MP_PRIME_FERMAT_C) # define MP_CLEAR_C # define MP_CMP_C @@ -848,6 +877,12 @@ #if defined(MP_SHRINK_C) #endif +#if defined(MP_SIEVE_CLEAR_C) +#endif + +#if defined(MP_SIEVE_INIT_C) +#endif + #if defined(MP_SIGNED_RSH_C) # define MP_ADD_D_C # define MP_DIV_2D_C @@ -1011,6 +1046,28 @@ # define MP_SUB_C #endif +#if defined(S_MP_ERATOSTHENES_C) +# define S_MP_ISQRT_C +# define S_MP_SIEVE_CLEAR_C +# define S_MP_SIEVE_GET_C +# define S_MP_SIEVE_SETALL_C +#endif + +#if defined(S_MP_ERATOSTHENES_INIT_C) +# define S_MP_ERATOSTHENES_C +#endif + +#if defined(S_MP_ERATOSTHENES_SEGMENT_C) +# define S_MP_ISQRT_C +# define S_MP_SIEVE_CLEAR_C +# define S_MP_SIEVE_NEXTSET_C +# define S_MP_SIEVE_SETALL_C +#endif + +#if defined(S_MP_ERATOSTHENES_SEGMENT_INIT_C) +# define S_MP_ERATOSTHENES_SEGMENT_C +#endif + #if defined(S_MP_EXPTMOD_C) # define MP_CLEAR_C # define MP_COPY_C @@ -1049,6 +1106,10 @@ #if defined(S_MP_GET_BIT_C) #endif +#if defined(S_MP_INIT_SINGLE_SEGMENT_WITH_START_C) +# define S_MP_ERATOSTHENES_SEGMENT_INIT_C +#endif + #if defined(S_MP_INVMOD_C) # define MP_ADD_C # define MP_CLEAR_MULTI_C @@ -1079,6 +1140,9 @@ # define MP_SUB_C #endif +#if defined(S_MP_ISQRT_C) +#endif + #if defined(S_MP_LOG_C) # define MP_CLEAR_MULTI_C # define MP_CMP_C @@ -1200,6 +1264,19 @@ #if defined(S_MP_RAND_PLATFORM_C) #endif +#if defined(S_MP_SIEVE_CLEAR_C) +#endif + +#if defined(S_MP_SIEVE_GET_C) +#endif + +#if defined(S_MP_SIEVE_NEXTSET_C) +# define S_MP_SIEVE_GET_C +#endif + +#if defined(S_MP_SIEVE_SETALL_C) +#endif + #if defined(S_MP_SQR_C) # define MP_CLAMP_C # define MP_CLEAR_C diff --git a/tommath_private.h b/tommath_private.h index 42920519c..8e3c1a26b 100644 --- a/tommath_private.h +++ b/tommath_private.h @@ -104,6 +104,23 @@ extern void *MP_CALLOC(size_t nmemb, size_t size); extern void MP_FREE(void *mem, size_t size); #endif +/* Size of the base sieve of mp_sieve*/ + +#define MP_SIEVE_PRIME_MAX 0xFFFFFFFFlu +#define MP_SIEVE_PRIME_MAX_SQRT 0xFFFFlu +#ifdef MP_64BIT +#define MP_SIEVE_FILL 0xFFFFFFFFFFFFFFFFllu +#else +#define MP_SIEVE_FILL 0xFFFFFFFFlu +#endif + +/* Size is in bits and must be of the form (2^n)-1! */ +#define MP_SIEVE_RANGE_A_B 0xFFF +#define MP_SIEVE_PRIME_NUM_BITS (sizeof(mp_sieve_prime)*CHAR_BIT) +#define MP_SIEVE_BASE_SIEVE_SIZE ((mp_sieve_prime)MP_SIEVE_PRIME_MAX_SQRT) + + + /* feature detection macro */ #ifdef _MSC_VER /* Prevent false positive: not enough arguments for function-like macro invocation */ @@ -208,6 +225,21 @@ MP_PRIVATE void s_mp_zero_buf(void *mem, size_t size); MP_PRIVATE void s_mp_zero_digs(mp_digit *d, int digits); MP_PRIVATE mp_err s_mp_radix_size_overestimate(const mp_int *a, const int radix, size_t *size); +MP_PRIVATE void s_mp_sieve_setall(mp_single_sieve *bst); +MP_PRIVATE void s_mp_sieve_clear(mp_single_sieve *bst, mp_sieve_prime n); +MP_PRIVATE mp_sieve_prime s_mp_sieve_get(mp_single_sieve *bst, mp_sieve_prime n) MP_WUR; +MP_PRIVATE mp_sieve_prime s_mp_sieve_nextset(mp_single_sieve *bst, mp_sieve_prime n) MP_WUR; +MP_PRIVATE mp_sieve_prime s_mp_isqrt(mp_sieve_prime n) MP_WUR; +MP_PRIVATE void s_mp_eratosthenes(mp_single_sieve *bst); +MP_PRIVATE mp_err s_mp_eratosthenes_init(mp_sieve_prime n, mp_single_sieve *bst) MP_WUR; +MP_PRIVATE void s_mp_eratosthenes_segment(mp_sieve_prime a, mp_sieve_prime b, mp_single_sieve *base, + mp_single_sieve *segment); +MP_PRIVATE mp_err s_mp_eratosthenes_segment_init(mp_sieve_prime a, mp_sieve_prime b, mp_single_sieve *base, + mp_single_sieve *segment) MP_WUR; +MP_PRIVATE mp_err s_mp_init_single_segment_with_start(mp_sieve_prime a, mp_single_sieve *base_sieve, + mp_single_sieve *single_segment, mp_sieve_prime *single_segment_a) MP_WUR; + + #define MP_RADIX_MAP_REVERSE_SIZE 80u extern MP_PRIVATE const char s_mp_radix_map[]; extern MP_PRIVATE const uint8_t s_mp_radix_map_reverse[];