《MCNP4B快速入门.doc》由会员分享，可在线阅读，更多相关《MCNP4B快速入门.doc（16页珍藏版）》请在三一办公上搜索。
1、MCNP4B快速入门(MCNP4B Quick Start Guide)2002年6月前 言本手册简单介绍了 Monte Carlo方法和 MCNP4B并附以实例说明 MCNP4B的使用。本文的目的在于告诉MCNP4B初学者如何快速入门，因此很多内容不曾涉及。手册包括3部分内容：Monte Carlo方法与MCNP4B简介；MCNP4B的使用；MCNP4B应用实例。由笔者水平所限，若有纰漏望批评指正。 严慧勇2002年6月第一章 Monte Carlo与 MCNP4B简介一、Monte Carlo方法简介1-6 Monte Carlo方法是诺斯阿拉莫斯实验室在总结其二战期间工作（曼哈顿计划）的
2、基础上提出来的。Monte Carlo的发明，主要归功于 Enrico Fermi、Von Neumann和Stanislaw Ulam等。自二战以来，Monte Carlo方法由于其在解决粒子输运问题上特有的优势而得到了迅速发展，并在核物理、辐射物理、数学、电子学等方面得到了广泛的应用。Monte Carlo的基本思想就是基于随机数选择的统计抽样，这和赌博中掷色子很类似，故取名 Monte Carlo。在Monte Carlo方法中，把要研究的对象称为母体、Monte Carlo方法就是从母体中抽取一些个体进行模拟，然后利用个体的特征来对母体的特征进行推断。从母体中抽取个体必须是随机的，每个
3、个体必须是独立的，抽样越多对母体的判断也就越准确。Monte Carlo方法非常适于解决粒子输运问题，从某种意义上讲，Monte Carlo方法可以部分的替代实验。Monte Carlo方法解决问题时，不从方程或者数学公式出发，非常直观而且不受几何条件的限制。这和确定性（deterministic）方法是截然不同的。确定性方法，如常用的分离坐标法，主要是求解表征大多数粒子行为的传输方程，而Monte Carlo并不求解某特定方程，而是模拟每个粒子的行为并记录它们的平均行为的某些方面（tally），然后利用中心极限定理便可推出粒子在物理系统中的平均行为。也就是说，这两种方法解决问题的途径是不同的
4、，确定性方法是通过求解传输方程得到问题的答案，而Monte Carlo方法则是通过模拟每个粒子的传输过程解决问题。在比较Monte Carlo方法和分离坐标法的异同时，人们常说，Monte Carlo方法是求解积分传输方程，而分离坐标法是求解积分一微分传输方程。这种说法很不准确，容易引起两种误解：首先，积分传输方程和积分一微分传输方程是同一方程的不同形式，如果一个求解了，另一个也就求解了。其次，Monte Carlo是通过模拟粒子的传输过程来解决问题的，并不求解传输方程Monte Carlo并不需要任何传输方程以解决传输问题。Monte Carlo方法非常适于解决复杂的三维问题，对于不能用确定
5、性方法解决的问题尤其有用，可以用来模拟核子与物质的相互作用。在粒子输运中，Monte Carlo技术就是跟踪来自源的每个粒子，从粒子产生开始，直到其消亡（吸收或逃逸等）。在跟踪过程中，利用有关传输数据经随机抽样来决定粒子每一步的结果。可以用图1-1来阐述Monte Carlo模拟过程。在图1所示的例子中，在事件1处入射中子发生碰撞后射向另外一个方向（散射）一中于散射方向是根据物理散射概率分布随机选择的；同时，一个光子产生并暂时存储起来。在事件2处，该人射中子发生裂变并消亡，同时产生两个新的中子和一个光子。裂变产生的光子和一个中子被存储起来，现在开始对另外一个中子进行跟踪。这个中子在事件3处被俘
6、获，这时取出那个刚才存储起来的中子，通过随机抽样，该中子 图1 中子入射到平板上与物质发生相互作用在事件4处逃逸出平板，MCNP结束对该中子的模拟，重新取出裂变产生那个光子进行跟踪。这个光子在事件5处发生碰撞并在事件6处逃逸出平板。现在取出最后一个粒子在事件1处产生的光子进行跟踪，该光子在事件7处被俘获吸收。至此入射中于的整个历史也就完成了。需要注意的是，MCNP在重新取回以前存储的粒子时，是按照这样的规则进行的：最后存储的粒子最先被取出，也就是“先进后出”或“后进先出”的堆栈操作。 由于 Monte Carlo在解决粒子输运问题方面的优越性，人们开发了很多种 Monte Carlo模拟程序，
7、如 EGS。MCNP、MORSE和 FLUKA等，其中，EGS4（Electron Gamma Shower，Version 4）和MCNP4B（Monte Carlo NParticle Transport Code，Version 4B）是应用最多并得到广泛验证的两种经典Monte Carlo模拟程序。EGS4是由美国斯坦福直线加速器中心。日本高能物理国家实验室和加拿大国家研究所联合推出的一套模拟电子和光子在物质中输运过程的通用Monte Carlo程序系统，是揭示电子和光子在物质中输运规律的有力和较为方便的模拟研究工具4。关于 EGS4和其他Monte Carlo程序系统，在文献5中有详
8、细的介绍，这里不再作过多叙述。MCNP是由美国诺斯阿拉莫斯国家实验室开发的通用目的、连续能量、依赖于时间的大型Monte Carlo程序，其开发工作量相当于 450多人年。MCNP可用于模拟中子、光子、电子或者中子一光子一电子的耦合输运，也能够计算重要系统的特征值。该程序能处理任意三维结构材料中的粒子输运问题。MCNP有几种模式：中子模式、光子模式、电子模式、电子光子和中子光子耦合模式，能模拟能量为 10-11MeV到 20MeV范围内的中子和 1keV到 1000MeV范围内的光子和电子。MCNP的最新版本MCNP4D还可以解决光于中子的耦合输运问题，能方便的解决光中子问题。需要注意的是，M
9、CNP在模拟电子输运的时候，使用的是连续慢化模型，考虑了正电子、韧致辐射和k层Xrays，但是没有考虑外部和电子自身引入的场的影响。二、MCNP4B与EGS4之简单比较 在本科毕业设计中，曾经使用过一段时间的EGS4。在这里据个人的理解，对这两种程序的异同作一个简单的比较。 MCNP4B使用起来更简单 与 EGS4相比,MCNP4B的使用是非常简单的，这表现在： (1) 运行EGS4之前，需要按照有关规范，建立一个Mortran或者Fortran。程序源文件，然后编译、连接、运行。也就是说需要编写程序，这对没有学过Fortran的人来说，可能就很麻烦了。而MCNP4B的运行则是非常简单，只需要
10、建立一个输入文件（文本文件），给出有关参数就可以了。 (2) 在EGS4中，用户必须用程序语言（Fortran或Mortran）直接描述粒子所需要经过的区域和几何空间，比如一个平板、一个锥台等等。当几何体很多而且不规则的时候，难度就可想而知了，就像如果我们用语言来清楚地描述一个屏蔽系统的具体形状和结构等，就不是那么容易的事情。MCNP4B是采用纯数学的方法来描述几何空间的。任何几何体，都是由若干个面围成的，而生活中常见物体的面一般都是比较简单的球面、锥面、柱面、平面等等，用数学公式描述这些面是非常轻松的事情。MCNP4B便是通过对“面”的描述来表达“体”的，而不是象EGS4那样直接对“体”进行
11、描述。 (3) 在EGS4中，用户需要编写代码来描述粒子如何从一个区域进入到另一个区域，而MCNP4B的输入文件中则不需要考虑这个问题。 （4）EGS4中对源的描述也是比较复杂的，比如描述一个均匀分布的源，就需要用户自己编写代码产生伪随机数。而在MCNP4B中，只需要填写有关的tally就可以了。另外，MCNP4B提供了很多的tally，用于结果处理，使用起来非常方便。EGS4的使用更灵活 前面曾经提到，EGS4需要用户自己编写程序代码，这给使用带来了一定困难。但是，也正是由于EGS4需要用户自己编写代码，才使得其应用显得更加灵活，也可以扩充EGS4的功能。而MCNP4B只能使用系统提供的那些
12、功能，有一定局限性。比如我们可以通过EGS4和MORSE结合实现光中子的模拟，这是MCNP4B做不到的。 MCNP4B和EGS4的应用有各自的侧重点 EGS4可以模拟用子和光子在物质中的传输，但不能模拟中子的输运，而MCNP4B则可以模拟电子、光子以及中子的输运过程。另外，EGS4能模拟粒子在电磁场中的输运，这是MCNP4B做不到的。三、MCNP4B网络资源 MCNP用户可以通过email列表服务器mcnp-forrumlanl.gov和mcnp-llanl.gov与其他使用者进行交流讨论，也可以通过访问http:/www-xdiv.lanl.gov/xtm/mcnp得到有关信息。用户要加入
13、MCNP email列表，应先和mcnplanl.gov联系。第二章MCNP4B的使用一、MCNP4B的安装与运行 MCNP4B的安装很容易，只需要将打包好的文件解压到某目录，比如F:MCNP4B，将生成如下两个目录： F:MCNP4BDlc189（数据文件和相应文档） F:MCNP4BDOS（应用程序和文档）注意：Dlc189目录下是有关的截面数据，要想使用这些数据，需保证这些数据文件都在MCNP4B当前工作目录下，如F:MCNP4BUser，敲入mcnp Inp=myinput xsdir=xsdir1便可运行输入又件了。 MCNP4B的运行分以下几步：a）建立输入文件，如test.inp
14、。在输入文件中定义几何区域、面、介质和源，并给出输出控制tally和输运模式等。b） 查看几何形状是否描述正确。test.inp中的几何描述是否和预期的一样呢？可 以通过敲入以下命令查看： mcnp ip itest.inp otest.out 这时，DOS进入图形模式，并显示提示符“Plot”，此时回车或键入“factor x” 将图形放大l/x 倍显示。敲入命令“help”可以得到有关帮助，结束后用命令“end”退出图形模式回到DOS系统。c）运行MCNP4B。敲以下命令开始运行： mcnp i=test.inp o=test.outd）查看结果test.oute）删除过时文件（run*,
15、com*） 需要注意的是，当运行MCNP4B时，如果程序检测到已经有同名输出文件存在，则会自动更改输出文件名而不会覆盖己有的文件。 二、MCNP4B输入文件的格式 在运行MCNP4B之前，首先要根据实际要解决的问题建立一个输入文件。输入文件需包括区域（cell）定义块、面（surface）定义块、物质定义块、粒子源定义块和输出 tally块等。输入文件格式如下：Message Block OptionalBlank Line DelimiterOne Line Problem Title CardCell CardsBlank Line DelimiterSurface CardsBlank
16、Line DelimiterData Cards (material and source definitions plus output tallies)Blank Line Terminator (optional)几何区域(cell)定义 在MCNP中，不直接对几何体进行描述，而是通过对围成该几何体的面进行描述来实现几何体或几何区域的定义。我们知道，任何一个几何区域都是由若干面围成的。区域（Cell）定义格式如下： 区域号 物质号 物质密度 区域定义（若干面围成） 比如，图2所示的区域1为一个长方体，由平面16围成，该区域为密度7.9g/cm3的物质1；区域2位于球面7以内但不包括区域1
17、，该区域为真空；区域3位于球面7以外且在球面8以内，由密度为8.9g/cm3的物质2填充。那么，在MCNP4B的输入文件中，图2所示几何区域可描述如下： 1 1 -7.9 1 -2 3 -4 5 -6 $区域1 2 0 -7 #1 $区域23 2 -8.9 7 -8 $区域3图3中，球面1和球面2的公共区（交集）为区域1，其余部分为区域2，那么，区域1和区域2（真空区域）可描述为： 1 0 -l -2 2 0 (-1 2):(-2 1) 图4图4中共3个区域，区域1由两个不相连的区域组成，面l5为平面，面69为圆柱面，设这3个区域的物质密度分别为 18.3g/cm3、8.9g/cm3和 11.
18、3g/cm3，物质号分别为2、1和4，则图4所示几何空间可描述为： 1 2 -18.3 (1 -2 -6):(4 -5 -9) 2 1 -8.9 2 -3 -7 3 4 -11.3 3 -4 -8 说明：物质号由物质定义部分决定；一般需要在物质密度前加负号；表示在某面以内或以下或以左（即小于某面）的区域，则需要在该面序号前加负号，反之，若在某面的以外或以上或以右（即大于某面）的区域，则应在该面序号前不加负号；“$”后面的文字是对前面的语句的注释；“c”表示该行为注释行；面序号由面定义部分给出；表示不在某区域（cell）以内，则在该区域号前加“#”；默认密度单位为g/cm3，默认长度单位为cm，
19、默认的粒子能量单位为Mev。面（surface）的定义 面定义格式如下： 面序号 面助记符（Mnemonic）入口参数（entries） 如图4所示的面l-9等面，可描述如下：l py -1.5 $平面，序号为 1，平面方程为 y=-1.5 cm2 py3 py 24 py 2.55 py 3.56 cy 1.5 $圆柱面，圆柱轴线和y轴重合，圆柱半径为1.5cm。7 cy 28 cy 2.59 cy 1 其他常用面的定义请看表1。三物质、粒子源和输出tally的定义 物质定义格式为： mn Z1A1.XXX物质 Z1百分含量 Z2A2.XXX物质 Z2百分含量 m为助记符，n为物质序号（与区
20、域定义块中的物质序号对应）Z为原子序数，A为质量数，当物质为单质时，A可取“000”，A为三位数；XXX表示反应截面数据库的序号与类型，有时候可以省略，在定义物质时该项可根据文献3第675715页决定。下面给出一个物质定义块的例子； m1 74000 1 $钨m2 79197.60c 1 $金m3 1001.50c 0.667 8016.60c 0.333 $水，物质序号为3 水中氢原子占2/3，氧占1/3。m4 29000 1 $铜源的定义相对比较复杂，可阅读文献3有关章节（第267278页）。下面仅简单说明一下源的定义方法：1）定义模式（Mode card） MCNP可以几种不同方式运行：
21、 Mode N neutron transport only（default） N P neutron and neutron-induced photon transport P photon transport only E electron transport only P E photon and electron transport N P E neutron, neutroninduced photon and electron transport 如果不定义Mode，则系统认为要模拟的粒子为中子。2）定义区域粒子重要性 格式如下： IMP：type keyl key2 key3.
22、keyn 从左到右，依次为助记符、粒子类型、粒子在区域1至区域n这n个区域中的重要性。如果key为0。则表示粒子一旦进入该区域便会被抛弃。举例说明如下： imp：e 1 1 1 0 $仅模拟电子在区域 1，区域 2和区域 3中的传输 imp：n 1 1 1 0 $模拟中子在区域1-3中的传输，中子进入区域4便结束对其模拟3）源的说明 定义源的特点，包括源的位置、能量、能谱、方向、类型、分布等等。下面举例说明简单的粒子源定义方法： mode e p sdef sur2 erg 15 pos 0 -5.5 0 rad d1 dir 1 vec 0 1 0 par e sil 1.5 以上定义了一个
23、电子源，入射电子在以（0，-5.5,0）为圆心、半径为 1.5cm的圆内均匀分布（该圆位于面2上），出射方向平行y轴。Sur（cell）给出粒子出发位置所在的面或区域号；PoS给出出发点坐标，默认为原点；erg给出初始能量或能谱，默认为 14MeV；Par给出源粒子类型（N for N,N P, N P E；P for P，P E；E for E）。si、sp与rad等配合，可给出源的空间分布；dir给出粒子初始运动的方向或方向分布(相对与参考矢量的夹角的余弦值）；vec给出方向参考矢量以原点与该点的连线方向作为参考方向。关于源的详细定义方法可阅读文献3的267278页。4）输出tally说明
24、 表2给出了常用tally。参考该表，如果我们要定义几个探测器，记录（0,1,1)和（10，23，20）和（-10，20，-30）三点的粒子通量（平均一个粒子入射时，在该点单位面积上通过的粒子数目）以及面2、面3和面4上的粒子通量（平均一个粒子人射时，在该面上的平均通量），则输出tally可以写作： F15:p 0 1 1 0 $定义点探测器1，记录该点上的光子通量 F25:n 10 23 20 0 $定义点探测器2 F35:e -10 20 -30 0 $定义点探测器3，探测器序号不可重复 F12:P 2 F22:n 3 $面探测2，记录面3上的平均中子通量 F32:e 4三、一个实例下面以
25、一个实例来说明MCNP4B的输入文件格式及MCNP4B的运行。Shielding System Designing for 15MeV CT，to shield Xrays $Title应从第一列开始1 2 -18.3 （1 -2 13 -16）:（2 -8 12 -14):（8 -9 -14） $W 定义区域。依次为区域号、物质号、密度、区域定义（由平面围成）2 1 -1 （1 -2 16 -17）：（2 -9 14 -15）：（9 -10 -15）$C2H43 2 -19.3 3 -4 -11 $Wu，定义复合靶区域4 3 -19.66 4 -5 -11 $AU5 4 -1.0 5 -6
26、-11 $H2O6 5 -8.9 6 -7 -11 $Cu7 6 -1.54 7 -8 -11 $Graphite8 0 -18 #l #2#3 #4 #5 #5 #6 #7 $表示区域8在面18以内，且不含1-7区域9 0 18 （区域定义完后，空一行，开始定义面。事实上，我们都是先定义面后定义区域的。只是将其放在区域后面而已）1 py -16 5 $y-16.5 cm,垂直于y的平面。面序号为 12 py -5.53 py -1.474 py -1.355 py -1.256 py -1.17 py -1.08 py 09 py 1110 py 3111 cy 0.5 $圆柱面,轴线和y轴
27、重合,半径为0.5 Cm12 cy 1.413 cy 6.014 cy 5.915 cy 20.916 cy 10.517 cy 25.5 $圆柱面，轴线和y轴重合，半径为 25.5 cm18 so 150 $球面，球心位于坐标原点，半径为150 cm 定义完面后 定义源，之间空一行）mode e p $模式,可以是e p n；e p；p n等imp:p 1 1 1 1 1 1 1 1 0 $在区域中的importance0，模拟；=0抛弃imp:e 1 1 1 1 1 1 1 1 0c sdef erg D1 pos 0 0 0 $注释,（默认）中子源定义, 位置，能谱, 各向同性c spl
28、 -3 1.025 2.926 $f=-3 Watt fission energy spectrum：p (E) =Cexp (-E/a) sinh (bE) (1/2)sdef sur=2 erg 15 pos 0 -5.5 0 rad dl dir 1 vec 0 1 0 ara 0.001767 par 3 $源位于平面0 -5.5 0 为圆心的半径 0.075，均匀分布，方向余弦为1，参考方向为（0,0,0）（0,1,0）sil 0.075 $半径为 0.075m1 6000.60c 0.317 1001.50c 0.633 5010.60c 0.05 $物质定义，含硼聚乙烯。物质序号
29、和前质序号要一致m2 74000 1 $钨m3 79197.60c 1 $金m4 1001.50c 0.667 8016.60c 0.333 $水m5 29000 1 $铜m6 6000 1 $碳cf5:P 0 100 0 0 $0 $输出结果控制，输出某点光子的通量率。Tally F5输出某点的剂量率。点由4个数组成，前3个数表示点的坐标，最后一个数0定义了一个点。fl5:p 0 99.62 8.716 0 $5f25:p 0 98.48 17.633 0 $10c f35:n 0 96.59 25.882 0 $15f45:p 0 93.97 34.202 0 $20c f55:n 0 9
30、0.63 42.262 0 $25f65:p 0 86.60 50.000 0 $30c f75:p 0 81.92 57.358 0 $35f85:p 0 76.60 64.278 0 $40c f95:p 0 70.71 70.711 0 $45f105:p 0 64.28 76.604 0 $50f115:p 0 50 86.603 0 $60f125:p 0 34.20 93.969 0 $70f135:p 0 17.63 98.480 0 $80f145:p 0 0 100.00 0 $90f155:p 0 -17.63 98.480 0 $100f165:p 0 -34.20 93
31、.969 0 $110f175:p 0 -50 86.603 0 $120f185:p 0 -64.28 76.604 0 $130f195:p 0 -76.60 64.278 0 $140f205:p 0 -86.60 50 0 $150f2l5:P 0 -93.97 34.202 0 $160f225:P 0 98.48 17.633 0 $170f235:P 0 100 0 0 $180f2: p 18f12:e 18 $面上的剂量通量de0 0.01 0.03 0.05 0.07 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0.55 0.6 0.65
32、 0.7 0.8 1.0 1.4 1.8 2.2 2.6 2.8 3.25 3.75 4.25 4.75 5.0 5.25 5.75 6.25 6.75 7.5 9.0 11 13 15 $能量插值点。log或者线性插值，de7与df0配合使用，df0给出对应的函数值df0 1.1E-08 1.617E-09 8.056E-10 7.166E-10 7.861E-l0 1.053E-09 1.392E-09 1.753E-09 2.108E-09 2.438E-09 2.736E-09 3.000E-09 3.250E-09 3.528E-09 3.778E-09 4.000E-09 4.22
33、2E-09 4.666E-09 5.500E-09 6.972E-09 8.306E-09 9.500E-09 1.061E-08 1.114E-08 1.225E-08 1.342E-08 1.453E-08 1.556E-08 1.611E-08 1.669E-08 1.769E-08 1.872E-08 1.975E-08 2.128E-08 2.436E-08 2.861E-08 3.278E-08 3.694E-08 $log interolation,dose rate is abtained.mSv/s, 给出的函数 值为photon fluxdose rate转换系数E0 0
34、141 15 $能量间隔。结果输出时。统计数据采用的能量间隔cnps 20000 $抽样数目cut:e j 0.611 $电子截至能量 输入文件建立好以后，存为S15.nip。需要注意的时，文件名长度不要超过8个字符。接下来查看一下该输入文件描述的几何空问和预期的是否一样，键如下命令： mcnp ip i=S15.inp oS15.out系统进入图形模式并显示提示符“plot”，此时敲“factor 0.5”将图形放大2倍显示，见附图。 图形显示无误，删除S15.out并运行mcnp4B： mcnp i=S15.inp oS15.out当程序运行完毕后，阅读S15.out，找出需要的信息。结果
35、文件开头部分为文件的内容，后面才是结果以及误差分析。最后删除文件run*和com*。四、Plot功能的使用 MCNP提供了较好的图形功能（参考文献3），为调试带来了很大的方便，这里仅简单说明一下 ploter的使用。敲入命令 mcnp ip i=输入文件o=输出文件，进入图形模式，显示图形提示符“plot”，此时敲如“help”可得到相关帮助而敲“end”命令则退出图形模式。下面列出另三个常用命令。 Origin x0 y0 z0以（x0，y0，z0）为原点显示图形，即将图形平移 Pz（py or px）z1（y1 or x1）沿平面zz1(yy1 or xx1）剖开几何体 Factor m将图形放大1/m倍显示 需要注意的是，MCNP中，水平方向为y轴。
