|
液体や気体中に浮遊する微粒子は、不規則な運動をくりかえしています。この現象は、1827年に植物学者ブラウンによって、水に浮く花粉の不規則な運動
として初めて発見されました。この現象は、酔っぱらいが、ちどり足で歩くのとにていることから、ランダムウォークあるいは酔歩の問題と呼ばれます。
1次元ランダムウォーク
まずはじめに、簡単な1次元のランダムウォークを考えてみます。図6.1.1のように1本の道の上で、酔っぱらいが、ふらふらといったりきたりしているとします。
図6.1.1 1次元ランダムウォーク |
|
酔っぱらいが、ある時刻で x の場所にいたとして、次の時刻に前へ行くか、後ろへ行くかは、まったくランダムで、予測が着きません。このような場合の、 x について、シミュレーションしてみましょう。
プログラムのポイント
プログラムS0610.javaは、1次元ランダムウォークのシミュレーション・プログラムで、各クロックごとに、(1)において確率1/2で、前へ進むか、後ろに進むかをきめています。その結果を図6.1.2に示します。横軸時刻
t 、縦軸 x のグラフとして表現しています。このランダムウォークは、3.1節でのべたベルヌーイ試行を時刻 t ごとにおこなっているものです。
またこれは、勝ったり負けたりのゲームで、持ち金の金額の変化をあらわしています。
図6.1.2 1次元ランダムウォークのシミュレーション結果 |
|
S0610.java 1次元ランダムウォーク | Download | 実行 | |
/* S0610.java
* 1次元ランダムウォーク
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0610 extends Applet implements ActionListener {
Button button0;
public void init() {
button0 = new Button(" 再実行 ");
add(button0);
button0.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
repaint();
}
public void paint(Graphics g) {
Window w ;
w = new Window();
int SPACE = 30;
int HIGHT = 400;
int WIDTH = 640;
long t = 0;
long x = 0;
long T_END = 1000; /* 終りの時刻 */
double P = 0.5; /* 確率1/2 */
/*グラフィックの準備*/
w.setWindow(0, 0.0,-100.0,T_END,100.0,
SPACE,HIGHT-SPACE,WIDTH-SPACE,SPACE);
w.axis(0, "t", T_END/10, "x", 10, g);
w.moveTo(0, 0.0, 0.0, g);
/*メイン*/
for (t = 0; t < T_END; t ++) {
if (P > Math.random()) {
//(1) 確率1/2で前へ進むか、後ろに進むか
x = x + 1;
} else {
x = x - 1;
}
g.setColor(Color.green);
w.lineTo(0, (double)(t), (double)(x), g);
}
stop();
}
}
|
2次元ランダムウォーク
平面上をブラウン運動する粒子を考えます。周囲の流体の分子から乱雑な衝撃を受け、ブラウン粒子はランダムな方向に運動をしますが、その方向は、 0
から 2π
の間に一様に分布すると考えられます。衝撃を受けたブラウン粒子は、そのたびごとにさまざまな距離を進むと考えられますが、ここでは、図6.2.1のよう
に一定の距離を進むものとします。このような、2次元のランダムウォークをシミュレーションしてみましょう。
図6.2.1 2次元ランダムウォーク |
|
プログラムのポイント
プログラムS0620.javaでは、(1)において 0 から 2π の乱数を発生させ、(2)のように x, y から、次の時点の x, y へ一歩進ませています。(3)でブラウン粒子の軌跡をディスプレイします。
シミュレーション結果を図6.2.2に示します。
図6.2.2 2次元ランダムウォークのシミュレーション結果 |
|
S0620.java 2次元ランダムウォーク | Download | 実行 | |
/* S0620.java
* 2次元ランダムウォーク
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0620 extends Applet implements ActionListener {
Button button0;
public void init() {
button0 = new Button(" 再実行 ");
add(button0);
button0.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
repaint();
}
public void paint(Graphics g) {
Window w ;
w = new Window();
int SPACE = 30;
int HIGHT = 400;
int WIDTH = 640;
double X_MAX = 80.0;
double Y_MAX = 50.0;
double x,y; /* 粒子の座標 */
double th; /* θ */
int t; /*時刻*/
long T_END = 10000; /* 終りの時刻 */
/*グラフィックの準備*/
w.setWindow(0, -X_MAX, -Y_MAX, X_MAX, Y_MAX,
SPACE,HIGHT-SPACE,WIDTH-SPACE,SPACE);
w.axis(0, "", 10,"", 10, g);
w.moveTo(0, 0, 0, g);
/*メイン*/
x = 0.0;
y = 0.0;
for (t = 1; t <= T_END; t ++) {
th = 2.0 * Math.PI * Math.random(); //(1)
x = x + Math.cos(th); //(2)
y = y + Math.sin(th); //(2)
g.setColor(Color.green);
w.lineTo(0, x, y, g); //(3)
}
stop();
}
}
|
拡散過程 (文献(10)参照)
たくさんの粒子を一ケ所に集めて、ブラウン運動をさせると、粒子はどんどん広がって、ついには、一様にひろがってしまうでしょう。このような過程を拡散過程といいます。本節では、この拡散過程をシミュレーションします。
ところで、単位時間に1度ずつ、衝撃を受けるブラウン運動で、 t=0 で原点にいた粒子が、 t=n にはどこにいる確率が高いでしょうか。図6.3.1のように、一回の衝撃で動く距離を1とすると、ブラウン粒子が、 n 回目の衝撃をうけると、原点からの距離 dn は、つぎのように求められます。
|
・・・・式6.3.1 |
となり、 n を大きくしていくと、平方根の中の2項目は n にくらべて無視できるようになります。したがって
|
・・・・・ |
式6.3.2 |
となり、原点からの平均距離は衝撃の回収がふえると、その回数の平方根に比例して遠くなることがわかります。
たくさんの粒子を一ケ所に集めて、ブラウン運動をさせると、粒子はどんどん広がって、ついには一様に分布してしまいます。このことを、実際にシミュレーションで確かめてみましょう。
図6.3.1 粒子の位置 |
|
プログラムのポイント
プログラムS0630.javaにおいては、2つのウインドウを用います。ウインドウ0では、粒子が運動する平面、ウインドウ1は、原点からの粒子の平均距離を表示することとします。
時刻 0 で、1000個の粒子は原点にあるとし、クロックを1進めるたびに、前節と同じように、 0〜2π の乱数を発生させ、ランダムウォークをさせます。1000個の粒子に対しこのことを実行し、原点からの粒子の平均距離 d を計算したあと、時刻を進めていきます。
実行結果
このプログラムを実行してみましょう。図6.3.2のように、はじめ原点にあった粒子は、あたかも生き物のように動きまわり、全体はしだいに拡散してい
く様子が表示されます。また、たしかに原点からの平均距離は、ウインドウ1に示すように、時刻の平方根に比例していることも示されます。
図6.3.2 拡散のシミュレーション結果(途中の図 長時間CPU占有注意) |
|
S0630.java 拡散のシミュレーション | Download | 実行 (長時間CPU占有注意)| |
/* S0630.java
* 拡散のシミュレーション
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0630 extends Applet implements ActionListener {
Button button0;
Button button1;
int sw =0;
public void init() {
button0 = new Button(" 実行 ");
button1 = new Button("クリア");
add(button0);
add(button1);
button0.addActionListener(this);
button1.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
if (label.equals(" 実行 ")) sw = 1; else sw = 0;
repaint();
}
public void paint(Graphics g) {
Window w ;
w = new Window();
int SPACE = 30;
int HIGHT = 400;
int WIDTH = 640;
int PARTICLE = 10000; /* 粒子の数 */
double MAX = 50.0;
long T_END = 1000;
double x[] = new double[PARTICLE]; /* すべての粒子は原点にある */
double y[] = new double[PARTICLE];
int i; /* 粒子のカウンタ */
long t = 0; /* 時刻 */
double th; /* θ */
double d1; /* 原点からの距離の和 */
double d; /* 原点からの距離の平均値 */
/*グラフィックの準備*/
w.setWindow(0, -MAX,-MAX,MAX,MAX,
SPACE,HIGHT-SPACE,HIGHT-SPACE,SPACE);
w.axis(0, "", MAX, "", MAX, g);
w.setWindow(1, 0,0,T_END,MAX,
HIGHT,HIGHT/2,WIDTH-SPACE,SPACE);
w.axis(1, "時間", T_END / 5, "原点からの平均距離", MAX / 5, g);
w.moveTo(1, 0.0, 0.0, g);
/*メイン*/
if ( sw == 1) {
while (t < T_END) {
d1 = 0.0;
for (i = 0; i < PARTICLE; i ++) {
th = 2.0 * Math.PI * Math.random();
if (x[i] <= MAX && x[i] >= -MAX &&
y[i] <= MAX && y[i] >= -MAX) {
g.setColor(Color.white);
/* 前の粒子の位置を消す */
w.putPixel(0, x[i], y[i], g);}
x[i] = x[i] + Math.cos(th);
y[i] = y[i] + Math.sin(th);
if (x[i] <= MAX && x[i] >= -MAX &&
y[i] <= MAX && y[i] >= -MAX) {
g.setColor(Color.blue);
/* 次の粒子の位置を表示 */
w.putPixel(0, x[i], y[i], g); }
d1 = d1 + Math.sqrt(x[i] * x[i] + y[i] * y[i]);
}
/* 消された軸を再度書く */
w.axis(0, "",MAX,"",MAX, g);
d = (double)(d1 / PARTICLE);
t = t + 1;
w.lineTo(1, t, d, g);
}
stop();
}
}
}
|
周辺部では遅く動く
前節のプログラムを長時間まわしていると、実は変なことが起こります。中心部にあった粒子が徐々に拡散していき、ついに MAX という枠をはみだしてしまいます。そこで、今回は枠まで粒子が到着したら、そこに壁があり、その壁ではねかえるようにします。
また、前節では、粒子が1クロックで動く長さを一定としていましたが、中心部では速く、周辺部では遅く動くようにしてみましょう。どんなことがおこるでしょうか 。
プログラムのポイント
プログラムS0640.javaにおいては、S0630.javaと同様に1000個の粒子をブラウン運動させます。ただし、初期状態として一様に分布させることとし、また粒子の動き方は、1クロックで進む長さを l とすれば、
|
・・・・・ |
式6.4.1 |
とします( MAX は壁までの距離、 xi, yi は粒子 i の位置)。すなわち、中心部で最も速く、
|
・・・・・ |
式6.4.2 |
でうごき、中心部から最も遠い隅の
|
・・・・・ |
式6.4.3 |
の時、停止状態となるものとします。
また、次のクロックにおける位置を計算した結果、その値が MAX をこえた場合、図6.4.1のように、壁のところで、完全反射するものとします。すなわち
if (x[i] >= MAX) { x[i] = MAX - (x[i] - MAX) ; }
のように MAX からはみだした分、ひき返すように計算します。 y 座標も同様です。
図6.4.1 粒子の動き方 |
|
実行結果
このプログラムを実行してみましょう。部屋のほこりは、中心部では風が吹いたり、人の動きがあるため、速く動きまわり、次第に風の弱い部屋の隅に集まってきます。このような様子が図6.4.2のように、シミュレーションによって確かめることができます。
図6.4.2 ほこりはすみにたまるシミュレーション結果(長時間CPU占有注意) |
|
S0640.java ホコリのシミュレーション | Download | 実行(長時間CPU占有注意) | |
/* S0640.java
* ホコリのシミュレーション
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0640 extends Applet implements ActionListener {
Button button0;
Button button1;
int sw =0;
public void init() {
button0 = new Button(" 実行 ");
button1 = new Button("クリア");
add(button0);
add(button1);
button0.addActionListener(this);
button1.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
if (label.equals(" 実行 ")) sw = 1; else sw = 0;
repaint();
}
public void paint(Graphics g) {
Window w ;
w = new Window();
int SPACE = 30;
int HIGHT = 400;
int WIDTH = 640;
int PARTICLE = 1000; /* 粒子の数 */
double MAX = 50.0;
long T_END = 10000;
double x[] = new double[PARTICLE];
double y[] = new double[PARTICLE];
int i; /* 粒子のカウンタ */
long t = 0;
double th;
double l;
/*グラフィックの準備*/
w.setWindow(0, -MAX,-MAX,MAX,MAX,
SPACE,HIGHT-SPACE,HIGHT-SPACE,SPACE);
w.axis(0, "", MAX, "", MAX, g);
/*開始*/
if ( sw == 1) {
for (i = 0; i < PARTICLE; i ++) {
/* 始めは一様に分布 */
x[i] = MAX * (2.0 * Math.random() - 1.0);
y[i] = MAX * (2.0 * Math.random() - 1.0);
g.setColor(Color.green);
w.putPixel(0, x[i], y[i], g);
}
/*メイン*/
while (t < T_END) {
for (i = 0; i < PARTICLE; i ++) {
g.setColor(Color.white);
w.putPixel(0, x[i], y[i], g);
th = 2.0 * Math.PI * Math.random();
l = (Math.sqrt(2.0) * MAX -
Math.sqrt(x[i] * x[i] + y[i] * y[i])) / MAX;
x[i] = x[i] + l * Math.cos(th);
y[i] = y[i] + l * Math.sin(th);
/* 反射 右の壁 */
if (x[i] >= MAX) { x[i] = 2*MAX - x[i];}
/* 反射 左の壁 */
if (x[i] <= -MAX) { x[i] = -2*MAX - x[i];}
/* 反射 上の壁 */
if (y[i] >= MAX) { y[i] = 2*MAX - y[i];}
/* 反射 下の壁 */
if (y[i] <= -MAX) { y[i] = -2*MAX - y[i];}
g.setColor(Color.blue);
w.putPixel(0, x[i], y[i], g);
}
w.axis(0, "", MAX, "", MAX ,g);
t = t + 1;
}
stop();
}
}
}
|
1. プログラムS0620.javaでは、歩幅は一定であったが、歩く度に歩幅がランダムに変わる場合をシミュレーションしなさい。(時間間隔一定、歩幅不定のランダムウォーク)
(解答) A0610.java | 実行 |
2. p が小さい場合のベルヌーイ試行をおこない、事象が起ったときに再度確率1/2で、右に行くか左に行くかきめるとします。このときの、軌跡を図示するプログラムを作りなさい。(歩幅一定、あるく時間が不定のランダムウォーク)
(解答) A0620.java | 実行 |
|
|