[課題21] べき剰余

課題

2以上の整数m，整数x(0≦x＜m)，0以上の整数nが与えられたときに，「xⁿをmで割ったときの余り」で定義される「べき剰余」を，以下で説明する要領で再帰的に計算するメソッドを作成してください．なお0⁰=1とします．メソッドはテンプレートに組み込んで作成して，完成したプログラムを適宜テストした上で提出して下さい．

べき剰余はインターネットでの暗号通信において一般的に利用されているRSA暗号で基礎となっている処理です．

べき剰余を計算するメソッドの定義形式
テストプログラム(テンプレート)
べき剰余の計算手順
技術要素
サンプルプログラム
[参考] RSA暗号
[参考] 逐次的な繰り返し処理による「べき剰余」計算
[参考] 再帰呼び出しの問題
[参考] 2進展開に基づく繰り返し処理による「べき剰余」計算

べき剰余を計算するメソッドの定義形式

「xⁿをmで割ったときの余り」をべき剰余といいます．べき剰余を計算するメソッドは次の形式で定義します．


  # xⁿをmで割ったときの余りを求める
  def modexp(x,n,m)


  end

次に注意して下さい．

メソッドでは「最後に評価する式の値」が，そのメソッドの値になる(別途指定しない限り)
if式では「最後に評価する式の値」が「if式の値」になる

たとえば次のメソッドf(x,y)ではif式がメソッドの最後の式です．そのif式の値がメソッドの値になります．


  def f(x,y)
    z = x+3*y
    if x.even? # xは偶数?
      z**2
    else
      2*z
    end
  end

  u = f(5,2) # u == 2*11 == 22
  v = f(2,3) # v == 11**2 == 121 (**はべき乗の演算)

「べき剰余」メソッドは次のように使うことになります．


  m = 7
  x = 5
  n = 2

  # メソッドmodexp(x,n,m)の値(5²を7で割った余り)をyに代入する
  # 「modexp(x,n,m)の値」はメソッドで最後に評価される式の値であるように定義してあるとする
  y = modexp(x,n,m) # y ==> 4

テストプログラム(テンプレート)

今回のプログラムは次のテンプレートにしたがって作成してください．

modexp_template.rb

テンプレートは次のように構成されています．

定数定義
modexp(x,n,m)のメソッド定義 [課題]
小さな数での予備テストの定義
大きな数でのテストの定義
テストの実行と結果判定

メソッドが定義できたらプログラムを実行してみてください．最初は小さな数で予備テストを行うようになっています．プログラムを実行すると計算の結果を順に表示していきます．すべてのテストにパスすると「CLEARED」と表示します．途中で問題が見つかった場合にはその時点で実行を中断して「FAILED」と表示します．



  $ ruby modexp.rb 
  2^0 % 22 = 1
  2^1 % 22 = 2
  2^2 % 22 = 4
    :
    :
  21^10 % 22 = 1
  21^11 % 22 = 21
  CLEARED

予備テストがうまくいったら，次のようにプログラム内の定数の定義(26行目)を変更して下さい．


  (変更前) FUNC=:preliminary
  (変更後) FUNC=:trial

変更後は大きな数でのテストを実行します．この場合もテストにパスすると「CLEARED」と表示します．途中で問題が見つかった場合にはその時点で実行を中断して「FAILED」と表示します．

modexpの1回のテストは，ごく短い時間で終了するはずです．modexpの1回の計算にかかる時間が長すぎる場合は，プログラムに問題があると考えられます．


  $ ruby modexp.rb 
  428014040==>329905401==>428014040 OK
  332757652==>585163767==>332757652 OK
   :
   :
  374370147==>456294692==>374370147 OK
  CLEARED

このテストでは指定された回数だけ整数x(2≦x＜N)をランダムに選んで，それぞれについて次の計算を行います．


  ## 以下の計算でx == zとなることが期待される
  x = rand(2...N) # 2,3,...,N-1のいずれか
  y = modexp(x,E,N)
  z = modexp(y,D,N)

2≦x＜Nを満たすxについて，上のような計算の結果として「x==z」が成り立つことが確かめられるはずです(x=0,1についても，もちろん成立します)．なおこの関係は勝手なE,D,Nで成立するわけではなく，E,D,Nは適切に選択する必要があります(参考: RSA暗号)．

[参考] べき剰余に基づく暗号

上の計算ではE,D,Nが適切に選択されているとき，「x」が(E,N)によって「y」に変換されて，また変換後の「y」に対して(D,N)を適用すると，元の「x」が得られるようになっています．

そこでこの仕組みを「暗号」として利用できます． 2つの「べき剰余」の計算過程を暗号化と復号の処理として用いるわけです．

[暗号化]	秘密データ(x)	→	暗号データ(y = x^E % N)	(E,Nによってxをyに暗号化，E,Nは公開してよい)
[復号]	暗号データ(y)	→	秘密データ(x = y^D % N)	(D,Nによってyからxを復号，Dは秘匿する)

「E,N」は公開して「D」は秘密にしておきます．第三者に知られたくない秘密のデータ「x」を相手から自分に送ってもらうことを考えます．その場合「E,N」によって「x」から「y」を計算して「y」を送ってもらいます．「E,N」は公開されているため，「y」を計算することは問題ありません(E,Nを知っていれば誰でも暗号化できます)．この暗号データ「y」を元の秘密データ「x」に戻すためには「D,N」を使います．すでに見た通り，「D」を知っていれば「y」から「x」が計算できます(暗号化と全く同様です)．ところが「D」を知らなければ(「E,N」を知っていたとしても)，「y」から「x」を計算することは極めて難しいと考えられています(簡単に解読できる方法が知られていません)．この仕組みがインターネットでも一般的に使われているRSA暗号に利用されています．

べき剰余の計算手順

コンピュータでテストプログラムのように大きな値についてべき剰余を計算することは単純には実現できません．計算に現れる値の大きさに注意する必要があります．以下ではべき剰余を計算する手順について説明します．

さてxⁿをmで割ったときの余りは，数学的には単純にxⁿを計算してからmで割れば求められます．


 x = 3
 n = 4
 m = 5

 # 「x**n」でxⁿを計算する
 (x**n) % m # ==> 81 % 5 == 1

ところが(xがそれほど大きくなくても)nが大きな値のときにはxⁿの値が大きくなりすぎて，コンピュータでは単純には扱えなくなくなります．実際にテストプログラムで定義されている次のEの値を使って「2^E」をRubyで計算させてみると次のようになります．


  E = 153888981

  # 2^Eを計算したい
  2**E # ==> Infinity

  warning: in a**b, b may be too big

このようなときにはどうすればよいのでしょうか．今回計算したいのはxⁿの値そのものではなく，それをmで割った値(0以上m未満の値)であることに注目します．またべき乗が再帰的に定義できることを利用します．そうすることで，(そのままでは扱えないような)大きな値の出現を避けることができます．

べき乗の再帰的定義

nを0以上の整数として，xⁿを次のように再帰的に定義できます．




xⁿ = x^k・x^k・x^b  n > 0, n=2k+b, b∈{0,1}


xⁿ = 1 n = 0 (0⁰=1とします)

ここでk,bはnを2で割ったときの商と余りです．

便宜上0⁰=1と定義します(Rubyでも「0**0==1」です)．


  [備考]
  n=2kのとき    xⁿ = x^k・x^k
  n=2k+1のとき  xⁿ = x^k・x^k・x

  これらをまとめると
  xⁿ = x^k・x^k・x^b, n=2k+b,b∈{0,1}

剰余計算の方法

さて一般に乗算してから剰余を計算した値は，剰余を計算してから乗算を行って，さらに剰余をとった値と等しくなります．


  # 「u % v」はuをvで割った余り
  # (※ ここでは「＝」を等式の意味で使っている)
  (a・b) % m ＝ ((a % m)・(b % m)) % m

「乗算→剰余」という順序で計算をする代わりに，「剰余→乗算」の順で計算して最後にもう一度剰余をとればよいというわけです．

剰余を先にとることで計算は複雑になるものの，これを利用すれば(コンピュータで扱うには)大きすぎる値を避けうることが分かります．

具体的には上に示したようにn＞0のときn = 2k+b(b∈{0,1})として，「xⁿ = x^b・x^k・x^k」と再帰的に定義できることから，次のように「xⁿ % m」(= modexp(x,n,m))を「x^k % m」(= modexp(x,k,m))と「x^b」(= xあるいは1)の値から計算できることが分かります．


　xⁿ = x^k・x^k・x^b (n=2k+b,b∈{0,1})

 プログラムでは
   xⁿ % m → modexp(x,n,m) 
   x^k % m → modexp(x,k,m) 
 → modexp(x,n,m)は，modexp(x,k,m)と「x**b」から計算できる

こうすることで計算の過程ではせいぜいm³までの値を扱うだけで済むようにできます(0≦x＜mとしていることに注意)．

技術要素

プログラムの実行中断

実行したプログラムがいつまで経っても終了しなさそうな場合，次のように端末に入力することで，プログラムの実行を中断して，強制終了させることができます．
```
 [Ctrl] + [c]
```
[Ctrl]キーを押しつつ，[c]を押します．

剰余演算(%)
a % bでaをbで割った余りを計算できます．
```
  a = 10
  b = 3
  c = a % b # c == 1
```

商と剰余の計算(divmod)
n.divmod(m)でnをmで割ったときの商と剰余が同時に得られます．
```
  a = 10
  b = 3
  q,r = a.divmod(b) # q == 3, r == 1
```

サンプルプログラム

サンプルプログラムを示します．ブラウザの画面で開いたときに文字化けしてしまう場合には，ダウンロードしてEmacs等で開いてみて下さい．

combination_rec.rb
$ ruby combination_rec.rb 6 2
factorial_rec.rb
$ ruby factorial_rec.rb 10

[参考] RSA暗号

インターネットで一般的に用いられているRSA暗号は，(適切に選ばれた)大きな素数の積については「素因数分解」が本質的に難しいであろうことを根拠として，「暗号化」する鍵を一般公開しても安全にデータをやりとりできる「公開鍵暗号」の一種です．

「公開鍵暗号」が発明される前の暗号(共通鍵暗号)では，暗号化のための鍵と復号(元に戻す)のための鍵が同一です．したがって，暗号化の鍵を公開してしまうと，第三者に暗号文を解読されることになります．一般にインターネットでの通信は傍受される危険性があるため，鍵をネット上で安全にやりとりすることはできません．そのため共通鍵暗号は(単独では)インターネットで安全に使うことが困難です(※公開鍵暗号の技術と共通鍵暗号を組み合わせて使うのが一般的です)．

しかしRSA暗号では暗号化のための鍵は公開しても構いません．これは復号するときに暗号化の鍵(公開鍵)とは別の秘密鍵を用いるためです．公開鍵と秘密鍵はペアで定義されます．このようにペアになる2つの鍵を使う暗号システムを一般に公開鍵暗号と呼びます．

さて今回のテストプログラムでは次のような計算をしています．


  y = x^E % N  # x --> y
  z = y^D % N  # y --> z == x

「x」を守るべき秘密のデータとして，「y = x^E % N」という計算によって，「EとN」という「鍵」で秘密のデータ「x」を「y」に暗号化したものと考えることができます．すでに見た通り今回使ったE,D,Nについては「z = y^D % N」のように計算するとzがxと一致します．この計算により「DとN」という「鍵」で暗号データ「y」を復号して，元の秘密データ「x」を復元したものとみなせます．

実際にこの仕組みを使うには，まずデータ受信者が(E,D,N)を適切に決めて，送信者に(E,N)を渡します． (E,N)は第三者も知りうる公開鍵となります．送信者は秘密のデータ「x」を公開鍵(E,N)で暗号化した結果の「y」を暗号データとして受信者に送ります．秘密の「D」を知っている受信者は，暗号データ「y」から「z = y^D% N」と計算することで簡単に秘密データ「x」(== z)を得ることができます．しかしたとえ第三者が公開鍵(E,N)を知ったとしても， (E,D,N)を適切に決めた場合， (E,N)の情報のみでは復号の鍵「D」を計算することは極めて困難であると考えられています．このとき暗号データ「y」と公開鍵(E,N)を知られても，秘密鍵「D」を秘匿することで秘密データ「x」が守られることになります．

P,Qを相異なる素数として，E，D，Nを次の性質を満たすように定めれば，それらを「鍵」として使えることが知られています．これがRSA暗号の原理です．


   N = P・Q  (P≠Q) P,Qは素数
   E s.t. gcd(E,φ(N)) == 1 # gcdは最大公約数(つまりEとφ(N)は互いに素)
   D s.t. (E・D) % φ(N) == 1 
   φ(N)=(P-1)(Q-1)

十分大きなP,Qを適切に選ぶことで，(E,N)を公開しても，(D,P,Q)を秘密にしておけば，暗号文yからもとの秘密のデータxを得ることは極めて難しいと信じられています(すでに見ている通りP,Qは計算には使いません)．一方秘密のDを知っていれば，今回のプログラムで実現したように元の秘密データを計算するのはごく簡単です．

公開されているNを素因数分解してP,Qを得ればRSA暗号を解読することができるのですが，十分に大きな素数P,Qを適切に選んだ場合(大きなP,Qであれば何でもよいというわけではありません)，たとえN=P・Qであると分かっていても，現実的な時間でNをPとQに素因数分解する方法は知られていません．原理的には素数Pの候補を2,3,5,7,...,と順に試していけば，いつかはNを割り切るPが見つかるはずなのですが，P,Qが大きい場合，この素朴な手法では絶望的に長い時間が必要となります．たとえば次のNを素因数分解することができるでしょうか．


 N = 114381625757888867669235779976146612010218296721242362562561842935706935245733897830597123563958705058989075147599290026879543541

じつはこのNはすでに次のように素因数分解されている「とても小さな」数です(129桁)．


 N = P・Q
 P = 3490529510847650949147849619903898133417764638493387843990820577
 Q = 32769132993266709549961988190834461413177642967992942539798288533

現在使われている暗号では，Nとして600桁程度(2048bit)の数を使います．

なお素因数分解の効率的な方法は見つかっていませんが，素因数分解が本当に難しいかどうかは証明されていません．また素因数分解以外の方法ではRSA暗号を解読できないかということも分かっていません．

[参考] 逐次的な繰り返し処理による「べき剰余」計算

次のようにして逐次的に乗算を行うことで，再帰呼び出しを使うことなく「べき剰余」を計算できます．


  def modexp_linear_iter(x,n,m)
    y = 1
    n.times do 
      y = (y*x) % m
    end
    y
  end

この方法はとても単純でわかりやすく，もちろん計算も実行可能です．しかしこの方法では乗算をn回繰り返します．一方で今回作成する「再帰的に乗数を半分ずつにしていく」計算方法では乗算の回数はlog₂(n)程度です． n=153888981(テストプログラムのE)のとき，log₂(n)は約27です．このことから想像できるように，逐次的な繰り返し処理では桁違いに長い時間が必要になります．実際，再帰ではほぼ一瞬で計算が終わるのに対して，逐次的な繰り返し処理では何十秒も必要となります．

[参考] 再帰呼び出しの問題

再帰呼び出しメソッドは問題の構造を素直に記述することで実現できます．これはプログラムを作成する上では大きなメリットです．しかし残念ながら再帰呼び出しの利用には制限があります．

次に簡単な例(factorial_rec.rb)を示します．


  def fact(n)
    puts "fact(#{n})"
    if n == 0
      1
    else
      n*fact(n-1)
    end
  end

  if ARGV.size < 1
    STDERR.puts "#{File.basename($0)} n"
    exit(1)
  end

  n = ARGV.shift.to_i
  y = fact(n)
  puts "#{n}! = #{y}"


  $ ruby factorial_rec.rb 10
  fact(10)
  fact(9)
  fact(8)
  fact(7)
  fact(6)
  fact(5)
  fact(4)
  fact(3)
  fact(2)
  fact(1)
  fact(0)
  10! = 3628800 # OK

  $ ruby factorial_rec.rb 10000
  fact(10000)
  fact(9999)
  fact(9998)
   :
   :
   :
  fact(1816)
  fact(1815)
  factorial_rec.rb:6: stack level too deep  # ERROR!!

このような現象が起きる原因について次のページで説明しています．

メソッドの再帰呼び出し

なお再帰呼出しを使わなければ「10000!」でも問題なく計算できます(factorial_iter.rb)．

[参考] 2進展開に基づく繰り返し処理による「べき剰余」計算

「べき剰余」の計算は再帰呼び出しを使わず，また逐次的な繰り返し処理も使わずに実現できます．

具体的には，乗数nをn=n₀・1+n₁・2+n₂・4+n₃・8+n₄・16+...のように2のべき乗の和に分解しながら(n_i∈{0,1})，つまりnを2進表現に展開しながら， x，x²%m，x⁴%m，...を順次求めることで「べき剰余」を繰り返し処理で計算できます．

[ プログラミング演習(Ruby) > プログラミング演習(Ruby) 課題(2024) > [課題21] べき剰余 ]

xⁿ =	x^k・x^k・x^b	n > 0, n=2k+b, b∈{0,1}
xⁿ =	1	n = 0 (0⁰=1とします)