蒙特卡洛方法 (MC) 方法

MCC 方法
蒙特卡罗方法的诞生
蒙特卡罗方法的产生可追溯到Buffon 投针实验。

法国数学家Buffon 用此实验来估算π值，它的原理是这样子的：在桌面上划一组间距为d 的平行线，然后向桌面上随意抛掷长度为L 的细针，从针与平行线相交的概率就可以得到π值。

其中 [0,)A d ∈ [0,)
x π∈ 由积分性质可得投针置于平行线上的概率为sin 1
2l d l p dAdx d π
θ
π
π
==
⎰
⎰
假如在N 次投针实验中，有M 次与平行线相交，则有2l M P d N
π=
=
图3.2
Buffon 的投针实验
图3.3 投针位置分析
1930年，费米利用蒙特卡罗方法研究了中子的扩散，并设计了一个蒙特卡罗机械装置，用于计算核反应堆的临界状态。

冯.诺依曼是蒙特卡罗方法的正式奠基者,他与Stanislaw Ulam 合作建立了概率密度函数、反累积分布函数的数学基础，以及伪随机数产生器，从而使得蒙特卡罗方法得以推广，成为科学领域一种常用的模拟方法。

蒙特卡罗方法的基本思想
对某一个待解决的物理问题（当这个物理问题可以抽象为数学问题时）建立一个概率模型，即确定某个随机事件X ，使得待求问题的解等于随机事件X 出现的概率或随机变量的数学期望值。

然后进行模拟实验，重复多次地模拟随机事件X 。

最后对随机实验结果进行统计平均，求出X 出现的频数作为问题的近似解。

这就是蒙特卡罗方法的基本思想。

具体来说：
假设所要求的量x 是随机变量的数学期望
，那么近似确定x 的方法是
对进行N
次重复抽样，产生相互独立的值的序列、、……、，并计算其算术平均值：
1
1
N
N
n
n N
ξξ
==∑
根据大数定理有
P (l i m )
N N x ξ→∞
==
因此，当N 充分大时，下式 ()N E x ξξ≈=
成立的概率为1,亦即可以用
作为所求量x 的估计值。

用蒙特卡罗方法求解时，最简单的情况是模拟一个发生概率为P 的随机事件A 。

考虑一个随机变量，若在一次试验中事件A 出现，则取值为1；若事件A 不出现，
则
取值为0。

令q=1-p ，那么随机变
量
的数学期
望
，此即一次试验中事件Ａ出现的概率。

的方差。

假设在Ｎ次试验中事件Ａ出现v次，那么观察频数v也是一个随机变量，其数学期望，方差。

令，表示观察频率，那么按照
加强大数定理，当Ｎ充分大时，式=p成立的概率等于1。

因此由上述模型得到的频率近似地等于所求量p。

这就说明了频率收敛于概率，而且可以用样本方差
2
(1) ()
1
p p p
N
σ
-=
-
作为理论方差2()p
σ的估计值[13]。

碰撞问题
等离子体中存在大量运动着的电子、离子、中性粒子，它们之间不断发生着各种类型的碰撞。

一般的蒙特卡罗碰撞是采用碰撞时间随机的方法。

而PIC方法中的MCC模型不同于一般的蒙特卡罗碰撞。

在PIC方法中，粒子推进和场推进的时间步长是固定的，采用在一个时间步长内随机决定粒子之间是否发生碰撞来实现蒙特卡罗碰撞[10]。

首先给定初始电磁场和初始粒子，在电场和磁场的作用下按照牛顿力学及洛伦滋方程处理碰撞粒子对的位置和速度，得到一个时间步长后的粒子的新位置和速度。

然后根据粒子运动前后的位置和速度，在空间网格上分配电荷、电流密度，最后利用Maxwell方程组求解新的电场和磁场。

再在新的电场和磁场下更新粒子位置和速度，如此循环下去，模拟出等离子体的动态物理过程。

并采用蒙特卡罗模型得到碰撞后的位置和速度。

图3.4 碰撞截面示意图（摘自[12]）
两粒子发生碰撞的碰撞截面定义如上图所示，可得碰撞截面
2(12)r r σπ=+ 每个时间步长内发生碰撞的概率：
,1exp(())c i j i P n
l σε=--∆
程序结构
图程序流程图。