3.7粒子輸运问题的蒙特卡洛模拟 棋拟法 真换棋拟油是基于粒子贛运过程的隨机鴕计特性的考建,认 为物理上的可巩测量就是大量粒子的行为失同货的统计结。 因此,该方法就是考一个一个粒子的传翁,棋拟它们在物质中 随机运动的历史。记录其在运动中对感兴越的物理棋拟量的货 就。在对单个粒子运动历史进行大量的量复模拟之后,就可以对 物理模拟量进行计平均,;到所曇的物狸結果。 中子与物作用后,一部分会被吸收。另一部分过多次散 射后会穿還物孜层射出去。中子与物孜层作用后可能会产生次 級粒子,仉不考虑这些次级粒子的迁移。如是中子和物质材料 中第m种原子核作用的全微面为 "(E)=o"(E)+aa(E) dm(E)和am(E)分别衰示中子被一个原子核m的散射和吸收微面。 如果单位体积第m种原子核的数量记为pn,则中子作用在单位体
3. 7 粒子输运问题的蒙特卡洛模拟 一、直接模拟法 直接模拟法是基于粒子输运过程的随机统计特性的考虑，认 为物理上的可观测量就是大量粒子的行为共同贡献的统计结果。 因此，该方法就是考虑一个一个粒子的传输，模拟它们在物质中 随机运动的历史，记录其在运动中对感兴趣的物理模拟量的贡 献。在对单个粒子运动历史进行大量的重复模拟之后，就可以对 物理模拟量进行统计平均，得到所需要的物理结果。 中子与物质作用后，一部分会被吸收，另一部分经过多次散 射后会穿透物质层透射出去。中子与物质层作用后可能会产生次 级粒子，我们不考虑这些次级粒子的迁移。如果中子和物质材料 中第m 种原子核作用的全截面为 ( ) E (E) (E) m a m s m σ t = σ +σ (E) m σ s 和 分别表示中子被一个原子核 的散射和吸收截面。 如果单位体积第 种原子核的数量记为 (E) m σ a m m ρ m，则中子作用在单位体
积内第m种元上的总氰面为 ∑"(E)=pno"(E) 假如村料中有多种元景,该中子与村料作用的总微面为 ∑,(E)=∑Em 假定中子与某一原子被散射后的分布表示为dG(E)4,当射 角分布对方位角是各向同性时,方角口可以被积萍,脊到微 分散射薇面do(E)/d(os0)。无论微分戴窗do(E)或者 dσ(E)d(os),我们部可以得到相蠃的理论公式。 设在O点有一个能量为E的中子兽直入射到物质层中。我仉 记录这时该中子的状世形为50=(x0=0, Eo. cos e=1),经过第一次 碰撞后散射到状庵世形S1=(x1,E1,cos1),再经过第二次磁撞后射 到状形s2=(x2,E2,cosB2),…,如此我们依次记下该中子在物 败层中动历史上的形点的孰遼 S0→)S1→S2…→SM 直到在米庵。该中子被物质屈吸收。射背射出来。或者在 s处该中子的能量E氐于某一,则程序就停止跟腺。 親在我仉讨论程序如何具体模拟暇麝这个运动历翟 物始状忞位形如前所逃已給出,假定为 s0=(x0=0,E0,cos0=1) 现在要由位形31=(x-1,E1osO-)确突下一个状庵位形 s,=(x,E;,cos)。 我们采用如下步骤来确定状乏的各个參数:
积内第m 种元素上的总截面为 (cos Ms G = 0 Ei , (E) (E m m t m Σt = ρ σ ). 假如材料中有多种元素，该中子与材料作用的总截面为: Σ ( ) = ∑Σ ( ) m m t E t E . 假定中子与某一原子核散射后的角分布表示为dσ (E)/ dΩ，当散射 角分布对方位角ϕ 是各向同性时，方位角ϕ 可以被积掉，得到微 分散射截面 dσ ( ) E / d(cosθ ) 。无论微分截面 dσ ( ) E / dΩ 或 者 dσ ( ) E / d θ )，我们都可以得到相应的理论公式。 设在O点有一个能量为 的中子垂直入射到物质层中。我们 记录这时该中子的状态位形为 E0 ( 0, ,cos 1) s0 = x0 = E0 θ 0 = G ，经过第一次 碰撞后散射到状态位形 ( ) 1 1 ,cos 1 1 s = x , E θ G 2 ，再经过第二次碰撞后散射 到状态位形 ( 2 2 2 s = x , E ,cosθ ) G ，…，如此我们依次记下该中子在物 质层中运动历史上的位形点的轨迹： M s s s s G G G G 0 → 1 → 2 ⋅⋅⋅ → . 直到在 状态，该中子被物质层吸收，透射或背射出来，或者在 sM 处该中子的能量 低于某一阈值，则程序就停止跟踪。 G EM 现在我们讨论程序如何具体模拟跟踪这个运动历程。 初始状态位形如前所述已经给出，假定为 s0 = ( ) x0 , E0 ,cosθ 0 = 1 G 。 现在要由位形 ( ) 1 1 1 1 , ,cos i− = i− Ei− i− s x θ G 确定下一个状态位形 ( ) i i i s = x , cosθ G 。 我们采用如下步骤来确定状态 i s G 的各个参数：
(1)首先确定坐标参数x。中子到达状痕以前,经历过第 i-1次碰贛后敵訇遠直織运动,其运动的自由翟y演足分布度 函数f(y)=∑,(E-)exp,(E)。我们可以采用直接抽样法得到自 由程y的抽样值 则x由下式给出 x=x+ycos6-1=x压cose (2)确定碰撞的原子核种类。中子与物质层中m原子核碰 的几率为 pm=m(E)、(E) 则可以由高散型分布随机变量的直接抽禅法很容易确定发生磁 撞的是何种原子核。 (3)确定碰撞的性质是吸收还是散射。中子与第m原子核发 生觉射的几事为 pm)=om(E-)σ"(E) 同样可以采用高散型随机变量的直接法抽取。若抽绪是为吸 临,则停止跟隙回到状,,开始下一个中子的跟腺;若抽桦绪 最为散封,则进入第(4)乡。 (4)确定中子散射角日和能量E。由于理论上一般给出的是质 心系中的散射微分微面公式d(E)dcos-,因此我们要首先 换照质心系的微分截面抽取散射角汆弦c0s,cosb足的分布
（1） 首先确定坐标参数 xi。中子到达 i s G 状态点以前，经历过第 次碰撞后做匀速直线运动，其运动的自由程 满足分布密度 函数 i −1 y ( ) [ ( ) 1 1 ( ) = Σt i− ⋅ exp − Σt f y E y y ] Ei− 。我们可以采用直接抽样法得到自 由程 的抽样值， ( )lnξ 1 Σ −1 = − t Ei y . 则xi由下式给出 ( ) 1 1 1 1 1 cos ln cos − − − − − Σ = + = − i t i i i i i E x x y x θ ξ θ （2） 确定碰撞的原子核种类。中子与物质层中第m 原子核碰撞 的几率为 ( ) ( ) 1 1 ( 1) − − − = Σ i Σt i m t i pm E E . 则可以由离散型分布随机变量的直接抽样法，很容易确定发生碰 撞的是何种原子核。 （3） 确定碰撞的性质是吸收还是散射。中子与第 原子核发 生散射的几率为 m ( ) ( ) 1 1 ( 1) , − − − = i m i t m s i pm s σ E σ E . 同样可以采用离散型随机变量的直接法抽取。若抽样结果为吸 收，则停止跟踪回到 状态，开始下一个中子的跟踪；若抽样结 果为散射，则进入第（4）步。 0 s G （4） 确定中子散射角θ i和能量 。由于理论上一般给出的是质 心系中的散射微分截面公式 Ei ( ) 1 cos dσ Ei−1 d θ i− ，因此我们需要首先 按照质心系的微分截面抽取散射角余弦 θ c cos ， θ c cos 满足的分布
粤度函数为 d cos e d 0c/J-1 d cos 0 理论上。射后的中子能量E由下弌计扩得到 E[(1+r)+(1-r)cos] 2 其中r=1-,A是原子核质量与中子质量之比。质心系散射 A+1 0可以用下面公式换算为对应的实验宜系的散射角日2 cosL=(1+ Acose\n+A+2Acos8c 开根据下面的歌面三角公式,通过验蜜系散射角来确定日。 b 6. cos 6, + sin 6, sin 6, cos 其中φ为方位角。由守们考的中子散射过程是各向同性的, 方肩φ是通过抽样φ=2n确定抽祥。 按照上面的计算步暻,我们就完成了从S到§状虍的取瞭。 量复上述中宁眾计算过程,真则中宁在物质暴中运动历程的终 点。如果们一失棋拟了N个中子的动过程,接下来,我们就 要对感兴慜的棋拟物颯量选行统计平均。 以计算投射率來计犷可观测的物理量 我们定义在拟过程中邛n个中子对造射率的贡就n为 ≥ x≤a,或被吸收 其中下标M为该中子在物败层中碰的次数。我们得到穿造物质 层的中子数N为
密度函数为 ( ) ( ) ( ) ∫− − − = 1 1 1 1 cos cos cos cos c c i c i c d d d E d d E f θ θ σ θ σ θ . 理论上，散射后的中子能量Ei由下式计算得到 [ ] ( ) ( ) i i c E E 1 r 1 r cosθ 2 1 = −1 + + − . 其中 2 1 1 + − = A A r ，A是原子核质量与中子质量之比。质心系散射角 θ c可以用下面公式换算为对应的实验室系的散射角θ L。 ( ) θ L A θ c A A θ c cos 1 cos 1 2 cos 2 = + + + . 再根据下面的球面三角公式，通过实验室系散射角θ L来确定θ i。 cosθ i = cosθ i−1 cosθ L + sinθ i−1 sinθ L cosϕ 其中ϕ 为方位角。由于我们考虑的中子散射过程是各向同性的， 方位角ϕ 是通过抽样ϕ = 2πξ 确定抽样值。 按照上面的计算步骤，我们就完成了从 i−1 s G 到 i s G 状态的跟踪。 重复上述中子跟踪计算过程，直到中子在物质层中运动历程的终 点。如果我们一共模拟了 个中子的运动过程，接下来，我们就 要对感兴趣的模拟物理量进行统计平均。 N 以计算投射率来计算可观测的物理量： 我们定义在模拟过程中第n个中子对透射率的贡献η n为 ≤ ≥ = 0 ,或被吸收 1 x a x a M M η n 其中下标M 为该中子在物质层中碰撞的次数。我们得到穿透物质 层的中子数N1为
由此孙到射亭的一个佑计值为 P n 在1-a量信水平下,P的误鑾佔计为 <t./N σn的均方。由于n是一个二项式分布的隨机变量,所以 P(-P)≈P(-P) 过对大量中子运动过程的跟腺,我们也很容易求出适射中子的 能量和肩分布。只要将能量E和极肩日分成普干个小区间,如: E0>E1>E2>…>Em, 0=b0<61<62<…<bn=/2 将适射中子的能量E和极扇0记入国中对应区间,计落入吝个 能量区间或肩度区间的中子数,并画出直方图。这桦觉们就到 相应的射中子的能量分布域角分布图形。 〖接拟法的优痕:模拟过程置了物理过程的机侧,模拟 思丸朴录、蒿阜。 直接模拟法的峡痕:当物赝层毓厚时,适射事会很小,导政 误大。 二、杈量法 实称上在模拟大量粒子穿过校尾物放层时,只有很少教的
∑ . = = N n N n 1 1 η 由此得到透射率的一个估计值为 ∑= = = N n n N N N P 1 1 1 η 在1−α 置信水平下，P 的误差估计为 P − P < tασ η / N . σ η是η n的均方差。由于η n是一个二项式分布的随机变量，所以 = P( ) 1− P ≈ P(1− P) 2 σ η 通过对大量中子运动过程的跟踪，我们也很容易求出透射中子的 能量和角分布。只要将能量E 和极角θ 分成若干个小区间，如： E0 E1 E2 > Emin > > > ⋅⋅⋅ ， 0 = θ 0 < θ1 < θ 2 < ⋅⋅⋅ < θ M = π / 2 将透射中子的能量E 和极角θ 记入图中对应区间，统计落入各个 能量区间或角度区间的中子数，并画出直方图。这样我们就得到 相应的散射中子的能量分布或角分布图形。 直接模拟法的优点：模拟过程重复了物理过程的机制，模拟 思想朴素、简单。 直接模拟法的缺点：当物质层较厚时，透射率会很小，导致 误差较大。 二、权重法 实际上，在模拟大量粒子穿过较厚物质层时，只有很少数的