ES-BASIC #7

  • (by K, 2019.09.21)

(9-1) ベンチマーク#1

  • ES-BASICが他の言語と比べてどのくらい速いのか、もしくはどのくらい遅いのかを比べてみようと思って、どんなプログラムで比較したらいいのかを考えてみました。
  • それで世間の人はどんなプログラムで比較しているのかなと思って調べてみたら、以下のようなものがありました。
    • [1]再帰を使いまくってfibo(42)を計算する。
    • [2]あえて収束の悪い級数を使って(arctan(1)の級数展開)、円周率を計算する。
    • [3]単純な10億回ループ。
    • [4]自然数を与えたときに、その平方根までの奇数で割り算して余りを出す方法で、素数かどうかを判定する関数を書き、その関数を使ってかなり大きな数まで素数かどうかを判定する。
  • これらのベンチマークは以下の問題があると思いました。
    • [1]フィボナッチ数を計算するのなら、普通は再帰なんか使わずに漸化式を使って計算する。そのほうが何億倍も速い。あえて極端に遅いアルゴリズムを選んで、それで優劣を決めるのは実際の性能を反映しているだろうか?フィボナッチ数を再帰で計算しようとすれば、それは相当な回数の関数呼び出しが発生するのだが、それはもはや関数呼び出し命令のベンチマークにしかなっていない気がする。普通のプログラムはそこまで関数呼び出しばかりするわけではなく、他の演算の性能のほうがずっと重要だと思う。
    • [2]普通、円周率を計算するのには、もっと収束の速い級数を使うので、このベンチマークも意味がない。プログラミング言語の作者は、アホなアルゴリズムで性能が出るかどうかは正直どうでもいいだろう。それよりも実際的なアルゴリズムでの性能を考えて言語設計をしたに違いない。力士を集めて優劣をつける時に、なぜ将棋で競わなければいけないのか。相撲でやればいいではないか。だから、あえてダメなアルゴリズムを持ってきて、それでベンチマークを取るのは有意義ではないと思う。
    • [3]単純10億回ループは、まあ最速なアルゴリズムではあると思うが、しかしループの中で何もしないのなら、そもそもループする必要がないと思うので、やはり不自然だ。10億回やらなければいけない処理があって、それで10億回ループしているのなら何の問題もないのだが・・・。
    • [4]単に素数探しをしたいだけなら、エラトステネスのふるいのアルゴリズムを使うのが普通だし、そのほうが何倍も速い。だからこれもベンチマークとしては適切ではないと思う。
  • ということで、まともなアルゴリズムを使っているのにそれなりに処理時間がかかるような処理を探さなければいけないわけです。・・・ということで、マンデルブロ集合の計算かなと思いました。せっかく計算してもグラフィック描画ができないと表示できないのですが、表示部分があるとグラフィックができない言語ではテストできなくなってしまうので、プログラムを改造して、画像の表示はしないことにして、でも演算結果の集計値を計算してprintすることにします。そうでないと、ちゃんと計算できているのかどうかがわからないので。
  • 私の手元でベンチマークしたらこんな感じです:
    言語処理時間備考
    Python 3.7.4324.276秒コードは(9-5)
    Ruby 2.6.498.800秒コードは(9-6)
    ES-BASIC (32bit)2.906秒コードは(9-3)、JITコンパイルに要した時間も含む
    gcc(C言語)(32bit)1.954秒コードは(9-2)、もちろん最適化レベル最大、コンパイル時間含まず
    ES-BASIC (32bit、レジスタ変数利用)1.820秒コードは(9-7)、JITコンパイルに要した時間も含む
    gcc(C言語)(64bit)1.195秒コードは(9-4)、もちろん最適化レベル最大、コンパイル時間含まず
  • ちなみにこのベンチマークプログラムは、以下の画像の各ピクセルの色を決定するための計算をしているだけです。
    • http://k.osask.jp/files/esb20190827b.png

(9-2) C言語版(32bit向け)

  • もちろんC++もこれで行けます。
    #include <stdio.h>
    #include <time.h>
    
    #define MAX 447
    #define STEP (14.0 / 0x100000)
    #define DN 4
    
    int main()
    {
        int x, y, n, sx, sy, sn, h = 0;
        double zx, zy, cx, cy, xx, yy;
        for (y = 0; y < 384; y++) {
            for (x = 0; x < 512; x++) {
                sn = 0;
                for (sx = 0; sx < DN; sx++) {
                    cx = (double) 0x4750 / 0x10000 + (x * DN + sx) * (STEP / DN);
                    for (sy = 0; sy < DN; sy++) {
                        cy = (double) 0x1e8 / 0x10000 + (y * DN + sy) * (STEP / DN);
                        zx = zy = 0;
                        for (n = 0; n < MAX; n++) {
                            xx = zx * zx;
                            yy = zy * zy;
                            if (xx + yy > 4) break;
                            zy = zx * zy * 2 - cy;
                            zx = xx - yy + cx;
                        }
                        sn += n;
                    }
                }
            //  n = sn / (DN * DN);
                n = sn >> 4;
                h += n; // 表示の代わりに集計
            }
        }
        printf("%d\n", h);
        printf("(%.3f[sec])\n", clock() / (double) CLOCKS_PER_SEC);
        return 0;
    }
  • うまく実行できると 28654738 が表示されると思いますが、浮動小数点の実装の微妙な違いで、もしかすると多少は違う値になるかもしれません。でもそれでもこの計算結果を画像化した場合は大差ないと思われますので問題ないと私は考えます。1%以上違う値になるのであれば、それはこのベンチマークを正しく再現できていない可能性が高いです。

(9-3) ES-BASIC版(32bit向け)

  • ES-BASICは現状のバージョンでは浮動小数点演算をサポートしていません(情けない)。整数演算のみです。・・・だから符号部1ビット、整数部7ビット、小数部24ビットの固定小数点を使って演算しています。
  • ES-BASICでも定数に名前を付ける機能(ALIAS命令)を持っているので、もう少しきれいに書くこともできるのですが、このプログラムはその機能を作る前から作っていたこともあり、面倒だったので清書していません。
       1 H=0
       2 FOR Y=0,383
       3   FOR X=0,511
       4     SN=0
       5     FOR SX=0,3
       6       CX=(X*4+SX)*56+4673536
       7       FOR SY=0,3
       8         CY=(Y*4+SY)*56+124928
       9         ZX=0
      10         ZY=0
      11         FOR N=0,446
      12           XX=ZX*>>ZX,24; // これはZXとZXをIMULで計算して64bitの結果を生成した後で、SHRD命令で24ビット右シフトして下位32bitを結果とする命令
      13           YY=ZY*>>ZY,24
      14           T=XX+YY
      15           IF T>0X4000000 GOTO SKIP
      16           ZY=ZY*>>ZX,23
      17           ZX=XX-YY+CX
      18           ZY=ZY-CY
      19         NEXT
      20         LABEL SKIP
      21         SN=SN+N
      22       NEXT
      23     NEXT
      24     N=SN>>4; // N=SN/16 (N=0...447)
      25     H=H+N; // 表示の代わりに集計
      26   NEXT
      27 NEXT
      28 PRINT H
  • これは実行すると 28656617 という結果になります。オリジナルの(9-2)と比べると0.007%ほど大きいことになります。これは固定小数点演算化したことによる誤差だと思われます。問題ないという認識です。

(9-4) C言語版(64bit向け)

  • もちろんC++もこれで行けます。
  • ES-BASICでは「浮動小数点演算が使えないから」という理由で整数演算だけで頑張りましたが、x64ならSHRD命令が使えなくても普通のシフト命令で代用できるので、C言語でも整数命令だけでこのベンチマークをやることができます。そうすると、(9-2)より速くなります。
    #include <stdio.h>
    #include <time.h>
    
    #define MAX 447
    #define DN 4
    
    int main()
    {
        long long x, y, sx, sy, zx, zy, cx, cy, xx, yy;
        int n, sn, h = 0;
        for (y = 0; y < 384; y++) {
            for (x = 0; x < 512; x++) {
                sn = 0;
                for (sx = 0; sx < DN; sx++) {
                    cx = (x * 4 + sx) * 56 + 4673536;
                    for (sy = 0; sy < DN; sy++) {
                        cy = (y * 4 + sy) * 56 + 124928;
                        zx = zy = 0;
                        for (n = 0; n < MAX; n++) {
                            xx = (zx * zx) >> 24;
                            yy = (zy * zy) >> 24;
                            if (xx + yy > 0x4000000) break;
                            zy = (zx * zy) >> 23;
                            zx = xx - yy + cx;
                            zy = zy - cy;
                        }
                        sn += n;
                    }
                }
            //  n = sn / (DN * DN);
                n = sn >> 4;
                h += n; // 表示の代わりに集計
            }
        }
        printf("%d\n", h);
        printf("(%.3f[sec])\n", clock() / (double) CLOCKS_PER_SEC);
        return 0;
    }
  • 結果としては当然のことながら(9-3)と同じく28656617 が表示されます。
  • 32bitのC/C++では、いろいろ工夫してもdoubleを使った場合よりも速くすることはできませんでした(インラインアセンブラでIMUL+SHRDを使えば速くできますが)。

(9-5) Python版

  • 私はPython初心者なので、もっといい書き方があるかもしれません。もっといい書き方を知っている人がいたら、ぜひ教えてほしいです。
    import time
    start = time.time()
    h = 0
    CX0 = 0x4750 / 0x10000
    CY0 = 0x1e8 / 0x10000
    STP = 14.0 / 0x100000 / 4
    for y in range(384):
      for x in range(512):
        sn = 0
        for sx in range(4):
          cx = CX0 + (x * 4 + sx) * STP
          for sy in range(4):
            cy = CY0 + (y * 4 + sy) * STP
            zx = zy = 0
            n = 0
            while n < 447:
              xx = zx * zx
              yy = zy * zy
              if xx + yy > 4: break
              zy = zx * zy * 2 - cy
              zx = xx - yy + cx
              n = n + 1
            sn += n
        n = sn >> 4
        h += n
    print(h)
    print(time.time() - start)
  • これを Python 3.7.4 で実行したところ、28654739になりました。(9-2)とはわずかに違いますが、まあ誤差の範囲だと思います。

(9-6) Ruby版

  • 私はRuby初心者なので、もっといい書き方があるかもしれません。もっといい書き方を知っている人がいたら、ぜひ教えてほしいです。
    start_time = Time.now
    h = 0
    CX0 = 0x4750 / 65536.0
    CY0 = 0x1e8 / 65536.0
    STP = 14.0 / 0x100000 / 4
    y = 0
    while y < 384 do
      x = 0
      while x < 512 do
        sn = 0
        sx = 0
        while sx < 4 do
          cx = CX0 + (x * 4 + sx) * STP
          sy = 0
          while sy < 4 do
            cy = CY0 + (y * 4 + sy) * STP
            zx = zy = 0
            n = 0
            while n < 447 do
              xx = zx * zx
              yy = zy * zy
              if xx + yy > 4 then
                break
              end
              zy = zx * zy * 2 - cy
              zx = xx - yy + cx
              n = n + 1
            end
            sn += n
            sy = sy + 1
          end
          sx = sx + 1
        end
        n = sn >> 4
        h += n
        x = x + 1
      end
      y = y + 1
    end
    p h
    p Time.now - start_time
  • これを ruby 2.6.4 で実行したところ、28654739になりました。(9-5)と同じになりました。

(9-7) ES-BASIC版(32bit向け、レジスタ変数利用)

  • ES-BASICはいくつかの変数をレジスタに割り当てることで、大幅に高速化することができます。これを使うとどのくらい速くなるのか、どのくらい可読性が落ちるのかを確認したかったので、書いてみました。
      1 H=0
      2 $EDI = EDI; // EDIレジスタを自由に使いたいときの前処理
      3 FOR Y=0,383
      4   FOR X=0,511
      5     SN=0
      6     FOR SX=0,3
      7       CX=(X*4+SX)*56+4673536
      8       FOR SY=0,3
      9         CY=(Y*4+SY)*56+124928
     10         ALIAS ZX:ECX, ZY:EBX, XX:ESI, YY:EDI
     11         ZX=0
     12         ZY=0
     13         FOR N=0,446
     14           XX=ZX*>>ZX,24
     15           YY=ZY*>>ZY,24; // 結果はEAXにも入っている.
     16           EAX=EAX+XX
     17           IF EAX>0X4000000 GOTO SKIP
     18           ZY=ZX*>>ZY,23
     19           XX=XX-YY+CX
     20           ZX=XX
     21           ZY=ZY-CY
     22         NEXT
     23         LABEL SKIP
     24         SN=SN+N
     25       NEXT
     26     NEXT
     27     N=SN>>4; // N=SN/16 (N=0...447)
     28     H=H+N; // 表示の代わりに集計
     29   NEXT
     30 NEXT
     31 PRINT H
  • これは実行すると 28656617 という結果になります。(9-3)と同じです。
  • 3行ほど増えていますが、これだけ速くなってくれるのであれば、私なら十分に許容範囲です。

(9-8) ES-BASIC版(64bit向け、レジスタ変数利用)

  • ES-BASICの64bit版ができたので試してみました。
      1 $DI = RDI; // RDIレジスタを自由に使いたいときの前処理
      2 ALIAS ZX:RCX, ZY:RDX, XX:RBX, YY:RSI, N:RDI, CX:R08, CY:R09
      3 ALIAS SX:R10, SY:R11, SN:R12, H:R13, X:R14, Y:R15 
      4 H=0
      5 FOR Y=0,383
      6   FOR X=0,511
      7     SN=0
      8     FOR SX=0,3
      9       CX=(X*4+SX)*56+4673536
     10       FOR SY=0,3
     11         CY=(Y*4+SY)*56+124928
     12         ZX=0
     13         ZY=0
     14         FOR N=0,446
     15           XX=(ZX*ZX)>>24
     16           YY=(ZY*ZY)>>24
     17           RAX=XX+YY
     18           IF RAX>0X4000000 GOTO SKIP
     19           ZY=(ZY*ZX)>>23
     20           XX=XX-YY
     21           ZX=XX+CX
     22           ZY=ZY-CY
     23         NEXT
     24         LABEL SKIP
     25         SN=SN+N
     26       NEXT
     27     NEXT
     28     N=SN>>4; // N=SN/16 (N=0...447)
     29     H=H+N; // 表示の代わりに集計
     30   NEXT
     31 NEXT
     32 PRINT H
  • これは実行すると 28656617 という結果になります。(9-3)と同じです。
  • 3行ほど増えていますが、これだけ速くなってくれるのであれば、私なら十分に許容範囲です。

こめんと欄


コメントお名前NameLink

トップ   編集 凍結 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2019-11-08 (金) 22:26:49 (28d)