计算物理学练习题及参考解答

计算物理学练习题及参考解答1．计算物理学的英文表示：computatioal physics 或者computer physics2．什么是计算物理学？它与理论物理、实验物理有什么区别和联系？答：计算物理是指以计算机及计算机技术为工具和手段，运用计算数学的方法解决复杂物理问题的一门应用科学。

计算物理方法是除理论方法和实验方法之外的第三种研究手段，计算物理现已成为物理学研究的三大支柱之一，它与实验物理和理论物理的关系如下图：3．计算物理学是物理学、数学、计算机科学三者结合的产物，它也是物理学的一个分支，与理论物理、实验物理有着密切的联系。

4．计算机在物理学中有哪些应用？答：计算机数值分析、计算机符号处理、计算机模拟、计算机实时控制5．计算机技术有各种各样的算法，可以概括为最基本的两类：串行计算和并行计算。

6．理论物理在实际计算中遇到许多困难：非线性问题求解和非对称问题的求解；自变量较多问题求解；非规则界面问题求解等。

7．计算物理的优点有：省时省钱；具有更大的自由度和灵活性；能够模拟极端条件下的实验。

8．第一原理方法是基于量子力学基本原理建立起来的；分子动力学方法是基于经典力学基本原理建立起来的；蒙特卡罗方法是基于统计力学基本原理建立来的。

9．计算机模拟一般有哪两种类型？答：随机模拟和确定性模拟，比如蒙特卡罗模拟和分子动力学模拟。

10．什么是蒙特卡罗模拟？它的应用一般有哪三种形式？答：通过不断产生随机数序列来模拟过程。

直接蒙特卡罗模拟、蒙特卡罗积分、Metropolis蒙特卡罗模拟。

11．蒙特卡洛方法的理论依据答：（1）大数法则：人们发现，在一个随机事件中，随着试验次数的增加，事件发生的频率趋于一个稳定值；人们同时也发现，在对物理量的测量实践中，测定值的算术平均也具有稳定性。

大数法则反映了大量随机数之和的性质。

（2）中心极限定理：中心极限定理，是概率论中讨论随机变量和的分布以正态分布为极限的一组定理。

这组定理是数理统计学和误差分析的理论基础，指出了大量随机变量近似服从正态分布的条件。

中心极限定理告诉我们：在有足够大，但又有限的抽样数n的情况下，蒙特卡洛估计值是如何分布的。

12．请简述著名的巴夫昂（Buffon）投针实验。

并写出用Matlab实现的代码。

（给出方程、算法框图、程序）l 的细针，那么从针与平行线相答：在平滑桌面上画一组相距为s的平行线，向此桌面随意地投掷长度s交的概率就可以得到π的数值。

NMp =，p 2M 2N =≈π clearS=1; %平行线间距 L=1; %针长N=1000000; %总投针次数 M=0; for i=1:Nx=rand*S/2; %针到最近平行线的距离 a=rand*pi/2; %偏角 if(x<sin(a)*L/2)M=M+1; %统计相交次数 end endtestpi=2*N/M %pi 的实验值 13．在考虑蒙特卡罗模拟的精确度时，不能只是简单地减少方差和增加模拟次数，还要同时兼顾计算费用，即机时耗费。

14．假定我们研究连续的随机变量，由随机变量的分布可以得到它取某给定值的概率，即][)(du u u u P du u g +<'<=)(u g 称为u 的概率密度分布函数，它表示随机变量u '取u 到du u +之间值的概率。

而⎰∞-=udx x g u G )()( 则称为u 的分布函数。

15．高斯分布可以由给定的期望值μ和方差2σ完全确定下来，通常用),N(2σμ来表示]2/)(ex p[21),N(222σμπσσμ--=x比如期望值为1，方差为1的高斯分布表达式为]2/)1(ex p[21)1,1N(2--=x π16．对物理问题的计算机模拟所需要的伪随机数应当满足什么样的标准？有哪些统计检验方法？ 答：良好的统计分布特性；高效率；循环周期长；产生程序可以移植性好；可以重复产生。

统计检验有：均匀性检验；独立性检验；组合规律检验；无连贯性检验；参数检验等等。

17．在蒙特卡洛方法应用中减小方差的基本技术：重要抽样法，分层抽样法，控制变量法和对偶变量法。

然而，单独使用这四种减小方差的技巧仍然有其局限性。

人们发展了一些复合蒙特卡洛计算技术，如适应性蒙特卡洛方法和多道蒙特卡洛抽样方法等。

这些蒙特卡洛技巧对于被积函数在积分范围内具有多个尖峰的情况，特别具有实用价值。

18．真随机数列答：真随机数列是不可预测的，因而也是不能重复产生的数据序列。

19．伪随机数列答：通过某些数学公式计算而产生的随机数列 20．同余产生器及程序代码答：一般通过如下的线性同余关系式来产生数列))(m od (x 1n m c ax n +=+m x n n /=ξfloat Random(int n, int m, float seed, float a, float b) {int i; float r; r = seed;for (i = 1; i <= n; i++) r = (a * r + b) % m; return r / m; }21．均匀性检验答：均匀性检验又称频率检验，它用来检验用随机数（样本值）确定的经验频率和均匀分布频率是否有显著性差异。

常用的统计检验方法有2x 检验和累积频率检验（K-S 检验）。

22．随机变量抽样答：指的就是由给定分布函数产生随机数的方法。

首先，在[0,1]区间抽取均匀分布的伪随机数列，再从中抽取满足给定分布密度函数的简单子样，并且各个伪随机数相互独立。

23．连续分布的随机变量抽样一般有哪些方法？答：直接抽样法；变换抽样法；舍选抽样法；复合抽样法；近似抽样方法 24．试述离散型分布的随机变量的直接抽样答：对于离散序列数 ,,,21i x x x 给定每个数的取值概率 ,,,21i p p p ，则我们可定义其分布函数F(x) 如下：∑≤=xx i i p x F )(。

在区间[0,1]上取均匀分布的随机数ξ，判断满足下式的j 值：)()(1j j x F x F <≤-ξ则抽样值η为j x ，η分布符合分布函数F(x)的要求为。

25、试述连续分布的随机变量的变换抽样法。

答：设连续型随机变量η的分布密度函数为)(x f 。

要对满足分布密度函数f(x)的随机变量η 抽样较难时可考虑通过其它已知函数的抽样来得到。

考虑变换)(x h y =，)(y g x =选择)(y φ，使得)())(()()(x h x h dxdyy x f '==φφ 则可对)(y φ抽样得到δ，通过变换)(δηg =，得到满足分布密度函数)(x f 的随机抽样。

更为一般的情况，设连续型随机抽样),(δη的分布密度函数为),(y x f 。

考虑变换),(1y x h u =，),(2y x h v =，J y x h y x h g J v u g y x f )),(),,((),(),(21==，yv x v y u x u J ∂∂∂∂∂∂∂∂=////，这样就可以通过抽取满足分布密度函数),(v u g 随机抽样),(δη''得到待求的满足分布密度函数),(y x f 的随机抽样),(δη，其中),(1δηηh ='，),(2δηδh ='。

26．试给出一个用随机数计算π的Matlab 程序。

（10分） 解：物理模型：如图第一项限中单位正方形内投点在圆内的概率即为单位圆面积的四分之一。

数学方程：⎰⎰--⋅=11222121)1(4x x dx dx θπ算法框图：产生随机点（ξ，η）M 个；统计其中满足条件122≤+ηξ的点的个数N ；计算π值M N /4=π。

Matlab 程序：P=4/100000*length(find(sum(rand(2,100000).^2)<1))27．对指数分布⎩⎨⎧>>=其他0,λx λe f(x)-λλ直接抽样答：第一步：求分布函数⎰⎰--∞-===λx λx x-e 1dx λe f(x)dx F(x)第二步：抽样]1,0[~ξ 第三步：计算分布函数的反函数ληηξ-e -1)F(=＝)1ln(1ξλη--=28．梅氏游走法计算氦原子中两电子间库伦作用的平均值。

（给出动力学方程、数值解方程、算法框图、程序）答：! 梅氏游走法计算氦原子中两电子间库伦作用的平均值 program main dimension x(6)!输入抽样参数Nt,Ng,Nf,Ns,dx Nt=1000 !热化步数Ng=100 !样本组数Nf=10 !抽样间隔Ns=10000 !每组抽样个数dx=0.1 !游走步长!初始化：随机设置初始值xdo 10 i=1,6call random(RND)10 x(i)=0.01*(rnd-0.5) !两电子在核附近!热化：消除初始化影响，趋于平衡分布do 20 j=1,Nt20 call walk(RND,dx,x)!分组间隔抽样，计算平均值和误差su=0su2=0sdu=0do 40 ig=1,Ng !样本分成Ng个组ug=0ug2=0do 30 k=1,Ns !Ns间隔Nf抽样，Ns个样本为一组call walk(RND,dx,x)if(mod(k,Nf).ne.0) goto 30x12=x(1)-x(4)y12=x(2)-x(5)z12=x(3)-x(6)r12=dist(x12,y12,z12)u=1/r12ug=ug+uug2=ug2+u*u !组内求和30 continueug=ug/Nssigmag2=ug2/Ns-ug*ugdug=sqrt(sigmag2/Ns) !组内平均、方差和误差su=su+ugsu2=su2+ug*ug40 sdu=sdu+dug !组间求和avu=su/Ngsigma2=su2/Ng-avu*avudu1=sqrt(sigma2/Ng)du2=sdu/Ngdel=du1-du2 !组间平均、方差和误差!输出avu,du1,du2,del100 open(12,file='out.dat')write(12,1000) Nt,Ng,Nf,Ns,dx,avu,du1,du2,del close(12)1000 format(4i10,5f15.4)end! 计算距离的函数子程序 function dist(x,y,z) dist=sqrt(x*x+y*y+z*z) return end! 计算权重的函数子程序 subroutine weight(x,f) dimension x(6) r1=dist(x(1),x(2),x(3)) r2=dist(x(4),x(5),x(6)) f=exp(-3.375*(r1+r2)) return end! 梅氏游动一步的子程序 subroutine walk(RND,dx,x) dimension x(6),x0(6) call weight(x,f0) do 10 i=1,6 x0(i)=x(i) call random(RND) ! 存旧10 x(i)=x(i)+dx*(RND-0.5) ! 生新 call weight(x,f) call random(RND) if(f.ge.f0*RND) goto 30 !游动 do 20 i=1,620 x(i)=x0(i) !不动 30 return End29．有限差分法答：微分方程和积分微分方程数值解的方法。