Hi,

I'm finally getting around collecting some of the crypto-related Ruby code
I have lying around into one package. While doing so I've given some
thought to the fact that Ruby does not have a native power modulo operator
(which for example Python has). Its a bit of a shame since its heavily
used in crypto algs which has some importance in modern-day computer use.
Maybe we should add it to bignum.c?

Anyway I currently use the simple, recursive and pretty fast

  # Calculate ((b**p) % m) assuming that b and m are large integers.
  def Math.power_modulo(b, p, m)
    if p == 1
      b % m
    elsif (p & 0x1) == 0 # p.even?
      t = power_modulo(b, p >> 1, m)
      (t * t) % m
    else
      (b * power_modulo(b, p-1, m)) % m
    end
  end

but since this is so speed-critical it can (in this abnormal case) be
wise to spend time optimizing it. Since many heads are wiser than one I
thought I'd pose this to you as a challenge:

  Given that p, b and m are integers, p is typically less than 2**18 while
  b and m are around 1000 or more bits what is the fastest pure-Ruby
  implementation of powmod, ie. calculating ((b**p) % m)?

Maybe its not really a "kata" since the (best) code will not be thrown
away ;).

If you want a test case here is one:

class TestPowMod < Test::Unit::TestCase
  def test_01_large_powmod
    expected = 24656295158149552129522791057189025935153769756028704310461488183165007781190622
66158248517786824586416789690271181089924166576876358272852153508275457030547406
50790047459930747867644527843993869970230088426805592858713543180070230688245872
57820097544818526648507725023614743971311566767966295145235315729887013427381259
85437819072366890714347945762217274178746870467684532926146890122994632184298568
39868247687863706332185080258360241983212661935661577867566186261592507929101422
20369343724700074136993030921203032483179240961705092122797953462669310767801367
5181607550763200080096436693225368019004915688931329
    assert_equal(expected, Math.power_modulo(1+2**1000, 65537, 2**2048))
  end
end

Oh, and don't run the testcase above with

Math.power_modulo(b, p, m)
  ((b**p) % m)
end

unless you have lots of time to wait... ;)

Regards,

Robert Feldt