unofficial mirror of guile-devel@gnu.org 
 help / color / mirror / Atom feed
From: Eric Hanchrow <offby1@blarg.net>
Subject: Re: Forwarded patch for modular exponentiation support (GMP powm)
Date: Tue, 10 Feb 2004 23:00:50 -0800	[thread overview]
Message-ID: <87ad3qm76l.fsf@offby1.atm01.sea.blarg.net> (raw)
In-Reply-To: 87hdxyctsa.fsf@zip.com.au

[-- Attachment #1: Type: text/plain, Size: 136 bytes --]

This isn't exactly a patch -- rather, it's a new file followed by a
patch.  Anyhow I've incorporated most of Kevin Ryde's suggestions.


[-- Attachment #2: semi-patch --]
[-- Type: application/octet-stream, Size: 4819 bytes --]

;;;; mexp.test --- test suite for Guile's modular exponentiation functions    -*- scheme -*-
;;;; Eric Hanchrow <offby1@blarg.net> --- January 2004
;;;;
;;;; Copyright (C) 2004 Free Software Foundation, Inc.
;;;; 
;;;; This program is free software; you can redistribute it and/or modify
;;;; it under the terms of the GNU General Public License as published by
;;;; the Free Software Foundation; either version 2, or (at your option)
;;;; any later version.
;;;; 
;;;; This program is distributed in the hope that it will be useful,
;;;; but WITHOUT ANY WARRANTY; without even the implied warranty of
;;;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;;;; GNU General Public License for more details.
;;;; 
;;;; You should have received a copy of the GNU General Public License
;;;; along with this software; see the file COPYING.  If not, write to
;;;; the Free Software Foundation, Inc., 59 Temple Place, Suite 330,
;;;; Boston, MA 02111-1307 USA

;; You can easily double-check these numbers by using perl with the Math::BigInt::GMP library.
(with-test-prefix "modulo-expt"
   (pass-if (= 1 (modulo-expt 17 23 47)))
   
   (pass-if (= 1 (modulo-expt 17 -23 47)))
   
   (pass-if (= 17 (modulo-expt 17 -22 47)))
   
   (pass-if (= 36 (modulo-expt 17 22 47)))
   
   (pass-if (= 183658794479969134816674175082294846241553725240 (modulo-expt 111122223333444455556666 111122223333444455556666 1153478690012629968439432872520758982731022934717)))

   (pass-if-exception
    "Proper exception with 0 modulus"
    (cons 'numerical-overflow "")
    (modulo-expt 17 23 0))

   (pass-if-exception
    "Proper exception when result not invertible"
    (cons 'numerical-overflow "")
    (modulo-expt 10 -1 48))
)
\f
Index: guile-core/libguile/numbers.c
===================================================================
RCS file: /cvsroot/guile/guile/guile-core/libguile/numbers.c,v
retrieving revision 1.221
diff -w -u -r1.221 numbers.c
--- guile-core/libguile/numbers.c	6 Jan 2004 21:55:29 -0000	1.221
+++ guile-core/libguile/numbers.c	11 Feb 2004 06:56:11 -0000
@@ -1519,6 +1519,84 @@
 }
 #undef FUNC_NAME
 
+static void
+coerce_to_big(SCM in, mpz_t out)
+{
+  if (SCM_BIGP(in))
+    mpz_set(out, SCM_I_BIG_MPZ(in));
+  else if (SCM_INUMP(in))
+    {
+      mpz_set_si(out, SCM_INUM(in));
+    }
+  else
+    scm_wrong_type_arg("modulo-expt", 1, in);
+}
+
+SCM_DEFINE(scm_modular_expt, "modulo-expt", 3, 0, 0,
+           (SCM n, SCM k, SCM m),
+            "Return @var{n} raised to the non-negative integer exponent\n"
+	    "@var{k}, modulo @var{m}.\n"
+	    "\n"
+	    "@lisp\n"
+	    "(modulo-expt 2 3 5)\n"
+	    "   @result{} 3\n"
+	    "@end lisp")
+#define FUNC_NAME s_scm_modular_expt
+{
+  SCM result = scm_i_mkbig();
+  mpz_t n_tmp; 
+  mpz_t k_tmp; 
+  mpz_t m_tmp; 
+  int needs_inverting = 0;
+
+  if (SCM_NFALSEP(scm_zero_p(m))) {
+    scm_num_overflow(FUNC_NAME);
+  }
+  
+  mpz_init(n_tmp);
+  mpz_init(k_tmp);
+  mpz_init(m_tmp);
+
+  coerce_to_big(n, n_tmp);
+  coerce_to_big(k, k_tmp);
+  coerce_to_big(m, m_tmp);
+  
+  /* if the exponent K is negative, and we simply call mpz_powm, we
+     might well get a divide-by-zero exception.  Since those are hard
+     to handle, we'll do the inversion ourselves -- because that way
+     we get a simple failure code, which is easy to handle. */
+  
+  if (-1 == mpz_sgn(k_tmp))
+    {
+      needs_inverting = 1;
+      mpz_neg(k_tmp, k_tmp);
+    }
+
+  mpz_powm (SCM_I_BIG_MPZ (result),
+            n_tmp,
+            k_tmp,
+            m_tmp);
+
+  if (needs_inverting)
+    {
+      mpz_t inverted;
+      mpz_init(inverted);
+      if (!mpz_invert (inverted, SCM_I_BIG_MPZ(result), m_tmp))
+        {
+          mpz_clear(inverted);
+          scm_num_overflow(FUNC_NAME);
+        }
+      mpz_set(SCM_I_BIG_MPZ(result), inverted);
+      mpz_clear(inverted);
+    }
+  mpz_clear(m_tmp);
+  mpz_clear(k_tmp);
+  mpz_clear(n_tmp);
+
+  return scm_i_normbig (result);
+}
+#undef FUNC_NAME
+
 SCM_DEFINE (scm_integer_expt, "integer-expt", 2, 0, 0,
             (SCM n, SCM k),
 	    "Return @var{n} raised to the non-negative integer exponent\n"
Index: guile-core/libguile/numbers.h
===================================================================
RCS file: /cvsroot/guile/guile/guile-core/libguile/numbers.h,v
retrieving revision 1.77
diff -w -u -r1.77 numbers.h
--- guile-core/libguile/numbers.h	18 Nov 2003 19:59:51 -0000	1.77
+++ guile-core/libguile/numbers.h	11 Feb 2004 06:56:11 -0000
@@ -201,6 +201,7 @@
 SCM_API SCM scm_logtest (SCM n1, SCM n2);
 SCM_API SCM scm_logbit_p (SCM n1, SCM n2);
 SCM_API SCM scm_lognot (SCM n);
+SCM_API SCM scm_modular_expt (SCM n, SCM k, SCM m);
 SCM_API SCM scm_integer_expt (SCM z1, SCM z2);
 SCM_API SCM scm_ash (SCM n, SCM cnt);
 SCM_API SCM scm_bit_extract (SCM n, SCM start, SCM end);

[-- Attachment #3: Type: text/plain, Size: 102 bytes --]


-- 
Hamburgers!  The cornerstone of any nutritious breakfast.
        -- Jules {From "Pulp Fiction"}

[-- Attachment #4: Type: text/plain, Size: 142 bytes --]

_______________________________________________
Guile-devel mailing list
Guile-devel@gnu.org
http://mail.gnu.org/mailman/listinfo/guile-devel

  reply	other threads:[~2004-02-11  7:00 UTC|newest]

Thread overview: 17+ messages / expand[flat|nested]  mbox.gz  Atom feed  top
2004-01-28 17:22 Forwarded patch for modular exponentiation support (GMP powm) Rob Browning
2004-02-09 23:15 ` Kevin Ryde
2004-02-09 23:31   ` Kevin Ryde
     [not found]   ` <87llnah2hf.fsf@offby1.atm01.sea.blarg.net>
2004-02-11  1:02     ` Kevin Ryde
2004-02-11  7:00       ` Eric Hanchrow [this message]
2004-02-11 23:44         ` Kevin Ryde
2004-02-12  4:56           ` Eric Hanchrow
2004-02-14  0:23             ` Kevin Ryde
2004-02-15  0:04               ` Eric Hanchrow
2004-02-15 22:08                 ` Kevin Ryde
2004-03-20 21:20             ` Marius Vollmer
2004-03-20 21:32               ` Kevin Ryde
2004-03-21  2:28                 ` Marius Vollmer
2004-03-21 22:18                   ` Kevin Ryde
2004-03-22  0:01                     ` Eric Hanchrow
2004-03-22  0:37                       ` Kevin Ryde
2004-03-25  0:02                       ` Kevin Ryde

Reply instructions:

You may reply publicly to this message via plain-text email
using any one of the following methods:

* Save the following mbox file, import it into your mail client,
  and reply-to-all from there: mbox

  Avoid top-posting and favor interleaved quoting:
  https://en.wikipedia.org/wiki/Posting_style#Interleaved_style

  List information: https://www.gnu.org/software/guile/

* Reply using the --to, --cc, and --in-reply-to
  switches of git-send-email(1):

  git send-email \
    --in-reply-to=87ad3qm76l.fsf@offby1.atm01.sea.blarg.net \
    --to=offby1@blarg.net \
    /path/to/YOUR_REPLY

  https://kernel.org/pub/software/scm/git/docs/git-send-email.html

* If your mail client supports setting the In-Reply-To header
  via mailto: links, try the mailto: link
Be sure your reply has a Subject: header at the top and a blank line before the message body.
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).