esbasic0012
の編集
https://essen.osask.jp/?esbasic0012
[
トップ
] [
編集
|
差分
|
バックアップ
|
添付
|
リロード
] [
新規
|
一覧
|
単語検索
|
最終更新
|
ヘルプ
]
-- 雛形とするページ --
BracketName
EssenRev4
FormattingRules
FrontPage
Help
InterWiki
InterWikiName
InterWikiSandBox
K
MenuBar
PHP
PukiWiki
PukiWiki/1.4
PukiWiki/1.4/Manual
PukiWiki/1.4/Manual/Plugin
PukiWiki/1.4/Manual/Plugin/A-D
PukiWiki/1.4/Manual/Plugin/E-G
PukiWiki/1.4/Manual/Plugin/H-K
PukiWiki/1.4/Manual/Plugin/L-N
PukiWiki/1.4/Manual/Plugin/O-R
PukiWiki/1.4/Manual/Plugin/S-U
PukiWiki/1.4/Manual/Plugin/V-Z
RecentDeleted
SDL2_01
SandBox
WikiEngines
WikiName
WikiWikiWeb
YukiWiki
a21
a21_acl01
a21_bbs01
a21_challengers
a21_count
a21_edu01
a21_edu02
a21_edu03
a21_edu04
a21_edu05
a21_edu06
a21_edu07
a21_edu08
a21_edu09
a21_edu10
a21_edu11
a21_hlx000
a21_hlx001
a21_hlx001_1
a21_hlx001_2
a21_hlx001_3
a21_hlx002
a21_hlx002_1
a21_hlx003
a21_hlx003_1
a21_hlx004_1
a21_memo01
a21_opt
a21_opt02
a21_opt03
a21_p01
a21_special
a21_tl9a
a21_todo
a21_txt01
a21_txt01_10
a21_txt01_1a
a21_txt01_2
a21_txt01_2a
a21_txt01_2b
a21_txt01_3
a21_txt01_4
a21_txt01_5
a21_txt01_6
a21_txt01_6a
a21_txt01_7
a21_txt01_8
a21_txt01_8a
a21_txt01_9
a21_txt01_9a
a21_txt02
a21_txt02_10
a21_txt02_10a
a21_txt02_10b
a21_txt02_11
a21_txt02_11a
a21_txt02_12
a21_txt02_12a
a21_txt02_12b
a21_txt02_1a
a21_txt02_1b
a21_txt02_2
a21_txt02_2a
a21_txt02_3
a21_txt02_3a
a21_txt02_4
a21_txt02_4a
a21_txt02_5
a21_txt02_5a
a21_txt02_6
a21_txt02_6a
a21_txt02_6b
a21_txt02_6b_rev0
a21_txt02_6x
a21_txt02_7
a21_txt02_7a
a21_txt02_8
a21_txt02_8a
a21_txt02_9
a21_txt02_9a
a22_acl2_01
a22_acl2_02
a22_edu12
a22_intro01
a22_intro02
a22_intro03
a22_memman01
a22_memman02
a22_memman03
a22_memman04
a22_memman05
a22_memman06
a22_memman07
a22_memo01
a22_mingw_debug
a22_txt03
a22_txt03_1a
a22_txt03_1b
a22_txt03_2
a22_txt03_2a
a22_ufcs01
a23_bbs
a23_ec001
a23_ec002
a23_intro00
a23_intro000
a23_intro01
a23_intro02
a23_intro03
a23_intro04
a23_intro05
a23_intro06
a23_intro07
a23_intro08
a23_intro09
a23_intro10
a23_intro10wk1
a23_intro10wk2
a23_intro10wk3
a23_intro11
a23_intro12
a23_intro13
a23_intro13wk1
a23_intro14
a23_intro15
a23_intro16
a23_intro90
a23_intro91
a23_neopixel1
a23_os01
a23_useSelfMade
a23_usm001
a23_usm002
a23_usm003
a23_usm004
a23_usm005
a23_usm006
a23_usm007
a23_usm008
a23_usm009
a24_memo01
a24_osc20240310
a24_osc20241026
a24_raspberrypi01
a24_useSelfMade
aclib00
aclib01
aclib02
aclib03
aclib04
aclib05
aclib06
aclib07
aclib08
aclib09
aclib10
aclib11
aclib12
aclib13
aclib14
aclib15
aclib16
aclib17
aclib18
aclib19
aclib20
aclib21
aclib22
aclib23
aclib24
aclib25
aclib_bbs
arm64_01
avm0001
edu0001
edu0002
edu0003
esb02b_hrb
esb_dbg
esbasic0001
esbasic0002
esbasic0003
esbasic0004
esbasic0005
esbasic0006
esbasic0007
esbasic0008
esbasic0009
esbasic0010
esbasic0011
esbasic0012
esbasic0013
esbasic0014
esbasic0015
esbasic0016
esbasic0017
esbasic02a
esc0001
escm0001
essen_hist
esvm0001
esvm0002
esvm0003
esvm0004
esvm0005
esvm0006
esvm_i0
hh4a
idea0001
idea0002
idea0003
impressions
jck_0000
jck_0001
kawai
kbcl0_0000
kbcl0_0001
kbcl0_0002
kbcl0_0003
kbcl0_0004
kbcl0_0005
kbcl0_0006
kbcl0_0007
kclib1_0000
kclib1_0001
kclib1_0002
kclib1_0003
kclib1_0004
kclib1_0005
kclib1_0006
kclib1_0007
kclib1_0008
kclib1_0009
kclib1_0010
kpap0001
members
memo0001
osask4g
osask4g_r2
p20200311a
p20200610a
p20200610b
p20200624a
p20200711a
p20200716a
page0001
page0002
page0003
page0004
page0005
page0006
page0007
page0008
page0009
page0010
page0011
page0012
page0013
page0014
page0015
page0016
page0017
page0018
page0019
page0020
page0021
page0022
page0023
populars
seccamp
seccamp2019
sechack
sechack2019
seclang01
sh3_2020
sh3_2020_kw
sh3_2020_nk
sh3_2021_kw
sh3_2021_nk
sh3_2022_kw
sh3_2023_kw
sh3_2024_kw
sh3_kw_hist
termux001
termux002
text0001
text0001a
text0002
text0002a
text0003
text0004
text0005
text0006
text0006a
text0007
text0008
text0010
text0011
text0012
text0013
text0014
text0015
text0016
text0017
text0018
text0019
text0020
text0021
tl1c
tl2c
tl3c
tl3d
* ES-BASIC #12 -(by [[K]], 2019.11.21) ** はじめに -このページは[[esbasic0007]]の内容が古くなってしまったので、全面的に書き直したものです。 ** (19-1) ベンチマーク ver.1.1 -ES-BASICが他の言語と比べてどのくらい速いのか、もしくはどのくらい遅いのかを比べてみようと思って、どんなプログラムで比較したらいいのかを考えてみました。 -それで世間の人はどんなプログラムで比較しているのかなと思って調べてみたら、以下のようなものがありました。 --[1]再帰を使いまくってfibo(42)を計算する。 --[2]あえて収束の悪い級数を使って(arctan(1)の級数展開)、円周率を計算する。 --[3]単純な10億回ループ。 --[4]自然数を与えたときに、その平方根までの奇数で割り算して余りを出す方法で、素数かどうかを判定する関数を書き、その関数を使ってかなり大きな数まで素数かどうかを判定する。 -これらのベンチマークは以下の''「問題」''があると思いました。 --[1]フィボナッチ数を計算するのなら、普通は再帰なんか使わずに漸化式を使って計算する。そのほうが何億倍も速い。あえて極端に遅いアルゴリズムを選んで、それで優劣を決めるのは実際の性能を反映しているだろうか?フィボナッチ数を再帰で計算しようとすれば、それは相当な回数の関数呼び出しが発生するのだが、それはもはや関数呼び出し命令のベンチマークにしかなっていない気がする。普通のプログラムはそこまで関数呼び出しばかりするわけではなく、他の演算の性能のほうがずっと重要だと思う。 --[2]普通、円周率を計算するのには、もっと収束の速い級数を使うので、このベンチマークも意味がない。プログラミング言語の作者は、アホなアルゴリズムで性能が出るかどうかは正直どうでもいいだろう。それよりも実際的なアルゴリズムでの性能を考えて言語設計をしたに違いない。力士を集めて優劣をつける時に、なぜ将棋で競わなければいけないのか。相撲でやればいいではないか。だから、あえてダメなアルゴリズムを持ってきて、それでベンチマークを取るのは有意義ではないと思う。 --[3]単純10億回ループは、まあ最速なアルゴリズムではあると思うが、しかしループの中で何もしないのなら、そもそもループする必要がないと思うので、やはり不自然だ。10億回やらなければいけない処理があって、それで10億回ループしているのなら何の問題もないのだが・・・。 --[4]単に素数探しをしたいだけなら、エラトステネスのふるいのアルゴリズムを使うのが普通だし、そのほうが何倍も速い。だからこれもベンチマークとしては適切ではないと思う。 -ということで、''「まともなアルゴリズムを使っているのにそれなりに処理時間がかかるような処理」を探さなければいけない''わけです。 ---- -・・・それで思いついたのが「マンデルブロ集合」です。これは具体的には --http://k.osask.jp/files/esb20190827b.png -この絵を計算して描画する処理です。 -プログラムは以下のようになります(blaライブラリを使っています)。 // "19-1a.c" #include <stdio.h> #include <time.h> #define STATIC static #include "bla.c" #define MAX 447 #define STEP (14.0 / 0x100000) #define DN 4 int bla_main(int argc, const char **argv) { int x, y, sx, sy, n, sn, c; float cx, cy, zx, zy, xx, yy, tx, ty; bla_Window *win = bla_openWin(512, 384, "mandel", 1); for (y = 0; y < 384; y++) { for (x = 0; x < 512; x++) { sn = 0; for (sx = 0; sx < DN; sx++) { cx = (float) 0x4750 / 0x10000 + (x * DN + sx) * (STEP / DN); for (sy = 0; sy < DN; sy++) { cy = (float) -0x1e8 / 0x10000 - (y * DN + sy) * (STEP / DN); zx = zy = 0.0; for (n = 0; n < MAX; n++) { xx = zx * zx; yy = zy * zy; if (xx + yy > 4.0) break; tx = xx - yy; ty = zx * zy * 2.0; zx = tx + cx; zy = ty + cy; } sn += n; } } n = sn / (DN * DN); c = bla_rgb(0, 0, 0); if (n < 256) c = bla_rgb(n, 0, 0); else if (n < MAX) c = bla_rgb(255, n - 255, 0); bla_setPix(win, x, y, c); } bla_leapFlushAll(win, 300); } printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); bla_end(); return 0; // dummy } -このプログラムを手元の環境で実行したら、1.815秒でした。 -しかしこのプログラムには問題があります。 --[1]グラフィック描画をサポートしていない言語には移植できない。 --[2]グラフィック描画処理が重い環境では、言語の純粋な計算力の比較ではなくなるのではないか? -ということで、上記のプログラムからグラフィック処理を全部抜き取って、各ピクセルごとのnの値の合計を計算することにしました。 --本当は合計を計算する必然性はなかったのですが、しかし各ピクセルごとに計算した値を使わないと、C言語などでは最適化によってnの計算そのものを省略することがあり得て、それだとさすがに困ってしまうので、合計を計算することにしたのです。 ---- -ということでできあがったのが以下のプログラムです。 // "19-1b.c" #include <stdio.h> #include <time.h> #define MAX 447 #define STEP (14.0 / 0x100000) #define DN 4 int main() { int x, y, sx, sy, n, sn, h = 0; float cx, cy, zx, zy, xx, yy, tx, ty; for (y = 0; y < 384; y++) { for (x = 0; x < 512; x++) { sn = 0; for (sx = 0; sx < DN; sx++) { cx = (float) 0x4750 / 0x10000 + (x * DN + sx) * (STEP / DN); for (sy = 0; sy < DN; sy++) { cy = (float) -0x1e8 / 0x10000 - (y * DN + sy) * (STEP / DN); zx = zy = 0.0; for (n = 0; n < MAX; n++) { xx = zx * zx; yy = zy * zy; if (xx + yy > 4.0) break; tx = xx - yy; ty = zx * zy * 2; zx = tx + cx; zy = ty + cy; } sn += n; } } n = sn / (DN * DN); h += n; } } printf("h=%d\n", h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28656157 time=1.807[sec] ---- -私は当初、この「19-1b.c」で満足していました。しかし、このプログラムを実行するには、float型(つまり浮動小数点演算)を言語がサポートしている必要があって、まあたいていの言語はサポートしているわけですが、しかし自作言語では最初は整数型だけになるのが普通で、そうなると言語の完成度がかなり上がらないとこのベンチマークは試せないということになります。 -それは私には不満でした。できるだけ早期の段階で性能テストできる方が、基本部分の設計を早い段階で練ることができるわけで、後になったらもう直したり試行錯誤したりするのが大変になるかもしれないからです。 -ということで、floatの利用をやめて、cx~tyは小数部が24ビットの固定小数点にしてみました。 -これで足し算や引き算はうまく行くのですが、掛け算は少し問題があります。というのは乗算してから24ビットシフトする必要があるのですが、乗算結果は一時的に32ビットを超えてしまうので、64bit整数にキャストしなければいけません。しかしC言語では、掛け算ごとにキャストして頑張るよりも、そもそもcx~tyを最初から64bit整数にしておく方がずっと高速でした。ということで、こんなプログラムになりました。 // "19-1c.c" #include <stdio.h> #include <time.h> int main() { int x, y, sx, sy, n, sn, h = 0; long long cx, cy, zx, zy, xx, yy, tx, ty; for (y = 0; y < 384; y++) { for (x = 0; x < 512; x++) { sn = 0; for (sx = 0; sx < 4; sx++) { cx = (x * 4 + sx) * 56 + 4673536; for (sy = 0; sy < 4; sy++) { cy = (y * 4 + sy) * -56 - 124928; zx = zy = 0; for (n = 0; n < 447; n++) { xx = (zx * zx) >> 24; yy = (zy * zy) >> 24; if (xx + yy > 0x4000000) break; tx = xx - yy; ty = (zx * zy) >> 23; zx = tx + cx; zy = ty + cy; } sn += n; } } n = sn / 16; h += n; } } printf("h=%d\n", (int) h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28656617 time=1.177[sec] -「19-1b.c」の結果と比べると、hが少しだけ増えています。これはfloatと24ビット固定小数点での精度の差に由来するもので、本質的なものではありません。実際この計算結果で画面描画してみても、違いを感じることはできませんでした。 -そして1.807秒と比較すると1.5倍以上も速くなっていたのです!・・・つまり、C言語のように浮動小数点演算を十分にサポートできている場合であっても、整数演算だけで代用できるのであれば、そのほうがずっと高速だということです。なるほど勉強になります。 -そういう代用を検討する場合は、64ビット整数が使えると非常に重宝するので、64ビット化というのはこういう風に活用すればいいんだなと悟りました(私はそれまで、整数なんて32ビットもあれば十分だ、64ビット整数なんてあっても使いきれない、メモリの無駄だ!とか思っていたのです)。 ---- -この「19-1c.c」での大幅な高速化が気に入って、私はもっと高速化できる方法はないかなと考えてみることにしました。すると、まず、txやtyがいらないと思いました。 zy = ((zy * zx) >> 23) + cy; zx = xx - yy + cx; -とでもすればそれで済む話だだからです。 -さらに、xやyを工夫すれば、sxやsyもなくせると思いました。 -さらに加えて、zx=zy=0で回り始めるのなら、n=0のときのif文でbreakすることなんてありえないわけで、すると zx=cx; zy=cy; になってn=1になるので、この初期条件で回してもいいはずです。その方が速いはずです。 -さらに、cxが変化しなければ、最初のzx*zxの値は変化しないはずです。だからcyのループに入る前にxx0として計算して置いてこれを流用するようにします。 -さらに、 n=1~447 でループするのではなく、n= -446~0 でループするようにします。これならループ中は0との比較をするだけですむのでわずかに高速化できるかもしれません。このつじつまを合わせるために、snの初期値は0ではなく、447*16=7152にします。 -ええとなぜこのような努力をするのかというと、「まともなアルゴリズムを使う」ことがベンチマークでは何よりも大事だと私は思っているからです。良いアルゴリズムが使えるのに、それを使わずに計算速度の優劣を競っても、有意ではないと思っているのです。 -これらすべての変更を盛り込んだ究極形が以下の「19-1d.c」です。 // "19-1d.c" #include <stdio.h> #include <time.h> int main() { int n, sn, h = 0; long long x, y, cx, cy, zx, zy, xx, yy, xx0; for (y = -125152; y > -211168; y -= 224) { for (x = 4673760; x < 4788448; x += 224) { sn = 7152; for (cx = x - 224; cx < x; cx += 56) { xx0 = (cx * cx) >> 24; for (cy = y + 224; cy > y; cy -= 56) { zy = cy; zx = cx; xx = xx0; yy = (cy * cy) >> 24; n = -446; do { zy = (zy * zx) >> 23; zx = xx + cx - yy; zy = zy + cy; n++; if (n >= 0) goto skip; xx = (zx * zx) >> 24; yy = (zy * zy) >> 24; } while (xx + yy <= 0x4000000); sn += n; skip: ; } } h += sn >> 4; } } printf("h=%d\n", h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28656617 time=1.115[sec] -hの値が全く変わっていないので、正しく計算できていることがよくわかります。「19-1c」と比較すると1.055倍くらい高速化できたことになります。 ---- -参考までに、floatを使いつつ精一杯高速化するとこうなります。これはC言語においては全然速くはありませんが、一時的に64ビットになるような乗算が不得意なら、こちらのアルゴリズムの方が速いかもしれません。 // "19-1e.c" #include <stdio.h> #include <time.h> int main() { int x, y, ix, iy, n, sn, h = 0; float cx, cy, zx, zy, xx, yy, xx0; for (y = -125152; y > -211168; y -= 224) { for (x = 4673760; x < 4788448; x += 224) { sn = 7152; for (ix = x - 224; ix < x; ix += 56) { cx = ix * (1.0 / 16777216.0); xx0 = cx * cx; for (iy = y + 224; iy > y; iy -= 56) { cy = iy * (1.0 / 16777216.0); zy = cy; zx = cx; xx = xx0; yy = cy * cy; n = -446; do { zy = zy * zx * 2.0; zx = xx + cx - yy; zy = zy + cy; n++; if (n >= 0) goto skip; xx = zx * zx; yy = zy * zy; } while (xx + yy <= 4.0); sn += n; skip: ; } } h += sn >> 4; } } printf("h=%d\n", h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28655896 time=1.761[sec] -hの値がまた少し変わっていますが、これもおそらく精度に由来する誤差であり、気にしなくてよいと思われます。 ** (19-2) gccの最適化について -今回、19-1b → 19-1c → 19-1d と改良してきたわけですが、この改良のつど、処理時間は確実に短くなっており、つまり人間のセンスによる最適化は、コンパイラによる自動的な最適化を十分に上回っていたことになります。もし上回ることができていなければ、コンパイラは最適化の結果として同じバイナリにたどり着いて同じ処理時間になっていたはずなのです。そうならなかったのは、コンパイラの最適化は、まだまだ人間には及ばないという証拠なのです。 -だから私はES-BASICにおいて言語が最適化を頑張ることを放棄して、人間による最適化を支援する方向で言語を設計しているわけです。 ** つづく -[[esbasic0013]]へつづく * こめんと欄 #comment
タイムスタンプを変更しない
* ES-BASIC #12 -(by [[K]], 2019.11.21) ** はじめに -このページは[[esbasic0007]]の内容が古くなってしまったので、全面的に書き直したものです。 ** (19-1) ベンチマーク ver.1.1 -ES-BASICが他の言語と比べてどのくらい速いのか、もしくはどのくらい遅いのかを比べてみようと思って、どんなプログラムで比較したらいいのかを考えてみました。 -それで世間の人はどんなプログラムで比較しているのかなと思って調べてみたら、以下のようなものがありました。 --[1]再帰を使いまくってfibo(42)を計算する。 --[2]あえて収束の悪い級数を使って(arctan(1)の級数展開)、円周率を計算する。 --[3]単純な10億回ループ。 --[4]自然数を与えたときに、その平方根までの奇数で割り算して余りを出す方法で、素数かどうかを判定する関数を書き、その関数を使ってかなり大きな数まで素数かどうかを判定する。 -これらのベンチマークは以下の''「問題」''があると思いました。 --[1]フィボナッチ数を計算するのなら、普通は再帰なんか使わずに漸化式を使って計算する。そのほうが何億倍も速い。あえて極端に遅いアルゴリズムを選んで、それで優劣を決めるのは実際の性能を反映しているだろうか?フィボナッチ数を再帰で計算しようとすれば、それは相当な回数の関数呼び出しが発生するのだが、それはもはや関数呼び出し命令のベンチマークにしかなっていない気がする。普通のプログラムはそこまで関数呼び出しばかりするわけではなく、他の演算の性能のほうがずっと重要だと思う。 --[2]普通、円周率を計算するのには、もっと収束の速い級数を使うので、このベンチマークも意味がない。プログラミング言語の作者は、アホなアルゴリズムで性能が出るかどうかは正直どうでもいいだろう。それよりも実際的なアルゴリズムでの性能を考えて言語設計をしたに違いない。力士を集めて優劣をつける時に、なぜ将棋で競わなければいけないのか。相撲でやればいいではないか。だから、あえてダメなアルゴリズムを持ってきて、それでベンチマークを取るのは有意義ではないと思う。 --[3]単純10億回ループは、まあ最速なアルゴリズムではあると思うが、しかしループの中で何もしないのなら、そもそもループする必要がないと思うので、やはり不自然だ。10億回やらなければいけない処理があって、それで10億回ループしているのなら何の問題もないのだが・・・。 --[4]単に素数探しをしたいだけなら、エラトステネスのふるいのアルゴリズムを使うのが普通だし、そのほうが何倍も速い。だからこれもベンチマークとしては適切ではないと思う。 -ということで、''「まともなアルゴリズムを使っているのにそれなりに処理時間がかかるような処理」を探さなければいけない''わけです。 ---- -・・・それで思いついたのが「マンデルブロ集合」です。これは具体的には --http://k.osask.jp/files/esb20190827b.png -この絵を計算して描画する処理です。 -プログラムは以下のようになります(blaライブラリを使っています)。 // "19-1a.c" #include <stdio.h> #include <time.h> #define STATIC static #include "bla.c" #define MAX 447 #define STEP (14.0 / 0x100000) #define DN 4 int bla_main(int argc, const char **argv) { int x, y, sx, sy, n, sn, c; float cx, cy, zx, zy, xx, yy, tx, ty; bla_Window *win = bla_openWin(512, 384, "mandel", 1); for (y = 0; y < 384; y++) { for (x = 0; x < 512; x++) { sn = 0; for (sx = 0; sx < DN; sx++) { cx = (float) 0x4750 / 0x10000 + (x * DN + sx) * (STEP / DN); for (sy = 0; sy < DN; sy++) { cy = (float) -0x1e8 / 0x10000 - (y * DN + sy) * (STEP / DN); zx = zy = 0.0; for (n = 0; n < MAX; n++) { xx = zx * zx; yy = zy * zy; if (xx + yy > 4.0) break; tx = xx - yy; ty = zx * zy * 2.0; zx = tx + cx; zy = ty + cy; } sn += n; } } n = sn / (DN * DN); c = bla_rgb(0, 0, 0); if (n < 256) c = bla_rgb(n, 0, 0); else if (n < MAX) c = bla_rgb(255, n - 255, 0); bla_setPix(win, x, y, c); } bla_leapFlushAll(win, 300); } printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); bla_end(); return 0; // dummy } -このプログラムを手元の環境で実行したら、1.815秒でした。 -しかしこのプログラムには問題があります。 --[1]グラフィック描画をサポートしていない言語には移植できない。 --[2]グラフィック描画処理が重い環境では、言語の純粋な計算力の比較ではなくなるのではないか? -ということで、上記のプログラムからグラフィック処理を全部抜き取って、各ピクセルごとのnの値の合計を計算することにしました。 --本当は合計を計算する必然性はなかったのですが、しかし各ピクセルごとに計算した値を使わないと、C言語などでは最適化によってnの計算そのものを省略することがあり得て、それだとさすがに困ってしまうので、合計を計算することにしたのです。 ---- -ということでできあがったのが以下のプログラムです。 // "19-1b.c" #include <stdio.h> #include <time.h> #define MAX 447 #define STEP (14.0 / 0x100000) #define DN 4 int main() { int x, y, sx, sy, n, sn, h = 0; float cx, cy, zx, zy, xx, yy, tx, ty; for (y = 0; y < 384; y++) { for (x = 0; x < 512; x++) { sn = 0; for (sx = 0; sx < DN; sx++) { cx = (float) 0x4750 / 0x10000 + (x * DN + sx) * (STEP / DN); for (sy = 0; sy < DN; sy++) { cy = (float) -0x1e8 / 0x10000 - (y * DN + sy) * (STEP / DN); zx = zy = 0.0; for (n = 0; n < MAX; n++) { xx = zx * zx; yy = zy * zy; if (xx + yy > 4.0) break; tx = xx - yy; ty = zx * zy * 2; zx = tx + cx; zy = ty + cy; } sn += n; } } n = sn / (DN * DN); h += n; } } printf("h=%d\n", h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28656157 time=1.807[sec] ---- -私は当初、この「19-1b.c」で満足していました。しかし、このプログラムを実行するには、float型(つまり浮動小数点演算)を言語がサポートしている必要があって、まあたいていの言語はサポートしているわけですが、しかし自作言語では最初は整数型だけになるのが普通で、そうなると言語の完成度がかなり上がらないとこのベンチマークは試せないということになります。 -それは私には不満でした。できるだけ早期の段階で性能テストできる方が、基本部分の設計を早い段階で練ることができるわけで、後になったらもう直したり試行錯誤したりするのが大変になるかもしれないからです。 -ということで、floatの利用をやめて、cx~tyは小数部が24ビットの固定小数点にしてみました。 -これで足し算や引き算はうまく行くのですが、掛け算は少し問題があります。というのは乗算してから24ビットシフトする必要があるのですが、乗算結果は一時的に32ビットを超えてしまうので、64bit整数にキャストしなければいけません。しかしC言語では、掛け算ごとにキャストして頑張るよりも、そもそもcx~tyを最初から64bit整数にしておく方がずっと高速でした。ということで、こんなプログラムになりました。 // "19-1c.c" #include <stdio.h> #include <time.h> int main() { int x, y, sx, sy, n, sn, h = 0; long long cx, cy, zx, zy, xx, yy, tx, ty; for (y = 0; y < 384; y++) { for (x = 0; x < 512; x++) { sn = 0; for (sx = 0; sx < 4; sx++) { cx = (x * 4 + sx) * 56 + 4673536; for (sy = 0; sy < 4; sy++) { cy = (y * 4 + sy) * -56 - 124928; zx = zy = 0; for (n = 0; n < 447; n++) { xx = (zx * zx) >> 24; yy = (zy * zy) >> 24; if (xx + yy > 0x4000000) break; tx = xx - yy; ty = (zx * zy) >> 23; zx = tx + cx; zy = ty + cy; } sn += n; } } n = sn / 16; h += n; } } printf("h=%d\n", (int) h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28656617 time=1.177[sec] -「19-1b.c」の結果と比べると、hが少しだけ増えています。これはfloatと24ビット固定小数点での精度の差に由来するもので、本質的なものではありません。実際この計算結果で画面描画してみても、違いを感じることはできませんでした。 -そして1.807秒と比較すると1.5倍以上も速くなっていたのです!・・・つまり、C言語のように浮動小数点演算を十分にサポートできている場合であっても、整数演算だけで代用できるのであれば、そのほうがずっと高速だということです。なるほど勉強になります。 -そういう代用を検討する場合は、64ビット整数が使えると非常に重宝するので、64ビット化というのはこういう風に活用すればいいんだなと悟りました(私はそれまで、整数なんて32ビットもあれば十分だ、64ビット整数なんてあっても使いきれない、メモリの無駄だ!とか思っていたのです)。 ---- -この「19-1c.c」での大幅な高速化が気に入って、私はもっと高速化できる方法はないかなと考えてみることにしました。すると、まず、txやtyがいらないと思いました。 zy = ((zy * zx) >> 23) + cy; zx = xx - yy + cx; -とでもすればそれで済む話だだからです。 -さらに、xやyを工夫すれば、sxやsyもなくせると思いました。 -さらに加えて、zx=zy=0で回り始めるのなら、n=0のときのif文でbreakすることなんてありえないわけで、すると zx=cx; zy=cy; になってn=1になるので、この初期条件で回してもいいはずです。その方が速いはずです。 -さらに、cxが変化しなければ、最初のzx*zxの値は変化しないはずです。だからcyのループに入る前にxx0として計算して置いてこれを流用するようにします。 -さらに、 n=1~447 でループするのではなく、n= -446~0 でループするようにします。これならループ中は0との比較をするだけですむのでわずかに高速化できるかもしれません。このつじつまを合わせるために、snの初期値は0ではなく、447*16=7152にします。 -ええとなぜこのような努力をするのかというと、「まともなアルゴリズムを使う」ことがベンチマークでは何よりも大事だと私は思っているからです。良いアルゴリズムが使えるのに、それを使わずに計算速度の優劣を競っても、有意ではないと思っているのです。 -これらすべての変更を盛り込んだ究極形が以下の「19-1d.c」です。 // "19-1d.c" #include <stdio.h> #include <time.h> int main() { int n, sn, h = 0; long long x, y, cx, cy, zx, zy, xx, yy, xx0; for (y = -125152; y > -211168; y -= 224) { for (x = 4673760; x < 4788448; x += 224) { sn = 7152; for (cx = x - 224; cx < x; cx += 56) { xx0 = (cx * cx) >> 24; for (cy = y + 224; cy > y; cy -= 56) { zy = cy; zx = cx; xx = xx0; yy = (cy * cy) >> 24; n = -446; do { zy = (zy * zx) >> 23; zx = xx + cx - yy; zy = zy + cy; n++; if (n >= 0) goto skip; xx = (zx * zx) >> 24; yy = (zy * zy) >> 24; } while (xx + yy <= 0x4000000); sn += n; skip: ; } } h += sn >> 4; } } printf("h=%d\n", h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28656617 time=1.115[sec] -hの値が全く変わっていないので、正しく計算できていることがよくわかります。「19-1c」と比較すると1.055倍くらい高速化できたことになります。 ---- -参考までに、floatを使いつつ精一杯高速化するとこうなります。これはC言語においては全然速くはありませんが、一時的に64ビットになるような乗算が不得意なら、こちらのアルゴリズムの方が速いかもしれません。 // "19-1e.c" #include <stdio.h> #include <time.h> int main() { int x, y, ix, iy, n, sn, h = 0; float cx, cy, zx, zy, xx, yy, xx0; for (y = -125152; y > -211168; y -= 224) { for (x = 4673760; x < 4788448; x += 224) { sn = 7152; for (ix = x - 224; ix < x; ix += 56) { cx = ix * (1.0 / 16777216.0); xx0 = cx * cx; for (iy = y + 224; iy > y; iy -= 56) { cy = iy * (1.0 / 16777216.0); zy = cy; zx = cx; xx = xx0; yy = cy * cy; n = -446; do { zy = zy * zx * 2.0; zx = xx + cx - yy; zy = zy + cy; n++; if (n >= 0) goto skip; xx = zx * zx; yy = zy * zy; } while (xx + yy <= 4.0); sn += n; skip: ; } } h += sn >> 4; } } printf("h=%d\n", h); printf("time=%.3f[sec]\n", (float) clock() / CLOCKS_PER_SEC); return 0; // dummy } -これを実行すると、手元の環境では以下のようになりました。 h=28655896 time=1.761[sec] -hの値がまた少し変わっていますが、これもおそらく精度に由来する誤差であり、気にしなくてよいと思われます。 ** (19-2) gccの最適化について -今回、19-1b → 19-1c → 19-1d と改良してきたわけですが、この改良のつど、処理時間は確実に短くなっており、つまり人間のセンスによる最適化は、コンパイラによる自動的な最適化を十分に上回っていたことになります。もし上回ることができていなければ、コンパイラは最適化の結果として同じバイナリにたどり着いて同じ処理時間になっていたはずなのです。そうならなかったのは、コンパイラの最適化は、まだまだ人間には及ばないという証拠なのです。 -だから私はES-BASICにおいて言語が最適化を頑張ることを放棄して、人間による最適化を支援する方向で言語を設計しているわけです。 ** つづく -[[esbasic0013]]へつづく * こめんと欄 #comment
テキスト整形のルールを表示する