This is an automated email from the git hooks/post-receive script. ppm-guest pushed a commit to annotated tag v0.14 in repository libmath-prime-util-perl.
commit 73b3a6e6141e8fbd9a881d9f0bb9a7f11532862f Author: Dana Jacobsen <[email protected]> Date: Thu Nov 22 17:56:40 2012 -0800 Move mulmod, powmod, etc. to separate file. A few changes. --- MANIFEST | 1 + factor.c | 101 +-------------------------------------------------------- mulmod.h | 110 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 112 insertions(+), 100 deletions(-) diff --git a/MANIFEST b/MANIFEST index 66357cb..4cf35c5 100644 --- a/MANIFEST +++ b/MANIFEST @@ -10,6 +10,7 @@ README TODO XS.xs ptypes.h +mulmod.h aks.h aks.c cache.h diff --git a/factor.c b/factor.c index b89ca99..deedbae 100644 --- a/factor.c +++ b/factor.c @@ -7,6 +7,7 @@ #include "factor.h" #include "util.h" #include "sieve.h" +#include "mulmod.h" /* * You need to remember to use UV for unsigned and IV for signed types that @@ -69,106 +70,6 @@ int trial_factor(UV n, UV *factors, UV maxtrial) } -#if (BITS_PER_WORD == 32) && HAVE_STD_U64 - - /* We have 64-bit available, but UV is 32-bit. Do the math in 64-bit. - * Even if it is emulated, it should be as fast or faster than us doing it. - */ - #define addmod(n,a,m) (UV)(((uint64_t)(n)+(uint64_t)(a)) % ((uint64_t)(m))) - #define mulmod(a,b,m) (UV)(((uint64_t)(a)*(uint64_t)(b)) % ((uint64_t)(m))) - #define sqrmod(n,m) (UV)(((uint64_t)(n)*(uint64_t)(n)) % ((uint64_t)(m))) - -#elif defined(__GNUC__) && defined(__x86_64__) - - /* GCC on a 64-bit Intel x86 */ - static UV _mulmod(UV a, UV b, UV c) { - UV d; /* to hold the result of a*b mod c */ - /* calculates a*b mod c, stores result in d */ - asm ("mov %1, %%rax;" /* put a into rax */ - "mul %2;" /* mul a*b -> rdx:rax */ - "div %3;" /* (a*b)/c -> quot in rax remainder in rdx */ - "mov %%rdx, %0;" /* store result in d */ - :"=r"(d) /* output */ - :"r"(a), "r"(b), "r"(c) /* input */ - :"%rax", "%rdx" /* clobbered registers */ - ); - return d; - } - #define mulmod(a,b,m) _mulmod(a,b,m) - #define sqrmod(n,m) _mulmod(n,n,m) - -#else - - /* UV is the largest integral type available (that we know of). */ - - /* Do it by hand */ - static UV _mulmod(UV a, UV b, UV m) { - UV r = 0; - while (b > 0) { - if (b & 1) { - if (r == 0) { - r = a; - } else { - r = m - r; - r = (a >= r) ? a-r : m-r+a; - } - } - a = (a > (m-a)) ? (a-m)+a : a+a; - b >>= 1; - } - return r; - } - - /* if n is smaller than this, you can multiply without overflow */ - #define HALF_WORD (UVCONST(1) << (BITS_PER_WORD/2)) - #define mulmod(a,b,m) (((a)|(b)) < HALF_WORD) ? ((a)*(b))%(m):_mulmod(a,b,m) - #define sqrmod(n,m) ((n) < HALF_WORD) ? ((n)*(n))%(m):_mulmod(n,n,m) - -#endif - -#ifndef addmod - #define addmod(n,a,m) ((((m)-(n)) > (a)) ? ((n)+(a)) : ((n)+(a)-(m))) -#endif - -/* n^power mod m */ -#ifndef HALF_WORD - static UV powmod(UV n, UV power, UV m) { - UV t = 1; - n %= m; - while (power) { - if (power & 1) t = mulmod(t, n, m); - power >>= 1; - if (power) n = sqrmod(n, m); - } - return t; - } -#else - static UV powmod(UV n, UV power, UV m) { - UV t = 1; - n %= m; - if (m < HALF_WORD) { - while (power) { - if (power & 1) t = (t*n)%m; - power >>= 1; - if (power) n = (n*n)%m; - } - } else { - while (power) { - if (power & 1) t = mulmod(t, n, m); - power >>= 1; - if (power) n = sqrmod(n,m); - } - } - return t; - } -#endif - -/* n^power + a mod m */ -#define powaddmod(n, p, a, m) addmod(powmod(n,p,m),a,m) - -/* n^2 + a mod m */ -#define sqraddmod(n, a, m) addmod(sqrmod(n,m), a,m) - /* Return 0 if n is not a perfect square. Set sqrtn to int(sqrt(n)) if so. * * Some simple solutions: diff --git a/mulmod.h b/mulmod.h new file mode 100644 index 0000000..bce8eb7 --- /dev/null +++ b/mulmod.h @@ -0,0 +1,110 @@ +#ifndef MPU_MULMOD_H +#define MPU_MULMOD_H + +#include "ptypes.h" + +#if defined(__GNUC__) || defined(_MSC_VER) + #define INLINE inline +#else + #define INLINE +#endif + +/* if n is smaller than this, you can multiply without overflow */ +#define HALF_WORD (UVCONST(1) << (BITS_PER_WORD/2)) + +#if (BITS_PER_WORD == 32) && HAVE_STD_U64 + + /* We have 64-bit available, but UV is 32-bit. Do the math in 64-bit. + * Even if it is emulated, it should be as fast or faster than us doing it. + */ + #define addmod(n,a,m) (UV)(((uint64_t)(n)+(uint64_t)(a)) % ((uint64_t)(m))) + #define mulmod(a,b,m) (UV)(((uint64_t)(a)*(uint64_t)(b)) % ((uint64_t)(m))) + #define sqrmod(n,m) (UV)(((uint64_t)(n)*(uint64_t)(n)) % ((uint64_t)(m))) + +#elif defined(__GNUC__) && defined(__x86_64__) + + /* GCC on a 64-bit Intel x86 */ + static INLINE UV _mulmod(UV a, UV b, UV c) { + UV d; /* to hold the result of a*b mod c */ + /* calculates a*b mod c, stores result in d */ + asm ("mov %1, %%rax;" /* put a into rax */ + "mul %2;" /* mul a*b -> rdx:rax */ + "div %3;" /* (a*b)/c -> quot in rax remainder in rdx */ + "mov %%rdx, %0;" /* store result in d */ + :"=r"(d) /* output */ + :"r"(a), "r"(b), "r"(c) /* input */ + :"%rax", "%rdx" /* clobbered registers */ + ); + return d; + } + #define mulmod(a,b,m) _mulmod(a,b,m) + #define sqrmod(n,m) _mulmod(n,n,m) + +#else + + /* UV is the largest integral type available (that we know of). */ + + /* Do it by hand */ + static INLINE UV _mulmod(UV a, UV b, UV m) { + UV r = 0; + a %= m; /* These are wasteful given that careful attention from the */ + b %= m; /* caller should make them unnecessary. */ + while (b > 0) { + if (b & 1) r = ((m-r) > a) ? r+a : r+a-m; /* r = (r + a) % m */ + b >>= 1; + if (b) a = ((m-a) > a) ? a+a : a+a-m; /* a = (a + a) % m */ + } + return r; + } + + #define mulmod(a,b,m) (((a)|(b)) < HALF_WORD) ? ((a)*(b))%(m):_mulmod(a,b,m) + #define sqrmod(n,m) ((n) < HALF_WORD) ? ((n)*(n))%(m):_mulmod(n,n,m) + +#endif + +#ifndef addmod + #define addmod(n,a,m) ((((m)-(n)) > (a)) ? ((n)+(a)) : ((n)+(a)-(m))) +#endif + +/* n^2 + a mod m */ +#define sqraddmod(n, a, m) addmod(sqrmod(n,m), a,m) +/* i*j + a mod m */ +#define muladdmod(i, j, a, m) addmod(mulmod(i,j,m), a, m) + +/* n^power mod m */ +#ifndef HALF_WORD + static INLINE UV powmod(UV n, UV power, UV m) { + UV t = 1; + n %= m; + while (power) { + if (power & 1) t = mulmod(t, n, m); + power >>= 1; + if (power) n = sqrmod(n, m); + } + return t; + } +#else + static INLINE UV powmod(UV n, UV power, UV m) { + UV t = 1; + n %= m; + if (m < HALF_WORD) { + while (power) { + if (power & 1) t = (t*n)%m; + power >>= 1; + if (power) n = (n*n)%m; + } + } else { + while (power) { + if (power & 1) t = mulmod(t, n, m); + power >>= 1; + if (power) n = sqrmod(n,m); + } + } + return t; + } +#endif + +/* n^power + a mod m */ +#define powaddmod(n, p, a, m) addmod(powmod(n,p,m),a,m) + +#endif -- Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git _______________________________________________ Pkg-perl-cvs-commits mailing list [email protected] http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-perl-cvs-commits
