第10回 アカデミックスキルII C言語(3) 乱数の発生と gnuplot の初歩

計算機はいつも決まった処理を行うように設計されているが， ゲームを作ったり，確率的シミュレーションを行なったりする場合には， サイコロのような確率的な値(乱数 と呼ばれる)が必要になることもある．

C言語で乱数を扱う関数や定数は <stdlib.h> ライブラリに定義されている． まずは, ~/cpp ディレクトリ内に以下のサンプルプログラム sample3-1.cpp を作り，コンパイル・実行してみよう．

• Emacs でファイルを開くには C-x C-f (file-file)
  • ~/cpp/sample3-1.cpp というファイルを開く/新しく作るには:

    C-x C-f ~/cpp/sample3-1.cpp RET
    

• Emacs 上にペーストするには C-y (yank)
• Emacs でファイルを保存するには C-x C-s (save-buffer)
• Emacs 上でコンパイル/実行するには M-! (shell-command)
  • sample3-1.cpp をコンパイルして sample3-1.o という実行ファイルを作るには:

    g++ sample3-1.cpp -o sample3-1.o
    

  • sample3-1.o を実行するには

    ./sample3-1.o
    

sample3-1.o を実行すると，以下のような表示されるはずだ：

r1: 16807
r2: 282475249
r3: 1622650073
RAND_MAX: 2147483647

r1, r2, r3 は，いずれも rand() という関数で初期化されているのに， 全く違うランダムな値が格納されていることを確認して欲しい．

rand() は, ある初期値(「種」と呼ばれる)から ある規則(線形合同法)に従う数列を順に返す． この「種」を変えることで，異なる疑似乱数を生成できる． 以下のプログラム sample3-2-seed.cpp を作成し，コンパイル・実行してみよう．

• このプログラムのように，ユーザー値を入力するのを「待つ」ようなプログラムは Emacs の M-! (shell-command) では実行できない．「端末」を使うか， M-x eshell として Eshell(Emacs 上のシェル)を使おう．

このプログラムでは，生成される数列のパターンを判別し易いように, あえて5つの乱数を表示させている．

このプログラムを何度か実行して，以下を確認してみよう．

1. 「種」として 1 を与えてみよう. sample3-1.o の結果と比較してみよう．
2. 「種」として 1 を与えた場合に生成される2番目の乱数(r2)を「種」として与えてみよう．
3. 「種」として 1, 10, 11, 100, 111, 1000, 10001 を与えた時の乱数列を比較してみよう
4. 何度も同じプログラムをコンパイルしたり実行したりするには「端末」や Emacs が 備える「履歴」機能を使うとよい．具体的には，端末や Emacs でコマンドが入力できる ところで カーソルの上下 を押すと，過去に入力したコマンドが順に表示されてくる．

• 10行目

  unsigned int seed;            // 乱数の「種」を格納する変数
  

  では，乱数の「種」を格納するために符号無し整数(unsigned int)型の変数 seed を定義している． 符号無し整数とは，全てのビットを数値の表現に使う型で 0 以上の整数の格納に用いられる．

• 12行目

  srand(seed);               // 乱数の「種」として seed を与える
  

  では, srand() という関数を用いることで，乱数の「種」を与えている． 具体的には, 関数 srand() を seed を 引数 として与えて呼び出すことで， 乱数の「種」の値に seed の値を代入している．

sample3-2-seed.cpp では「種」を入力させることで， 毎回異なる疑似乱数を発生させることに成功した． しかし, 実行のたびにユーザーに「種」を入力させるのは面倒だし， 同じ「種」からは毎回同じ乱数が生成されてしまう(予測できてしまう)． そこで，この「種」を自動的に，かつ，予測不能な形で生成する古典的な方法として， プログラムが呼び出された際の「時間」を用いる方法を紹介する． 以下のプログラム sample3-2-time.cpp を作成し，コンパイル・実行してみよう．

このプログラムを何度か実行してみよう．

1. 間隔を1秒以上あけて実行してみよう．
2. 間隔をなるべく短くして実行してみよう(カーソルの上下で過去のコマンドが入力できることを活用しよう)
3. 1つ目の乱数 r1 と，2つ目以降の乱数 r2, r3, r4, r5 の変り方を比べてみよう

4行目:

#include <time.h>               // 時間に関するライブラリを読み込む

で，時間に関するライブラリ <time.h> を読み込んでいる．

10行目:

unsigned int seed = (unsigned int) time(NULL); // 乱数の「種」を現在時刻で初期化

では, time() 関数を用いて現在の時刻(世界標準時の 1970年1月1日 0:00 からの経過秒数)を取得し， それを用いて seed を初期化している． それだけの処理をするだけなのに, time の前に (unsigned int) がついていたり, 引数に見慣れない変数 NULL が与えられていたり，とごちゃごちゃしている．

本来, 関数 time() は, time_t 型(時間値を表す型で <time.h> で定義されている)の ポインタ を引数とし， time_t 型の値を返す関数である． まず, (unsigned int) を前につけることで, time の戻り値(time_t 型)を srand が引数として受け取れる符号無し整数(unsigned int 型)に変換している． 次に，ここでは引数に何も指定する必要が無いのだが， time を引数無しで呼び出すとエラーになるため, NULL という空のポインタを与えている．

なお，変数のポインタというのはその変数の「仮想的な番地」に相当するもので， これを用いることで，メモリを無駄に使うことなく関数やオブジェクト間で変数の値を共有できる． その詳細については本講義では解説しないが， C/C++ において極めて重要な役割を果たす仕様であることだけ覚えておいて欲しい．

sample3-2-time.o を何度か繰り返して実行してみると, r1 だけは大きく変わらないことが判る． これは, srand に与える乱数の「種」(≒時刻)が殆ど変わらないため， その直後に呼び出された rand() が返す値が似てしまうのが原因である(sample3-2-seed.o で乱数の種に 10 と 11, 100 と 101, 1000 と 1001 を入れた場合の最初の変数が似てしまうのと同じ)． この問題は srand() の直後に rand() を何度か「空打ち」することで解決できる． 以下のプログラム sample3-2.cpp を作成し，コンパイル・実行してみよう．

• sample3-2.o と sample3-2-seed.o を何度か実行し，結果を比較してみよう．

rand() を使うことで「0以上 RAND_MAX 以下のランダムな整数」が生成できるのだが， これをそのままゲームやシミュレーションに使うのは難しい． RAND_MAX の値に意味が無い上，計算機環境によってその値が異なるからだ．

そこで, 古典的な方法として, rand() で生成された乱数を RAND_MAX+1 で割ることで， 「0以上1未満のランダムな実数」(と思われるもの)を作る方法を紹介しよう． 次のサンプルプログラム sample3-3.cpp を作り，コンパイル・実行してみよう．

15〜17行目

double r1 = (double) rand() / ( (double)RAND_MAX + 1.0 );
double r2 = (double) rand() / ( (double)RAND_MAX + 1.0 );
double r3 = (double) rand() / ( (double)RAND_MAX + 1.0 );

の各行では, rand() を RAND_MAX + 1 で割ったもので実数(double 型)変数 r1, r2, r3 を初期化している． ここで, rand() の出力や RAND_MAX の前に (double) をつけることで， もともと整数であるこれらの値を実数(double)型に型変換している(1.0 は最初から実数なので型変換は不要)．

なお，これらの分母を (double) RAND_MAX に変えれば 「0以上1以下」 の実数乱数が生成できる．

~/cpp ディレクトリの下に「0以上1以下の」実数乱数 x, y を発生させる 次のサンプルプログラム sample3-4.cpp を作り，コンパイル・実行してみよう．

sample3-4.o の実行結果は大変地味で， 2つの実数乱数を空白で区切ったものが3セット表示されるだけだ． これを gnuplot と呼ばれるグラフツールでプロットしてみることにしよう．

「端末」を起動し，以下のように入力しよう

cd ~/cpp
gnuplot

すると,

G N U P L O T
Version 4.6 patchlevel 5    last modified February 2014
Build System: Darwin x86_64

Copyright (C) 1986-1993, 1998, 2004, 2007-2014
Thomas Williams, Colin Kelley and many others

gnuplot home:     http://www.gnuplot.info

のような文字が表示された後，プロンプト(入力待ち状態での表示)が

gnuplot>

と変わるはずだ．これは，端末上で gnuplot が起動されていることを表している．

以下を順に入力してみよう(gnuplot> はプロンプトなので Del キーなどで消そうとしなくてもよい):

set xrange [0:1]
set yrange [0:1]
set size square
plot "<./sample3-4.o"

新しいウィンドウが開いて，3つの点がプロットされているだろうか？

何度か gnuplot> プロンプトに対して

plot "<./sample3-4.o"

を繰り返し入力してみよう(ここでもカーソルキーの上下で過去に入力したコマンドを表示させられる)．

もし点のサイズが小さいと感じるようなら

plot "<./sample3-4.o" pointsize 2

としてみよう. pointsize は ps と省略してもよい．また，2より大きな値を入れてもよい．

gnuplot 用のウィンドウに点が表示されている状態で, gnuplot プロンプトに対して

replot - x + 1 linetype 3 title ""

と入力してみよう. linetype は lt と省略してもよい．

gnuplot で円をプロットするには，少し工夫が必要である．

set parametric
plot [ 0: 0.5*pi ] cos(t), sin(t)

とした後,

plot "<./sample3-4.o" pointsize 2
replot cos(t), sin(t) lt 3 t ""

とすると，四分円とランダムな点をプロットできる．

gnuplot は非常に高機能なグラフ描画アプリケーションである． Excel のようにボタン1つでグラフの形状を変えたりはできないが， 3次元等高線や2次元ベクトル図など科学技術計算結果を簡単に表示させられる． 下記のようなサイトが参考になるだろう．

gnuplot homepage

http:www.gnuplot.info

gnuplot の初歩

http:graph.pc-physics.com

gnuplot コマンド集

http:www.gnuplot-cmd.com