蒙特卡罗方法与MCNP程序入门

更新时间:2023-04-27 16:54:01 阅读量: 实用文档 文档下载

说明:文章内容仅供预览,部分内容可能不全。下载后的文档,内容与下面显示的完全一致。下载之前请确认下面内容是否您想要的,是否完整无缺。

蒙特卡罗方法与MCNP程序入门

目录

第1章蒙特卡罗方法简介 (1)

§1.1蒙特卡罗方法的基本原理 (1)

§1.2 蒙特卡罗方法的解题手续和特点 (6)

§1.3 用蒙特卡罗方法模拟粒子输运 (7)

第2章MCNP程序入门 (9)

§2.1 MCNP简介 (9)

§2.2 MCNP程序的组成及特点 (11)

§2.3 MCNP4C程序安装、运行与源程序编译 (13)

§2.4 MCNP输入文件 (16)

第3章MCNP程序中的几何构建 (23)

§3.1 基础知识 (23)

§3.2 几何描述卡 (26)

§3.3 有效地构建几何 (30)

第4章MCNP程序的数理基础 (33)

§4.1 物理 (33)

§4.2 记数 (41)

§4.3 减小方差技巧 (47)

第5章MCNP程序中的数据卡 (56)

§5.1 问题类型卡 (56)

§5.2 栅元参数和曲面参数卡 (56)

§5.3 源的描述 (65)

§5.4 记数方式的指定 (74)

§5.5 材料的指定 (88)

§5.6 能量和热处理方式的指定 (90)

§5.7 问题截断条件 (93)

§5.8 外围卡 (95)

§5.9 MCNP输入文件综述 (98)

第6章经验 (102)

§6.1 一般应用步骤 (102)

§6.2 需注意的问题 (102)

第7章应用实例 (104)

§7.1 医学物理中的应用 (104)

§7.2 反应堆物理计算中的应用 (129)

附录 (157)

连续能量中子截面库ENDL851数据目录 (157)

中子热截面库BMCC1数据目录 (158)

离散中子截面库D91数据目录 (159)

光子截面库MCPLIB1数据目录 (160)

特殊材料S(α,β)热截面库TMCC1数据目录 (162)

参考文献 (164)

- i -

蒙特卡罗方法与MCNP 程序入门

- 1 - 第1章 蒙特卡罗方法简介

蒙特卡罗方法,又称随机抽样方法,是一种与一般数值计算方法有本质区别的计算方法,属于试验数学的一个分支,起源于早期的用几率近似概率的数学思想,它利用随机数进行统计试验,以求得的统计特征值(如均值、概率等) 作为待解问题的数值解。 随着现代计算机技术的飞速发展,蒙特卡罗方法已经在原子弹工程的科学研究中发挥了极其重要的作用,并正在日益广泛地应用于物理工程的各个方面,如气体放电中的粒子输运过程等。

§1.1 蒙特卡罗方法的基本原理

就数学特性而言,蒙特卡罗方法的发展可以追溯到18 世纪著名的蒲丰问题。 1777年,法国科学家蒲丰(Buffon) 提出用投针试验计算圆周率π值的问题。 这里我们用蒲丰问题来初步说明蒙特卡罗方法的基本原理和解决问题的基本手续。

蒲丰问题是这样一个古典概率问题: 在平面上有彼此相距为2a 的平行线,向此平面任意投一长度为2l 的针,假定l

等概率落入,即M l 距平行线的距离x 均匀分布在区间[0,a]之内; (2) 针与线的夹角θ均匀分布在区间[,]22

ππ-之内;(3) x 与θ互相独立。 如图1.1 所示,建立与平行线垂直且原点在某一条平行线上的x 轴,不失一般性,假定针的中心处于图示中的x 轴上。由于对称性,我们只需分析针中心处在(0,)x a ∈范围的情况即可。令探针中心的坐标值为x,显然,只有x ≤l 时才可能发生相交的事件。我们来分析在条件x ≤l 满足时, 针与线相交的概率: 只有当0cos x ar l

θθ≤=时才能相交,且相交的概率为 12cos x P arc l

π= (1.1) 下面再来分析针中心位置在轴上的分布, 显然, 这是一个均匀分布, 即针中心处于区间(x,x+dx)内的概率为

蒙特卡罗方法与MCNP 程序入门

- 2 -

2dx

dP a

=

(1.2) 这样,一次投掷,针中心落入(x,x+dx)且与线相交的概率为

122arccos x

dP PP dx a l π==

(1.3) 则一次投掷,针与线相交的总概率为

022arccos l

x l

P dP dx a l a

ππ===??

(1.4) 即:

2l

Pa

π=

(1.5) 从(5)式可见, 可利用投针试验计算π值:设投针N 次,其中n 次针与线相交,则可用频率值n/N 作为概率P 的估计值,从而求得π的估计值为

2l N

a n

π≈

(1.6) 这就是早期的用频率值作为概率近似值的方法的应用实例,表1是在历史上一些有名的用投针试验计算π值的结果,其中针长以a 为单位。

表1.1 投针试验计算π值的结果

实验者 时间(年份) 针长 投针次数 相交次数 π的估值 Wolf 1850 0.8 5,000 2,532 3.1596 Smith 1855 0.6 3,204 1,218.5 3.1554 De

1860

1.0

600

382. 5

3.137

图 1.1

蒲丰问题的概率分

蒙特卡罗方法与MCNP 程序入门

- 3 - Fox

1884 0.75 1,030 489 3.1595 Lezzerini

1901 0.83 3,408 1,808 3.1415929 Reina 1925 0.5419 2,520 859 3.1795

需要指出的是,上述由投针试验求得π的近似值的方法,是进行真正的试验,并统计试验结果,要使获得的频率值与概率值偏差小,就要进行大量的试验,这在实际中,往往难以做到。可以设想,对蒲丰问题这样一个简单的概率问题,若要进行10万次投针试验,以每次投针、作出是否相交判断并累加相交次数用时5 秒钟计算,则需用时50万秒,即大约139个小时。 那么,可以设想,对于象上述确定条件下的核裂变、直流气体放电中粒子的输运过程及粒子输运的总效应,若要用多次掷骰子的方法近似求出就是不可能的了。所以,在现代计算机技术出现之前,用频率近似概率的方法——抑或称为雏形时代的蒙特卡罗方法——并没有得到实质上的应用。

若用数值模拟方法代替上述的真正的投针试验,是利用均匀分布于(0,1) 之间的随机数序列,并构造出随机投针的数学模型,然后进行大量的随机统计并求得π的近似值。

如图1.2建立坐标系,平面上一根针的位置可以用针中心M l 的坐标x 和针与平行线的夹角

θ来决定, 在y 方向上的位置不影响相交性质。

任意投针,意味着x 与θ都是任意取的。但θ的范围可限于[0,π], x 的范围可限于[0,a]。在这种情况下, 针与平行线相交的数学条件是

sin ,0x l x a θ≤≤≤ (1.7)

其次,怎样模拟投针呢? 亦即如何产生任意的[x,θ]。x 在[0,a ]任意取值,意味着x 在

图1.2 用数值模拟方法计算蒲丰问

蒙特卡罗方法与MCNP 程序入门

- 4 - [0,a]上取哪一点的概率都一样,即x 的概率密度函数为

11,0x a ()0,x f x a ?≤≤?=???当时当为其他值时

(1.8) 类似的,θ的概率密度函数为

21,0()0,f θπθπθ?≤≤?=???当时当为其他值时

(1.9) 由此,产生任意(x,θ)的过程就变为由f 1(x)抽样x,由f 2(θ)抽样θ的过程。 容易得到 1

2x a ξθπξ== (1.10)

式中,ξ1,ξ

2 均为(0,1)上均匀分布的随机数。只要随机数的均匀性和独立性良好, 如此构造的数值模型就很好地模拟了实际试验中的一次投针, 并用下式判断是否相交且记录统计结果:

1,sin (,)0,i i i i x l i s x x θθ≤?=??

当时当为其他值时 如果投针N 次,那么

1

1(,)N

N i i i s s x N θ==∑ (1.11) 是相交几率P 的估计值。这样就实现了用数值方法模拟真正的投针试验。用此方法计算的π的近似值的情况如表1.2所示。

表1.2 用蒙特卡罗方法计算的π的近似值

表2中的计算结果表明,随着模拟投针次数的增大,所计算的π的近似值越来越接近于其真值,而要进行这样的数值模拟,就需要很大的计算量,只有利用计算机才能实现。

从蒲丰问题可以看出,用蒙特卡罗方法求解问题时,应建立一个概率模型,使待解问题与此概率模型相联系,然后通过随机试验求得某些统计特征值作为待解问题的近似解。与此相似,在一些物理问题,如核裂变、直流气体放电等过程中,粒子的输运过程及粒子输运的总效应,也是可以与某些概率过程联系起来,例如,电子与原子、分子、离子的碰撞过程,实际上就是与碰撞截面有关的概率过程,这样,从数学物理特征来说,类似于用随机投针方

蒙特卡罗方法与MCNP 程序入门

- 5 - 法计算π的近似值,确定条件下的核裂变、直流气体放电中粒子的输运过程及粒子输运的总效应可以用多次掷骰子的方法近似求出。

随着现代计算机技术的出现和飞速发展,用计算机模拟概率过程,实现多次模拟试验并统计计算结果,进而可获得所求问题的近似结果. 计算机的大存储量、高运算速度使得在短时间内,获得精度极高且内容丰富的模拟结果。在历史上,也正是原子弹工程研究初期阶段的工作,为模拟裂变物质的中子随机扩散,提出了运用大存储量、高运算速度计算机的要求,这也成为当时推动计算机技术发展的重要动力,也就是在第二次世界大战期间,冯.诺依曼和乌拉姆两人把他们所从事的与研制原子弹有关的秘密工作——对裂变物质的中子随机扩散进行直接模拟——以摩纳哥国的世界闻名赌城蒙特卡罗(Monte Carlo) 作为秘密代号来称呼。用赌城名比喻随机模拟,风趣又贴切,很快得到广泛接受,此后,人们便把这种计算机随机模拟方法称为蒙特卡罗方法。需要指出的是,正是由于广泛领域的物理问题中存在着大量的随机过程,如粒子间的碰撞等,使得蒙特卡罗方法在计算机物理和物理工程中得到日益广泛的应用,并成为沟通理论与实验研究的一个桥梁。

需要指出的是,蒙特卡罗方法不仅在处理具有概率性质的问题方面获得广泛的应用,对于具有确定性问题的计算也因其程序简单等优点获得了广泛的应用。这里以定积分的计算简要说明其处理确定性问题的手续。

对于定积分

()b

a s f x dx =? 通过变量替换,可以转换为下面的形式

1

0()s k g x dx =? (1.12) 其中g(x)∈(0,1),当x ∈(0,1)时即转换为求积分1

0()g x dx ?,亦即求边长为1的正方形中一个曲边梯形的面积的问题,如图1.3所示。

蒙特卡罗方法与MCNP 程序入门

- 6 - 我们可以设想这样一种随机投点求定积分10

()g x dx ?的方法:在一个边长为1的正方形上并以其两边分别为坐标轴画出曲线g(x),实际上就是图3,然后随机地正方形投掷小球,那么,小球击中g(x)曲线下部分的概率就等于所要求的积分1

0()g x dx ?, 这样就将确定性的定积分问题转化为一个概率问题,同样可以通过数值模拟方法——蒙特卡罗方法求得其近似解。用此方法,我们计算了积分1

20x dx ?, 当投球数为1万次时, 得到的积分近似值为0.332800 , 与其真值13

极为接近。 §1.2 蒙特卡罗方法的解题手续和特点

在用蒙特卡罗方法解算问题时,一般需要这样几个过程:构造或描述概率过程,对于本身就具有随机性质的问题,如粒子输运问题,主要是正确地描述和模拟这个概率过程。对于本来不是随机性质的确定性问题,比如计算定积分、解线性方程组、偏微分方程边值问题等,要用求蒙特卡罗方法求解,就必须事先构造一个人为的概率过程,它的某些参量正好是所要求问题的解。即要将不具有随机

性质的问题,转化为随机性质的问题。这构成了蒙特卡罗方法研究与应用上的重要问题之一,然后建立各种估计量,使其期望值是所要求解问题的解,最后根据所构造的概率模型编制计算程序并进行计算,获得计算结果。

图1.3 用蒙特卡罗方法求定积分

蒙特卡罗方法与MCNP程序入门

与其他的数值计算方法相比,蒙特卡罗方法有这样几个优点:

(1) 收敛速度与问题维数无关,换句话说,要达到同一精度,用蒙特卡罗方法选取的点数与维数无关;计算时间仅与维数成比例。但一般数值方法,比如在计算多重积分时,达到同样的误差,点数与维数的幂次成比例,即计算量要随维数的幂次方而增加。这一特性,决定了对多维问题的适用性;

(2) 受问题的条件限制的影响小;

(3) 程序结构简单,在电子计算机上实现蒙特卡罗计算时,程序结构清晰简单,便于编制和调试;

(4) 对于模拟象粒子输运等物理问题具有其他数值计算方法不能替代的作用。

蒙特卡罗方法的弱点是收敛速度慢,误差大的概率性质,这一情况在解粒子输运问题中仍然存在,除此之外,经验证明,只有当系统的大小与粒子的平均自由程可以相比较时,一般在10个平均自由程左右,这方法算出的结果较为满意。而对于大系统深穿透问题,算出的结果往往偏低。对于大系统,其他数值方法往往很适应,能算出较好的结果。因此,已有人将数值方法与蒙特卡罗方法联合起来使用,克服这种局限性,取得了一定的效果。

自二战以来,随着计算机技术的飞速发展,蒙特卡罗方法以其独特的优点广泛应用于计算物理和物理工程领域,特别是在核物理、辐射物理、数学、电子学等方面得到了广泛的应用。对蒙特卡罗方法的推广必将使其对物理学学科的发展发挥更大的作用。

§1.3 用蒙特卡罗方法模拟粒子输运

中子和光子在物质中输运的宏观表现是大量粒子与原子核微观作用的平均结果,蒙特卡罗方法通过逐一模拟和记录单个粒子的历程来求解输运问题。要得到比较合理的平均结果需要跟踪大量的粒子,至于单个粒子在其生命中的某一阶段如何度过,可以在已知统计分布规律的前提下通过抽取随机数来决定。

- 7 -

蒙特卡罗方法与MCNP程序入门

图1.4 单个中子随机历程示意

图1.4 显示了模拟中一个中子射入物质后的随机历程。首先根据中子与物质作用的物理规律(分布函数),选取一个随机数决定中子在何处与原子核碰撞,本例中在事件1处碰撞;然后再用抽取随机数的方法决定中子与原子核发生了哪种反应,这里抽出的是非弹性散射反应;散射中子的能量和向哪个方向飞行也是用抽取随机数的方法从已知分布函数中决定的;碰撞过程中是否产生光子以及光子的能量、飞行方向等参数还是要通过抽取随机数从已知分布中决定,这里产生了一个光子,在MCNP程序中该光子暂时被存储起来,散射后的中子在事件2处与原子核发生(n,2n)裂变反应,同时产生两个新的中子和一个光子,裂变产生的光子和一个中子被存储起来,现在开始对另外一个中子进行跟踪。这个中子在事件3处被俘获,这时取出那个刚才存储起来的中子,通过随机抽样,该中子在事件4处逃逸出平板,进入探测器,MCNP 结束对该中子的模拟,重新取出裂变产生那个光子进行跟踪。这个光子在事件5处与原子核发生了一次散射,在事件6处逃逸出平板。现在取出最后一个粒子——在事件1处产生的光子进行跟踪,该光子在事件7处与原子核碰撞并被吸收。至此,入射中子的整个历史也就完成了, 其有一个中子到达了探测器,感兴趣的结果被记录下来。跟踪越来越多的入射粒子历程后,平均结果就能反映出宏观效果。通过以上描述,读者不难领略蒙特卡罗方法如何通过跟踪粒子历程的方法计算问题,也了解了随机数在蒙特卡罗计算中的独特作用。需要注意的是,MCNP 程序在重新取回以前存储的粒子时,是按照这样的规则进行的:最后存储的粒子最先被取出,也就是“先进后出”或“后进先出”的堆栈操作。

- 8 -

蒙特卡罗方法与MCNP程序入门

第2章MCNP程序入门

由于蒙特卡罗(Monte Carlo)方法在解决粒子输运问题方面的优越性,人们开发了很多种Monte Carlo模拟程序,如EGS、MCNP、MORSE和FLUKA等,其中,EGS4(Electron Gamma Shower,Version 4)和MCNP4C(Monte Carlo N—Particle Transport Code,Version 4C)是应用最多并得到广泛验证的两种经典Monte Carlo模拟程序。EGS4是山美国斯坦福直线加速器中心、日本高能物理国家实验室(KEK)和加拿大国家研究所(NRC)联合推出的——套模拟电子和光子在物质中输运过程的通用Monte Carlo程序系统,是揭示电子和光子在物质中输运规律的有力和较为方便的模拟研究工具,目前KEK发布了EGS5版本。

§2.1 MCNP简介

MCNP全名为Monte Carlo Neutron and Photo Transport Code (蒙特卡罗中子-光子输运程序),它是由美国Los Alamos国家实验室应用理论物理部(X部)的Monte Carlo小组(X-6X小组)经过数十年的研究开发的一个基于蒙特卡罗方法的大型的多功能Monte Carlo粒子输运程序。从1977年开始产生到现在历经十几个版本,解决了核能领域很多关键性问题,功能也越来越强大。

现在MCNP可在微机的的UNIX、LINUX、DOS、WINDOWS 98、Windows XP等操作系统下工作。现其最新版本的MCNP-5 1.30具有如下功能:

(1)非带电粒子成相技术。在用户指定的栅格中,MCNP-5使用多个点探测器来确定某个象素区域的粒子流量。用户可以根据需求设置尽可能多的探测点以便生成尽可能平滑的图象。

(2)随机几何能力。该能力可用来分析颗粒燃料,还可用来研究燃料核在石墨矩阵中的随机位置。

(3)可处理复杂三维几何系统的输运问题,几何界面除任意平面和二阶曲面外,也可包括四阶椭环面。

(4)粒子输运方式可以是中子输运、光子输运、电子输运、中子-光子联合输运、光子-电子联合输运、电子-光子联合输运、中子-光子-电子联合输运。既可用于求解通常的

- 9 -

蒙特卡罗方法与MCNP程序入门

输运方程,也可解多群共轭输运方程。MCNP-5已经能够处理低能光子相互作用的不连续散射问题。

(5)既可计算穿透问题,也可计算临界特征值问题。对临界特征值的计算,给出了KEFF 、预期寿命和生存时间的计算方法,还可计算各种记数关于介质成分、密度或截面数据的一阶、二阶微扰量。

(6)配备的截面数据覆盖了所有常用的核素和同位素,并可选用点截面方式或多群截面方式。可处理的中子能量范围为10E-11至20MeV,光子和电子为0.001至1000MeV。

(7)有多种物理量的计算选择,包括点通量、界面通量、任意独立栅格的粒子流及通量、几何体上的通量及能量沉积,可给出按空间、时间、能量的谱(分布)和联合分布,粒子流还可增加角度分布。另外还有脉冲高度谱的计算。所有量的计算结果都同时给出了统计误差。

(8)有十余种降低方差技巧可以选用。

(9)具有较好的图形输出功能。可作几何结构任意方位的二维剖面图,计算结果的二维X-Y曲线图、等值线图、三维曲面图、截面数据图及粒子径迹图。既可交互方式作图,也可用批命令方式作图,也可在中途暂挂计算进行作图。MCNP-5可绘制最多有64种颜色的彩图。

(10)既可在串行环境下运行,也可在并行环境下运行。

另外,根据具体应用的需要,MCNP还有一些专门用途的版本。如MCNP-X和MCNP-PoliMi。MCNP-X(Monte Carlo N-particle code)是由美国Los Alamos 实验室开发的用于多粒子传输模拟计算的蒙特卡罗程序。它是由MCNP 4B版本发展而来, 可以跟踪的粒子数达到40余种,可跟踪的中子最高能量达到GeV以上。在辐射防护研究中也得到了广泛的应用。MCNP-PoliMi (Politecnico di Milano)是由意大利米兰工艺学院核工程学部根据MCNP 4C 开发用于(时间)相关性测量的蒙特卡罗软件。它可以很详细的记录粒子-核子的相互作用信息,比如每个粒子关于每一次碰撞的类型、目标核、能量沉积、位置等。MCNP-PoliMi在处理中子输运问题时,首先抽取样品中子碰装类型, 然后产生光子,这与标准 MCNP采用的技术相反。一个与MCNP-PoliMi配套的后加工处理程序也已经开发,用它能够定制模拟特定特征的探测器。该软件可用于模拟塑料闪烁体对快中子或伽吗射线的探测问题。

MCNP程序的应用范围十分广阔,主要包括:反应堆设计、核临界安全、辐射屏蔽和核防护、探测器的设计与分析、核测井、个人剂量与物理保健、加速器靶的设计、医学物理与放射性治疗、国家防御、废物处理、射线探伤等。

- 10 -

蒙特卡罗方法与MCNP程序入门

§2.2 MCNP程序的组成及特点

A、组成

MCNP是主程序模块,它根据正在运行的问题之需要分别调用IMCN、 MCRUN、XACT和 PLOT 等主要模块。各模块之间的关系,如图1所示。

图2.1 MCNP程序模块构成简图

IMCN模块总是要被调用的,其任务是读入输入文件INP,并对输入数据进行分析处理。它根据输入信息定义动态分配存储空间的大小,对相应的动态分配数组建立相对地址索引表,随后将调用所需要的模块如RDPROB、TRFMAT、ITALLY 和VOLUM等。RDPROB模块主要是依序读人输入文件INP的卡片映像,通过分析处理送入相应的存储位置。TRFMAT模块对非标准形式输入的曲面描述卡,完成坐标变换,给出曲面标准形式的方程系数;对输入的源分布处理成便于抽样的形式。ITALLY模块对输入的记数及材料说明信息进行加工处理。VOLUME模块则主要是计算栅元的体积及曲面的面积,但限于具有对称轴的情况。MCRUN模块是输运计算的实体,是MCNP的核心部分。它将根据需要调用两个模块TRNSPT和OUTPUT。TRNSPT模块按用户指定的样本数NPS或时间限CTME完成对粒子历史的模拟。OUTPUT模块则实现对计算结果的编辑输出。XACT模块的任务是从截面数据文件中读出问题中所用核素的截面数据,并根据用户给出的信息删去所关心能量范围以外的中子截面数据。PLOT模块实现在各种图像设备上绘制或显示问题的几何图形。它分别调用SETUP模块和WRAPUP模块。SETUP模块根据用户的终端键盘信息绘制几何图形,产生图形文件PIX。 WRAPUP模块在不同的图像设备上实现PIX文件的图形输出。

显然,上述各模块并不是运行任何问题都必须全部调用的,MCNP从分析用户在输入文件INP上指定的信息,决定所需要调用的模块。

- 11 -

蒙特卡罗方法与MCNP程序入门

B、特点

MCNP程序中物理量的单位规定为:

表2.2 MCNP程序中物理量的单位规定

物理量单位

长度cm

能量MeV

时间shake, =10-8sec

温度MeV (kT)

原子密度atoms/barn-cm

质量密度g/cm3

截面barns (10-24cm2)

加热量MeV/collision

此外,原子质量按中子质量为1.0计算,这种单位下阿佛伽德罗常数是0.59703109;程序运行时间以分钟为单位。

MCNP程序具有很强的通用性,主要体现在:

1.可以处理任意三维几何结构的问题。在输入文件INP中,空间被曲面(surface)分割成相互邻接的区域,称为栅元(cell),可以给栅元填充各种物质。栅元的界面可以是平面、二阶曲面或某些四阶曲面(如椭圆环状面)。

2.可以模拟中子输运、光子输运和二者联合输运。

3.用户可以非常方便地在任何位置指定体源、面源、线源或点源,设置源粒子位置、能量、时间、飞行方向等参数的分布。

4.程序提供多种记录模拟结果方法,包括通过某一界面的粒子流量或通量、进入某一栅元的通量、沉积能量和点通量。模拟结果在MCNP中称为记数(tally),可以按位置、能量、时间、粒子来向和粒子种类记数。

5.程序包携带了大量核反应数据库文件,包括连续和离散的中子截面库、光子点截面库、热中子点截面库等,几乎可对所有天然物质进行计算。程序能比较精细地模拟中子和光子输运过程,并对一些特定的物理过程允许用户选择使用哪种方式进行处理,如对热中子处理可选用自由气体模型或S(α,β)模型,对低能光子处理可以考虑或忽略相干散射等。

6.为了提高计算时效,给用户提供了许多可选用的减小方差(variance)技巧,主要包括:重要抽样、权重截断和轮盘赌、时间和能量截断、模拟俘获、指数变换、强迫碰撞、能量分

- 12 -

蒙特卡罗方法与MCNP程序入门

裂和轮盘赌、源的偏倚、点探测器记数、确定论输运、权窗等。

7.用户可通过设置源粒子数或运行时间来通知程序何时终止运行,还可以在原有计算结果的基础上接续运行程序。一些结果不会因计算的意外中断而丢失。

8.在输出文件OUTP中给用户提供丰富的信息,包括输入列表、使用的截面表、粒子生成和丢失表、栅元中的粒子活动情况、中子诱发光子表、记数和记数涨落表等,还可以根据用户要求给出其它信息。

9.提供了简单的问题调试工具。

§2.3 MCNP4C程序安装、运行与源程序编译

A、安装

通常我们能拿到的MCNP原始文件有两种:一种是可以直接运行的,这个无所谓安装,只要你把它放在本地硬盘的一个合适位置就可以了,多数人都是采用这种方式;另一种是有源代码的,里面含有安装文件C700DOS2.EXE。对于这种具体安装步骤如下:

1.找到并双击C700dos

2.exe;

2.把文件解压到目标文件夹下,如d:\mcnp(或者其它);

3.在MCNP目录下创建一个新的文件夹work;

4.在d:\mcnp\exe路径下,复制mcnp_dvf.exe并且改名为mcnp.exe ;

5.把mcnp.exe从d:\mcnp\exe移动到d:\mcnp\work下;

6.在d:\mcnp\exe路径下,复制xsdir2_dvf并且改名为xsdir;

7.修改xsdir文件第一行:把DATAPATH=\mcnp\xs修改为DATAPATH=..\xs ;

8.把xsdir从d:\mcnp\exe移动到d:\mcnp\work下;

9.运行时把输入文件也放在work目录下(就是说输入文件要和mcnp.exe、xsdir放在相同的目录下)。

B、运行

下面以前面安装的工作为基础简单介绍在XP下的运行方法。若用户不修改源程序,MCNP 的输入文件包括截面数据库文件、截面库目文件XSDIR、问题输入文件INP等。INP文件是用户要填写的主要输入文件,一般把该文件特指为输入文件。OUTP是MCNP的主要输出文件

- 13 -

蒙特卡罗方法与MCNP程序入门

(文本格式),其它输出文件还有转储文件RUNTPE(二进制格式)、运行信息文件OUTPUT等。所有文件的文件名不能超过7个字符。

运行的方法是:点击“开始”——“运行”,输入cmd。在命令提示符下输入d:,回车;输入cd mcnp\work,回车;输入cd work,回车,进入到运行目录了。

(1)如果输入卡文件名是INP,即使用缺省文件名,则运行的命令可以是:

d:\mcnp\work> mcnp

(2)如果不是使用以上缺省文件名,则在命令行中重新指定。假设输入卡文件名是mcin,那么运行命令可以是:

mcnp inp=mcin outp=mcout runtpe=mcruntpe

或:

mcnp i=mcin o=mcout ru=mcrntpe

如果在当前文件夹中存在与mcnp将要产生的文件名相同,则mcnp自动按字母顺序将文件名最后一个字符改变为字母表中下一字符。

(3)如果想让在一次运行中所产生的文件,有相似的文件名,使用如下方式:

mcnp name=job1

运行后产生的OUTP 文件将为JOB1O,RUNTPE 为JOB1R。如果文件已经存在,MCNP将不会覆盖它们,而是给出一条消息并中断运行。

(5)如果要看几何绘图,先删除work目录下刚生成的job1o和job1r两个文件,再输入下面命令并回车:

mcnp ip n=job1

如果能绘图的话会弹出plot提示符,在此例如你可以输入ex=100,就可以看到自己的几何描述了,详细的参数看手册。

C、源程序的编译

MCNP采用适应多环境能力的设计思想,即通过选择不同参数对程序文件进行加工处理,可得到能分别在UNIX(SUN、HP、IBM RS/6000、SGI、DEC等机型)、CRAY(UNICOS、COS、CTSS等系统)、DEC、VAX等多种不同环境下编译、运行的源程序。为了能对程序作诸如计数修改、添加有关子程序、截面参数调用程序编制和移植开发等方面的一些工作,用户首先必须根据自己的计算机系统,从MCNP安装包中分离出适应用户机器环境的可编译源程序,然后才能进行相关工作。

- 14 -

蒙特卡罗方法与MCNP程序入门

(1)系统配置及相关文件

编译处理的MCNP程序版本为4C版本,需具备表2.3中所列出的文件。

表2.3 MCNP4C源程序相关文件

源程序文件说明

MCSETUP.ID Fortran程序,用于制作编译配置文件

PRPR.ID Fortran程序,用于对MCNP程序预处理

MAKXS.ID Fortran程序,用于处理截面数据库

MCNPC.ID C程序,MCNP的C程序部分

MCNPF.ID Fortran程序,MCNP的Fortran程序部分

出于维护方便MCNP源程序主体是用标准Fortran 77 写成的,非标部分比如:系统相关特性,时间相关,图形,内存动态分配部分是用C语言写成的,为了能对这两种语言进行编译,对于4C版本WINDOWS XP下可安装Digital Visual Fortran ,版本6.0;Microsoft Visual C++ ,版本6.0。假设它们都被安装在E:盘Program Files下。

(2)编译器环境设置

建议先创建一个目录,如F:\MCC,以下操作均可在该目录下进行。为使编译过程简化,可在命令提示符下使用行编译命令,为此需要将E:\Program Files\Microsoft Visual Studio\DF98\BIN下的DFVARS.BAT,及E:\Program Files\Microsoft Visual Studio\VC98\Bin下的VCVARS32.BAT复制至F:\MCC目录下。在命令提示符窗口中,设定F:\MCC为当前工作目录,执行上述两个批处理文件,这一步非常重要,它将直接影响下面操作的顺利进行。

(3)产生编译配置文件

将安装盘中的MCSETUP.ID、PRPR.ID、MAKXS.ID、MCNPC.ID、MCNPF.ID复制到MCC目录下,并将MCSETUP.ID改名为MCSETUP.for,在命令提示符窗口中执行下面的编译命令,F:\MCC>f90 mcsetup.for

生成可执行文件mcsetup.exe,运行该程序,进入MCNP编译设置主菜单界面,选择适应目标机器环境的设置,具体设置如下:

- 15 -

蒙特卡罗方法与MCNP 程序入门

- 16 - ? 输入“1.1”,在出现的计算机系统描述菜单中选择“2”(PC DVF WINDOWS );

? 输入“3.3”,在出现的绘图选项菜单中选择“5”(DVF QuickWin );

? 输入“4.1”,在出现的截面数据路径菜单中选择“1”,出现路径输入提示行,假设机器上截面数据被存放在D:\MCNP\XS 下,应该输入D:\MCNP 。

以上设置完成后输入“P ”(Process ),回车,让程序完成最后处理工作,最终将在当前目录F:\MCC 下生成文件install.ans 、makemcnp.bat 、patchc 、patchf 四个文件。其中install.ans 记录了各设置选项,makemcnp.bat 是MCNP 程序编译批处理文件,patchc 与patchf 是源程序预处理标识文件。

(4) 编译源程序

在命令提示符窗口中,执行批处理makemcnp.bat ,系统将自动完成预处理及编译任务,生成可执行程序prpr.exe 、makxsf.exe 、mcnp.exe ,同时生成目录flib 与olib ,在flib 目录下存放的是MCNP 源程序的各个模块文件,olib 目录下存放的是各个模块编译后的目标文件。用户可以根据实际需要修改对应源程序,比如使用自定义源或记数箱,之后只需对修改的模块重新编译成目标文件,然后与其它模块的目标文件重新连接即可。

§2.4 MCNP 输入文件

A. 输入文件的基本形式

INP 文件有初始运行和接续运行两种形式,用户须在INP 中描述问题的几何、材料、源、记数和其它要求。INP 由一些被空行分隔的输入块组成,主要的输入块是信息块、标题和栅元块、曲面块和数据块等。输入块又由一些被称为卡的输入行组成。初始运行输入文件的格式如下:

信息卡 } 可选

空行分隔符

蒙特卡罗方法与MCNP程序入门

标题卡

栅元卡

.

.

.

空行分隔符

曲面卡

.

.

.

空行分隔符

数据卡

.

.

.

空行结束符

其它

信息卡的1~8列应填写MESSAGE:,后面跟着用空格分隔的参数项。可用A=B参数项更改输出文件名,如OUTP = MYOUT。信息块是可选的。

在信息块之后的第一行是问题的标题卡,它仅限于一行,占用1~80列,可以是任何信息,将作为OUTP文件中各个输出表的标题被复制。

用户在栅元块和曲面块中描述问题的几何。栅元由栅元卡描述。空间必须由彼此相邻的栅元填满,栅元之间不能重叠,也不能出现无栅元的空区,否则会出现错误。构建栅元的曲面由曲面卡定义,曲面卡在曲面块中给出。

曲面块之后是数据块,在数据块中用户描述源、记数方式、材料等。

无论数据块后是否有空行结束符,MCNP都能运行。用户可以把希望保留的一些附加信息写在数据块的空行之后,MCNP会将它们复制到OUTP文件末尾。

- 17 -

蒙特卡罗方法与MCNP 程序入门

- 18 - B. 一个简单的例子

为说明如何填写INP 文件,这里例举一个简单问题。如图2.2所示,在一个边长10cm 的石墨立方体3中有两个半径0.5cm 的球形空间,球1中

充满氧气,球2是铁球。在球1中置一14MeV 各向同

性中子点源,计算球2外表面与能量相关的中子通

量。建立的INP 文件如下:

SAMPLE PROBLEM INPUT DECK

1

1 –0.0014 -7 2

2 –7.86 -8 3

3 –1.60 1 –2 –3

4 –

5

6

7

8 4 0 -1:2:3:-4:5:-6

←空行

1 PZ -5

2 PZ 5

3 PY 5

4 PY -5

5 PX 5

6 PX -5

7 S 0 –4 –2.5 .5

8 S 0 4 4 .5

←空行

MODE N

IMP:N 1 1 1 0

SDEF POS=0 –4 –2.5 ERG=14

F2 8

E0 1E-5 1E-4 1E-3 .01 .1 1 2 3 4 5 6 7 8 9 10 11 12 13

M1 8016 1

M2 26000 1

M3 6000 1

4

321874321y 图2.2 例子的几何示意图

蒙特卡罗方法与MCNP 程序入门

- 19 - NPS 10000

空行

本例中没有信息块,第一行是标题卡,之后至空格前为栅元块。栅元卡上依次填写栅元号、材料号、密度和构成栅元界面的曲面号(带正负号),这里定义了4个栅元:栅元1由球面7围成,里面填充材料1(16O 2气体),密度是0.0014g/cm 3;栅元2由球面8围成,填充材料2(铁),

密度7.86g/cm 3;栅元3由平面1、2、3、4、5、6围成,不包括球面7、8以内的空间,填充材料3(石墨),密度1.6g/cm 3;栅元4是栅元3以外的空间,为真空。

曲面卡上需要填写曲面号、曲面类型和曲面参数,本例中定义了8个曲面,前6个为与原点距离5cm 垂直于各坐标轴的平面,后两个是半径0.5cm 的球面,球心分别在(0,-4,-2.5)和(0,4,4)。

数据块中指定了问题类型、源、记数方式、材料和运行粒子数,各卡数据项的意义如下:

MODE 卡 问题类型是中子输运

IMP 卡 4个栅元的中子重要性分别是 1 1 1 0 SDEF 卡 位于(0,-4,-2.5)、能量14MeV 的各向同性点源

F2卡

在曲面8上做中子通量记数 E0卡 对记数能量分区,1~14MeV 之间间隔为1MeV ,1MeV~10-5MeV

之间间隔为一个数量级

M1卡

材料1是16O 核素 M2卡

材料2是Fe 元素 M3卡

材料3是C 元素 NPS 卡

运行源粒子数10000

以上例子仅用于说明INP 文件格式,有关各输入卡的详细内容,读者在此不必深究。 C. 接续运行的输入文件

为了继续一个早先被终止的计算,或想对早先终止的计算进行重新编辑打印(指向文件输

蒙特卡罗方法与MCNP程序入门

出,即PRINT语句,不是向打印机输出),可以接续运行MCNP。

接续运行INP文件中必须有信息块,信息卡中必须有C或C m参数,其中m表示从第m次转储开始接续运行,m缺省表示从最后一次开始。接续运行输入文件的格式为:MESSAGE: C m

空行分隔符

CONTINUE

数据卡

.

.

空行结束符

其它

此时INP中不能有标题卡、栅元块和曲面块。数据块第一行必须是CONTINUE,允许使用的数据卡只是NPS、CTME、PRINT、FQ、DD、PRDMP、LOST和DBCN卡。NPS卡参数为累积运行的总源粒子数,CTME卡参数相对于本次运行的开始时间。如果给NPS置负参数,表示仅从转储文件中提取数据重新打印。

接续运行的另一个必要条件是上次运行的转储文件RUNTPE还没被删除。因每次运行都生成新的转储文件,若想对一次运行结果多次重新打印,则应把RUNTPE文件做备份。各次转储时运行的源粒子数和时间可从OUTPUT文件中查到。

D. 卡片格式

INP文件的每一行都限于使用1~80列。卡片都可以按行填写,数据卡也允许按列填写。

a) 行输入格式

通常卡片的1~5列用于填写栅元号、曲面号或数据卡的助记名,6~72列填写卡片参数,73~80列为注释,$符号之后也为注释。序号或卡片助记名可以写在1~5列的任何地方。带有粒

- 20 -

本文来源:https://www.bwwdw.com/article/erpq.html

Top