esbasic0012
をテンプレートにして作成
[
トップ
] [
新規
|
一覧
|
単語検索
|
最終更新
|
ヘルプ
]
開始行:
* 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億回ループは、まあ最速なアルゴリズムではあると...
--[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 ...
for (sy = 0; sy < DN; sy++) {
cy = (float) -0x1e8 / 0x10000 - (y *...
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_...
bla_end();
return 0; // dummy
}
-このプログラムを手元の環境で実行したら、1.815秒でした。
-しかしこのプログラムには問題があります。
--[1]グラフィック描画をサポートしていない言語には移植でき...
--[2]グラフィック描画処理が重い環境では、言語の純粋な計算...
-ということで、上記のプログラムからグラフィック処理を全部...
--本当は合計を計算する必然性はなかったのですが、しかし各...
----
-ということでできあがったのが以下のプログラムです。
// "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 ...
for (sy = 0; sy < DN; sy++) {
cy = (float) -0x1e8 / 0x10000 - (y *...
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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28656157
time=1.807[sec]
----
-私は当初、この「19-1b.c」で満足していました。しかし、こ...
-それは私には不満でした。できるだけ早期の段階で性能テスト...
-ということで、floatの利用をやめて、cx~tyは小数部が24ビ...
-これで足し算や引き算はうまく行くのですが、掛け算は少し問...
// "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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28656617
time=1.177[sec]
-「19-1b.c」の結果と比べると、hが少しだけ増えています。こ...
-そして1.807秒と比較すると1.5倍以上も速くなっていたのです...
-そういう代用を検討する場合は、64ビット整数が使えると非常...
----
-この「19-1c.c」での大幅な高速化が気に入って、私はもっと...
zy = ((zy * zx) >> 23) + cy;
zx = xx - yy + cx;
-とでもすればそれで済む話だだからです。
-さらに、xやyを工夫すれば、sxやsyもなくせると思いました。
-さらに加えて、zx=zy=0で回り始めるのなら、n=0のときのif文...
-さらに、cxが変化しなければ、最初のzx*zxの値は変化しない...
-さらに、 n=1~447 でループするのではなく、n= -446~0 で...
-ええとなぜこのような努力をするのかというと、「まともなア...
-これらすべての変更を盛り込んだ究極形が以下の「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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28656617
time=1.115[sec]
-hの値が全く変わっていないので、正しく計算できていること...
----
-参考までに、floatを使いつつ精一杯高速化するとこうなりま...
// "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_...
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億回ループは、まあ最速なアルゴリズムではあると...
--[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 ...
for (sy = 0; sy < DN; sy++) {
cy = (float) -0x1e8 / 0x10000 - (y *...
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_...
bla_end();
return 0; // dummy
}
-このプログラムを手元の環境で実行したら、1.815秒でした。
-しかしこのプログラムには問題があります。
--[1]グラフィック描画をサポートしていない言語には移植でき...
--[2]グラフィック描画処理が重い環境では、言語の純粋な計算...
-ということで、上記のプログラムからグラフィック処理を全部...
--本当は合計を計算する必然性はなかったのですが、しかし各...
----
-ということでできあがったのが以下のプログラムです。
// "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 ...
for (sy = 0; sy < DN; sy++) {
cy = (float) -0x1e8 / 0x10000 - (y *...
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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28656157
time=1.807[sec]
----
-私は当初、この「19-1b.c」で満足していました。しかし、こ...
-それは私には不満でした。できるだけ早期の段階で性能テスト...
-ということで、floatの利用をやめて、cx~tyは小数部が24ビ...
-これで足し算や引き算はうまく行くのですが、掛け算は少し問...
// "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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28656617
time=1.177[sec]
-「19-1b.c」の結果と比べると、hが少しだけ増えています。こ...
-そして1.807秒と比較すると1.5倍以上も速くなっていたのです...
-そういう代用を検討する場合は、64ビット整数が使えると非常...
----
-この「19-1c.c」での大幅な高速化が気に入って、私はもっと...
zy = ((zy * zx) >> 23) + cy;
zx = xx - yy + cx;
-とでもすればそれで済む話だだからです。
-さらに、xやyを工夫すれば、sxやsyもなくせると思いました。
-さらに加えて、zx=zy=0で回り始めるのなら、n=0のときのif文...
-さらに、cxが変化しなければ、最初のzx*zxの値は変化しない...
-さらに、 n=1~447 でループするのではなく、n= -446~0 で...
-ええとなぜこのような努力をするのかというと、「まともなア...
-これらすべての変更を盛り込んだ究極形が以下の「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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28656617
time=1.115[sec]
-hの値が全く変わっていないので、正しく計算できていること...
----
-参考までに、floatを使いつつ精一杯高速化するとこうなりま...
// "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_...
return 0; // dummy
}
-これを実行すると、手元の環境では以下のようになりました。
h=28655896
time=1.761[sec]
-hの値がまた少し変わっていますが、これもおそらく精度に由...
** (19-2) gccの最適化について
-今回、19-1b → 19-1c → 19-1d と改良してきたわけですが、こ...
-だから私はES-BASICにおいて言語が最適化を頑張ることを放棄...
** つづく
-[[esbasic0013]]へつづく
* こめんと欄
#comment
ページ名: