Ruby - 多倍長浮動小数点数の加減算！

前回は、C++ による多倍長浮動小数点数同士の加減算について紹介しました。

• C++ - 多倍長浮動小数点数の加減算！

今回は、同じことを Ruby で試してみました。

0. 前提条件

• Linux Mint 14 Nadia (64bit) での作業を想定。
• Ruby 2.0.0-p0 を使用。

1. 考え方

考え方は C++ 版と同じなので、「C++ - 多倍長浮動小数点数の加減算！」を参照のこと。

2. Ruby スクリプト作成

• 今回、 A - C < 0 (A < C) になることは想定していない。

File: add_big_float.rb

#! /usr/local/bin/ruby
#*********************************************
# 浮動小数点数加減算
#*********************************************
#
class AddBigFloat
  D_MAX = 1001  # 有効桁数(小数点以下＋べき指数)

  def initialize
    # 使用する被加減数・加減数を設定（テストなので適当な乱数を使用）
    #   ( A + B, A - C で A > C となることが条件だが
    #     稀に条件に合致しない値が生成されるかもしれない )
    @a, @b, @c = Array.new, Array.new, Array.new
    @a << rand(D_MAX - 1) + 1      # Aのべき指数部
    @b << rand(D_MAX - 1) + 1      # Bのべき指数部
    @c << @a[0] - rand(@a[0]) - 1  # Cのべき指数部
    (D_MAX - 1).times do |_|
      @a << rand(10)               # Aの小数部
      @b << rand(10)               # Bの小数部
      @c << rand(10)               # Cの小数部
    end
  end

  # 計算
  def calc
    # 計算する被加減数・加数・減数を出力
    puts "A ="; fdisp(@a)
    puts "B ="; fdisp(@b)
    puts "C ="; fdisp(@c)
    # 計算(Z = A + B) ・結果出力
    puts "A + B ="
    fdisp(fadd(@a, @b, 1))
    # 計算(Z = A - C) ・結果出力
    puts "A - C ="
    fdisp(fadd(@a, @c, -1))
  rescue => e
    raise
  end

  private

  # ========================================================
  # 加減算
  # ( 固定長整数同士の加減算, 正規化含む )
  #
  #   [引数]   a   : 被加減数配列
  #            b   : 加減数配列
  #            id  : 1: 加算(A + B(, -1: 減算(A - B)
  #   [返り値] z   : 計算結果配列
  # ========================================================
  def add(a, b, id)
    z = Array.new  # 結果格納配列

    begin
      # 計算
      if id >= 0  # 加算
        a.size.times { |i| z << a[i] + b[i] }
      else        # 減算
        a.size.times { |i| z << a[i] - b[i] }
      end
      # 正規化
      z = norm(z)
      return z
    rescue => e
      raise
    end
  end

  # ========================================================
  # 加減算
  # ( 浮動小数点数同士の加減算, 正規化含む )
  #
  #   [引数]   a  : 被加減数配列
  #            b  : 加減数配列
  #            id : 1: 加算(A + B), -1: 減算(A - B)
  #   [返り値] z  : 出力データ配列(Z = A + B or Z = A - B)
  # ========================================================
  def fadd(a, b, id)
    z = Array.new  # 結果格納配列

    begin
      # 計算
      if a[0] >= b[0]  # A のべき指数 >= B のべき指数
        k = a[0] - b[0]
        z = a[0...(k + 1)]
        z += add(a[(k + 1)..-1], b[1..(a.size - 1 + k)], id)
      else             # A のべき指数 <  B のべき指数
        z << b[0]
        k = b[0] - a[0]
        # 位の異なる部分
        if id >= 0  # 加算
          z += b[1..k]
        else        # 減算
          1.upto(k) { |i| z << -b[i] }
        end
        # 位の同じ部分の加減算
        z += add(a[1..(-k - 1)], b[(k + 1)..-1], id)
      end
      # 正規化
      z = [z[0]] + norm(z[1..-1])  # 固定長整数
      z = fnorm(z)                 # 浮動小数点数
      return z
    rescue => e
      raise
    end
  end

  # ========================================================
  # 正規化
  # ( Fixed Normalize(A) by Base-Value )
  #
  #   [引数]   a : 正規化前データ配列
  #   [返り値] a : 正規化後データ配列
  # ========================================================
  def norm(a)
    cr = 0         # 繰り上がり

    begin
      (a.size - 1).downto(1) do |i|
        cr        = a[i] / 10
        a[i]     -= cr * 10
        a[i - 1] += cr
      end
      (a.size - 1).downto(1) do |i|
        if a[i] < 0
          a[i - 1] -= 1
          a[i]     += 10
        end
      end
      return a
    rescue => e
      raise
    end
  end

  # ========================================================
  # 正規化
  # ( Floating Normalize(A) )
  #
  #   [引数]   a : 正規化前データ配列
  #   [返り値] a : 正規化後データ配列
  # ========================================================
  def fnorm(a)
    k = 0  # 上位の 0 の個数

    begin
      # 小数点以下第１位が桁あふれする場合、右に１桁シフト＆べき指数加算
      if a[1] >= 10
        (a.size - 2).downto(2) { |i| a[i + 1] = a[i] }
        a[2] = a[1] % 10
        a[1] = a[1] / 10
        a[0] += 1
      end
      # 正規化前のべき指数を退避
      exp = a[0]
      # 上位の 0 の個数をカウント
      1.upto(a.size) do |i|
        unless a[i] == 0
          k = i - 1
          break
        end
      end
      # 上位の 0 の部分に以降の数字をシフト
      1.upto(a.size - k) { |i| a[i] = a[i + k] }
      # シフト後の下位に 0 をセット
      (a.size - k).upto(a.size - 1) { |i| a[i] = 0 }
      # べき指数セット
      a[0] = exp - k
      return a
    rescue => e
      raise
    end
  end

  # ========================================================
  # 結果出力 (指数表示)
  #
  #   [引数]   s    : 結果出力データ配列
  #   [返り値] なし
  # ========================================================
  def fdisp(s)
    if s[1] < 0
      puts "ANS. < 0\nPlease retry!"
      exit
    end
    print "0."
    # 1行に50桁出力
    1.upto(s.size - 1) do |i|
      print s[i]
      print " "    if i % 10 == 0
      print "\n  " if i % 50 == 0
    end
    # べき指数
    puts s[0] < 0 ? " * 10^(#{s[0]})" : " * 10^#{s[0]}"
    puts
  rescue => e
    raise
  end
end

if __FILE__ == $0
  begin
    # 計算クラスインスタンス化
    obj = AddBigFloat.new
    # 計算
    obj.calc
  rescue => e
    $stderr.puts "[#{e.class}] #{e.message}"
    e.backtrace.each{ |tr| $stderr.puts "\t#{tr}" }
  end
end

• Gist - Ruby script to compute big-digit floats.

3. 実行

まず、実行権限を付与。

$ chmod +x add_big_float.rb

そして、実行。

$ ./add_big_float.rb
A =
0.2783737238 1303344608 9347169684 3282329196 5307669513 
  7488984027 9152724850 2953114516 5239024752 8881375559 
  4776030762 3200933909 2347805946 9027550863 4980016846 
  6437587376 3724863387 8289052079 2546386916 9463209606 
  9430104727 1417817305 6061340308 5886569658 3934807534 
  9022223169 9237629594 8589623932 8572577452 2995044468 
  1434355236 5967148315 0728551183 9326565879 0903552919 
  7193617021 7185307099 3793766010 3818257051 2043669198 
  7660839436 6124646218 0098359769 6767128550 1117072199 
  4079639090 8053259429 4657834539 9199989746 9142523390 
  0927673070 1972138156 9454265642 0547769513 2273878983 
  5658314613 3225502634 8727654863 4352393649 2731691913 
  2857300925 7959917752 5874362423 1196284003 9762780233 
  7490294640 8727767922 5976319116 4839032889 7603116162 
  9036273941 1385721188 1044449806 1580096886 0196620642 
  0616850330 5270598331 7963643161 0714343433 0403405273 
  4531243862 8516049185 0402843189 1629731926 2095630644 
  3475105637 9206505869 5169996557 5174741190 6613444696 
  4361713529 1213966077 1722382995 0672113765 2961341066 
  4279110855 7882795131 3033069262 8216293523 9587644749 
   * 10^473

B =
0.0022590812 3710169960 1089793981 3064983884 8976399063 
  0174741110 7085151138 3069839274 5570373653 4920990302 
  1324232041 2352590657 5877070452 1072953500 2407669857 
  6132879125 2529649102 6269263549 9717346260 0662788373 
  1945326059 4405971259 6353011882 7202514913 1548805895 
  5874073241 7889888186 6570490556 1949198176 0614375744 
  8623094975 0675069190 3405049410 1785796064 3289120677 
  3095849091 2012628735 6174677007 9165521273 3298469487 
  1592269693 2551540460 8420790624 5236765502 2036639189 
  8255945043 2041013035 9744058240 4809179002 9289582357 
  8719417926 8618216275 8595304962 6226968526 4466104488 
  7167567081 5603149454 0633372385 3552837529 8976482766 
  0420057267 1270892710 4581766110 4824554239 7046866365 
  5911003836 6570217932 5684673486 0072372723 3905035754 
  8771140381 0586414934 2960402420 2668236561 2089882061 
  6747293260 9461701944 1168686412 8630504241 5542791687 
  4284354351 7828275583 6839070337 5711975924 4086873321 
  2040262977 9734719019 5731178504 8276668426 2971662580 
  5066388788 2366427242 1882430337 6182753573 8337974385 
  8914788122 5122202068 9952880775 7626255357 6530388836 
   * 10^349

C =
0.0223515925 3368391270 4463511258 6838833759 5344753262 
  9298194326 6849626597 4092467072 1883731160 2936293254 
  9591015466 2862138558 7921550146 9830278274 7052198091 
  3162472250 5165091329 3265422066 8905699592 7889242878 
  1589146326 4595601809 5825637063 8874329409 5548438178 
  7824023837 6842833089 0563265810 5163775672 4359175871 
  7523234617 0939926181 4347155887 0876924642 6704926060 
  4752210118 4714967153 2454327713 5988714622 6605249346 
  9447198350 6644704623 1670956736 6329244574 4852326315 
  1782508575 2386060636 4910666220 6204822996 5681950263 
  0784586598 4507427723 1537926607 8576097458 7197469594 
  4129706598 5813054378 3653792810 3124948815 8397367015 
  3030694408 6362597428 1053575333 6750544320 4213368161 
  0932897812 6921253272 7919690724 6679955380 2214087007 
  7801741483 1836299808 9203585085 2352469586 0418839310 
  4977156521 7471613380 8878610577 3104178078 3844178350 
  8862935375 6050332904 3349661461 0753487181 5552636745 
  6878809641 7832467427 3930317837 8647324176 7295023847 
  8564697338 5348993215 3433521128 8695055944 3107985729 
  3300541411 4102409842 5813750909 6607397015 7712412612 
   * 10^448

A + B =
0.2783737238 1303344608 9347169684 3282329196 5307669513 
  7488984027 9152724850 2953114516 5239024752 8881375559 
  4776030762 3200933909 2347808205 9839921880 4940125826 
  0418893874 7609761027 7352069553 3657095432 0601516590 
  8704661764 5071309404 6363472731 7927804917 4592395241 
  9474330465 2737870361 8447237220 7697830417 2097671394 
  4984326971 2227214593 9101745716 5386006476 2163188220 
  9076337273 2098461979 9689353417 7060046040 0230326247 
  8217034356 4300707655 5843222079 1742196057 0307412704 
  3489817670 4117588341 5335144124 8291191009 7878140857 
  7935589622 3245468003 8941424869 0241024667 2734721062 
  6282838289 8727706298 7917480457 9395597750 5767666319 
  1097781843 6962846710 8232234364 9123145825 6038639764 
  2452917337 7254214533 0465035873 1920593204 7057179500 
  1421629224 8915618836 3810491811 8847223975 2907078818 
  6727332785 9510303018 4329234261 4551000454 8335973740 
  8017251100 1239439688 6157720303 2010790567 7029926684 
  5895372461 5767714857 7231671286 8435687360 8557561565 
  0774576579 5455520356 3409811430 5023896592 8545024973 
  4616682053 3807203818 6354273289 1194266995 8607217866 
   * 10^473

A - C =
0.2783737238 1303344608 9347167449 1689795512 6180624878 
  6363115639 5776771402 7626821534 5806356256 6221634634 
  8068811925 0084904546 3022310036 7480922242 1124137631 
  1422889073 5897392865 8479920454 5321335266 0330276952 
  7223415670 1458538413 1773524417 1253923702 3753849278 
  5315834426 6296674110 4771745692 6188809023 9686138835 
  4853303598 8399904723 3141375951 5864856479 8285409448 
  1604908252 4721040050 1187718488 2806409901 5328344655 
  4889479549 4662380165 5163665297 6932062103 0654755489 
  8405975798 3595810906 2026316714 8342465886 3078874283 
  4305611021 9672481337 4427957796 1887924438 9501563604 
  2997528852 3479630660 1768213566 3692535518 7293855375 
  3576269676 3078333778 9172832116 1755420378 0019969697 
  9956927135 4295725788 9160209787 5057763677 2275836965 
  9963807141 5847699047 2343671788 7431778523 0215728606 
  2108326805 8311994143 4032593389 5062168716 9065316487 
  3473512821 0708210743 2567754559 8092171422 8805197147 
  7328998103 0488350343 1495427769 4210562865 9870705393 
  2577927055 8796293126 9337597348 0938260275 3639806731 
  2166223905 2288364051 4460136257 4075152499 8603386612 
   * 10^473

桁を合わせて加減算してみると、結果が合っていることが分かる。

---

簡単な内容ですが、多倍長演算を行う際に必要になることもあるでしょう。

以上。