Numerical Modeling of Gamma Radiation from Galaxy Clusters

SAGCM-APD增益影响因子分析与优化胡大鹏，郭方敏，朱晟伟，熊大元【摘要】应用二维漂移扩散模型研究具有分立吸收层、渐变层、电荷层和倍增层结构（SAGCM）的InGaAsP-InP雪崩光电探测器（APD），仿真分析了不同电荷层、倍增层厚度和掺杂浓度对电场分布、电流响应及击穿电压的影响，特别是参数变量对增益计算模型的影响，载流子传输过程的时间依赖关系和倍增层中所处位置的影响。

仿真结果表明：较高掺杂浓度和较薄电荷层结构可以改变器件内部的电场分布，进而提高增益值。

当入射光波长为1.55 μm，光功率为500 W/m2时，光电流响应量级在10-2A；阈值电压降低到10 V以下，击穿电压为42.6 V时，器件倍增增益值大于100。

【期刊名称】红外与激光工程【年(卷),期】2010(039)002【总页数】6【关键词】SAGCM-APD；雪崩光电探测器；倍增增益；电流响应；仿真分析0 引言雪崩光电探测器（APD）主要应用于高比特的光波传输系统中，相比传统的PIN光电探测器，这种具有内部增益的APD器件是长距离光纤通信系统中的重要器件[1]。

近年来，具有分立吸收层、渐变层、电荷层和倍增层结构（SAGCM）的APD由于其结构优势而被广泛研究[2-8]。

图1是典型的SAGCM-APD结构，把倍增层和吸收层分离是为了得到较低的偏压和低噪声，因此，吸收作用只发生在低带隙材料层中，而应用宽带隙材料作为倍增层则有利于发生倍增效应[3]。

加入高掺杂浓度的电荷层可以使倍增层得到较高增益带宽积所需要的高电场分布[4]。

此外，在吸收层和倍增层之间加入渐变层结构可以减小空穴在异质结界面传输时发生的捕获效应，增加进入倍增层中的载流子数，提高电流响应。

虽然很多文献对SAGCM-APD进行理论研究，但是很少讨论电荷层和倍增层的厚度、掺杂浓度对电流响应及击穿电压的影响。

文中应用二维漂移扩散模型讨论了不同电荷层、倍增层厚度和掺杂浓度对SAGCM-APD电流响应及击穿电压的影响，研究了器件结构参量改变对倍增增益的影响。

翻译技巧第九节被动语态的译法英汉对比研究表明，英语中被动语态的使用极为普遍，尤其是以科技英语为主。

这是英语区别于汉语的显著特点之一。

在英语中，大凡为了强调受事，以突出其鲜明位置，而无需提及主动者、无意点明主动者、无从说出主动者，为了上下文的衔接与连贯等等，或是出于礼貌措辞圆通等方面的考虑不愿说出动作的执行者是谁的，往往都采用被动语态。

例如：(1)When will we be interviewed?我们什么时候来参加面试呢?(2)New British and American films are often shown to English major juniors.常给英语专业的三年级的同学放映新的英美电影。

(3)The happy man cannot be harried.吉人自有天相。

(4)You are cordially invited to a dinner party to be given at the Workers’ Club in Xidan at 7 P. m.Nov. 23.谨定于11月23日晚上7时在西单工人俱乐部举行晚宴，敬请光临。

(5)The plan was especially supported by those who wished to have more chance to speak French. 这个计划特别受到愿多有机会说法语的人的支持。

在英语中，被动结构能给人以间接、客观的印象，因而在科技文献、政论文和新闻报道等文体中使用尤为频繁。

相比之下，被动语句在汉语中的使用范围小得多。

主要原因在于，在历史上，汉语被动语句主要用来表达对主语来说是如意或不希望发生的事情。

现代汉语的这种限制虽然在西方语言的影响下有所动，被动句式有时也用来表达中性甚至希冀之事，但是汉语采用的主要是语义型句式，大多数被动意义并不一定非得通过被动语句来表达，而可以通过形式上主动、语义上被动的句式予以体现。

第50 卷第 9 期2023年9 月Vol.50，No.9Sept. 2023湖南大学学报（自然科学版）Journal of Hunan University（Natural Sciences）遍布节理模型在层状岩体模拟中的适用性研究黄娟1，周世杰1，贾朝军1†，宋银涛1，张建2（1.中南大学土木工程学院，湖南长沙 410075；2.中铁五局集团有限公司，湖南衡阳 420002）摘要：为弥补遍布节理模型未考虑节理长度、间距及节理刚度的不足，利用三轴压缩数值试验和参数校准准则，对有限差分软件FLAC3D中的遍布节理模型进行参数校准. 通过圆形洞室开挖的算例，对比分析了遍布节理模型与3DEC块体离散元模型的计算结果在位移、塑性区以及最大主应力上的差异. 依托具有典型层状围岩的新华山隧道工程，采用校准的遍布节理模型和离散元方法分析隧道开挖和初期支护后的力学响应. 最后探讨了层理角度对围岩变形和塑性区的影响，进一步验证校准后的遍布节理模型在工程中的适用性. 研究表明，经过校准的遍布节理模型能够较好地描述层状岩体的各向异性行为，可应用于类似工程之中.关键词：岩石力学；各向异性；遍布节理模型；隧道开挖中图分类号：TU45 文献标志码：AApplicability of Ubiquitous-Joint Model in Layered Rocks Simulation HUANG Juan1，ZHOU Shijie1，JIA Chaojun1†，SONG Yintao1，ZHANG Jian2（1.School of Civil Engineering， Central South University， Changsha 410075， China；2.China Railway No.5 Engineering Group Co.， Ltd.， Hengyang 420002， China）Abstract：To address the limitations of the Ubiquitous-Joint model which does not consider the effects of joint length， joint spacing， and joint stiffness， the parameters of the Ubiquitous-Joint model in FLAC3D were calibrated using triaxial compression numerical tests and parameter calibration criteria. The distinctions of modeling results between the Ubiquitous-Joint model and the 3DEC model including the deformation，the plasticity zone，and the maximum principle stress were compared and analyzed by an example of circular tunnel excavation. Based on the Xinhua Mountain Tunnel project with typical layered rock mass，the deformation and failure mode were analyzed with the calibrated Subiquitous model and discrete element method after the tunnel was excavated and primary support finished. Finally，the deformation and plastic zone of surrounding rock influenced by bedding angle was discussed，which further verified the applicability of the calibrated Subiquitous model in engineering. The results confirm that the calibrated Subiquitous model is capable to well describe the anisotropic behavior of layered rock，which can be applied to similar engineering projects.Key words：roke mechanics；anisotropy；Ubiquitous-Joint model；tunnel excavation∗收稿日期：2022-08-02基金项目：国家自然科学基金资助项目（U1934211），National Natural Science Foundation of China（U1934211）作者简介：黄娟（1977—），女，湖北荆州人，中南大学副教授，博士† 通信联系人，E-mail：******************.cn文章编号：1674-2974（2023）09-0131-11DOI：10.16339/ki.hdxbzkb.2023109湖南大学学报（自然科学版）2023 年层状岩体是岩土与地下工程建设中经常遇到的一类岩体，在自然界中广泛分布. 长期地质构造作用下所形成的层理面使岩体在强度和变形等方面都表现出明显的各向异性，这对隧道开挖时岩体锚固［1］、衬砌开裂、仰拱隆起［2］等工程问题有着显著的影响.因此，层状岩体的力学行为与响应机制研究具有重要意义与研究价值.近年来，随着材料本构不断完善以及计算机技术的更新迭代，越来越多的数值模拟技术用于岩体力学特性的研究，为室内试验或现场测试的局限性提供了补充和解决办法. 王培涛等［3］应用颗粒流软件PFC2D研究了不同层理角度的黑云变粒岩的强度特性. Singh等［4］通过UDEC探究了节理岩体在单轴加载条件下产生高侧向应变的原因. 刘爱华等［5］采用有限元软件ANSYS模拟了不同层面倾角的岩体抗拉、抗压试验下的破坏形态. 此外一部分学者还将有限元法［6-7］、有限差分法［8］、离散元法［9］、有限-离散元法［10］、真实破裂过程分析方法［11］等数值方法应用于模拟层状围岩地下洞室的变形和破坏机理等方面.虽然数值模拟方法繁多，但相比之下，采用离散元法能够最有效地描述层状岩体等不连续材料的力学性能［12］. 然而，考虑到离散元法计算的效率，若要模拟全部的节理或层理构造以进行某些大型地下工程的开挖掘进是不太可取的［13］. 近年来，有学者研发了高效颗粒离散元软件MatDEM［14］，但颗粒离散元软件很大程度上依赖于本构参数的准确标定，且该软件暂未广泛应用于层状岩体模拟之中. 为了避免这些限制，通常可以采用FLAC3D中的遍布节理模型来表示一些层状各向异性岩体. 例如，蒋青青等［15］采用FLAC3D内置的Ubiquitous-Joint模型分析了层状岩质边坡开挖过程中层理倾角和倾向与安全系数之间的关系. 朱泽奇等［16］、周鹏发等［17］采用改进的Ubiquitous-Joint模型模拟了层状围岩地下洞室开挖时的变形和破坏. Sainsbury等［18］针对岩体中普遍存在的随机节理，通过建立与主节理或层理方向正交的节理集，并提出遍布节理模型参数修正准则，较好地描述了自然界中各向异性岩体强度和变形特性.然而，目前遍布节理模型中参数的物理意义不够明确，不能由试验结果直接获取. 由于遍布节理模型没有考虑节理间距和节理刚度，如果直接将其材料参数与离散元模型的岩块和结构面参数赋值一致，模拟结果不能真实地反映实际工程中的岩体强度或变形. 因此，需要对遍布节理模型的参数进行校准修正，以便为工程设计或施工提供有意义的参考.本文通过总结部分学者对层状各向异性岩体的研究，在Sainsbury研究的基础上，分别采用FLAC3D 中的Ubiquitous-Joint模型和Softening-ubiquitous模型（考虑应变软化的Ubiquitous-Joint模型，以下简称Subiquitous模型）以及块体离散元软件3DEC对层状岩体的力学特征进行模拟并作对比分析. 基于校准后的Subiquitous模型，通过分析新华山隧道开挖和支护过程，揭示层理对围岩变形和破坏特征的影响，验证遍布节理模型的适用性.1 （应变软化）遍布节理模型在FLAC3D中，遍布节理模型有Ubiquitous-Joint 和Subiquitous模型两种. Ubiquitous-Joint模型［19］对应于摩尔-库仑模型，即在摩尔-库仑体中加入节理面，该节理面也服从摩尔-库仑屈服准则，使材料表现出强度各向异性. Ubiquitous-Joint模型同时考虑了岩石基质和节理的物理力学属性，必须在模型的指定区域内同时赋予基质和节理的参数.节理面的破坏包括拉伸和剪切破坏，如图1所示，其中剪切破坏包络线AB表示为f s=0：f s=τ+σ3′3′tanϕj-c j=0.（1）拉伸破坏包络线BC表示为f t=0：f t=σ3′3′-σt j=0.（2）式中：ϕj、c j和σt j分别为节理面的内摩擦角、黏聚力和抗拉强度；σ3′3′为节理面上的正应力.该模型的计算公式中未涉及节理的间距、长度以及岩层的弯曲刚度等. 如果不对相应的力学参数进行校准，可能会得到错误的岩体强度和变形响应. Subiquitous模型［19］是广义的Ubiquitous-Joint模型，该模型中基质和节理强度符合双线性摩尔-库仑准则，且允许材料基质和节理的强度发生硬化或软图1 节理面破坏准则Fig.1 Joint failure criterion132第 9 期黄娟等：遍布节理模型在层状岩体模拟中的适用性研究化. Subiquitous 模型和Ubiquitous-Joint 模型都是先根据摩尔-库仑准则检测基质的屈服，并进行相应的塑性修正，然后分析在新的应力状态下节理面上的破坏，在材料未达到极限强度前力学行为一致.在遍布节理模型中，弱面对岩体强度的影响通常与Jaeger 提出的单弱面理论［20］进行比较. 单弱面理论指出，当1-tan ϕtan β>0时，若满足式（3）则会发生结构面的剪切破坏.σ1≥σ3+2()c +σ3tan ϕ(1-tan ϕtan β)tan β.（3）式中：c 、ϕ分别为结构面的黏聚力和内摩擦角；β为结构面的倾角. 当1-tan ϕtan β<0时，岩体不会沿结构面破坏，只会发生基质的破坏. 故该理论只允许出现沿结构面的剪切滑移破坏和基质的破坏两种破坏模式.图2为Ubiquitous-Joint 模型［19］和Jaeger 单弱面理论的岩体承载强度与结构面倾角的关系的对比，可以看出两者紧密匹配.图2中ϕw 为结构的内摩擦角，其中曲线为带有“肩部”的“U ”形曲线. 当β<ϕ或β=90°时，岩体强度与弱面无关.图3为部分已有的层状岩体三轴压缩试验研究［21-24］，由图3可知，岩体的强度随着层面倾角连续变化，这一特征也得到了许多研究人员的验证. 而单弱面理论不能充分描述自然存在的层状岩体的各向异性. 遍布节理模型也存在同样的局限性，故需要进一步探讨其适用性.2 遍布节理模型与离散元模型的对比为了探讨遍布节理模型对层状岩体模拟的有效性，针对已有的层状页岩三轴压缩试验结果，采用FLAC3D 建立与试样同等规模的数值模型，用其内置的Ubiquitous-Joint 模型和Subiquitous 模型进行分析计算，并与3DEC 的模拟结果作比较.2.1 三轴压缩试验模拟中的比较2.1.1 块体离散元方法和Ubiquitous-Joint 模型为了研究层状岩体的强度和变形特性，参考页岩［22］的三轴压缩试验数据（如图4所示），使用3DEC 建立了直径50 mm 、高100 mm 的标准圆柱体模型. 层理倾角分别设置为0°、15°、30°、45°、60°、75°和图2 Ubiquitous-Joint 模型三轴抗压强度值与Jaeger 解析解的比较Fig.2 Comparison of triaxial compressive strengthvalues-Ubiquitous-Joint model versus analytical solution （a ）层状砂岩［23］（b ）层状页岩［22， 24］（c ）层状片岩［21］图3 层状岩体三轴压缩强度随倾角变化特性Fig.3 Variation of triaxial compressive strengthof layered rock mass with bedding angle133湖南大学学报（自然科学版）2023 年90°，层厚5 mm ，岩体参数标定结果见表1. 同时基于Ubiquitous-Joint 本构模型建立了类似的FLAC3D 模型，将表1中的岩体参数直接用作模型中岩石基质和节理的参数输入，3DEC 模型和Ubiquitous-Joint 模型的强度响应如图4所示.正如预期，离散元模型随着β角的增大而遵循连续变化的强度曲线. 其与室内试验不同倾角下的峰值强度相对误差小于8%，结果基本吻合. Ubiquitous-Joint 模型在β角小于15°时其强度不受节理的影响，与室内试验结果相差超过20%，这种“U ”形强度曲线上的肩部清楚地表明了模型的局限性.提取较为典型的层理倾角为60°时岩石破坏模式的试验结果与模拟结果，如图5所示.可知此时岩石表现为沿层理面的滑移破坏，其中从离散元模型结果可以看到层理面的错动，与试验结果一致. 而Ubiquitous-Joint 模型显示大量的节理剪切破坏，但无法得知具体的破裂面位置和破裂形态.图6表示了不同倾角下离散元模型和Ubiquitous-Joint 模型的弹性模量和应力-应变曲线.由图6可知，Ubiquitous-Joint 模型没有体现出峰后的应变软化行为. 当直接在模型中采用3DEC 岩石块体的刚度参数时，所得到的弹性模量明显高于3DEC 的模拟结果，这是Ubiquitous-Joint 模型未考虑节理刚度和节理间距导致的，在实际工程中要特别注意这一点.（a ）弹性模量变化曲线（b ）应力-应变曲线（0°~45°）（c ）应力-应变曲线（60°~90°）图6 两种模型的弹性模量和应力-应变曲线Fig.6 Elastic modulus and stress-strain curve of discontinuumand Ubiquitous-Joint model图4 离散元和Ubiquitous-Joint 模型强度各向异性曲线Fig.4 Anisotropic strength curves of discontinuum andUbiquitous-Joint model表1 岩石和节理力学参数Tab.1 Rock and joint mechanical properties层间结构面黏聚力/MPa22（°）19（GPa•m -1）20（GPa•m -1）10 （a ）试验结果 （b ）3DEC 模型 （c ）UB-Joint 模型图5 层理倾角60°时岩石破坏模式Fig.5 Failure mode of rock with bedding angle of 60°134第 9 期黄娟等：遍布节理模型在层状岩体模拟中的适用性研究因此，建议不要将3DEC 中的岩石块体和节理参数直接用作Ubiquitous-Joint 模型的参数，为使其产生有意义的结果，需要对岩石基质和节理参数进行校准，以匹配离散元模型的结果. 以下将对此进行探讨.2.1.2 考虑应变软化的Subiquitous 模型参数校准与Ubiquitous-Joint 模型相比，Subiquitous 模型在校准岩石基质和节理参数方面提供了更大的灵活性. 通过双线性软化关系，可以更好地表示层状岩体的强度和变形特性. 其参数校准准则如下［18］：1）将离散元模拟结果视为实际层状岩体的各向异性行为.2）节理黏聚力和内摩擦角的初始值不变，岩体达到峰值后，节理黏聚力与岩体基质黏聚力以相同的速率软化至0.3）校准岩石基质的强度和变形响应，以补偿节理刚度和节理间距参数的缺失.β在0°和90°的情况下，试样的峰值强度取决于岩石基质的黏聚力和内摩擦角，这些参数对应于β=0°时的离散元模型的强度响应进行校准. 岩体基质和节理黏聚力的软化速率参考离散元模型的结果.在整个校准过程中，强度和刚度参数以及试样的破坏过程都得以考虑. 比较离散元模型和Subiqui‑tous 模型的破坏模式，将其分为劈裂张拉破坏（β为90°时）、剪切滑移破坏（β为60°时）和复合破坏（β为30°时）. 通过监测加载过程中基质的屈服和节理的滑移剪切，可以揭示试样的破坏机制.前文中的三轴压缩数值试验已用Subiquitous 模型重建，采用经过校准的参数，具体取值见表2.离散元模型和Subiquitous 模型的各向异性“U ”形曲线如图7所示，并与开始的Ubiquitous-Joint 模型的结果进行了比较. 经过校准后的Subiquitous 模型随着β角的增大同样遵循连续变化的强度曲线，与离散元模型的结果更加贴切.图8显示了Subiqui‑表2 校准的Subiquitous 模型力学参数Tab.2 Calibrated mechanical properties of Subiquitousmodel弹性模量/GPa 32泊松比0.25黏聚力/MPa 47.5内摩擦角/（°）29节理黏聚力/MPa 22节理内摩擦角/（°）19图7 离散元和Subiquitous 模型强度各向异性曲线Fig.7 Anisotropic strength curves of discontinuumand Subiquitous model（a ）弹性模量变化曲线（b ）应力-应变曲线（0°~45°）（c ）应力-应变曲线（60°~90°）图8 弹性模量变化曲线及不同角度下的应力-应变曲线Fig.8 Elastic modulus and stress-strain response ofdiscontinuum and Subiquitous model135湖南大学学报（自然科学版）2023 年tous模型在不同层理倾角下的应力-应变曲线和弹性模量的变化，都与离散元模型更紧密地匹配.2.2 二维圆形隧道开挖分析为了验证2.1节中开发的校准后的Subiquitous模型在工程中的应用效果，建立了一个圆形隧道模型，研究隧道开挖后的力学响应，该模型是在不考虑重力加速度的各向同性应力场中模拟的. 为了比较模拟效果，建立了岩层厚度0.5 m的3DEC模型和等效的Ubiquitous-Joint模型. 模型参数取值见表3和表4. 图9比较了隧道开挖完成时每个模型的塑性区、位移和最大主应力.3DEC模型中显示隧道侧壁中有少量岩石块体的拉伸破坏，层理的剪切滑移破坏主要在洞顶和底部沿垂直层理方向延伸约2 m. 位移场分布明显受到了层理的影响，最大位移约为20 mm. 图10显示的是3DEC模型放大20倍的变形状况，在隧道顶部和底表3 模型材料参数Tab.3 Details of model material层间岩体层间结构面弹性模量/GPa2.69黏聚力/MPa0.08泊松比0.28内摩擦角/（°）20黏聚力/MPa1.64法向刚度/（GPa•m-1）3内摩擦角/（°）45.5剪切刚度/（GPa•m-1）1.15表4 校准的Subiquitous模型力学参数Tab.4 Calibrated mechanical propertiesof Subiquitous model弹性模量/GPa1.8泊松比0.34黏聚力/MPa1.5内摩擦角/（°）47节理黏聚力/MPa0.08节理内摩擦角/（°）20图9 5 m直径的圆形隧道开挖后的塑性区、位移和最大主应力Fig.9 Plastic zones， displacement and major principal stress around 5 m-diameter tunnel136第 9 期黄娟等：遍布节理模型在层状岩体模拟中的适用性研究部可以清楚地看到岩层的弯曲.Ubiquitous-Joint 模型（直接对岩石基质和节理赋予和3DEC 中块体和节理相同的参数）中没有显示出基质的屈服破坏，而主要为节理的滑移和拉伸破坏，在隧道顶部和底部延伸约4 m ，模型的最大位移明显小于3DEC 模型. 节理拉伸破坏导致区域周围出现显著的应力重分布，其破坏机制是因为遍布节理模型的公式中没有考虑岩层的弯曲刚度. 经过校准的Subiquitous 模型的隧道侧壁上也有少量的基质拉伸破坏，使得隧道周围出现更具有代表性的应力重分布，其位移场也更接近3DEC 模型.3 工程验证为了更好地研究Subiquitous 模型在实际工程中的使用性能，以新华山隧道为例，探讨隧道开挖以及在支护结构作用下围岩的变形和破坏特性，并通过现场实测数据验证模型的可靠性.3.1 工程概况和工程地质新华山隧道位于湖南省张家界市和湖南省湘西州永顺县境内. 该隧道为单洞双线隧道，起止里程为DK26+104.00-DK32+034.49，全长5 930.49 m ，最大埋深约为383 m ，开挖高度和宽度分别为12.64 m 和14.96 m.新华山隧道所处地貌为剥蚀低山地貌，地势较起伏，山坡自然坡度一般为30°～70°. 隧道穿越地层受区域构造影响严重，节理裂隙发育、岩体破碎. 本文以新华山隧道进口段DK26+490断面附近为研究对象. 根据前期地质勘查资料，新华山隧道围岩主要为层状特征较明显的炭质页岩，由于其所具有的各向异性和开挖后风化较快等特殊工程特性，使得隧道的开挖引起软弱围岩向洞内发生不均匀对称的变形.3.2 模型建立根据纵断面图可以发现，所模拟区段的埋深约110 m ，运用FLAC3D 建立如图11所示模型.为降低模型中的单元数量，仅在模型中创建部分上覆岩体，并通过在地层上表面施加荷载模拟其余上覆岩体的自重应力. 设定模型x 、y 、z 三个方向上的尺寸分别为100 m 、50 m 和100 m ，采用位移边界条件，除上表面外，其余5个边界面约束法向位移. 模型中，岩体层理倾角采用现场调查得到的层理倾角，即为75°. 隧道采用三台阶法开挖，模拟区段并未施做二次衬砌，故模型中支护体系仅包括锚杆和初期支护，相关力学计算参数根据支护设计方案确定（见表5）. 采用3DEC 建立同等规模的模型，根据现场测试以及《铁路隧道设计规范》（TB 10003―2016）取得如表6所示参数. 其中节理刚度参数参考文献［25］，并执行2.1节的校准程序取得Subiquitous 模型的参数，如表7所示.图10 3DEC 中显示的岩层弯曲变形（放大20倍）Fig.10 Bending deformation of bedding rock sown in 3DEC（magnified 20 times）图11 数值模型及细部构造（单位：m ）Fig.11 Numerical model and detailed construction （unit ： m ）表5 支护结构计算参数Tab.5 Parameters for the support system锚杆初衬截面积/mm 2153厚度/cm 28弹性模量/GPa 200密度/（kg•m -3）2 400砂浆刚度/MPa 50弹性模量/GPa 28砂浆黏聚力/kPa 400泊松比0.2砂浆摩擦角/（°）60137湖南大学学报（自然科学版）2023 年3.3 数值模拟结果与分析根据上述参数和模型，计算得到隧道中部横截面处（Y=25 m）开挖并施加初期支护后的围岩变形和块体塑性区以及节理塑性区情况如图12所示，从中可以看出：1）两种模拟方法的围岩变形和塑性区响应非常接近，说明经过校准后的Subiquitous模型能够较好地体现层状岩体的力学特性.2）隧道开挖完成后实测拱顶沉降和水平收敛分别为259.9 mm和173.5 mm，而3DEC中拱顶沉降和水平收敛分别为243.9 mm和160.2 mm，与实测值的差异分别为-5.4%和-7.6%，FLAC3D中分别为282.5 mm和189.2 mm，与实测值的差异分别为8.7%和9.0%，差异性较小，表明建立的模型能够较好地反映新华山隧道开挖及初期支护后的变形情况. 3）受层理的影响，拱顶和拱底都朝着层理倾角方向产生位移梯度，围岩位移场呈现出显著的不对称性，这也与现场观察到的非对称变形情况相符合. FLAC3D中边墙附近围岩位移比3DEC稍大，是因为Subiquitous模型无法表示完整岩层的屈曲变形，而岩层的厚度对围岩位移有重要影响.4）隧道开挖扰动作用下，围岩产生了节理剪切破坏、节理张拉破坏、岩石基质剪切破坏和岩石基质张拉破坏4种类型的破坏，主要处于节理和基质的剪切破坏状态，且大部分围岩体同时出现了多种破坏模式. 围岩塑性区分布也显示为极不对称性，围岩深部的塑性区主要集中在左拱脚和右拱肩. 3DEC中少量的节理张拉破坏沿洞周分布，FLAC3D中节理张拉破坏更少，这也与Subiquitous模型无法解释岩层间距有关. 基质的张拉破坏只出现在隧道底部，拱顶的塑性区范围都很小.图13给出了3DEC和FLAC3D模型（与实际掘进过程一致）Y=25 m断面处的拱顶沉降监测曲线与现场监测数据的比较，可以发现，3DEC与实测数据更为接近，而FLAC3D中采用Subiquitous模型的计算结果也能较好地吻合.综合以上分析，校准的Subiq‑uitous模型在工程中有较好的实用性，能为相应工程表6 岩石和层理面力学参数Tab.6 Rock and bedding plane mechanical properties层间岩体层理面弹性模量/MPa250黏聚力/kPa60泊松比0.37内摩擦角/（°）15黏聚力/kPa150法向刚度/（GPa•m-1）1内摩擦角/（°）23剪切刚度/（GPa•m-1）0.5表 7 校准的Subiquitous模型力学参数Tab.7 Calibrated mechanical propertiesof Subiquitous model弹性模量/MPa190泊松比0.40黏聚力/kPa140内摩擦角/（°）25节理黏聚力/kPa60节理内摩擦角/（°）15图12 离散元和校准的Subiquitous模型位移和塑性区对比Fig.12 Comparison of displacement and plastic zones of discontinuum and calibrated Subiquitous model138第 9 期黄娟等：遍布节理模型在层状岩体模拟中的适用性研究提供参考，且层理的存在对围岩的变形和破坏有重要影响. 此外，就计算效率而言，两者计算时长相差30~40倍.3.4 层理倾角对隧道开挖的影响当隧道施工穿越炭质页岩地层时，围岩和支护结构的变形很可能因围岩层理倾角的变化产生显性差异. 为分析层理倾角对围岩和支护结构变形的影响，用FLAC3D 依次建立层理倾角为0°、30°、45°、60°和90°等5种工况的仿真模型. 采用校准的Subiquitous 模型，除倾角外其余参数保持不变. 计算得到岩体围岩和支护结构的变形以及围岩塑性区分布，如图14所示，提取各个角度下拱顶的沉降得图15所示曲线. 由图14、图15可知：1）围岩和支护结构的变形显著受到层理倾角的影响. 层理倾角从0°到90°变化过程中，初期支护拱顶沉降呈现出倒“V ”形变化，即先增大后减小，45°时达到最大值. 位移场分布随着倾角改变，只有0°和90°时存在对称性.2）隧道开挖引起的塑性区形状和范围与层理倾角密切相关. 0°时塑性区范围最小，当层理倾角小于30°时，围岩深部塑性区沿垂直于层理方向发展；而当倾角为75°~90°时，深部塑性区主要沿层理方向发展；层理倾角为45°~60°时，塑性区呈现出“X ”形状，且范围较大，与前文所述岩体在45°~60°时强度较低相对应，表明该倾角范围内易使隧道围岩产生破坏.4 结论本文通过对比分析遍布节理模型与离散元模型在层状岩体三轴压缩以及层状围岩隧道开挖应用中的模拟效果，探讨采用等效参数的遍布节理模型在层状岩体模拟中的适用性，得出以下结论：1）离散元模型能够更好地体现层状岩体的变形图13 实测和模拟的拱顶沉降（Y =25 m ）Fig.13 Measured and simulated vault settlement （Y =25 m）图14 不同层理倾角下围岩变形和塑性区Fig.14 Deformation and plastic zone of adjacent rock mass at various bedding angles拱顶沉降/m m层理倾角/(°)图15 拱顶沉降随层理倾角的变化Fig.15 Vault settlement varies with bedding angles139湖南大学学报（自然科学版）2023 年和破坏特性，但若考虑计算效率，更适合于描述中小尺度层状岩体力学性质；而遍布节理模型由于其本身对节理裂隙考虑的不足，在模拟层状岩体时，需要对部分参数（弹性模量、泊松比以及岩石基质的黏聚力和内摩擦角）进行修正才能用于工程分析，且更适用于大尺度工程岩体的力学行为研究.2）对于新华山隧道工程而言，两种模型在网格单元划分接近的情况下，计算效率相差30~40倍，而校准的遍布节理模型得到的围岩位移与实测结果分别相差8.7%和9.0%，差异性较小，表明该模型兼顾效率的情况下准确度良好.3）层理弱面的抗剪强度和抗拉强度较低，故层状岩体在工程扰动的情况下，容易造成层理剪切滑移破坏以及张拉破坏，在工程中要重点关注.参考文献［1］GAO M，LIANG Z Z，JIA S P，et al．An equivalent anchoring method for anisotropic rock masses in underground tunnelling［J］．Tunnelling and Underground Space Technology，2019，85：294-306．［2］CHEN Z Q，HE C，XU G W，et al．A case study on the asymmetric deformation characteristics and mechanical behaviorof deep-buried tunnel in phyllite［J］．Rock Mechanics and RockEngineering，2019，52（11）：4527-4545．［3］王培涛，杨天鸿，于庆磊，等．含层理构造黑云变粒岩单轴压缩试验及数值模拟［J］．东北大学学报（自然科学版），2015，36（11）：1633-1637．WANG P T，YANG T H，YU Q L，et al．Uniaxial compression testand numerical simulation of stratified biotite granulite［J］．Journal of Northeastern University （Natural Science），2015，36（11）：1633-1637．（in Chinese）［4］SINGH M，SINGH B．High lateral strain ratio in jointed rock masses［J］．Engineering Geology，2008，98（3/4）：75-85．［5］刘爱华，董蕾，董陇军．节理岩体强度参数的数值模拟及工程应用［J］．中南大学学报（自然科学版），2011，42（1）：177-183．LIU A H，DONG L，DONG L J．Numerical simulation andengineering application of strength parameters of jointed rock mass［J］．Journal of Central South University （Science andTechnology），2011，42（1）：177-183．（in Chinese）［6］DO N A，DIAS D，DINH V D，et al．Behavior of noncircular tunnels excavated in stratified rock masses - Case of undergroundcoal mines［J］．Journal of Rock Mechanics and GeotechnicalEngineering，2019，11（1）：99-110．［7］赵勐，肖明，陈俊涛，等．地震动斜入射下层状岩体隧洞接触响应分析［J］．湖南大学学报（自然科学版），2021，48（5）：129-139．ZHAO M，XIAO M，CHEN J T，et al．Analysis on contactresponse of tunnel in layered rock mass subjected to obliquelyincidence earthquake［J］．Journal of Hunan University （NaturalSciences），2021，48（5）：129-139．（in Chinese）［8］左双英，叶明亮，唐晓玲，等．层状岩体地下洞室破坏模式数值模型及验证［J］．岩土力学，2013，34（S1）：458-465．ZUO S Y，YE M L，TANG X L，et al．Numerical model andvalidation of failure mode for underground Caverns in layered rockmass［J］．Rock and Soil Mechanics，2013，34（S1）：458-465．（inChinese）［9］SUN X M，ZHAO C W，ZHANG Y，et al．Physical model test and numerical simulation on the failure mechanism of the roadway inlayered soft rocks［J］．International Journal of Mining Scienceand Technology，2021，31（2）：291-302．［10］DENG P H，LIU Q S，HUANG X，et al．FDEM numerical modeling of failure mechanisms of anisotropic rock masses arounddeep tunnels［J］. Computers and Geotechnics，2022，142：104535．［11］WANG S Y，SLOAN S W，TANG C A，et al．Numerical simulation of the failure mechanism of circular tunnels intransversely isotropic rock masses［J］. Tunnelling andUnderground Space Technology，2012，32：231-244．［12］RIAHI A，CURRAN J H．Full 3D finite element Cosserat formulation with application in layered structures［J］．AppliedMathematical Modelling，2009，33（8）：3450-3464．［13］KARAMPINOS E，HADJIGEORGIOU J，HAZZARD J，et al．Discrete element modelling of the buckling phenomenon in deephard rock mines［J］. International Journal of Rock Mechanics andMining Sciences，2015，80：346-356．［14］刘春，乐天呈，施斌，等．颗粒离散元法工程应用的三大问题探讨［J］．岩石力学与工程学报，2020，39（6）：1142-1152．LIU C，LE T C，SHI B，et al．Discussion on three major problemsof engineering application of the particle discrete elementmethod［J］. Chinese Journal of Rock Mechanics and Engineering，2020，39（6）：1142-1152．（in Chinese）［15］蒋青青，胡毅夫，赖伟明．层状岩质边坡遍布节理模型的三维稳定性分析［J］．岩土力学，2009，30（3）：712-716．JIANG Q Q，HU Y F，LAI W M.Three-dimensional stabilityanalysis of stratified rock slope based on ubiquitous-jointmodel［J］．Rock and Soil Mechanics，2009，30（3）：712-716．（inChinese）［16］朱泽奇，盛谦，梅松华，等．改进的遍布节理模型及其在层状岩体地下工程中的应用［J］．岩土力学，2009，30（10）：3115-3121．ZHU Z Q，SHENG Q，MEI S H，et al．Improved ubiquitous-jointmodel and its application to underground engineering in layeredrock masses［J］．Rock and Soil Mechanics，2009，30（10）：3115-3121．（in Chinese）140。

·核科学与工程·放射性同位素伽玛源准直照射辐射场模拟研究*黄宇晨1,2， 钱易坤1,2， 冯　鹏2， 刘易鑫1， 张　颂1,2，何　鹏2， 魏　彪2， 毛本将1， 朱亚地1（1. 中国工程物理研究院 核物理与化学研究所，四川 绵阳 621900； 2. 重庆大学 光电技术及系统教育部重点实验室，重庆 400044）摘 要： 针对各向同性伽玛源参考辐射场尺寸关键技术问题，GB/T 12162系列标准虽然进行了相关规定，但是该规定并未对准直照射状态下照射室尺寸提出具体要求。

为减小用于辐射检测或监测类仪器仪表检定与量值校准时伽玛参考辐射场内散射影响，本文采用蒙特卡罗方法，研究了同位素放射源准直照射时，照射室尺寸变化对检验点处的剂量率值与能量分布的影响情况，获得了准直照射时伽玛辐射场照射室尺寸的边界条件，建立并完善了伽玛参考辐射场边界研究方法及相关标准细节，为准直照射状态下照射室尺寸设计提供了一种新方法或途径。

关键词： 核仪器仪表校准； 参考辐射； 准直照射； 蒙特卡罗； 边界条件 中图分类号： TL72 文献标志码： A doi : 10.11884/HPLPB202133.200294Simulation of radiation field from isotopic gamma source collimationHuang Yuchen 1,2， Qian Yikun 1,2， Feng Peng 2， Liu Yixin 1， Zhang Song 1,2，He Peng 2， Wei Biao 2， Mao Benjiang 1， Zhu Yadi 1（1. Institute of Nuclear Physics and Chemistry , CAEP , Mianyang 621900, China ;2. Key Laboratory of Optoelectronic Technology & Systems , Ministry of Education , Chongqing University , Chongqing 400044, China ）Abstract ： Aiming at the key technical issues of the reference radiation field size of isotropic sources, GB/T 12162 series of GB standards stipulates the size of the reference radiation field when using an isotropic source.However, there is no specific regulation for the size of the irradiation under collimation. To reduce the influence of scattering in the gamma reference radiation field when used for radiation detection or monitoring instrument verification and value calibration, Monte Carlo simulation was carried out to explore the effect of the size change of the irradiation chamber on the energy distribution and dose rate value during the radiation source collimation. The boundary conditions of the collimated gamma radiation irradiation chamber were obtained, and the details of the gamma reference radiation field boundary research method and related standards are established and improved. The study provides a new method or approach for the size design of the irradiation chamber under the collimated irradiation state.Key words ： nuclear instrument calibration ； reference radiation ； collimated irradiation ； Monte Carlo ；boundary conditions众所周知，核电站与核设施场所中，伽玛射线辐射剂量（率）仪表不仅是辐射防护工作必需的重要工具，更是监测和保障核应用场所与射线应用装置安全的必要条件。

a rXiv:n ucl-t h /441v28Ma y2The decay ρ0→π++π−+γand the coupling constant g ρσγA.Gokalp ∗and O.Yilmaz †Physics Department,Middle East Technical University,06531Ankara,Turkey(February 8,2008)Abstract The experimental branching ratio for the radiative decay ρ0→π++π−+γis used to estimate the coupling constant g ρσγfor a set of values of σ-meson parameters M σand Γσ.Our results are quite different than the values of this constant used in the literature.PACS numbers:12.20.Ds,13.40.HqTypeset using REVT E XThe radiative decay processρ0→π++π−+γhas been studied employing different approaches[1,5].There are two mechanisms that can contribute to this radiative decay: thefirst one is the internal bremsstrahlung where one of the charged pions from the decay ρ0→π++π−emits a photon,and the second one is the structural radiation which is caused by the internal transformation of theρ-meson quark structure.Since the bremsstrahlung is well described by quantum electrodynamics,different methods have been used to estimate the contribution of the structural radiation.Singer[1]calculated the amplitude for this decay by considering only the bremsstrahlung mechanism since the decayρ0→π++π−is the main decay mode ofρ0-meson.He also used the universality of the coupling of theρ-meson to pions and nucleons to determine the coupling constant gρππfrom the knowledge of the coupling constant gρter,Renard [3]studied this decay among other vector meson decays into2π+γfinal states in a gauge invariant way with current algebra,hard-pion and Ward-identities techniques.He,moreover, established the correspondence between these current algebra results and the structure of the amplitude calculated in the single particle approximation for the intermediate states.In corresponding Feynman diagrams the structural radiation proceeds through the intermediate states asρ0→S+γwhere the meson S subsequently decays into aπ+π−pair.He concluded that the leading term is the pion bremsstrahlung and that the largest contribution to the structural radiation amplitude results from the scalarσ-meson intermediate state.He used the rough estimate gρσγ≃1for the coupling constant gρσγwhich was obtained with the spin independence assumption in the quark model.The coupling constant gρππwas determined using the then available experimental decay rate ofρ-meson and also current algebra results as3.2≤gρππ≤4.9.On the other hand,the coupling constant gσππwas deduced from the assumed decay rateΓ≃100MeV for theσ-meson as gσππ=3.4with Mσ=400MeV. Furthermore,he observed that theσ-contribution modifies the shape of the photon spectrum for high momenta differently depending on the mass of theσ-meson.We like to note, however,that the nature of theσ-meson as a¯q q state in the naive quark model and therefore the estimation of the coupling constant gρσγin the quark model have been a subject ofcontroversy.Indeed,Jaffe[6,7]lately argued within the framework of lattice QCD calculation of pseudoscalar meson scattering amplitudes that the light scalar mesons are¯q2q2states rather than¯q q states.Recently,on the other hand,the coupling constant gρσγhas become an important input for the studies ofρ0-meson photoproduction on nucleons.The presently available data[8] on the photoproduction ofρ0-meson on proton targets near threshold can be described at low momentum transfers by a simple one-meson exchange model[9].Friman and Soyeur [9]showed that in this picture theρ0-meson photoproduction cross section on protons is given mainly byσ-exchange.They calculated theγσρ-vertex assuming Vector Dominance of the electromagnetic current,and their result when derived using an effective Lagrangian for theγσρ-vertex gives the value gρσγ≃2.71for this coupling ter,Titov et al.[10]in their study of the structure of theφ-meson photoproduction amplitude based on one-meson exchange and Pomeron-exchange mechanisms used the coupling constant gφσγwhich they calculated from the above value of gρσγinvoking unitary symmetry arguments as gφσγ≃0.047.They concluded that the data at low energies near threshold can accommodate either the second Pomeron or the scalar mesons exchange,and the differences between these competing mechanisms have profound effects on the cross sections and the polarization observables.It,therefore,appears of much interest to study the coupling constant gρσγthat plays an important role in scalar meson exchange mechanism from a different perspective other than Vector Meson Dominance as well.For this purpose we calculate the branching ratio for the radiative decayρ0→π++π−+γ,and using the experimental value0.0099±0.0016for this branching ratio[11],we estimate the coupling constant gρσγ.Our calculation is based on the Feynman diagrams shown in Fig.1.Thefirst two terms in thisfigure are not gauge invariant and they are supplemented by the direct term shown in Fig.1(c)to establish gauge invariance.Guided by Renard’s[3]current algebra results,we assume that the structural radiation amplitude is dominated byσ-meson intermediate state which is depicted in Fig. 1(d).We describe theρσγ-vertex by the effective LagrangianL int.ρσγ=e4πMρMρ)2 3/2.(3)The experimental value of the widthΓ=151MeV[11]then yields the value g2ρππ2gσππMσ π· πσ.(4) The decay width of theσ-meson that follows from this effective Lagrangian is given asΓσ≡Γ(σ→ππ)=g2σππ8 1−(2Mπ2iΓσ,whereΓσisgiven by Eq.(5).Since the experimental candidate forσ-meson f0(400-1200)has a width (600-1000)MeV[11],we obtain a set of values for the coupling constant gρσγby considering the ranges Mσ=400-1200MeV,Γσ=600-1000MeV for the parameters of theσ-meson.In terms of the invariant amplitude M(Eγ,E1),the differential decay probability for an unpolarizedρ0-meson at rest is given bydΓ(2π)31Γ= Eγ,max.Eγ,min.dEγ E1,max.E1,min.dE1dΓ[−2E2γMρ+3EγM2ρ−M3ρ2(2EγMρ−M2ρ)±Eγfunction ofβin Fig.5.This ratio is defined byΓβRβ=,Γtot.= Eγ,max.50dEγdΓdEγ≃constant.(10)ΓσM3σFurthermore,the values of the coupling constant gρσγresulting from our estimation are in general quite different than the values of this constant usually adopted for the one-meson exchange mechanism calculations existing in the literature.For example,Titov et al.[10] uses the value gρσγ=2.71which they obtain from Friman and Soyeur’s[9]analysis ofρ-meson photoproduction using Vector Meson Dominance.It is interesting to note that in their study of pion dynamics in Quantum Hadrodynamics II,which is a renormalizable model constructed using local gauge invariance based on SU(2)group,that has the sameLagrangian densities for the vertices we use,Serot and Walecka[14]come to the conclusion that in order to be consistent with the experimental result that s-waveπN-scattering length is anomalously small,in their tree-level calculation they have to choose gσππ=12.Since they use Mσ=520MeV this impliesΓσ≃1700MeV.If we use these values in our analysis,we then obtain gρσγ=11.91.Soyeur[12],on the other hand,uses quite arbitrarly the values Mσ=500 MeV,Γσ=250MeV,which in our calculation results in the coupling constant gρσγ=6.08.We like to note,however,that these values forσ-meson parameters are not consistent with the experimental data onσ-meson[11].Our analysis and estimation of the coupling constant gρσγusing the experimental value of the branching ratio of the radiative decayρ0→π++π−+γgive quite different values for this coupling constant than used in the literature.Furthermore,since we obtain this coupling constant as a function ofσ-meson parameters,it will be of interest to study the dependence of the observables of the reactions,such as for example the photoproduction of vector mesons on nucleonsγ+N→N+V where V is the neutral vector meson, analyzed using one-meson exchange mechanism on these parameters.AcknowledgmentsWe thank Prof.Dr.M.P.Rekalo for suggesting this problem to us and for his guidance during the course of our work.We also wish to thank Prof.Dr.T.M.Aliev for helpful discussions.REFERENCES[1]P.Singer,Phys.Rev.130(1963)2441;161(1967)1694.[2]V.N.Baier and V.A.Khoze,Sov.Phys.JETP21(1965)1145.[3]S.M.Renard,Nuovo Cim.62A(1969)475.[4]K.Huber and H.Neufeld,Phys.Lett.B357(1995)221.[5]E.Marko,S.Hirenzaki,E.Oset and H.Toki,Phys.Lett.B470(1999)20.[6]R.L.Jaffe,hep-ph/0001123.[7]M.Alford and R.L.Jaffe,hep-lat/0001023.[8]Aachen-Berlin-Bonn-Hamburg-Heidelberg-Munchen Collaboration,Phys.Rev.175(1968)1669.[9]B.Friman and M.Soyeur,Nucl.Phys.A600(1996)477.[10]A.I.Titov,T.-S.H.Lee,H.Toki and O.Streltrova,Phys.Rev.C60(1999)035205.[11]Review of Particle Physics,Eur.Phys.J.C3(1998)1.[12]M.Soyeur,nucl-th/0003047.[13]S.I.Dolinsky,et al,Phys.Rep.202(1991)99.[14]B.D.Serot and J.D.Walecka,in Advances in Nuclear Physics,edited by J.W.Negeleand E.Vogt,Vol.16(1986).TABLESTABLE I.The calculated coupling constant gρσγfor differentσ-meson parametersΓσ(MeV)gρσγ500 6.97-6.00±1.58 8008.45±1.77600 6.16-6.68±1.85 80010.49±2.07800 5.18-9.11±2.64 90015.29±2.84900 4.85-10.65±3.14 90017.78±3.23Figure Captions:Figure1:Diagrams for the decayρ0→π++π−+γFigure2:The photon spectra for the decay width ofρ0→π++π−+γ.The contributions of different terms are indicated.Figure3:The pion energy spectra for the decay width ofρ0→π++π−+γ.The contri-butions of different terms are indicated.Figure4:The decay width ofρ0→π++π−+γas a function of minimum detected photon energy.Figure5:The ratio Rβ=Γβ。

为了便于一些同学阅读矿井通风与空调方面的英文参考资料和为以后撰写英文论文发表，下面给出了一些常见的矿井通风与空调中英文专业词汇。

Abandoned workings 废弃坑道Absolute pressure 绝对压力Acceptable accuracy 允许精度Active regulation 主动调节（增压调节）Actual characteristic curves实际特征曲线Adiabatic and isentropic processes等熵线绝热的过程Adiabatic saturation process 绝热饱和过程Aerofoils风板Aerosol particles 气溶胶粒子Air crossings 风桥Air mover 鼓风机Air power空气动力Air pressure management 风压管理Air quantity survey空气质量调查Air regulators 风窗Airborne pollutants空气污染物Airflow measurements 风流测定Airflow reversal反向风流Airlock 气闸Airlocks 风门Airway resistance curve风路阻力曲线Alpha, beta and gamma radiation阿尔法、贝塔和伽玛辐射Altimeters 高度计Angular velocity角速度Asbestos 石棉Atkinson equation阿特金森方程式Atmospheric conditions 大气状态Atmospheric pressure 大气压力Auxiliary ventilation 辅助通风Axial fan轴流风机Axial impeller轴向式叶轮Backfill material 充填材料Barometers 气压计Barometric pressure at inlet 入口气压Becquerel (Bq) 贝克勒尔Bernoulli's equation for ideal fluids 理想流体伯努力方程Biot number 比奥数Blackdamp 窒息气体Blast fume炮烟Booster fans 局扇Boreholes 钻孔Branch resistance分支阻力Branch tree分支树Brattice curtain 风帘Brattices 风帘Bronchioles 细支气管Brownian motion 布朗运动Buoyancy (natural draft) effect浮力作用Burying the fire掩埋火源Cage and skip 罐笼和箕斗Carbon dioxide produced生成二氧化碳Carbon dioxide 二氧化碳Carbon monoxide 一氧化碳Carcinogenic (cancer causing) dusts 致癌粉尘Carnot cycle 卡诺循环Centrifugal fan离心风机Centrifugal impeller离心叶轮Chemical absorption化学吸收Chézy-Darcy equation谢兹-达西方程Chilled water spray chamber 冷却液体喷雾室Choke effect瓶颈效应Circular airway 循环风路Closed loop闭环Closed path回路Coal workers' pneumoconiosis (CWP) 煤工尘肺病Coefficient of drag阻力系数Coefficient of dynamic viscosity动力粘度系数Coefficient of friction摩擦系数Compressed air-assisted sprays 压气助喷雾Compressible flow可压缩流Computational fluid dynamics计算流体力学Condenser cooling tower 凝气器降热塔Condenser 冷凝器Consolidation 固结Contaminants 污染物Continuity equation 连续方程Controlled partial recirculation 受控开路循环通风Controlled recirculation in headings 掘进面受控循环通风Convected energy 扩散能Convective heat transfer 对流换热Conveyance运输工具Copper orebody 铜矿体Cross section of a duct or airway 管道或风路断面Curie, Ci 居里Cylindrical cyclone重力旋流器Dealing with a spontaneous heating 处理自热Degrees Celsius 摄氏度Degrees Kelvin绝对温度Density of gases 气体密度Desorption kinetics解吸动力学Dew point hygrometers 露点毛发湿度计Diaphragm gauge 隔膜片仪表Diesel emissions柴油机排放物Diesel exhaust fume柴油机尾气Diesel particulate matter柴油机颗粒物质Differential pressure instruments 微压差计Dimensionless无量刚Disaster management 灾害管理District systems 分区通风系统Dose rates 剂量率Downcast shaft入风井Droplet diameter雾滴直径Duct system风管系统Dust suppression 降尘Dynamic behavior of molecules 分子运动特征Electrochemical methods电化学方法Electrostatic precipitators 电除尘器Emanation of radon 氡的辐射Empirical method 经验方法Energy recovery device 能量回收装置Enthalpy of moist air潮湿空气的焓Enthalpy 焓Entry and exit losses 入口和出口阻力损失Environmental engineering 环境工程Equivalent length当量长度Equivalent resistance等效风阻Equivalent resistance等效阻力Equivalent sand grain roughness相当砂粒粗糙度Escape way 逃生通道；安全通道Euler's equation欧拉方程Evaporator蒸发器Excavating the fire挖掘火源Exhausting air 抽出空气Exhausting system 抽出式通风系统Explosive dusts 爆炸粉尘Explosives炸药Fan characteristic curve风机特征曲线Fan maintenance 风机维护Fan performance 风机性能Fan static pressure风机静压Fan total pressure风机全压Fan velocity pressure风机速度压Fibrogenic dusts 矿渣粉尘Filament and catalytic oxidation (pellistor) detectors丝状催化氧化探测器Fire triangle 火三角Firedamp 甲烷Firefighting with water 以水灭火First law of thermodynamics 热力学第一定律Fixed point measurement固定点测量Fixed quantity branch固定风量分支Flame safety lamps灯具安全火焰Flexible tubing 柔性风筒Flooded orifice scrubber 水淹孔洗涤器Flooding and sealing off 溢出和密封作用Flow work 流动功Fluid mechanics 流体力学Fluid pressure 流体压力Fog 雾Fogged air 雾气Forcing air压入空气Forcing or blowing system 压入式通风系统Fourier number傅里叶数Fragmented rock 破碎岩石Free crystalline silica (quartz, sand stones, flint)游离硅晶体Friction factor摩擦系数Frictional flow 摩擦流动Frictional losses摩擦损失Frictional pressure drop摩擦压降Frictional resistance 摩擦阻力Frictionless manner 无摩擦状态Gas adsorbents 气体吸收剂Gas chromatography气相色谱Gas constants 气体常数Gas drainage 瓦斯抽放Gas laws 气体定律Geothermic gradient 地热梯度Gob drainage采空区抽放气体Grab samples 样品收集Gravitational field 重力场Gravitational settlement of particles 引力沉降颗粒Gravitational settlement 重力沉降Hair hygrometers 毛发湿度计Hardy-Cross technique哈代克劳斯技术Haulage airways 运输风路Haulage level 运输平巷Heat capacity 热容Heat cramps 中暑痉挛Heat diffusivity 热扩散系数Heat exchanger 换热器Heat exchange换热Heat exhaustion 热量消耗Heat fainting 热昏厥Heat flux 热通量Heat illness 中暑Heat rash 热疹Heat stroke 中暑Heat tolerance 耐热性Heat transfer coefficient 传热系数High expansion foam高倍数泡沫High pressure tapping高压测压孔Hoisting shaft提升竖井Hot wire anemometer热线风速仪Hydraulic radius水力半径Hydrogen sulfide硫化氢气Hydrolift system 水力提升系统Hydropower 水电Ice system 冷却系统Ideal gas 理想空气Ideal isothermal compression理想恒温压缩Immediate response 应急反应Induction 感应Industrial Hygienists 工业卫生学家Inhalation rate吸入速度Initiation of explosions引发爆炸Injection of inert gases注射惰性气体Inlet and outlet ducts入口和出口管In-situ measurement 现场测量Intake airway 进风风路Interception and electrostatic precipitation 截留和静电沉淀Interference factor干扰因素Interferometers干涉计Internal Energy 内能Ionization smoke detectors离子感烟探测器Iron pyrites黄铁矿Jet fan 射流风机Junction节点Kata thermometer 卡它计Kinetic energy 动能Kirchhoff's Laws 基尔霍夫定律Laminar and turbulent flow层流和紊流Laminar resistance 层流阻力Laminar sublayer层流次边界层Laser spectroscopy激光光谱学Latent (or hidden) heat of the air空气的潜热Layout of mine 矿井布置Leakage control漏风控制Legislation 法规Level workings阶段工作面Loading station装运站Longitudinal fittings纵向装备Longwall长壁开采法Machine mounted gas monitors悬挂式气体检测器Main fans 主扇Main haulage route主运输道Main return 主（总）回风道Manometers 压差计Mass flow 质量流量Mass spectrometers质谱仪Mean free path 平均自由程Mean velocity of air 平均风速Mesh selection网孔选择Mesh网Metabolic heat balance 代谢热量平衡Metabolic heat代谢热Metal mine fires金属矿井火灾Meteorology 气象Methane drainage瓦斯排放Methane 甲烷Method of mining 采矿方法Mine climate 矿井气候Mine resistance 矿井阻力Mine ventilation 矿井通风Mist eliminator 除雾器Mist 雾Moisture content (specific humidity) of air空气的含湿量Momentum 动量Monitoring systems 监测系统Moving traverses运动线路Natural ventilating effect自然通风影响Natural ventilation 自然通风Neutral skin temperature 中性表皮温度Nikuradse's curves 尼库拉则曲线Nondispersive infra-red gas analyzer非分散红外线气体分析仪Nuisance dusts 粉尘污染Numerical method数值方法Nusselt number努塞尔数Old workings老工作面One standard atmosphere 一个标准大气压Open and concealed fires 明火和隐蔽火灾Ore pass 放矿溜井Ore production矿石生产Orebody deposit 矿体Outbursts from roof and floor 顶板和底板瓦斯突出Overlap systems of auxiliary ventilation 混合式局部通风Oxides of nitrogen氧化氮Oxygen Consumption耗氧量Parallel network or circuit并联网络或回路Paramagnetic analyzer 顺磁分析仪Passive regulator 可调风窗Pellistor methanometers瓦斯检定器Peripheral velocity圆周速度Permanent environmental monitors 持久环境监控Permeability 渗透率Personal dosemeters 个人剂量计Personal respirators 个体呼吸器Phases of oxidation氧化反应阶段Photometric (light-scattering) methods 分光光度Physical adsorption物理吸附Physical thermodynamics 物理热力学Pick face flushing and jet-assisted cutting 锯齿面冲洗与喷气助推器切割Piezoelectric instruments 压电仪器Pitot-static tube皮托静压管Polyvinyl chloride (PVC)聚氯乙烯Potential energy 势能Prandtl number 普兰特尔数Precautions against spontaneous combustion自燃预防Pressure energy 压能Pressure head 压头Pressure surveys压力调查Pressure transducers 压力传感器Pressure-volume surveys压力容积测量Profilometer轮廓仪Psychrometric chart 温湿图Psychrometric measurements 干湿度测量Push-pull system 压-抽混合式通风系统Radial velocity径向速度Radiation 辐射Radiative heat transfer 辐射传热Radioactive decay and half-life放射衰变和半衰期Radon daughters氡子体Radon decay constant 氡的衰变常数Radon, Rn氡气Ramp 斜坡道Rates of heat production 生产率Rates of oxygen consumption 氧气消耗率Refrigerant fluid 制冷液Refrigeration cycle制冷循环Refrigeration systems 制冷系统Refuge chambers避难洞室Regulator 调节器Relative humidity and percentage humidity相对湿度和湿度率Removal of dust from air 气体除尘Re-opening a sealed area重开封闭区Respirable dust呼吸性粉尘Respiratory system 呼吸系统Return airway 回风巷Reynolds Number雷诺数Room and pillar房柱式Rotating vane anemometer旋转叶片风速表Rough pipes 粗管Roughness粗糙度Safety and Health 安全卫生Saturation vapor pressure 饱和蒸汽压Sealants密封剂Seals 密闭Second law of thermodynamics 热力学第二定律Self-heating temperature (SHT) 自热温度Self-rescuers 自救器Sensible heat of the air 空气的显热Series network or circuit串连网络或回路Shaft fittings 井筒装备Shaft wall井壁Shear stress 剪切应力Shock loss factor冲击损耗系数Shock losses 冲击损失Short-Term Exposure Limit (STEL) 短时间接触阈限值SI system of units 国际标准单位体系Sigma heat 西格玛热Smoke tube烟筒Smoking and flame safety lamps 烟火安全灯Smooth concrete lined光滑混凝土内衬Specific heat (thermal capacity)比热容Specific heats 比热Spontaneous combustion of sulfide ores硫化矿自燃Spontaneous combustion自燃Spontaneous heating 自热Spot cooler 现场冷却器Spray fan 喷雾风机Steady flow energy equation稳流能量方程Steady flow physical thermodynamics稳流物理热力学Steady-flow thermodynamics 稳定流热力学Stokes' diameter斯托克斯粒径Stoping areas 回采区Stoppings 密闭Subsurface openings 地下空间Subsurface ventilation 地下通风Sulfide dust explosions 硫化矿粉尘爆炸Sulfur dioxide二氧化硫Sulfuric acid vapor硫酸雾Swinging vane anemometer摆动叶片风速表Tangential velocity at outlet出口切向速度Temperature-entropy diagram温熵图Temporary stopping暂时停止Terminal velocities 自由沉降速度The square law平方定律Thermal conductivity of insulation 绝缘导电温度Thermal conductivity导热系数Thermal equilibrium 热平衡Thermodynamic state 热力学状态Thermoluminescent dosemeters (TLD) 热释光剂量计Thermoregulation 体温调节Threshold limit values (TLV) 阈限值Through-flow ventilation 贯穿通风Time-Weighted Average (TWA)时间加权平均Total energy balance 总能量守恒Total shaft resistance 井筒总阻力Tube bundle systems 管束系统Turbulent resistance紊流阻力U tube manometers U型压差计U tube U型管Uncontrolled recirculation 无控循环通风Underground ventilation system 地下通风系统Unloading station卸载站Upcast shaft出风井Uranium mines 铀矿Vasodilation 血管舒张Velocity contour等流速线Velocity limit速度限值Velocity pressures 动压Velometer速度计Ventilation circuit 通风回路Ventilation door 风门Ventilation engineers 通风工程师Ventilation network analysis通风网络分析Ventilation planning 通风设计Ventilation raise 通风天井Ventilation survey team 通风测量术语Ventilation survey通风测量Venturi scrubber文丘里洗涤器Vertex顶点Viscosity 粘度Viscous drag粘性阻力Volume flow 体积流量Volumetric efficiency 容积效率Vortex-shedding anemometer漩涡式风速表Water gauge pressure 水柱压力Water infusion 注水（水封孔）Water mass flowrate 水质量流量Water vapor content 水蒸气含量Wet and dry bulb hygrometers (psychrometers)干湿球温度表Wet bulb thermometer 湿球温度计Wet Kata thermometer湿球卡他温度表Wet scrubbers湿式除尘器Wetting agents 润湿剂Worked-out area采空区Working face工作面Working level month, WLM 工作水平月Working Level 开采水平Zinc blende闪锌矿——上述词汇摘录自：吴超主编。

Ž.Journal of Hazardous Materials 591998211–233Mathematical modelling of accidental gas releases Helena Montiel ),Juan A.Vılchez,Joaquim Casal,Josep Arnaldos ´()Centre d’Estudis del Risc Tecnologic CERTEC ,Department of Chemical Engineering,Uni Õersitat `Politecnica de Catalunya-Institut d’Estudis Catalans Diagonal 647,08028Barcelona,Catalonia,Spain `Received 16April 1997;accepted 26October 1997AbstractThe increasing use of natural gas entails large networks,located mostly in highly populated zones.Therefore,any accidental releases of gas that might occur through damaged pipes imply a risk which must be controlled.In such cases,the prediction of release flow-rate and time of duration of the emergency is very important.However,the mathematical models currently available for performing this prediction show a major gap in the range of conditions over which they can be applied with a reasonable degree of accuracy.Furthermore,they do not take into consideration the unsteady state existing when a safety device is closed after a certain time of release.In this paper,a new model is developed as a combination of the classical ‘hole’and ‘pipe’models,for the calculation of gas releases in distribution systems at medium and low pressures.Ž.This model can be applied to these two cases small hole in a tank or a pipe and full bore holes X Ž.Ž2.Abbreviations:A ,Parameter defined by Eq.47;A ,Area of the cross section of the pipe m ;A ,c or Ž2.Ž.Ž.Hole area m ;B ,Parameter defined by Eq.47;C ,Velocity of the sound m r s ;for ideal gases ŽŽ..1r 2X Ž.Ž.C s k P M r R P T ;C ,Parameter defined by Eq.47;CPR,critical pressure ratio -;D ,Pipe diameter Ž.Ž.Ž.Ž.Ž.m or mm ;d ,Hole diameter mm ;F ,Function defined by Eq.45-;f ,Friction factor -;G ,Mass flux or Ž2.Ž.Ž.Ž.kg r m s ;k ,Heat capacity ratio -;K ,Pressure drop coefficient for each fitting -;L ,Pipe length m ;L ,i e Ž.Ž.Ž.Equivalent pipe length m ;M ,Molecular weight kg r kmol ;m ,Mass of gas contained in the pipe kg ;Ma,Ž.Ž.Ž.ŽMach number,u r C -;N ,number of fittings -;P ,Pressure Pa or bar abs ;P ,critical pressure Pa or bar i c .Ž3.Ž.abs ;Q ,Mass discharge rate kg r s or Nm r h ;R ,Ideal gas constant J r kmol K ;Re,Reynolds number,m Ž.Ž.Ž.Ž.Ž.Du r r m -;T ,Temperature K ;t ,Time s ;t ,Critical time s ;u ,Velocity of the gas m r s ;V ,Pipe c P Ž3.Ž.Ž.Ž.volume m ;Y ,Parameter defined by Eq.14-;Greek letters:a ,Parameter defined by Eq.36;b ,Ž.Ž.Ž3.ŽParameter defined by Eq.39;F ,Parameter defined by Eq.47;r ,Density kg r m ;m ,Viscosity kg r m .s ;Subscripts:0,Steady-state;0,Steady state or initial values in subsonic flow;1,At the beginning of the sub pipe;2,Inside the pipe and on a level with the hole;3,At the hole;a,In the surroundings,at atmospheric pressure;E,Regulator inlet conditions;i,Initial conditions;j ,Any point )Corresponding author.0304-3894r 98r $19.00q 1998Elsevier Science B.V.All rights reserved.Ž.PII S 0304-38949700149-0()212H.Montiel et al.r Journal of Hazardous Materials591998211–233and also to the intermediate situation,a large hole in a pipe.Its application to various accident scenarios is discussed.q1998Elsevier Science B.V.Keywords:Mathematical modelling;Gas releases;Accidental1.IntroductionThe consumption of natural gas in the industrialized countries has increased steadily in recent years.Its clean-combustion characteristics and its ease of distribution at lowŽpressures justify its wide use both for the generation of electricity and as a domestic .fuel,which has doubled from1990to1994,resulting in a worldwide consumption of 2170=1012m3in1994.This use has meant installing and maintaining complex piping systems to transport and distribute gas,a great deal of them located in highly populated zones.Due to these facts,accidents caused by the loss of containment of natural gas can involve substantial economic losses and even victims amongst the population.Polyethylene tubing has largely replaced the steel or cast iron gas pipes which were used traditionally.Polyethylene pipe has several advantages:it is easier to install,has good mechanical properties and it is not corroded by damp soils.Thousands of kilometers of this pipe have been installed and its use is still increasing.Gas pipes,which are mostly installed underground,can be damaged by various activities:underground work,bull-dozers,etc.If there is a partial or complete breaking of the pipe,the result will always be a loss of containment.Such a situation can be no more than an incident if the leak is small or is rapidly detected,or can lead to an accident depending on the circumstances.w xA recent survey of185accidents involving natural gas1showed that,of the total, 131were caused during transportation,either by road,railway,ship or pipeline.The analysis of these data clearly shows the relatively high frequency of accidents in pipes: 127of them occurred in piping systems.The most frequent causes of the accidents were mechanical failure,impact failure,human error and external events.Amongst theŽ.accidents arising from impact failure39,the most frequent specific cause wasŽ.Ž.Ž. excavating machinery21accidents,followed by vehicles5and heavy objects5.Ž. Other specific causes due to external events were ground subsidence4cases and Ž.sabotage r vandalism4cases.In the analysis of such accidents,or when forecasting the consequences of this type of hypothetical situation,the first stage is to estimate the flow-rate at which the gas isw xbeing released through the damaged pipe2.This is a complex problem,which in fact begins with the definition of the‘hole diameter’:a particular opening size must be assumed;it can be a small hole,a large one or even the whole cross-section of the pipe if it is completely broken.Then,the flow-rate through this hole must be calculated,bothŽfor the case of constant flow-rate and for the case of decreasing velocity if an.emergency valve is automatically closed,for example.In this paper,a mathematical model is described which allows the prediction of the flow-rate of gas through a damaged pipe.()H.Montiel et al.r Journal of Hazardous Materials591998211–2332132.Accidental releaseIn the flow of compressible fluids there can be significant changes in fluid density; this implies significant variation as well of pressure and temperature.Therefore,the analysis of such systems involves four equations:the equation of state and those of continuity,momentum and energy.This makes the analysis rather complicated.To simplify it,it is usually assumed that the flow is reversible and adiabatic,these conditions implying isentropic flow.Furthermore,it is often assumed that the fluid is aŽ.perfect gas with constant specific heat average value.Ž.In the case of gas release the existing models describe two situations:a gas flow through a hole,the pipe being considered to be like a tank:usually called‘hole models’;Ž.b gas flow through a hole which corresponds to the complete breaking of the pipe:‘pipe models’.w x Hole models have been widely treated in the literature.Woodward and Mudan3 also developed a model for the calculation of liquid and gas discharge rates through holes in process vessels,which takes into account the decrease in pressure as a function w xof time.Levenspiel4developed an interesting model which was later adapted to thew xestimation of accidental releases by Crowl and Louvar5.It considers the release of a gas contained in a tank through a hole;two assumptions are made:pressure inside the tank is constant and gas expansion is isentropic.This model can be applied to the loss of containment of a gas from a pipe if the hole is small;however,for large holes the model overpredicts the flow-rate.The pressure inside the pipe is considered to be constant,and therefore unsteady state is not taken into account.Neither does the model take into account either changes in gas exit velocity Ž.thus implying constant outlet temperature or pressure drop in the pipe.All these aspects make this model adequate for the prediction of release through a hole in a tank, but not for the description of release through a hole in a pipe.w xThe same authors4,5also describe a model for the case of complete rupture of theŽpipe,assuming adiabatic flow i.e.without any heat loss or heating from the surround-.ings.This model assumes isentropic release;a constant pressure is assumed at an initial point in the pipe,and pressure drop along it is taken into account.The model gives good predictions for the case of complete breaking of the pipe,but it cannot be applied to the flow through holes with a diameter smaller than the pipe diameter.Numerical solutionsw x w x w x w x¨have been treated by Gyori6,Cochran7and Farina8.Giraldo9applied the pipe model to a non-adiabatic case considering heat transfer between the pipe and the surroundings.These two models leave a major gap:all the cases ranging from a relatively smallŽhole up to a large orifice with a diameter smaller than that of the pipe hole models can. be applied to small orifices,and pipe models to the complete breaking of the pipe. Nevertheless,this range of hole diameters corresponds to a situation which can be found in those accidents caused by accidental releases from pipes.Furthermore,the aforemen-tioned models do not take into account the unsteady-state regime found,for example, when a safety device is closed and pressure and gas flow-rate decrease progressively.ŽAlthough the hole models and pipe models are often sufficient for hazard analysis case()214H.Montiel et al.r Journal of Hazardous Materials591998211–233.of an hypothetical incident,a more accurate model can be useful,especially for incidentŽ. investigation purposes when more precise information on the release is available.A new model bridging this gap would therefore be very useful for the case of accidental gas releases.In the following sections a model is described which takes into account all these aspects,thus covering all the possible hole diameters—assuming diverse hypothesis described in Section3—for the release of gas from a damaged pipe. This new model,which is a combination of the two aforementioned models,should leadŽto the same predictions as the hole models for the case of small holes as,in this case,.pressure inside the pipe will undergo little variation;and should likewise give the same prediction as the pipe model for those cases in which the hole diameter is equal to the pipe diameter.This model has been developed for gas distribution systems operating at medium and Ž.low pressures i.e.pressures at which ideal gas behaviour exists.It is therefore different from other models intended for dealing with accidental releases occurring in long,w xhigh-pressure gas pipelines10,11,in which the propagation of pressure disturbances within the pipe should be also taken into account.3.Hole–pipe generalized modelThe system analyzed here is shown schematically in Fig.1.From the distributionŽ.main,there is a length of pipe L after which there is a hole with a certain diameter through which the release takes place.Fig.1.Accident scenario.()H.Montiel et al.r Journal of Hazardous Materials 591998211–233215The conditions at the various points are considered:Point 1,at the beginning of the pipePoint 2,inside the pipe and on a level with the holePoint 3,at the holepoint in the surroundings,at atmospheric pressure.Ž.The following hypothesis are assumed:a the model is devised for gas distribution Žsystems operating at low and medium pressure,with relatively small pipe diameters up .to 350mm i.d.and lengths up to approximately 5km.This means that ideal gas Ž.properties are assumed with C s constant .The model can be applied then at pressures p Ž.at which ideal gas behaviour exists;b the model does not take into account the initial Ž.depressurization before the gas release due to the pressure drop undergone by the gas while flowing through the pipe.In fact,this effect is small and,in any case,neglecting it Ž.implies a certain margin of safety;c isentropic flow at the release point and adiabatic Ž.flow in the pipe are assumed;d a model of flow essentially 1D is assumed.The conditions downstream of the release point are relatively important only in the case of total rupture of the pipe.But even in this case,the contribution to the release of Žthis downstream length of pipe will be small and will last a very short time of the same .order of magnitude as that corresponding to the unsteady state emptying of the pipe .By applying the energy and momentum equations to the adiabatic flow through a pipe,the following equation is obtained:k q 1P T M P 2P 24fL 1221e ln q y q s 01Ž.2ž/ž/ž/k P T T T D RG 2121where L is the equivalent length of the pipe,which can be calculated using the e following expression:DL s L q ÝN P K P 2Ž.e i i ž/f In this expression,the Fanning friction factor can be calculated as follows:For Re -100000,with the classical Blasius equation:f s 0.079Re y 0.253Ž.For Re )100000,the Colebrook equation can be used.In this case,for smooth w x polyethylene pipes,this equation has been transformed into an explicit expression 12:f s 0.0232Re y 0.15074Ž.The mass discharge rate at the orifice can be calculated using the following expression,obtained from the continuity equation and the law of ideal gases for an isentropic expansion:2k q 12M kP P k k a a Q s A P P P y 5Ž.m or 2)ž/ž/R P T k y 1P P 222()H.Montiel et al.r Journal of Hazardous Materials 591998211–233216Fig.2.Different possibilities of release.An empirical discharge coefficient is often used in this expression;for situations of Ž.uncertainty for example,when the exact shape of the hole is unknown ,a conservative value of 1is usually recommended;this value has been taken here.The value of the flow-rate at the hole will depend on whether the flow is sonic or Ž.subsonic.This will be established by the critical pressure ratio CPR ,which indicates,for each case,the transition from sonic flow to subsonic flow:k P 2ak y 1CPR s s 6Ž.ž/P k q 12c P is the critical pressure at Point 2,above which there will always be sonic flow at 2c Ž.the outlet release into the environment,i.e.,towards atmospheric pressure .If the pressure inside the pipe becomes equal to or lower than P ,the flow at the outlet will 2c be subsonic.For the case of natural gas,with k s 1.28,and for a release towards atmospheric pressure,the value of the critical pressure is P s 1.82bar abs.2c The various possibilities are shown in Fig.2.If there is subsonic flow in the pipe,the flow at the orifice can be either sonic or subsonic,depending on the ratio between the Ž.diameter of the pipe and that of the hole.If there is sonic flow in the pipe at Point 2,there will always be sonic flow at the orifice.This leaves,therefore,three different situations to be analyzed.3.1.Subsonic flow in the pipe,sonic flow at the holeIf it is assumed that all the gas flowing through the pipe also flows through the hole,the following mass balance can be established:Q s A P u P r s A P u P r s A P u P r 7Ž.m or 33c 22c 11In this case,although the flow in the pipe is subsonic,the ratio between pipe diameter and hole diameter has a value such that gas release takes place with sonic velocity at the hole.For sonic flow at the hole,Q does not depend on P ,and taking into account Eq.m a Ž.Ž.6,Eq.5becomes:k q 1M 2k y 1Q s A P P k P 8Ž.)m or 2ž/RT k q 12()H.Montiel et al.r Journal of Hazardous Materials 591998211–233217By applying the continuity law to the flow of gas through the pipe and the orifice,and expressing the mass flux as a function of Mach number,the following expression is obtained for the mass flux through the pipe and the hole:k q 1A M 2kM kM ork y 1G s P k s Ma P s Ma P 9Ž.)21122((A RT k q 1RT RT c 212Ž.By introducing this value of G into Eq.1,we obtain:22k q 1P T M P P 4fL 1221e ln q y q s 0k q 12k P T T T D 2121A M 2ork y 12R P k 2ž/ž/A RT k q 1c 210Ž.Taking into account the following relationships:Y 1T sP T 11Ž.21Y 2Ma 11P s P P 12Ž.21(Ma Y 22Ma 12r sP r 13Ž.21(Ma Y 21where k y 12Y s 1q P Ma 14Ž.j j ž/2Ž.and substituting them into Eq.10,the following expression of the equation of energy is obtained which defines the present system:2A c 22ž/k q 1Ma Y Ma 4fL A 212e or ln q 1y q s 015Ž.2k q 12ž/2k Y D Ma Ma 2112k y 1k k q 1In these conditions,the following relationships apply:2Y 1P )P P Ma 16Ž.(211k q 1P a-CPR 17Ž.P 2()H.Montiel et al.r Journal of Hazardous Materials 591998211–233218The parameters corresponding to the point of release,Point 3,are defined for sonic flow as follows:k 2k y 1P sP P 18Ž.32ž/k q 12T s P T 19Ž.32ž/k q 112k y 1r s P r 20Ž.32ž/k q 1and the parameters defining the gas at atmospheric pressure,for the same situation,are finally:k y 1P k a T sP T 21Ž.a 3ž/P 31P k ar s P r 22Ž.a 3ž/P 33.2.Subsonic flow in the pipe and at the orificeIn this case,the following relationship applies:P aG CPR 23Ž.P 2Ž.The mass discharge rate is now given by Eq.5.By applying the continuity equation to the flow of gas through the pipe and through the hole,and expressing the mass flux as a function of Mach number,the following expression is obtained:2k q 1A M 2k P P kM k k ora a G s P y s Ma P 211()ž/ž/A RT k y 1P P RT c 2221kMs Ma P 24Ž.22(RT 2()H.Montiel et al.r Journal of Hazardous Materials 591998211–233219Ž.Ž.Taking into account again Eqs.11–13,and introducing the value of G given by Ž.Eq.24,the equation of energy can be expressed finally as:2A c 2ž/k q 1Ma Y k y 1A 21or ln q 22k q 12k Y 2k Ma 21P Ma Y P Ma Y k k a 22a 22y ((ž/ž/P Ma Y P Ma Y 111111Ma 24fL 2e =1y q s 025Ž.2ž/D Ma 1ŽŽ.This expression,together with the equation of continuity Eq.5,describes the phenomenon for this case.As the flow is subsonic,the conditions at Point 3and the atmospheric ones have the Ž.Ž.same values,given by Eqs.21and 22by replacing P ,T and r with the values 333corresponding to the conditions at Point 2.In the case of full bore hole,Point 2and Point 3are the same.Therefore,there is no Ž.Ž.isentropic expansion between these two points i.e.P s P s P ;Eq.5cannot be 23a ŽŽed and the mass balance Eq.24becomes:kM kM G s Ma P s Ma P 24b Ž.1122((RT RT 12Ž.Ž.Ž.This is the expression which must be used together with Eq.1and Eqs.11–13to solve this case.()3.3.Sonic flow in the pipe at Point 2and at the hole2Y 1When P -P P Ma the flow in the pipe becomes sonic at Point 2.In this case,(211k q 1the flow through an orifice can only be sonic.Now,the starting point is the expression which gives the mass discharge rate through Ž.an orifice sonic flow :k q 1M 2k y 1Q s A P Ma P P P k P 26Ž.)m or 11ž/RT k q 11Following the same path as in the previous cases,and applying the following relationships:P 22Y 212s Ma P 27Ž.12k q 1P 1P 2T 112s 28Ž.22P T Ma 2112P T 2Y 121s 29Ž.2P T k q 1Ma Ž.211()H.Montiel et al.r Journal of Hazardous Materials 591998211–233220the following equation is obtained for the energy balance:2A c ž/k q 12Y 14fL A 1e or ln q 1y q s 030Ž.2k q 12ž/2k D k q 1Ma Ma Ž.112k y 1k k q 1The conditions at Point 2are now given by the following relationships:2Y 1P s Ma P P 31Ž.(211k q 12Y 1T s P T 32Ž.21k q 1r s Ma P r 33Ž.211(2Y 1And the parameters at Point 3and for the atmospheric conditions are defined by Eqs.Ž.Ž.18–22.Ž.Ž.Ž.Eqs.15,25and 30must be solved using an iterative method.In this case,the secant method has been used.In the case of subsonic flow in the pipe and at the hole,w x the Wegstein method to accelerate convergence 13has been used to solve faster the system formed by the equation of continuity and the equation of energy.Table 1shows the sets of equations to be solved for each case in steady state,when Q -Q .m max 4.Examples of applicationTo check the validity of the proposed model,it has been applied to two different Žaccident scenarios in the first one,together with the aforementioned ‘hole’and ‘pipe’w x .models 5.In the following paragraphs the results are presented and discussed.Case A:Partial or complete breaking of a pipe in which natural gas is flowing,with the following conditions:Ž.pressure inside the pipe at the initial point ,P s 5bar abs;1temperature of the gas inside the pipe,T s 288K;1i.d.of the pipe,D s 163.6mm;distance between the point at which the pipe is damaged and the initial point,L s 1000m;molecular weight of natural gas,M s 17.4kg r kmol;density of natural gas at Point 1,r s 3.68kg r m 3;1viscosity of natural gas,m s 1.01=10y 5kg r ms;1Ž.Ž.fittings:it is supposed that there is a 908T K s 0.5every 75m i.e.13T’s ;Ž.it is supposed that the regulator located at the feeding point 1has no limit for the flow-rate,P thus always being constant.1()H.Montiel et al.r Journal of Hazardous Materials 591998211–233221T a b l e 1S e t o f e q u a t i o n s t o b e s o l v e d i n s t e a d y s t a t e a n d u n s t e a d y s t a t eŽ.C a s eS u b s o n i c f l o w i n t h e p i p e ,S u b s o n i c f l o w i n t h e p i p e ,S o n i c f l o w i n t h e p i p e P o i n t 2,abas o n i c f l o w a t t h e h o l e s u b s o n i c f l o w a t t h e h o l e s o n i c f l o w a t t h e h o l e S t e a d y s t a t e ,U n k n o w n v a r i a b l e s :u ,u ,P ,T U n k n o w n v a r i a b l e s :u ,u ,P ,T U n k n o w n v a r i a b l e s :u ,Q 122212221mŽ.Ž.w i t h o u t a n y e q u a t i o n s :E q u a t i o n s :E q .11E q u a t i o n s :E q .26Ž.l i m i t a t i o n o n t h e E q .30f l o w -r a t e o r Q -Q m m a xk q 1A M 2k M o r k y 1Ž.P k s M a P 9a (w x '222A R T k q 1R T c222k q 1k M k M A M 2k P P k M o r a a k kŽ.Ž.M a P s M a P 9b P y s M a P 24a '''Ž.Ž.1122222R T R T A R T k y 1P P R T (12c 2222Ž.Ž.E q .11E q .24b Ž.Ž.E q .15E q .25S t e a d y s t a t e ,U n k n o w n v a r i a b l e s :P ,u ,P ,T U n k n o w n v a r i a b l e s :P ,u ,P ,T U n k n o w n v a r i a b l e s :P ,Q 122212221mw i t h a l i m i t a t i o n o n t h e f l o w -r a t e Ž.Q G Q m m a x E q u a t i o n s :T h e s a m e s e t a s i n t h e p r e v i o u s c a s e E q u a t i o n s :T h e s a m e s e t a s i n t h e p r e v i o u s c a s eE q u a t i o n s :T h e s a m e s e t a s i n t h e p r e v i o u s c a s eU n s t e a d y s t a t eU n k n o w n v a r i a b l e s :P ,T ,u ,u ,P ,T U n k n o w n v a r i a b l e s :P ,T ,u ,u ,P ,T U n k n o w n v a r i a b l e s :P ,T ,u ,Q 111222111222111mE q u a t i o n s :T h e s a m e s e t a s i n t h e E q u a t i o n s :T h e s a m e s e t a s i n t h e E q u a t i o n s :T h e s a m e s e t a s i n t h e Ž.Ž.Ž.Ž.Ž.Ž.p r e v i o u s c a s e p l u s E q s .40a n d 41p r e v i o u s c a s e p l u s E q s .46a n d 48p r e v i o u s c a s e p l u s E q s .40a n d 41aŽ.I n t h e s e c a s e s ,i n u n s t e a d y s t a t e ,E q .35m u s t b e u s e d t o c a l c u l a t e t .c bŽ.Ž.Ž.Ž.I n t h i s c a s e ,f o r f u l l b o r e h o l e P i s k n o w n P s P ,E q .1m u s t b e u s e d i n s t e a d o f E q .25a n d E q .24a i s n o t r e q u i r e d .22a()222H.Montiel et al.r Journal of Hazardous Materials591998211–233The value of the release flow-rate,as calculated according to the three models,has been plotted as a function of the hole diameter in Fig.3.As can be observed,for small values of hole diameter,the new model gives results very close to those obtained with the‘hole model’;this is due to the existence of very small friction loss and a practically negligible expansion of the gas inside the pipe.However,for large diameter holes,the values given by this model approach those of the‘pipe model’.In the case of complete Ž.breaking of the pipe diameter of the hole s diameter of the pipe,Fig.3shows that the new model gives the same value than would be obtained from the‘pipe model’.Over the whole range of possible hole diameters,the new model gives results which lay between those of the‘pipe’and‘hole’models and which are much closer to the real values of release flow-rate,as it takes into consideration the change in the flow inside the pipe,and the variation in pressure and temperature.Figs.4and5show the prediction of the pressure and temperature of the gas at the point of release,as a function of hole diameter.The different trends existing for the zones corresponding to sonic and subsonic flow can be observed.Case B:In practice,the feeding points of a natural gas distribution net are controlledŽ.by flow-rate regulation valves see the Appendix,which impose a limit on theFig.3.Variation in release flow-rate as a function of hole diameter,according to hole model,pipe model and Ž.the new hole–pipe model.()H.Montiel et al.r Journal of Hazardous Materials591998211–233223Fig.4.Prediction of the pressure of the gas at the point of release,according to the pipe–hole model.maximum possible flow-rate in the network.Therefore,the release flow-rate cannot be higher than this maximum value;this means that if the release flow-rate reaches this maximum value,the pressure at the feeding point will no longer be constant and theŽpressure throughout the pipe will decrease although,as in the previous example,the.Ž. situation will correspond to a steady-state condition.In this case Q G Q,the setm maxof equations to be solved is the same than for the previous case for all the possiblesituations,but with a change in the boundary conditions.Now the value of u is known1Ž.and the value of P is unknown Table1.1The conditions will now be assumed to be the same as in Example1,but with the3Ž.flow-rate limited to a maximum value of9130Nm r h P s9bar abs by a TartariniEseries FL-BP DN50regulator.Fig.6shows the variation in the release flow-rate as a function of the diameter of thehole;the flow-rate reaches a maximum value of9130Nm3r h at a d of approximatelyorŽ65mm.This implies that from this value of d the pressure at Point1the feedingor.Ž.point will decrease if the hole diameter is larger Fig.7.It can be observed also in Fig.Ž7that the limitation in the flow-rate gives rise to a change in the flow at the hole from .sonic to subsonic.This change occurs at a hole diameter smaller than that correspond-()224H.Montiel et al.r Journal of Hazardous Materials591998211–233Fig.5.Prediction of the temperature of the gas at the point of release,according to the pipe–hole model.ing to the case in which the flow-rate has not reached the limiting value or there is noregulator.This is why the temperature at Point3and T reach the same value at aaŽ.smaller hole diameter Fig.8;these two temperatures and the temperature at Point2 have also the same value at a smaller hole diameter.Fig.8shows also that at large hole Ž.diameters but smaller than pipe diameter these three temperatures have already the same value.The variation of gas density at the different points as a function of the hole diameterŽ. has been plotted in Fig.9.As can be observed,for large but smaller than pipe diameter holes the densities at Point2and Point3are the same than that corresponding at atmospheric conditions.5.Unsteady stateThe above analysis is related to a stationary state,i.e.the release flow-rate has a constant value as a function of time.。

a r X i v :a s t r o -p h /9812370v 1 20 D e c 1998Constraining the IMF using TeV gamma rayabsorptionJ.S.Bullock a R.S.Somerville b D.MacMinn cand J.R.Primack aa PhysicsDepartment,University of California,Santa Cruz b RacahInstitute of Physics,Hebrew University,Jerusalemc Deceased 1Calculating the EBL and constraining the IMFSemi-analytic merging-tree (SAM)models of galaxy formation incorporate parameterized treatments of astrophysical processes such as gas cooling,star formation,supernovae feedback and dust absorption within the hierarchical structure formation ing the SAM models developed in Ref.[6],we can efficiently model the origin of the EBL in a variety of cosmological scenarios in a physical way.We can then use observations to contrain the nature of the IMF and the effects of dust (see [4,5]and [7],hereafter PBSM,for details).Fig.1.The attenuation factor forΛCDM models using a Salpeter(solid)and Scalo(dashed)IMF as a function of gamma ray energy Eγ,whereτ(Eγ)is the optical depth of the universe.The Salpeter IMF produces more high mass stars,thus more ultraviolet light and more∼100µm light reradiated by dust.For distant sources (redshift z s∼1)the Salpeter model’s increased UV light causes a noticeably larger attenuation of gamma rays.For nearby sources,the increase in reradiated light rel-ative to Scalo implies significantly increased gamma ray attenuation at Eγ>∼10TeV. In this contribution we focus on a model set in aflatΛCDM universe(Ωm= 0.4,ΩΛ=0.6),and investigate how the nature of the IMF affects the expected gamma-ray attenuation.We present results using two commonly used formsfor the IMF,Scalo[8]and Salpeter[9].Our results will primarily constrain the ratio of high mass to low mass stars produced.This quantity is of general interest in the context of understanding supernovae rates,metal production,and high redshift galaxies which are typically identified in the far UV.Note that there is some degeneracy between the effects of the IMF and the wave-length dependence of the dust extinction curve(see PBSM for a description of our treatment of dust).Figure1shows the gamma ray attenuation factor,exp(−τ),from ourΛCDM models using a Salpeter(solid)and Scalo(dashed)IMF.The optical depth isa function of gamma ray energy,τ(Eγ),and is most strongly influenced by the EBL at wavelengthsλEBL∼(Eγ/TeV)µm,where the cross section for pair production is maximized.The numbers next to each pair of curves indicatethe redshift of the source,z s.As discussed in detail in PBSM(seefigures2and4),the Salpeter IMF pro-duces more high mass stars,thus more ultraviolet light than does the ScaloIMF,and therefore a larger optical depth to<∼1TeV gamma rays.The increased UV light,in addition,produces more∼100µm light reradiated by dust,anda larger optical depth to>∼10TeV photons.Distant sources(z s>∼0.5)should provide an interesting probe of the IMF using<∼1TeV gamma ray telescopes. For nearby sources,the excess reradiated light relative to Scalo implies excess gamma ray attenuation at Eγ>∼10TeV,the relevant range for ground-based air-shower detectors.2ConclusionsGamma ray attenuation at Eγ<∼1TeV and>∼10TeV is significantly affected by the IMF,specifically the ratio of high to low mass stars,although there is some degeneracy associated with uncertainties in the modeling of dust absorption and reradiation.Observations of gamma ray absorption below∼1TeV for high redshift sources(z s>∼0.5)and above∼10TeV for nearby sources(z s∼0.03) will provide a useful probe of the nature of the IMF and galaxy formation.For more details,includingfitting functions of the optical depth of the universe as a function of redshift,cosmology,IMF,and gamma ray energy,see[4].References[1]R.J.Gould,G.P.Schr´e der,Phys.Rev.155(1967)1404.[2] F.W.Stecker,O.C.DeJager,M.H.Salamon,ApJ390(1992)L49;O.C.DeJager,F.W.Stecker,M.H.Salamon,Nature369(1994)294.[3] D.MacMinn,J.R.Primack Space Science Reviews75(1996)q413.[4]J.S.Bullock,R.S.Somerville, D.MacMinn,J.R.Primack in preparation(1998).[5]R.S.Somerville,J.S.Bullock,J.R.Primack,in preparation(1999).[6]R.S.Somerville,J.R.Primack MNRAS,accepted(1998).[7]J.R.Primack,J.S.Bullock,R.S.Somerville,D.MacMinn,these proceedings.[8]J.M.Scalo,Fundam.Cosmic Phys.11(1986)1.[9] E.E.Salpeter,ApJ121(1955)61.。