2011年11月29日火曜日

このブログのGitHubレポジトリを作りました


このブログのGitHubのレポジトリを作りました。

blog-codes@GitHub

これまでのブログエントリで公開していたコードを、GitHubのレポジトリに移しました。

粒子法のプログラム第1回(概要)
sph.cpp - GitHub

エントリはこちら。粒子法のプログラムを解説しています。論文や本の数式をプログラムに落とすところやパラメータの設定が分かるよう、ソースコードも合わせて公開しています。

Haskell、OCamlでSPH法
sph.hs - GitHub
sph.ml - GitHub
sph.cpp - GitHub

エントリはこちら。粒子法のコードを、Haskell、OCaml、C++で書いて比較しています。

【OCaml】スタンフォードバニーの読み込みライブラリ
PlyReader.ml - GitHub
PlyReader.mli - GitHub

エントリはこちら。OCamlで、スタンフォードバニーの読み込みライブラリを書きました。

粒子法シミュレーションを記述するための言語を設計中
Simulation.hs - GitHub

エントリはこちら。粒子法シミュレーションを記述するための言語を考えており、それを紹介しています。

CPUキャッシュによるプログラム高速化の例
code1.c - GitHub
code2.c - GitHub

エントリはこちら。CPUのキャッシュメモリを使うで、プログラムの実行性能をあげることができます。キャッシュメモリを活用することでどの程度性能が向上するのかを紹介しています。


--

2011年11月24日木曜日

Acceleration of programs using CPU caches effectively


Effective use of CPU caches make the performance of programs improved. This entry shows how much improvement is given with effective use of CPU caches.

Motivation
Computer data storages have hierarchy which consists of registers, cache, main memory and so on, distinguished by response time. Registers and cache memory are smaller but faster than main memory.

CPUs can directly operate registers only. So, between registers and main memory, there is about x100 gap in access speed, and to fill the gap, CPUs have cache memory inside. With effectively using CPU cache, we can improve the performance of programs.

However, programmers can not make CPUs directly to use the cache memory in the way they want, so they must write programs considering how the CPUs read/write the data on the memory to indicate CPUs indirectly how to use the cache.

* SPARC64 VIIIfx used in the K computer of RIKEN in Japan has "sector cache" functionality and programmers can directly indicate it how to use the cache. With sector cache, the cache in a CPU is divided in two sectors and programmers can instruct which sector they use to store data.

This entry demonstrates an example of performance improvement using CPU cache effectively.

Example
I show the optimization in efficiency of a program which computes the multiplication of 1000x1000 matrices using cache blocking technique.

The code before optimized (code1.c):
#include <stdio.h>
#include <stdlib.h>
 
#define N 1000
 
int main ()
{
  double A[N][N], B[N][N], C[N][N];
  int i, j, k;
  
  // initialize the matrices
  for ( i = 0; i < N; i++ )
    for ( j = 0; j < N; j++ )
      {
        A[i][j] = 1.0;
        B[i][j] = 1.0;
        C[i][j] = 0.0;
      }
  
  // matrix multiplication
  for ( i = 0; i < N; i++ )
    for ( j = 0; j < N; j++ )
      for (k = 0; k < N; k++ )
        C[i][j] += A[i][k] * B[k][j];
  
  printf ( "%f\n", C[0][0] );
  
  return 0;
}

I rewrite this code as below with cache blocking (code2.c). Of course, the result of multiplication is same as to that of code1.c.
$ cat code2.c
#include <stdio.h>
#include <stdlib.h>
 
#define N 1000
 
int main ()
{
  double A[N][N], B[N][N], C[N][N];
  int i, j, k;
  
  // initialize the matrices
  for ( i = 0; i < N; i++ )
    for ( j = 0; j < N; j++ )
      {
        A[i][j] = 1.0;
        B[i][j] = 1.0;
        C[i][j] = 0.0;
      }
  
  // matrix multiplication
  int ibl = 100;
  int ib, jb, kb;
  for (ib=0; ib<N; ib+=ibl)
    for (jb=0; jb<N; jb+=ibl)
      for (kb=0; kb<N; kb+=ibl)
        for (i=ib; i<ib+ibl; i++)
          for (j=jb; j<jb+ibl; j++)
            for (k=kb; k<kb+ibl; k++)
              C[i][j] += A[i][k] * B[k][j];
  
  printf ( "%f\n", C[0][0] );
  
  return 0;
}

Demonstration
The times taken to execute the unoptimized code and the optimized code are as below respectively:
$ gcc -Wl,-stack_size,1000000000 -O2 -o code1 code1.c
$ time ./code1 > /dev/null

real    0m7.595s
user    0m7.494s
sys     0m0.035s

$ gcc -Wl,-stack_size,1000000000 -O2 -o code2 code2.c
$ time ./code2 > /dev/null

real    0m2.509s
user    0m2.462s
sys     0m0.030s
* The -Wl option is to designate ld to allocate large arrays on the stack

The result represented as a chart:



Cache blocking improve the performance by about x3.

Comparing the performance as the size of cache blocking
Next, I show how the performance changes as the size of block. I changed the size of the block to 100, 200, 250, 500, then timed to calculation.

As a chart:


With the block size of 500, the performance degrades massively, which represents that the expansion of the block size leads to overflow of data from the cache.

The bytes of the 500x500 sized block is 2MB, with double typed array, so, considering other programs running in background, it agrees with the CPU cache size of my computer :)

--

2011年11月23日水曜日

CPUキャッシュによるプログラム高速化の例


CPUのキャッシュメモリを使うで、プログラムの実行性能をあげることができます。このエントリでは、キャッシュメモリを活用することでどの程度性能が向上するのか見てみます。

動機
コンピュータの記憶領域は、レジスタ、キャッシュメモリ、メインメモリなどからなる階層構造になっています。レジスタやキャッシュメモリはメインメモリより容量が小さいものの、非常に高速です。

CPUが直接触れるのはレジスタのみです。レジスタとメインメモリでは、データアクセスの速度に約100倍の差があり、その差を埋める目的で、CPU内部にキャッシュメモリがあります。キャッシュメモリを活用できれば、プログラムの実装性能を大きく高めることができます。

ただし、基本的に、CPUがどのようにキャッシュメモリを使うのかをプログラマが直接指示することはできません。なので、キャッシュメモリを活用するには、CPUがどのようにメモリ上のデータを使うのかを考慮したプログラムを書くことで間接的にキャッシュメモリを使うことになります。

※理化学研究所の京で使われているSPARC64 VIIIfxにはセクタキャッシュという機構があり、プログラマがキャッシュを制御できます。セクタキャッシュを使うと、キャッシュメモリを2つの領域に分割し、データごとにキャッシュに使う領域をプログラマが指定できます。

このエントリでは、キャッシュメモリを活用することで、実際にどの程度性能が向上するのか、一例を見てみます。


サンプルコード
1000×1000行列の積を求めるプログラムについて、キャッシュブロック化によって実行性能を最適化します。

以下が、最適化前のCのコードです。(code1.c)
#include <stdio.h>
#include <stdlib.h>
 
#define N 1000
 
int main ()
{
  double A[N][N], B[N][N], C[N][N];
  int i, j, k;
  
  // initialize the matrices
  for ( i = 0; i < N; i++ )
    for ( j = 0; j < N; j++ )
      {
        A[i][j] = 1.0;
        B[i][j] = 1.0;
        C[i][j] = 0.0;
      }
  
  // matrix multiplication
  for ( i = 0; i < N; i++ )
    for ( j = 0; j < N; j++ )
      for (k = 0; k < N; k++ )
        C[i][j] += A[i][k] * B[k][j];
  
  printf ( "%f\n", C[0][0] );
  
  return 0;
}

これを、キャッシュブロック化によって以下のコードに書き換えます(code2.c)。もちろん、計算の結果そのものは変わりません。
$ cat code2.c
#include <stdio.h>
#include <stdlib.h>
 
#define N 1000
 
int main ()
{
  double A[N][N], B[N][N], C[N][N];
  int i, j, k;
  
  // initialize the matrices
  for ( i = 0; i < N; i++ )
    for ( j = 0; j < N; j++ )
      {
        A[i][j] = 1.0;
        B[i][j] = 1.0;
        C[i][j] = 0.0;
      }
  
  // matrix multiplication
  int ibl = 100;
  int ib, jb, kb;
  for (ib=0; ib<N; ib+=ibl)
    for (jb=0; jb<N; jb+=ibl)
      for (kb=0; kb<N; kb+=ibl)
        for (i=ib; i<ib+ibl; i++)
          for (j=jb; j<jb+ibl; j++)
            for (k=kb; k<kb+ibl; k++)
              C[i][j] += A[i][k] * B[k][j];
  
  printf ( "%f\n", C[0][0] );
  
  return 0;
}


実行
最適化前、最適化後それぞれについて、実行にかかる時間は以下のようになりました。
$ gcc -Wl,-stack_size,1000000000 -O2 -o code1 code1.c
$ time ./code1 > /dev/null

real    0m7.595s
user    0m7.494s
sys     0m0.035s

$ gcc -Wl,-stack_size,1000000000 -O2 -o code2 code2.c
$ time ./code2 > /dev/null

real    0m2.509s
user    0m2.462s
sys     0m0.030s
※大きな配列をスタック上に確保するために、ldにオプションを渡してスタックサイズを大きくしています

グラフにすると以下のようになります。



キャッシュブロック化によって、約3倍の性能向上が得られました。


ブロックサイズによる性能比較
次に、ブロック化のサイズによって性能がどのように変わるかを調べてみます。ブロック化のサイズを100、200、250、500に変えて、実行時間を計ると、以下のようになりました。
100 : 2.311s
200 : 2.390s
250 : 2.509s
500 : 6.588s
※ブロック化のサイズの単位は、配列の要素数です。バイト数に直すと、double型の配列なので8を掛けたものになります

グラフにすると、以下のようになります。


ブロック化のサイズが500になると、急激に性能が悪化しています。このことから、ブロック化のサイズが大きくなりキャッシュにデータが収まらなくなった様子が分かります。

500×500のブロックのデータ量が2MBなので、このプログラム以外が動いていることも勘案すると、使っている環境のCPUの持つキャッシュメモリサイズと一致した結果といえます。


--

2011年11月22日火曜日

Common Lispでの、依存関係を考慮したファイルロード


Common Lispにおいて、asdfを使って、依存関係を考慮してプロジェクトの構成ファイルをロードする方法をまとめました。

動機
よっぽど小さなプログラムでない限り、プログラムの規模が大きくなってくると、ファイルを複数に分けることになります。

困ったことに、Common Lispには、依存関係をみて必要なファイルをロードする機能がありません。なので、ファイルを複数に分けると、必要なファイルを1つ1つロードするはめになります。

asdfを使うことで、ファイルの依存関係をみて必要なファイルを1度にロードすることができます。

asdfとは?
asdfは、"Another System Definition Facility"の略で、Common Lispプログラムを1つの「モジュール」として管理し、読み込むことができるツールです。現在はASDF2が広く使われています。

makeや、JavaのAntに相当するもので、PerlのCPANやRubyのRubyGemsの1つ下のレイヤに位置づけられます。CPANやRubyGemsにあたるものとしては、asdf-installがあります。

以下のような「sample」プロジェクトがあるとします。
sample
 - caller.lisp  ;; callee.lispで定義されている関数を使う
 - callee.lisp  ;; caller.lispから呼ばれる
 - sample.asd   ;; asdfモジュール定義ファイル

各ファイルの内容は以下のようになっています。
$ cat caller.lisp
(in-package :common-lisp)

(defpackage sample.caller
  (:use :common-lisp
        :sample.callee))
(in-package :sample.caller)

(hello-world)

$ cat callee.lisp 
(in-package :common-lisp)

(defpackage sample.callee
  (:use :common-lisp)
  (:export :hello-world))
(in-package :sample.callee)

(defun hello-world ()
  (print "Hello, World!"))

$ cat sample.asd 
(defsystem "sample"
  :components ((:file "callee")
               (:file "caller"))
  :serial t)

sample.asdがasdfから参照されるように設定し、common lispを起動します。

以下の関数を評価することで、sampleプロジェクトをロードすることができます。sample.asdの:serialオプションによって、ファイルの依存関係を考慮した形でファイルがロードされます。
(asdf:load-system 'sample)


--

2011年11月7日月曜日

TPPと混合診療の関係を簡単にまとめ


TPPと混合診療の関係を簡単に調べてみたので、そのまとめ。

重要なポイントとして、TPPによるアメリカの狙いの最も大きなものの1つが、サービス(金融)の自由化によって日本の個人資産を喰うこと。これは間違いない。わざわざサービス(金融)を後から捩じ込むだけのボリュームがあるターゲットは、日本だけだから。

アメリカ観点(混合診療を解禁しろ派)
混合診療解禁
→皆保険枠と混合診療枠ができる
→混合診療枠が医療保険のターゲット
→医療保険のノウハウを持っているのはアメリカの保険会社
→アメリカの保険会社はウマウマ

ref. http://haetarou.web.fc2.com/TPP/TPP_Kaisetu.html

※サービス(金融)はもともとTPPには入っておらず、アメリカが後から追加したもの。間違いなく、日本の個人資産を喰おうと狙っている

医師会観点(混合診療を解禁するな派)
混合診療解禁
→「政府は、財政難を理由に、保険の給付範囲を見直そうとしている」
→「混合診療を認めることによって、現在健康保険で見ている療養までも『保険外』とする」
→高額な診療が保険外になる
→高額な診療を受ける患者が減る
→医者の取り分が減る。困る

※医者は、国の医療費からたくさん診療費をもらいたい

ref. http://www.med.or.jp/nichikara/kongouqa/

厚生労働省観点(混合診療を解禁しろ派)
建前上、財政健全化したい
→医療費を減らしたい
→混合診療を解禁することによって、保険対象を減らすことで、医療費を減らす

※現状では、混合診療を認めていないため、その分手厚く保険対象にしてある

--