Chapel とは？
Chapel とは、米国 Cray 社の開発している HPC (High Performance Computing) 用の言語で、PGAS (Partitioned Global Address Space) 概念に基いています。
PGAS言語
PGAS は日本語では区分化大域アドレス空間と訳されています。素朴には分散メモリーに統一したアドレスを与えてアクセスしやすくしようとうアイデアです。ここで局所的なメモリーと遠方のメモリーを、プログラム側で意識的に分けて扱えないと演算性能が出ないことが分かっているので、そういう仕組みを含んだ概念だと思われます。
PGAS 言語には、古くには CAF (Co-Array Fortran) や UPC (Unified Parallel C) などがあります。やや新しくは、IBM の X10、Cray の Chapel、Sun の Fortress、Berkeley の Titanium などがあります。
言語的には、CAF と UPC はそれぞれ Fortran と C99 の拡張、X10 と Titunium は Java の拡張、Chapel と Fortress は最近流行りの軽量言語の特徴を入れた独自言語となっています。
このうち CAF は Fortran2008 の機能として ISO 規格に正式に取り入れられました。X10 と Chapel は細々と開発が続いています。Fortress は資金が打ち切られたため開発中止となりました。Titanium もあまり話を聞きません。UPC は生きているようです。
Chapel
Chapel のサイトは http://chapel.cray.com/ にあります。ここで言語の概略などを知ることができます。入門用の教材を豊富に用意するなど新規参入者への敷居を低くする努力をしています。
しかし笛吹けど踊らず、新奇なもの好きのネット界の言語マニアからも無視される状態が続いています。
興味
私の Chapel に対する興味は、そのデータ構造の取り方にあります。Chapel では Domain と称するデータ構造を記述する抽象を導入して柔軟なデータ構造を扱おうとしています。この仕組みがどのようなもので、どの程度有効なのかを知りたいと思います。
Fortran は配列処理に長けていますが、その特徴を生かせるデータ構造が矩形に限られます。丁度かつて "What's wrong with APL?" で APL が批判されたような状況になっています。また分散型のプログラムでは、縁（Halo）の扱いなどが常に問題になりますが、その扱いも煩雑です。
また Fortran では data parallel 用の機能は多く用意されていますが、task parallel 用の機能は基本的に存在せず、OpenMP の拡張があるのみです。これに対して Chapel は data parallel/ task parallel とも多く機能を持っています。
Chapel の知見がコンセプト的に Fortran に応用できるのか？何か参考にならないか？というあたりに興味があります。
MandelBrot 図形
以下では、いつもの如く Fortran program の移植で BMP ファイルを出力して図を書くことにします。
chapel のサイトではベンチマークがありますが、それによればベース言語としての chapel は Fortran や C の４倍くらい遅いようです。ここのプログラムは普通の時間で結果が返ってきます。
ここでは単純な矩形配列しかとりませんでしたが、Domain による配列宣言は確かに楽な気がしました。Mandelbrot 図形の計算では各座標での計算が独立になっています。このため Chapel の機能を試すには単純すぎるので、今後も色々試してみたい気がしました。
実行結果
プログラム改 (H26.10.1)
構造体の定義、バイナリ・ファイルの I/O 等が必要になります。OOP仕様にしてみました。
module m_mandel {
record t_bmp_file_header {
var bfType : int(16) = 19778; //"BM"
var bfSize : int(32);
var bfReserved1 : int(16) = 0;
var bfReserved2 : int(16) = 0;
var bfOffBits : int(32);
}
record t_bmp_info_header {
var biSize : int(32) = 40; // Z'28'
var biWidth : int(32);
var biHeight : int(32);
var biPlanes : int(16) = 1; // always 1
var biBitCount : int(16);
var biCompression : int(32) = 0; // 0:nocompression
var biSizeImage : int(32);
var biXPelsPerMeter : int(32) = 3780; // 96 dpi
var biYPelsPerMeter : int(32) = 3780; // 96 dpi
var biClrUsed : int(32) = 0;
var biClrImportant : int(32) = 0;
}
record t_rgb {
var b, g, r : uint(8);
}
class t_bmp {
var D: domain(2);
var file_header: t_bmp_file_header;
var info_header: t_bmp_info_header;
var rgb: [D] t_rgb;
proc wr(fn:string) {
var f = open(fn, iomode.cw);
var w = f.writer(kind=ionative);
var nx: int(32) = this.D.high(1):int(32);
var ny: int(32) = this.D.high(2):int(32);
this.file_header.bfSize = 14 + 40 + 0 + nx * ny * 3;
this.file_header.bfOffBits = 14 + 40;
this.info_header.biWidth = nx;
this.info_header.biHeight = ny;
this.info_header.biBitCount = 24;
this.info_header.biSizeImage = nx * ny * 3;
w.write(this.file_header.bfType );
w.write(this.file_header.bfSize );
w.write(this.file_header.bfReserved1);
w.write(this.file_header.bfReserved2);
w.write(this.file_header.bfOffBits );
w.write(this.info_header.biSize );
w.write(this.info_header.biWidth );
w.write(this.info_header.biHeight );
w.write(this.info_header.biPlanes );
w.write(this.info_header.biBitCount );
w.write(this.info_header.biCompression );
w.write(this.info_header.biSizeImage );
w.write(this.info_header.biXPelsPerMeter);
w.write(this.info_header.biYPelsPerMeter);
w.write(this.info_header.biClrUsed );
w.write(this.info_header.biClrImportant );
for j in D.dim(2) do
for i in D.dim(1) do
w.write(this.rgb(i, j).b, this.rgb(i, j).g, this.rgb(i, j).r);
w.close();
}
}
proc mandel(c:complex) {
const maxiter : int = 255;
var z: complex = c;
var i, k: int;
for i in 1..maxiter {
k = i;
if abs(z) > 2.0 then break;
z = z * z + c;
}
return k;
}
proc main() {
const nx: int = 1440, ny: int = 1080;
const x0: real = -2.0, y0: real = -1.5,
x1: real = 2.0, y1: real = 1.5;
var x: [1..nx] real, y: [1..ny] real;
forall ix in 1..nx do x(ix) = (x1 - x0) / (nx - 1) * (ix - 1) + x0;
forall iy in 1..ny do y(iy) = (y1 - y0) / (ny - 1) * (iy - 1) + y0;
const D = {1..nx, 1..ny};
var c: [D] complex;
forall (ix, iy) in D do c(ix, iy) = x(ix)+ y(iy) * 1.0i;
var imandel: [D] int;
imandel = mandel(c);
var bmp = new t_bmp(D);
bmp.rgb.r = (255 - imandel):uint(8);
bmp.rgb.g = (255 - imandel):uint(8);
bmp.rgb.b = (255 - imandel):uint(8);
bmp.wr("test.bmp");
}
}
プログラム
構造体の定義、バイナリ・ファイルの I/O 等が必要になります。不勉強でよく分かっていないので、OOP 式に構成するのは、やめにしました。
module m_mandel {
record t_bmp_file_header {
var bfType : int(16) = 19778; //"BM"
var bfSize : int(32);
var bfReserved1 : int(16) = 0;
var bfReserved2 : int(16) = 0;
var bfOffBits : int(32);
}
record t_bmp_info_header {
var biSize : int(32) = 40; // Z'28'
var biWidth : int(32);
var biHeight : int(32);
var biPlanes : int(16) = 1; // always 1
var biBitCount : int(16);
var biCompression : int(32) = 0; // 0:nocompression
var biSizeImage : int(32);
var biXPelsPerMeter : int(32) = 3780; // 96 dpi
var biYPelsPerMeter : int(32) = 3780; // 96 dpi
var biClrUsed : int(32) = 0;
var biClrImportant : int(32) = 0;
}
record t_rgb {
var b, g, r : uint(8);
}
proc wr_bmp(D, bmp:[D] t_rgb, fn:string, nx:int(32), ny:int(32)) {
//var D = {1..nx, 1..ny};
var f = open(fn, iomode.cw);
var w = f.writer(kind=ionative);
var bmp_file_header : t_bmp_file_header;
var bmp_info_header : t_bmp_info_header;
bmp_file_header.bfSize = 14 + 40 + 0 + nx * ny * 3;
bmp_file_header.bfOffBits = 14 + 40;
bmp_info_header.biWidth = nx;
bmp_info_header.biHeight = ny;
bmp_info_header.biBitCount = 24;
bmp_info_header.biSizeImage = nx * ny * 3;
w.write(bmp_file_header.bfType );
w.write(bmp_file_header.bfSize );
w.write(bmp_file_header.bfReserved1);
w.write(bmp_file_header.bfReserved2);
w.write(bmp_file_header.bfOffBits );
w.write(bmp_info_header.biSize );
w.write(bmp_info_header.biWidth );
w.write(bmp_info_header.biHeight );
w.write(bmp_info_header.biPlanes );
w.write(bmp_info_header.biBitCount );
w.write(bmp_info_header.biCompression );
w.write(bmp_info_header.biSizeImage );
w.write(bmp_info_header.biXPelsPerMeter);
w.write(bmp_info_header.biYPelsPerMeter);
w.write(bmp_info_header.biClrUsed );
w.write(bmp_info_header.biClrImportant );
for j in 1..ny do
for i in 1..nx do
w.write(bmp(i, j).b, bmp(i, j).g, bmp(i, j).r);
w.close();
}
proc mandel(c:complex) {
const maxiter : int = 255;
var z: complex = c;
var i, k: int;
for i in 1..maxiter {
k = i;
if abs(z) > 2.0 then break;
z = z * z + c;
}
return k;
}
proc main() {
const x0: real = -2.0; const x1: real = 2.0;
const y0: real = -1.5; const y1: real = 1.5;
const nx: int(32) = 1440;
const ny: int(32) = 1080;
var x: [1..nx] real;
var y: [1..ny] real;
const D = {1..nx, 1..ny};
var c: [D] complex;
var imandel: [D] int;
var ix, iy: int;
forall ix in 1..nx {
x(ix) = (x1 - x0) / (nx - 1) * (ix - 1) + x0;
};
forall iy in 1..ny {
y(iy) = (y1 - y0) / (ny - 1) * (iy - 1) + y0;
}
forall (ix, iy) in D {
c(ix, iy) = x(ix)+ y(iy) * (0 + 1i);
}
imandel = mandel(c);
var bmp : [D] t_rgb;
bmp.r = (255 - imandel):uint(8);
bmp.g = (255 - imandel):uint(8);
bmp.b = (255 - imandel):uint(8);
wr_bmp(D, bmp, "test.bmp", nx, ny);
}
}