Z80での高速な乗算処理 [Z80]
Twitter のタイムラインで見かけた
の等式を利用し、テーブルを使って二乗計算を速く行うことで乗算処理を高速にできる・・と言うことがずっと頭の片隅に残っていました。
時間もあるので少し検討してみたいと思います。以降上記の手法をプラマイ二乗法と記します。
★追記 2023/02/15
「Z80での高速な乗算処理（その２）最適化」の記事に最適化して更に高速化したことについて記載しました。
[TOP] [ 前へ ] 連載記事一覧 [ 次へ ]
|
(x+y)^2 - (x-y)^2 = 4*x*y ⇒ x*y = ((x+y)^2 - (x-y)^2) / 4 |
|---|
の等式を利用し、テーブルを使って二乗計算を速く行うことで乗算処理を高速にできる・・と言うことがずっと頭の片隅に残っていました。
時間もあるので少し検討してみたいと思います。以降上記の手法をプラマイ二乗法と記します。
- プラマイ二乗法の検討
8bit の二乗は最大で 16bit になるので二乗の結果を収めた 512 バイトのテーブルが必要になります。上位バイトと下位バイトを必ずペアでアクセスするので上位アドレスの使い回しができるように下位バイトのテーブルと上位バイトのテーブルに分けました。
コードを考える上で問題となるのが (x + y) が最大 9 ビットになることです。9 ビットになった場合、MSBとそれ以外の 8 ビットを分離して処理しようとすると
(2^8 + x)^2
= 2^16 + 2^8*2*x + x*x
= 2^16 + 2^9*x + x*x
となり、x の二乗を求めた後に 9 bit シフトして加算する必要があり、その後 4 で割るので 18 ビットの演算が必要になります。
場合分けしてもできそうですが、面倒になる根源は x+y という項があることです。そこで 加算の二乗を使わない方法を考えてみました。
- マイナス二乗法
- 各方式の比較
冒頭で紹介した Twitter メッセージには URL も書いてあったのでプラマイ二乗法の処理を確認して見ました。
二乗テーブルの持ち方は前項で書いたマイナス二乗法と同じですね。x+y が 9 bit になった場合の対処として偶数と奇数の場合に場合分けして処理していました。処理のステート数については「CPC Cycles: 136-172 (154 on average)」と書いてあります。
三つの方式のステート数の一覧が下表になります。
★変更 2023/02/13 Normal method のステート数を修正No. Method Max Min Average 1 Normal method 363 315 339.0 2 Minus square method 154 151 152.5 3 Plamai square method 172 136 154.0
今回作ったマイナス二乗法はノーマル処理の約２倍の速さと言うことになります。
ネットにあったプラマイ二乗法の処理を下に示しますが、ステート数を数えると Max:158, Min:129, Ave:143.5 になります。
下記のソースでは引数とリターン値を入れるレジスタのアサインを不自然に感じますが、本来の処理は引数とリターン値を入れるレジスタがリーズナブルな物であり、その分のステート数を加えたものが上記の表内の数値になるのかもしれません（マイナス二乗法でもリターン値をHLに揃えなくてもよければそれだけで数ステート短くなります）。
★追記 202/02/19 {
プラマイ二乗法が CPCWiKi に追加される際の議論を見つけました。この中でシフト命令等を変更してより高速化されています。最終ソースのステート数が異なるのはこの為かも知れません。
}
ネットで公開されていた高速乗算処理（Z80アセンブラ） ; Web version finded in Internet ; A,B <- data0,data1 ; A:E -> data0*data1 019A B8 IMul: cp b ; 4 019B 30 03 jr nc,l1 ; 12/7 019D 5F ld e,a ; 4 019E 78 ld a,b ; 4 019F 43 ld b,e ; 4 01A0 4F l1: ld c,a ; 4 01A1 90 sub b ; 4 01A2 1F rra ; 4 01A3 57 ld d,a ; 4 01A4 79 ld a,c ; 4 01A5 80 add a,b ; 4 01A6 1F rra ; 4 01A7 6F ld l,a ; 4 01A8 26 02 ld h,high(MULTBL) ; 7 01AA 7E ld a,(hl) ; 7 01AB 5D ld e,l ; 4 01AC 6A ld l,d ; 4 01AD 30 0F jr nc,l2 ; 12/7 01AF 96 sub (hl) ; 7 odd 01B0 6B ld l,e ; 4 01B1 5F ld e,a ; 4 01B2 24 inc h ; 4 loads high(sqrhi) 01B3 7E ld a,(hl) ; 7 01B4 6A ld l,d ; 4 01B5 9E sbc a,(hl) ; 7 01B6 57 ld d,a ; 4 01B7 7B ld a,e ; 4 01B8 80 add a,b ; 4 01B9 5F ld e,a ; 4 01BA 7A ld a,d ; 4 01BB CE 00 adc a,0 ; 7 01BD C9 ret ; 10 01BE 96 l2: sub (hl) ; 7 even 01BF 6B ld l,e ; 4 01C0 5F ld e,a ; 4 01C1 24 inc h ; 4 01C2 7E ld a,(hl) ; 7 01C3 6A ld l,d ; 4 01C4 9E sbc a,(hl) ; 7 01C5 C9 ret ; 10 ; Max 158 = 4+(7+4+4+4)+4+4+4+4+4+4+4+4+7+7+4+4+(7+7+4+4+4+7+4+7+4+4+4+4+4+7+10) ; Min 129 = 4+(12)+4+4+4+4+4+4+4+4+7+7+4+4+(12+7+4+4+4+7+4+7+10) ; Ave 143.5
- 応用編
8 ビットの高速な乗算処理ができたので「ＧＡＭＥ言語での ASCIIART の高速化」の記事で書いた ASCIIART に適用してみました。
先の記事では 16 ビットの固定小数点の乗算処理を作って高速化しましたが、今回は手抜きをして次の方法で 8 ビット乗算で 16 ビット乗算を作ってみました。
(256*x1 + x0)*(256*y1 + y0)
= 2^16*x1y1 + 2^8*(x0*y1 + x1*y0) + x0*y0
従来のアセンブラで書いた ASCIIART（asciiar2）を実行後に、今回の乗算処理を組み込んだ ASCIIART（fasciiar）を実行した画面が、下図になります。実行環境は Z80 が 12MHz で動作している環境です。
TeraTerm のマクロで実行時間を計測した結果も表示されていますが、それを見ると実行時間はほぼ同じ（寧ろ若干遅くなっている）ですね。16 ビット化するための上記の展開式の処理のオーバーヘッドが元々の 2 倍速を食いつぶしてしまったようです。
★追記 2023/02/17 {
上記で使用した固定小数点演算では 2byte x 2byte の結果の 4byte の中央の 2byte のみ必要なため、内部では 24 bit の計算をしているので通常のオーバーフロー無しの乗算処理（結果が 32 bit）よりは速いものです。
}
それにしても意外だったのは新乗算処理を組み込んだ後、一発で ASCIIART が表示されたことです。気分のいいところでそろそろ筆を置きたいと思います。
ASCIIART 実行結果
|
(x - y)^2 = x^2 -2*x*y + y^2 ⇒ x*y = (x^2 + y^2 - (x - y)^2) / 2 |
|---|
アセンブラで書いたソースが下記になります。それほど最適化していないのでまだのりしろがあると思います。
| 今回考案したマイナス二乗法での乗算処理（Z80アセンブラ） |
|
参考として、テーブルを使わない一般的な乗算処理は下記のようになるのではないかと思います。
| 一般的な unsigned 8 ビットの乗算処理（Z80アセンブラ） |
|
★追記 2023/02/15
「Z80での高速な乗算処理（その２）最適化」の記事に最適化して更に高速化したことについて記載しました。
[TOP] [ 前へ ] 連載記事一覧 [ 次へ ]
)
){kind=link}