公開日:2021/06/12
title.png

16th Season


#16-6 半精度浮動小数点演算ユニット(弐)

 前回基本設計を済ませたHFPU(Half precision FPU)のシミュレーションによる動作検証と、オリジナル16ビットCPUコアMoscoviumへの組み込みを済ませて完成させたいと思います。完成と言うからには、FPGAボード(Sipeed Lichee Tang Primer)に書き込んでの実機検証も済ませます。

 早速、検証方法から検討して行きます。

検証方法の検討

 どんな回路もそうですが、設計するなら検証方法についてもちゃんと検討しておかなければなりません。

 今回のHFPUをどうやって検証したら良いのか?

 PCの浮動小数点数値はfloat型やdouble型だし、half型の期待値を生成する方法を考えないとなりません。float型やdouble型だと、確かにビット数は異なるのだけれど、実行する浮動小数点演算自体は同じことをしている訳で、単に演算精度の問題に帰着するわけです。
 Verilog-HDLにも浮動小数点演算可能な型が有って、real型と呼んでいます。real型の内部は、どうやらdouble型(64ビット浮動小数点数)みたいです。real型を使って期待値を計算させ、HFPUを検証するのが良さそうです。

 検証方針を以下の通りとしました。

a166_verify.png

 ポイントとなる数値の演算をシミュレーション波形で細かくチェックした後、演算対象の数値Ah,Bhをランダム($random)に発生させて、スティミュラスで自動的に期待値比較させます(ロードテスト)。

half型→real型変換task

 以下、Verilog-HDLで作成したhalf型→real型変換taskのコードを示して、その下に解説を加えます。

  task fpu_half2real;
  input [15:0] hx;
  reg sx;
  reg signed [4:0] ex;
  reg [10:0] fx;
    begin
      // half hx[15:0] to real rx conversion


 half型からreal型への変換taskのヘッダ部分です。入力となるhalf型バイナリデータがhx[15:0]です。task内部で使用する符号(sx)や指数(ex[4:0])、仮数(fx[10:0])等のreg宣言があります。
 変換先のreal型変数は、グローバル宣言してあります(real rx;)。

      // conversion
      if (hx[14:10]==5'h0)
        begin
          // zero
          rx=0.0;
        end


 Ex’=0の場合は、ゼロまたは非正規化数ですが、HFPUでは非正規化数には非対応なので、ゼロ(rx=0.0;)を返します。
 なお、Verilog-HDLで数値に小数点(.)を付ける場合、小数点の上位側と下位側両方に数字を付加する必要があります(0.0等)。C言語だと、上位側または下位側のいずれかに数字が有れば良いため、"0."や".0"という記述が可能でした。

      else if (hx[14:10]==5'h1f)
        begin
          // inf or nan
          if (hx[9:0]==10'h0)
            rx=rinf;
          else
            rx=rnan;
        end


 Ex’=31の場合は、INF(無限大)またはNaN(非数)です。別なinitialステートメントの中で、real型のグローバル変数rinfとrnanに対して、「rinf=1.0/0.0;rnan=0.0/0.0;」としてINFとNaNを代入してあります。この辺り、非有限値関係の演算ルールは、今回のhalf型FPUでも踏襲しています。

      else
        begin
          // not zero, finite
          sx=hx[15];
          ex[4:0]=hx[14:10] - 5'd15;
          fx[10:0]={ 1'b1,hx[9:0] };
          rx=fx[10:0];
          if (ex<5'sd0)
            rx=rx/32'h400/(32'h1<<(-ex));
          else
            rx=rx/32'h400*(32'h1<<ex);
        end


 half型の定義に従って、real型の浮動小数点数値に変換しています。
 指数ex[4:0]はsigned宣言してあるので、符号付数(2の補数)として演算されます。ex[4:0]等のベクタ値は、単に"ex"と記載することもできるのですが、"ex"と書くと、ベクタ数なのかスケイラ数なのか?一見して解らないのが嫌いなので、普段は"ex[4:0]"と記述しています。
 ところが、"[4:0]"等のレンジ指定と言う記法をsigned宣言された識別子に使うと、unsignedとして取り扱われてしまうので、要注意です。上コードのif文条件の中で、exをゼロ(5'sd0、5ビット符号付0)と比べて正/負を判定していますが、この部分を「ex[4:0]<5'sd0」や「ex<5'd0(5ビット符号無0)」と書いてしまうと、符号無数値として比較が実行されるので期待通りに動きません。
 筆者だけの問題かもしれませんが……、一応注意喚起しておきます。

      if (hx[15])
        rx=-rx;
    end
  endtask


 最後に符号ビットを確認して、rxの符号を設定して終了です。

real型→half型変換task

 前節に引き続き、Verilog-HDLで作成したreal型→half型変換taskのコードを示して、その下に解説を加えます。

  task fpu_real2half;
  output [15:0] hx;
  reg sx;
  reg signed [5:0] ex;
  reg [10:0] fx;
  integer idx;
    begin
      // real rx to half hx[15:0] conversion


 real型からhalf型への変換taskのヘッダ部分です。入力となるreal型の変数rxは、グローバル宣言してあります。変換後のhalf型バイナリデータをhx[15:0]に出力します。task内部で使用する符号(sx)や指数(ex[5:0])、仮数(fx[10:0])等のreg宣言があります。

      sx=0;
      if (rx<0.0)
        sx=1;


 符号ビット(sx)をrxの正/負に従って決定しました。

      if (rx==0.0)
        hx[15:0]={ sx,15'h0 };


 rx=0.0だった場合のhalf型バイナリデータ=0x0000です。Verilog-HDL(iverilog)が−0.0を返さない・判定できないので、half型の−0=0x8000は出力できません。よって、0x8000との期待値比較ができません。ここだけ注意して使います。

      else if (rx==rinf || rx==-rinf)
        hx[15:0]={ sx,5'd16 + 5'd15,10'h0 };


 rxがINFだった場合のhalf型バイナリデータを返します。
 蛇足ですが、iverilogでは、INF値の比較(rx==rinf)が機能した様です。数学的には、INFとINFの比較はできないのですが……。実はこれが機能しなかった場合、この部分(↑)の記述が成立しなかったので、結果オーライです。

      else
        begin
          // real --> half conversion
          ex[5:0]=6'sd0;
          if (rx<0.0)
            rx=-rx;


 上で処理した以外の有限値だった場合の変換処理が始まります。指数ex[5:0]を初期化して、rxの絶対値を取り、準備します。

          while (rx<1.0 && ex>=-6'sd14)
            begin
              ex[5:0]=ex - 6'sd1;
              rx=rx*2;
            end


 rx<1.0の間、rxを左シフト(×2倍)しながら、exを減算します。正規化処理(の前半)です。

          while (rx>=2.0 && ex<=6'sd15)
            begin
              ex[5:0]=ex + 6'sd1;
              rx=rx/2;
            end


 rx≧2.0の間、rxを右シフト(×1/2倍)しながら、exを加算します。1つ上の処理と合わせて正規化処理が完了します。結果として1.0≦rx<2.0に調整され、指数exが求まりました。

          if (ex>6'sd15)
            begin
              // overflow -> INF
              ex[5:0]=6'sd16;
              fx[10:0]=11'h0;
            end


 ex>15となった場合、half型で表現可能な絶対値最大値を超えているので、INFとします。

          else if (ex<-6'sd14)
            begin
              // underflow -> 0
              ex[5:0]=-6'sd15;
              fx[10:0]=11'h0;
            end


 ex<−14となった場合、half型で表現可能な絶対値最小値を下回っているので、ゼロとします。

          else
            fx[10:0]=(rx - 0.5/1024)*11'h400;


 rxから仮数fx[10:0]を求めました。fx[10:0]は小数点上の"0x1."と合わせて11ビットですので、小数点以下10ビットが仮数Fx’となります。尚、Verilog-HDLだと、real型→整数型への代入では小数点以下を四捨五入してしまうので、HFPUの演算に合わせて、「小数点以下切り捨て」になる様な処理(−0.5/1,024)としています。

          if (fx[10:0]===11'hx)
            begin
              // nan
              sx=1'b1;
              ex[5:0]=6'sd16;
              fx[10:0]=11'h200;
            end


 rxがNaNだった場合、fx[10:0]が不確定値(11'hx)なので、これを判定して、half型のバイナリデータNaNを返します。Verilog-HDLでは、不確定値('hx)やハイインピーダンス値('hz)との比較をする場合、イコール(=)を3つ重ねた演算子(===)を用います。

          hx[15:0]={ sx,ex[4:0] + 5'd15,fx[9:0] };
        end
    end
  endtask


 個別に生成した符号(sx)・指数(ex)・仮数(fx)を組み合わせて、half型バイナリデータをhx[15:0]へ設定して終了します。

単体検証

 これで、real型⇔half型の双方向変換が可能となりました。ただし、Verilog-HDLが−0(負号付きゼロ)を返してくれないため、half型の−0との期待値比較だけは実行できませんでした。

  // check
  err=!( hc[15:0]===hx[15:0] ||
      (hc[15:0]===16'h8000 && hx[15:0]===16'h0000) );


 上コードでは、期待値Ehがhx[15:0]、HFPUによる演算結果Chがhc[15:0]となっています(こちらの図を参照)。Verilog-HDL(iverilog)では、期待値hx[15:0]として0x8000(−0)を返せないため、演算結果hc[15:0]=0x8000(−0)だった時、期待値hx[15:0]=0x0000(0)ならば、pass判定してあります。
 −0を返す演算に関しては、別途個別に検証してありますが、乱数を流し込むロードテストは実行できていません。


 以下、HFPU単体検証(ロードテスト)の状況です。

  検証状況:【 】内はCygwinシェルウインドウからの実行コマンド
  【iverilog -f test_ut_hfpu.lst ; if ( ! $status ) vvp a.out】
    :
  0x7b0d + 0x998d = 0x7b0c (0x7b0c): OK (57760.000000 -0.002710 57728.000000 57759.997290)
  0x8465 + 0x5212 = 0x5211 (0x5211): OK (-0.000067 48.562500 48.531250 48.562433)
  0xe301 + 0xcd0d = 0xe329 (0xe329): OK (-896.500000 -20.203125 -916.500000 -916.703125)
  0xf176 + 0xcd3d = 0xf178 (0xf178): OK (-11184.000000 -20.953125 -11200.000000 -11204.953125)
  0x57ed * 0xf78c = 0xfc00 (0xfc00): OK (126.812500 -30912.000000 -inf -3920028.000000)
  0xe9f9 * 0x24c6 = 0xd320 (0xd320): OK (-3058.000000 0.018646 -57.000000 -57.020203)
  0x84c5 * 0xd2aa = 0x1bf2 (0x1bf2): OK (-0.000073 -53.312500 0.003880 0.003880)
  0xf7e5 * 0x7277 = 0xfc00 (0xfc00): OK (-32336.000000 13240.000000 -inf -428128640.000000)
  0xd612 / 0xdb8f = 0x366c (0x366c): OK (-97.125000 -241.875000 0.401367 0.401550)
  0x69f2 / 0x96ce = 0xfc00 (0xfc00): OK (3044.000000 -0.001661 -inf -1832299.278990)
  0x7ae8 / 0x4ec5 = 0x6814 (0x6814): OK (56576.000000 27.078125 2088.000000 2089.361800)
  0x495c / 0x28bd = 0x5c86 (0x5c86): OK (10.718750 0.037018 289.500000 289.556472)
  total errors: 0 ( 12)


 加減算・乗算・除算をランダムに4回づつ実行した時のシミュレーションログです。
 左から、演算データ1、演算子(+,×,÷)、演算データ2、演算結果、期待値、期待値との照合結果(OK)、各half型バイナリデータをreal型に変換した値、と、並んでいます。inf等、非有限値データとの期待値照合も一致していることが解ります。
 この調子で数を増やして数10万パタン実行してみたところでは、照合エラー(ERR)は発生しませんでした。

 16ビットのhalf型入力データが2個と演算子が3種類の、計216+16×3通りの組み合わせがあります。総当たりでシミュレーションすると、実行できないことは無いかもしれないですね(ちょっと大変だけど)……。ただし、全パタンの検証は実施していませんので、ご承知おきください。

Moscoviumに追加する命令

 HFPUの単体検証まで終了しました。次は、MoscoviumにHFPU命令を追加します。追加する命令は以下の通りとしました。

 詳しい説明は、Moscovium命令マニュアルでご確認ください。
 実装自体は、コプロセッサHFPUに対してコマンドccmd[4:0]を発行して、結果を回収するだけなので、簡単です。特記事項はありません。


 アセンブラの記述を容易にするため、Moscoviumアセンブラ(asm_mcvm.olv)に以下の命令を追加します。

 こちらも、Moscovium命令マニュアルでご確認ください。

half型の大小比較方法

 half型浮動小数点数値同士で大小比較する場合の基本は、hsub命令です。引き算して見て、結果の正/ゼロ/負で判定します。

  例)half型汎用レジスタr0とr1の大小を比較する
  hsub r0,r1
  mov r0,r0  // hsub命令ではsrレジスタの演算フラグが変化しないので必要
  beq equal  // r0=r1の時、分岐
  bpl larger  //r0≧r1の時、分岐
  bmi smaller  // r0<r1の時、分岐
  bne not_equal  // r0≠r1の時、分岐


 mov命令でセットされるsrレジスタのフラグは、ネガティブNフラグとゼロZフラグだけなので、有効な分岐命令はbeq命令bne命令bpl命令bmi命令です。
 r0,r1がINF(無限大)やNaNなど非有限値の場合や、演算結果が−0(負号付きゼロ)だった場合の大小比較には注意して下さい。hsub命令(上のコード例)では大小比較不能なケースがあります。

 特に、負号付きゼロは「−10×0=−0」など、比較的頻繁に出現します。half型バイナリデータでは、half(−0)=0x8000であり、「mov r0,r0」によって、状態レジスタsrのネガティブNフラグが立って、同ゼロZフラグがネゲートされるため、beq命令では分岐せず、bmi命令で分岐します(即ち、扱いは「非ゼロの負数」)。

 この様な特別な値(INFやNaN等を含む)については、以下の様にuint16値として直接比較することを推奨します。

  例)half型汎用レジスタr0がINFかどうか?
  ldwi r1,0x7c00  // 0x7c00はhalf型では+INFの表現
  cmp r0,r1  // r0が+INFか?
  beq r0_is_inf


※half型のバイナリ値は、整数型のcmp命令によって一致・不一致比較はできますが、大小比較はできません。

 C言語などを使ってプログラミングしていると、この辺りの例外処理はC言語が処理してくれています。float型やdouble型だと、演算結果がオーバフロウしてINFが出力されるケースも少ないでしょう。half型では、絶対値が65,504を超えただけでINFが出力されるので、このような非有限値の例外処理を実行すべきケースが増加しています。
 しかも、Moscoviumプログラミングの基本はアセンブラなので、十分注意してプログラムする必要があります。以下、気になる例を示しておきます。

 いずれも、こちらに示した非有限値との加減算が基本となるので、よく理解した上でプログラミングしてください。


※「#16-10 それは比較の問題だ」にて、half型浮動小数点数の比較命令を追加しました。かなり便利になったと思います。

設計データの更新

 HFPU本体と、HFPU命令に対応した新バージョンのMoscovium CPUコア、HFPU搭載版マイコン(MCOC115FP0440)等のデザインデータを以下リンク先に更新しておきました。

btn_moscovium.png


 今回、RAMサイズを最大の40Kバイトまで拡張したのに伴って、コンパイルオプションでRAMサイズを8Kバイト単位で設定できるようにしておきました(↓)。

  // compile option
  //`define MCOC_RAM_16K
  //`define MCOC_RAM_24K
  //`define MCOC_RAM_32K
  `define MCOC_RAM_40K

  //`define MCVM_COPR_NOFPU


 最上位モジュールmcoc115fp0440.vの先頭付近に、上の様な`defineマクロが用意されているので、必要なRAMサイズの`defineを有効化してコンパイルしてください。複数の`defineマクロを有効化した場合、その中で最大のRAMサイズが選択されます。
 すべての`defineを無効化(コメント化)すると、デフォルトの8Kバイト版となります。

 また、HFPUに関しても、他のコプロセッサ(MULC16、DIVC32)と同様に、有効化/無効化できるようになっています。HFPUが不要な時には、MCVM_COPR_NOFPUの`defineを有効化してください。かなり論理規模が大きくなりますが、乗除算コプロセッサ+HFPUというリッチな構成で使用することもできます。MCOC115FP0440では、除算コプロセッサは無効化しておきました(用途が無かったため)。

 HFPU単体でのコンパイル結果もご参考までに記しておきます(↓)。

実機検証を兼ねた、実験

 ESP32のXtensa LX6コアに搭載されたfloat型FPUの実力については、以前実験済みで周知のとおりです(「#16-1 見せてもらおうか、Xtensa LX6コアの性能とやらを」参照)。
 上にリンクした記事では、ESP32(やArduino DUE)にマンデルブロ集合を計算させて、その計算時間を使って浮動小数点演算性能を調べました。

 ESP32では、有り余るRAM容量(520Kバイト)を使って240×240ドット×2バイト=112.5Kバイトの記憶領域を確保して、演算したデータをすべてRAM上に保存しておき、演算終了後、一気にシリアルモニタにデータを送信しました。演算しながらシリアルモニタに送信する方法では、シリアルモニタの通信速度で律速してしまうためです。

 これをHFPU搭載版MCOC115FP0440(RAM 40Kバイト版)でも実行したいのですが、最大40KバイトのRAMではメモリが足りません。そもそも、16ビットCPUであるMoscoviumにはメモリ空間が64Kバイトしかないので、112.5Kバイトのメモリは扱えません。
 止むを得ないので、140ドット×140ドット×2バイト≒38Kバイトの領域を確保して小さいサイズで計算し、計測した時間を240×240ドットに換算して比べてみようと思います。

 以下、Moscoviumアセンブラで作成したマンデルブロ集合計算プログラムをリンクしておきます。

  マンデルブロ集合計算プログラム:mandelbrot.asm


 ESP32やDUEで使ったC言語版(mandelbrot.ino)のCコードをコメントで残して、その下にアセンブラのコードを記述して行くスタイルで、C言語のコードをMoscoviumアセンブラに変換しました(手👐変換!以下、コード例)。

  // mdbr_tr=sq (mdbr_zr) - sq (mdbr_zi) + mdbr_cr;
  ldwi r7,mdbr_zr
  ldw r0,[r7]
  addi r7,2  // mdbr_zi
  ldw r1,[r7]
  ldwi r7,mdbr_cr
  ldw r2,[r7]
  addi r7,2  // mdbr_ci
  ldw r3,[r7]
  mov r4,r0
  hmul r4,r4
  mov r5,r1
  hmul r5,r5
  hsub r4,r5
  hadd r4,r2
  ldwi r7,mdbr_tr
  stw [r7],r4


 HFPUが無ければ、一々浮動小数点演算ルーチンをコールしなければならないですし、もっと大変なコーディングになった筈ですが、上コード例の様に比較的容易に、まるで普通の整数演算のようにコーディングすることが出来ます。


 では実験です。mandelbrot.asmと、MCOC115FP0440を書き込んだFPGAボードと、FTDI等のUSB⇔UART変換モジュールをこちらの様に接続します(下図参考)。

a158_expr.jpg


 PC側では、ESP32の時に使用したスクリプト(mandelbrot_smon.olv)をOlive+で起動して、開いたコマンド入力ウインドウから以下の様に操作します。

  実験サンプル:【 】内はmandelbrot_smon.olvからの送信データ
  【mdbr】
  OK: Mandelbrot (half) 0x0016be42[us]

a161_mandelbrot.jpg

 上の様な画像が表示されたら、正常に動作しています。


 計算に要した時間を以下表にまとめます。

  ESP32 @240MHz DUE @84MHz MCOC115FP @24MHz
実測 @1MHz 実測 @1MHz 実測 @1MHz
マンデルブロ集合half 4.38s 105s
マンデルブロ集合float 1.95s 468s 15.31s 1,286s
マンデルブロ集合double 4.16s 997s 17.84s 1,499s

※ESP32とDUEに関しては240×240ドットの演算に要した時間。MCOC115FPに関しては140×140ドットの演算に要した時間を240×240ドットに換算した(実測=0x0016be42μs≒1.490s @140×140ドット)。

 「#16-1 見せてもらおうか、Xtensa LX6コアの性能とやらを」で計測したESP32とDUEの結果も併せて示してあります。ただし、ESP32はfloat型FPUだし、MCOC115FPはhalf型FPUなので、条件の異なる二者を比較するのはフェアではありません。単なる参考値として比較してみた、と言う事をご理解ください。
 サイクル性能(@1MHz換算値)が……、イヤ、特にコメントを挟みません(∵参考値)ので、上表からご推察ください。


 前々回のロジアナでもハードウェアの実力をご覧頂けたかと思いますが、今回のHFPUにしても、もし、Moscoviumのソフトウェアで演算ルーチンを組んで使うとなったら、1演算当たり数100サイクルはかかります。今回もハードウェアの底力をご覧いただけたものと確信します。余談ですが筆者にとっては、ソフトウェア(アセンブラ)でプログラムを組むより、今回のようにハードウェアで設計してしまった方がずっと楽なのでした。


btn_arduino.png

btn_prev.png btn_next.png

home