《计算物理》第三章习题参考答案

第三章《蒙特卡罗方法的若干应用》习题参考答案

1. 解：一般地[1]，首先定义n-维欧几里得空间球坐标r，为此，引入1个径向坐标，n-1个角坐标。以xi表示笛卡尔坐标，则 x1 rcos( 1),x2 rsin( 1)cos( 2),x3 rsin( 1)sin( 2)cos( 2),

xn 1 rsin( 1) sin( n 2)cos( n 1),xn rsin( 1) sin( n 2)sin( n 1).反变换

tan( 1) tan( n 2) tan( n 1)

1

n 2xn

.xn 1

其中， n 1 [0,2 ], [0, 1, 2, , n 2 ].

n-维欧几里得空间球体积元由下列Jacobian变换得出 dVn det

(xi)

drd 1d 2 d n 1

(r, j)

=rn 1sinn 2( 1)sinn 3( 2) sin( n 2)drd 1d 2 d n 1而n-维欧几里得空间球得体积由下列积分给出

R

2

Vn

0

1 0

dVn

n 2 0n 1 0

解析地，n-维欧几里得空间球得体积为

Vn CnRn, Cn

n

2

n ( 1)2

.

4 3 248 25 36

特别地，V3 R,V4 R,V5 R,V6 R.

32156

对于用MC方法计算多维球体积，可以参照多维独立随机变量分布情形抽样。下面以六维为例说明之。n 6, R 1, 六维球体积元

dV6 sin4( 1)sin3( 2)sin2( 3)sin( 4)drd 1d 2d 3d 4d 5 其中，r [0,1], 1,2,3,4 [0, ], 5 [0,2 ]. 改写上述体积元

dV6 f2( 2) f3( 3) f4( 4) f5( 5)drd 1d 2d 3d 4d 5

f(r) 1, f4( 1) sin4( 1), f3( 2) sin3( 2),

其中， 2

1.f2( f5 5) 3) sin( 3),f1( 4) sin( 4), (

引入6个均匀分布的伪随机数 , 1,2,3,4,5 [0,1],使得

r , f(r) 1;

5 2 5, f5( 5) 1;

对于fn( ) sinn , n=1,2,3,4，由

fn( ) 1 cos

2

n 12

dcos

(cos ) cos (h( )) h ( ) d ,x [ 1,1];

(x) (1 x)

2

n 12

其中， cos 1x g(x),

g 1( ) cos h( ).

首先可以产生满足 (x)分布的伪随机数序列，为此，我们注意到

0 1 x 1, (x) (1 x)

2

2n 12

(y) yn 1,y [0,1]

用第一类舍选法抽样，可得

y max( 1, 2, , n),

x cos 1 2max( 1, 2, , n).由此得所求的伪随机数序列

cos 1 1 2max( 1, 2, , n) , fn( ) sinn . 具体地，

1 cos 1 1 2max( 1, 2, 3, 4) , f4( 1) sin4 1; 2 cos 1 1 2max( 1, 2, 3) , f3( 2) sin3 2;

1 2max( 1, 2) , f2( 3) sin 3;

4 cos 1 1 2 1 , f1( 4) sin 4.

3 cos

1

2

[1] see the reference from http://wapedia.mobi/en/N-sphere

2. 解： MC计算步骤：

f(x)

xedx f(x)dx xe xdx f (x)e xdx, f (x) x5/2.

e0000

i). [0,1], i ( i) 1, i 1, ,n;

5/2 x

ii). 首先对偏倚密度函数g(x) e x抽样：

x

F1(x) e xdx 1 e x,

set i F1( i) i ln i;

f(x)

iii). (求出f x) 在各抽样点的值：

g(x)

55

f( i)2

f( i) i ( ln i)2;

g( i)

5

1n

iiii). { ln i,i 1, ,n},I= ( ln i)2.

ni 1

3. 解：函数A(u)的期望值定义为

E{A} A(u)dG(u) A

方差为

V{A} E{(A A)2} E{A2 2AA A} =E{A} A A A

2

2

2

2

2

2

另一方面，设V(Ai) A,i 1, ,N, 则有

1 V{N} 2

N

2

1 A2

V{Ai} 2 N A NNi 1N

当N足够大时，有limV(N) V(A), 即

N

A A

2

2

2 A

N

1. N

2

4. 解：由于f( v) f(v)，故先讨论f(v) Cv2e v,v 0，并设C 0, 0 .

C

ifmax(v0) e 1;选择初始位置： 0 v0

ii). [0,1], 1 ( 1) 1,

def. step i 1 [0, );

iii). [0,1], 2 ( 2) 1, 引入过渡概率： (vi vi 1) min{1,

C(vi i)e

Cve

22

(vi i)2

2 vi2

i

v i (2vi i i2)

min{1, i} e

vi

judge: 2 (vi vi 1),

if it's true, i 1 vi 1 vi i, then

goto ii) and walk for the next step vi 1 vi 2; if not, goto ii) and walk for the step vi vi 1 again.

iiii). ={ 0, 1, , i, , N}, Mathematica plotting program:

ListLinePlot[{{ 0,f( 0)}, ,{ N,f( N)}}], 其中N 为投点数。

5. 解：Metropolis方法可以对无法归一化的分布密度函数进行抽样。

i). 0,选择初始位置： 0 x0 fmax(x0) A, A 0;

ii). [0,1], 1 ( 1) 1,

def. step i 1 [0, );

iii). [0,1], 2 ( 2) 1, 引入过渡概率： (vi vi 1) min{1,

e

(xi i)2

2 x2 i2

e

min{1, e

(2xi i i2)

2

}

judge: 2 (xi xi 1), if it's true, i 1 xi 1 xi i, then

goto ii) and walk for thenext step xi 1 xi 2; if not, goto ii) and walk for the step xi xi 1 again. iiii). ={f( i) Ae 0, 1, , i, , N},

i2

2

.