This is an automated email from the git hooks/post-receive script. ppm-guest pushed a commit to annotated tag v0.41 in repository libmath-prime-util-perl.
commit 1dc5d3c1b2b0f1eecf1b0a2dfea67e7b2f160898 Author: Dana Jacobsen <d...@acm.org> Date: Fri May 2 18:10:20 2014 -0700 Replace complicated invmod with simpler version that works properly for all n --- lib/Math/Prime/Util.pm | 6 +++--- util.c | 44 ++++++++++++-------------------------------- xt/pari-compare.pl | 11 +++++++++++ 3 files changed, 26 insertions(+), 35 deletions(-) diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm index 533e36a..51dce96 100644 --- a/lib/Math/Prime/Util.pm +++ b/lib/Math/Prime/Util.pm @@ -1165,9 +1165,9 @@ Objects can be passed to functions, and allow early loop exits. forcomposites { say } 2000,2020; Given a block and either an end number or a start and end pair, calls the -block for each composite in the inclusive range. The composites are the -numbers greater than 1 which are not prime: -C<4, 6, 8, 9, 10, 12, 14, 15, ...> +block for each composite in the inclusive range. The composites, +L<OEIS A002808|http://oeis.org/A002808>, are the numbers greater than 1 +which are not prime: C<4, 6, 8, 9, 10, 12, 14, 15, ...> =head2 fordivisors diff --git a/util.c b/util.c index 6286f44..3937473 100644 --- a/util.c +++ b/util.c @@ -1481,38 +1481,18 @@ UV znprimroot(UV n) { return 0; } -/* Calculate 1/a mod p. From William Hart. */ -UV modinverse(UV a, UV p) { - IV u1 = 1, u3 = a; - IV v1 = 0, v3 = p; - IV t1 = 0, t3 = 0; - IV quot; - while (v3) { - quot = u3 - v3; - if (u3 < (v3<<2)) { - if (quot < v3) { - if (quot < 0) { - t1 = u1; u1 = v1; v1 = t1; - t3 = u3; u3 = v3; v3 = t3; - } else { - t1 = u1 - v1; u1 = v1; v1 = t1; - t3 = u3 - v3; u3 = v3; v3 = t3; - } - } else if (quot < (v3<<1)) { - t1 = u1 - (v1<<1); u1 = v1; v1 = t1; - t3 = u3 - (v3<<1); u3 = v3; v3 = t3; - } else { - t1 = u1 - v1*3; u1 = v1; v1 = t1; - t3 = u3 - v3*3; u3 = v3; v3 = t3; - } - } else { - quot = u3 / v3; - t1 = u1 - v1*quot; u1 = v1; v1 = t1; - t3 = u3 - v3*quot; u3 = v3; v3 = t3; - } - } - if (u1 < 0) u1 += p; - return u1; +/* Calculate 1/a mod n. */ +UV modinverse(UV a, UV n) { + IV t = 0; UV nt = 1; + UV r = n; UV nr = a; + while (nr != 0) { + UV quot = r / nr; + { UV tmp = nt; nt = t - quot*nt; t = tmp; } + { UV tmp = nr; nr = r - quot*nr; r = tmp; } + } + if (r > 1) return 1; /* No inverse */ + if (t < 0) t += n; + return t; } UV divmod(UV a, UV b, UV n) { /* a / b mod n */ diff --git a/xt/pari-compare.pl b/xt/pari-compare.pl index 666cf99..5234c26 100755 --- a/xt/pari-compare.pl +++ b/xt/pari-compare.pl @@ -196,6 +196,17 @@ while (1) { # TODO: carmichael lambda? Pari doesn't have it. + # TODO: if we ever export out modular inverse, here's a test: + # use Math::Pari qw/PARI lift/; + # for (1..1000) { + # my $p = random_nbit_prime(64); + # for $i (2..100) { + # my $mpu = znorder($i,$p); + # my $pari = lift(PARI "Mod(1/$i,$p)"); + # die "$i $p" unless $pari == $mpu; + # } + # } + { my $m = int(rand($n-1)); my $mn = PARI "Mod($m,$n)"; my $order = znorder($m, $n); -- 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 Pkg-perl-cvs-commits@lists.alioth.debian.org http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-perl-cvs-commits