Modular exponentiation

Modular exponentiation is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Find the last 40 decimal digits of $a^{b}$, where

• $a=2988348162058574136915891421498819466320163312926952423791023078876139$
• $b=2351399303373464486466122544523690094744975233415544072992656881240319$

A computer is too slow to find the entire value of $a^{b}$. Instead, the program must use a fast algorithm for modular exponentiation: $a^{b}\mod m$.

The algorithm must work for any integers $a,b,m$ where $b\geq 0$ and $m\neq 0$.

C

Given numbers are too big for even 64 bit integers, so might as well take the lazy route and use GMP: <lang c>#include <gmp.h>

int main() { mpz_t a, b, m, r;

mpz_init_set_str(a, "2988348162058574136915891421498819466320" "163312926952423791023078876139", 0); mpz_init_set_str(b, "2351399303373464486466122544523690094744" "975233415544072992656881240319", 0); mpz_init(m); mpz_ui_pow_ui(m, 10, 40);

1. if 0 /* Don't actually need these two lines: hopefully GMP does it. */

mpz_mod(a, a, m); /* obvious */ mpz_mod(b, b, m); /* Fermat's little theorem */

1. endif

mpz_init(r); mpz_powm(r, a, b, m);

gmp_printf("%Zd\n", r); /* 16808958343740453059 */

mpz_clear(a); mpz_clear(b); mpz_clear(m); mpz_clear(r); return 0; }</lang>

Go

<lang go>package main

import (

   "fmt"
   "math/big"

)

func main() {

   a, _ := new(big.Int).SetString(
       "2988348162058574136915891421498819466320163312926952423791023078876139", 10)
   b, _ := new(big.Int).SetString(
       "2351399303373464486466122544523690094744975233415544072992656881240319", 10)
   m := big.NewInt(10)
   r := big.NewInt(40)
   m.Exp(m, r, nil)

   r.Exp(a, b, m)
   fmt.Println(r)

}</lang> Output:

1527229998585248450016808958343740453059

J

Solution:<lang j> m&|@^</lang> Example:<lang j> a =: 2988348162058574136915891421498819466320163312926952423791023078876139x

  b =: 2351399303373464486466122544523690094744975233415544072992656881240319x
  m =: 10^40x

  a m&|@^ b

1527229998585248450016808958343740453059</lang> Discussion: The phrase m&|@^ is the natural expression of a^b mod m in J, and is recognized by the interpreter as an opportunity for optimization, by avoiding the exponentiation.

Java

This example is in need of improvement: Tell algorithm it uses under the hood.

<lang java>import java.math.BigInteger;

public class PowMod {

   public static void main(String[] args){
       BigInteger a = new BigInteger(
     "2988348162058574136915891421498819466320163312926952423791023078876139");
       BigInteger b = new BigInteger(
     "2351399303373464486466122544523690094744975233415544072992656881240319");
       BigInteger m = new BigInteger("10000000000000000000000000000000000000000");
       
       System.out.println(a.modPow(b, m));
   }

}</lang> Output:

1527229998585248450016808958343740453059

Perl

<lang perl>use bigint;

my $a = 2988348162058574136915891421498819466320163312926952423791023078876139; my $b = 2351399303373464486466122544523690094744975233415544072992656881240319; my $m = 10 ** 40; print $a->bmodpow($b, $m), "\n";</lang> Output:

1527229998585248450016808958343740453059

PHP

<lang php><?php $a = '2988348162058574136915891421498819466320163312926952423791023078876139'; $b = '2351399303373464486466122544523690094744975233415544072992656881240319'; $m = '1' . str_repeat('0', 40); echo bcpowmod($a, $b, $m), "\n"; ?></lang> Output:

1527229998585248450016808958343740453059

Python

<lang python>a = 2988348162058574136915891421498819466320163312926952423791023078876139 b = 2351399303373464486466122544523690094744975233415544072992656881240319 m = 10 ** 40 print(pow(a, b, m))</lang> Output:

1527229998585248450016808958343740453059

Ruby

Ruby's core library has no modular exponentiation. OpenSSL, in Ruby's standard library, provides OpenSSL::BN#mod_exp. To reach this method, we call Integer#to_bn to convert a from Integer to OpenSSL::BN. The method implicitly converts b and m.

Library: OpenSSL

<lang ruby>require 'openssl'

a = 2988348162058574136915891421498819466320163312926952423791023078876139 b = 2351399303373464486466122544523690094744975233415544072992656881240319 m = 10 ** 40 puts a.to_bn.mod_exp(b, m)</lang>

Or we can implement a custom method, Integer#rosetta_mod_exp, to calculate the same result. This method does exponentiation by successive squaring, but replaces each intermediate product with a congruent value.

<lang ruby>class Integer

 def rosetta_mod_exp(exp, mod)
   exp < 0 and raise ArgumentError, "negative exponent"
   prod = 1
   base = self % mod
   until exp.zero?
     exp.odd? and prod = (prod * base) % mod
     exp >>= 1
     base = (base * base) % mod
   end
   prod
 end

end

a = 2988348162058574136915891421498819466320163312926952423791023078876139 b = 2351399303373464486466122544523690094744975233415544072992656881240319 m = 10 ** 40 puts a.rosetta_mod_exp(b, m)</lang>

Tcl

While Tcl does have arbitrary-precision arithmetic (from 8.5 onwards), it doesn't expose a modular exponentiation function. Thus we implement one ourselves.

Recursive

<lang tcl>package require Tcl 8.5

1. Algorithm from http://introcs.cs.princeton.edu/java/78crypto/ModExp.java.html
2. but Tcl has arbitrary-width integers and an exponentiation operator, which
3. helps simplify the code.

proc tcl::mathfunc::modexp {a b n} {

   if {$b == 0} {return 1}
   set c [expr {modexp($a, $b / 2, $n)**2 % $n}]
   if {$b & 1} {

set c [expr {($c * $a) % $n}]

   }
   return $c

}</lang> Demonstrating: <lang tcl>set a 2988348162058574136915891421498819466320163312926952423791023078876139 set b 2351399303373464486466122544523690094744975233415544072992656881240319 set n [expr {10**40}] puts [expr {modexp($a,$b,$n)}]</lang> Output:

1527229998585248450016808958343740453059

Iterative

<lang tcl>package require Tcl 8.5 proc modexp {a b n} {

   set c 1
   for {for {set t 1} {$t<$b} {set t [expr {$t<<1}]} {}} {$t} {set t [expr {$t>>1}]} {

if {$b & $t} { set c [expr {($c**2 * $a) % $n}] } else { set c [expr {$c**2 % $n}] }

   }
   return $c

}</lang> Demonstrating: <lang tcl>set a 2988348162058574136915891421498819466320163312926952423791023078876139 set b 2351399303373464486466122544523690094744975233415544072992656881240319 set n [expr {10**40}] puts [modexp $a $b $n]</lang> Output:

1527229998585248450016808958343740453059

TXR

<lang txr>@(bind result @(exptmod 2988348162058574136915891421498819466320163312926952423791023078876139

                       2351399303373464486466122544523690094744975233415544072992656881240319
                       (expt 10 40)))</lang>

$ ./txr rosetta/modexp.txr
result="1527229998585248450016808958343740453059"