Geant4单位系统全解析:从G4double赋值到G4BestUnit输出 某次帮同事排查一个模拟结果对不上的问题前后折腾了大半天。现象是探测器灵敏体积里沉积的能量比预期高了整整两个数量级几何、材料、物理列表都翻了个遍最后发现原因简单到让人想砸键盘几何尺寸的构造函数里直接填了个1.5作者脑子里想的是 1.5 cm而 Geant4 的单位系统把它当成 1.5 mm 处理探测器的厚度差了十倍立体角、立体体积、能量沉积一路全错。这个坑不算新鲜几乎每个刚开始写 Geant4 应用的人都会撞上一次撞完之后才会真正理解这个框架为什么要在所有接口上都强制铺一层单位约定。这篇就把 Geant4 的单位系统从头到尾拆一遍内部单位制为什么这么设计、核心单位常量长什么样、代码里从赋值到输出该怎么写、宏命令里的单位字符串怎么被解析、自定义单位怎么注册以及我这几年踩过的那些坑。不管你是刚跑通第一个 example 的新手还是已经在维护一套完整模拟程序的老手里面应该都有能直接拿去用的东西。1. 内部单位制到底在解决什么问题1.1 数字本身不带单位这是所有麻烦的根源Geant4 的绝大部分接口参数类型都是G4double。G4double是什么本质就是一个double一个纯粹的浮点数它自己不知道自己是长度、能量还是时间。你调用G4Box(det, 10, 10, 10)的时候框架收到的是三个 10它唯一能做的判断就是“这三个数得按我约定的标尺来理解”。也就是说单位系统在 Geant4 里不是一个类型系统不是编译器能帮你检查的东西它是一份全体开发者共同遵守的约定一份写死在头文件里的换算表。理解这一点很重要因为它决定了很多事情的做法。既然是约定而不是类型那就意味着编译器不会在你忘记乘单位时报错也不会在你乘错单位时警告。所有的错误都会在物理结果上体现出来而且往往体现得非常隐蔽几何尺寸错了会改变立体角密度错了会改变阻止本领能量阈值错了会改变次级粒子产生率最后汇总成一个“结果好像不太对”的模糊印象。所以唯一的防线就是纪律每一个表示物理量的字面量后面都必须跟着单位常量。这套约定的好处也很直接。框架内部所有的算法、所有的常量、所有的数据表都在同一个标尺下运行不需要在函数边界上做任何换算。你传进去的是内部单位的数值框架内部算完吐出来的还是内部单位的数值中间没有转换损耗也没有转换开销。代价就是你必须记住这套标尺是什么并且在所有对外交互的边界上读入用户参数、输出结果、和外部数据比对做显式转换。1.2 为什么标尺选 mm、ns、MeV 这一套Geant4 的内部基本单位是长度用毫米时间用纳秒能量用兆电子伏电荷用基本电荷角度用弧度温度用开尔文。这个组合不是随便挑的背后有几层考虑。第一层是数值范围。双精度浮点数的有效范围大约到 1e±308看着很宽但真正麻烦的是精度损失和可读性。如果用国际单位制能量得用焦耳1 MeV 等于 1.602e-13 J那么一个 1 GeV 的粒子能量写成焦耳就是 1.602e-10 J。数值本身没问题但你打开一个日志文件满屏都是 1e-10、1e-13 这种量级人眼根本没法快速判断哪个数大哪个数小。换成 MeV1 GeV 就是 1000一眼就能看出量级关系。探测器尺寸同理用 m 的话一个 5 mm 的硅层要写 5e-3用 mm 就是 5。第二层是物理常量的“顺手程度”这也是我个人认为最关键的一点。在 mm-ns-MeV 这套单位下光速恰好是 299.792458 mm/ns约等于 300。这个数字在心算估计的时候太好用了一个粒子飞行 300 mm 大约需要 1 ns飞行 1 m 需要 3.3 ns。如果换成国际单位制光速是 2.998e8 m/s你每次做时间-of-flight 的粗略估计都要在脑子里处理指数。同样约化普朗克常量在 MeV·ns 下是 6.582e-22在 MeV·mm 体系下配合光速能直接算出 197.327 MeV·fm也就是那条经典的 197 公式。第三层是和外部数据源的对接。粒子物理领域的数据表、截面数据、评价核数据库能量单位普遍是 eV 系列长度普遍用 cm 或 mm时间用 ns 或 s。内部单位选成 MeV 和 mm等于把最常见的两种外部单位之一直接做成了内部单位转录时的换算次数最少出错概率也最低。用生活里的类比说这就像装修时统一用毫米做图纸单位。你不会在图纸上同时出现“3 米 5 厘米”和“0.035 百米”所有数字都是同一个量纲下的纯数工人拿到任何一个数字都知道该怎么量。Geant4 做的事情一模一样只不过它的“图纸”覆盖了整个探测器的几何、材料和物理过程。1.3 G4 与 CLHEP 的分水岭版本差异必须搞清楚这里有一个历史包袱必须讲清楚因为它直接决定你会不会遇到一堆莫名其妙的编译错误。早期版本的 Geant4 并不自己维护单位常量而是直接使用 CLHEP 库提供的SystemOfUnits.h那些常量都放在CLHEP命名空间里。所以那个年代的老代码里你经常能看到文件开头一行using namespace CLHEP;然后下面所有单位直接裸写。从 Geant4 10.0 开始单位和物理常量都被搬进了 Geant4 自己的代码库里对应的头文件是G4SystemOfUnits.hh和G4PhysicalConstants.hh常量定义在全局作用域不再需要using namespace CLHEP。这个变化带来的直接影响是如果你的项目里同时存在新旧两套写法比如某个陈旧的工具类还留着using namespace CLHEP;而你新写的文件用的是 G4 版本的头文件那就可能出现同名常量的二义性编译器会告诉你mm有歧义。更麻烦的情况是两个头文件都被间接包含某些编译器下会报重定义。判断自己用的哪一套很简单打开对应的头文件看一眼如果里面写的是static constexpr G4double millimeter 1.;这样的形式说明是 G4 自维护版本如果看到的是static const double millimeter 1.;而且外面裹着namespace CLHEP {那就是老版本。新建项目一律用 G4 版本的头文件老项目迁移时最省事的办法是把所有#include CLHEP/Units/SystemOfUnits.h替换掉然后逐个删掉using namespace CLHEP;重新编译一次看哪里报错就改哪里。2. 单位常量的组织方式与换算关系2.1 基本单位与派生单位的分层结构整套单位表是分层的。底层是一批“基本单位”它们在数值上直接定义为 1其他所有单位都是通过乘除这些基本单位派生出来的。在 Geant4 沿用的这套体系里被当作基本单位处理的包括毫米、纳秒、兆电子伏、基本电荷、弧度、开尔文、摩尔等。你在代码里写的millimeter、nanosecond、MeV、eplus、radian、kelvin、mole它们的值全都是 1.0。这个设计的巧妙之处在于派生单位的定义直接写成了物理定义式。比如电流的安培定义成“每秒通过一个基本电荷”代码里就是ampere eplus/nanosecond。磁场的特斯拉定义成“每平方米每秒一伏特”代码里就是tesla volt*second/meter2而伏特又是焦耳每库仑。你把这条链一路展开最后所有东西都落到那几个基本单位上系数由常量表达式在编译期算好。这正是为什么你写1.0*tesla的时候完全不用关心它背后展开成什么——编译器帮你算完了。派生链长了之后某些量的数值会变得相当反直觉。比如压强的帕斯卡在这个体系里要展开成能量密度也就是 MeV/mm³ 量级的东西磁场的单位虽然是特斯拉但它在底层展开后系数是一长串。这些都不需要你操心你只要记住一件事永远用符号名永远不要自己手算系数。我见过有人在代码里写3.0e-4*tesla来表示 3 高斯理由是“1 高斯等于 1e-4 特斯拉我算过了”。这种写法不算错但它把一个可以从符号名读出来的语义变成了一个魔数后面维护的人看到 3.0e-4 根本不知道这是磁场也不知道是高斯换算来的。写成3.0*gauss就一目了然。2.2 常用单位换算速查表下面这张表是我自己整理过很多遍的日常写代码时基本够用。表里的“相对值”指的是相对于内部基本单位的数值比如长度类的相对值是相对于毫米能量类相对于 MeV。这些数值都可以从头文件里直接查证我建议你自己打开G4SystemOfUnits.hh对一遍比记表更靠谱。类别常用写法相对内部单位的值长度mm/millimeter1长度cm10长度m/meter1000长度km1e6长度um/micrometer1e-3长度nm/nanometer1e-6长度angstrom1e-7长度fermi1e-12时间ns/nanosecond1时间ps1e-3时间us1e3时间ms1e6时间s/second1e9能量MeV1能量eV1e-6能量keV1e-3能量GeV1e3能量TeV1e6能量joule6.241509074e12动量MeV/c形式的MeV1质量g/gram5.60958865e26质量kg/kilogram5.60958865e29角度rad/radian1角度deg/degree0.017453292519943295角度mrad1e-3电荷eplus1电荷coulomb6.241509074e18磁场tesla1磁场gauss1e-4磁场kilogauss0.1关于质量和动量那两行有个细节值得单独说一下。在这套单位体系里光速被当作 1 来处理严格说是被吸收进了单位定义所以质量和能量在量纲上是同一回事动量和能量也是同一回事。这就是为什么质量的单位可以直接写成 MeV动量的单位写成 MeV/c 但c在数值上是 1所以两者在代码里其实就是同一个常量。物理上当然不是一回事但在写代码的时候你确实可以直接写G4double mass 938.27*MeV;框架内部就是这么处理的。密度那一行我没有放进去因为它的写法和别的不同。密度不是单个单位常量而是一个表达式g/cm3在代码里是gram除以cm3cm3是立方厘米对应的常量。所以你写2.7*g/cm3是合法的写成2700*kg/m3也是合法的两者在数值上完全相等。我个人偏好前者因为材料手册上的密度基本都是 g/cm³直接抄不用换算。2.3 物理常量头文件千万别漏掉单位常量只解决“量纲怎么标”的问题物理常量解决的是“这些带量纲的数具体是多少”的问题。这两个头文件必须一起用G4SystemOfUnits.hh和G4PhysicalConstants.hh。很多人只包含前者然后在代码里手写G4double c 3e8;这就是典型的错误来源。在内部单位制下光速不是 3e8而是 299.792458 mm/ns差了将近千分之一。这个误差在做飞行时间计算时会直接变成一个系统偏差而且因为它太“合理”了你很难在结果里发现。常用的物理常量我列几个在下面这些值都是按内部单位预先换算好的直接拿来用就行。常量名内部单位下的数值单位含义c_light299.792458mm/nsc_square89875.51787mm²/ns²h_Planck4.135667696e-21MeV·nshbar_Planck6.582119569e-22MeV·nseplus1基本电荷electron_mass_c20.51099895000MeVproton_mass_c2938.27208816MeVneutron_mass_c2939.56542052MeVk_Boltzmann8.617333262e-11MeV/KAvogadro6.02214076e231/molefine_structure_const0.0072973525693无量纲classic_electr_radius2.8179403262e-12mmelectron_Compton_length2.42631023867e-12mmBohr_radius5.29177210903e-8mm注意最后三行长度类常量的单位是毫米所以数值都带一个很大的负指数。这个不是错误是内部单位决定的。你在代码里打印这些量的时候用G4BestUnit会自动显示成 fm 或者 Å比裸看2.8e-12舒服得多。3. 代码里从赋值到输出的完整写法3.1 赋值环节让每个数字后面都长着单位这是整个单位系统里最基础也最容易被忽视的一环。规则只有一条任何一个表示物理量的字面量都要乘上对应的单位常量。这条规则听起来简单但在实际项目里能严格执行的人并不多原因是有时候数字“看起来”已经是内部单位了乘不乘都一样比如G4double length 10*mm;和G4double length 10;在数值上完全相同。但这两者绝不是等价的。前者的可读性、可维护性、可审查性都远超后者。半年后有人改成厘米只需要把mm改成cm数值自动跟着变而后者得去查文档确认这个 10 到底是毫米还是厘米。更现实的问题是 review 的时候你看到一堆裸数字根本没法判断对错。我现在的习惯是如果某个文件里出现了表示物理量的裸字面量我会直接标出来要求改不管它当前的值对不对。#include G4SystemOfUnits.hh #include G4PhysicalConstants.hh // 几何尺寸全部显式带单位 G4double halfX 5.0*cm; G4double halfY 2.5*cm; G4double halfZ 100.0*um; G4Box* sensitiveLayer new G4Box(SensLayer, halfX, halfY, halfZ); // 材料密度g/cm3 是材料手册的常用单位 G4double density 2.3290*g/cm3; G4Material* silicon new G4Material(Silicon, 14., 28.0855*g/mole, density); // 能量阈值keV 是人眼最容易核对刻度的量级 G4double cutEnergy 1.0*keV; // 磁场别手算特斯拉 G4double fieldValue 1.5*tesla;这段代码里有一个地方特别值得说28.0855*g/mole。摩尔质量这个量的量纲是质量除以物质的量在内部单位下它是克每摩尔。写的时候一定要把g和mole都带上别图省事写一个裸的 28.0855因为那样写之后某天有人把材料换成铅207.2就会下意识写 207.2数值本身没错但可读性就没了。还有一个非常隐蔽的坑整数除法。1/2*cm的结果是 0不是 0.5 cm。因为1/2在 C 里先做整数除法得到 0再乘cm还是 0。这个错误在几何里表现为整个固体退化成一个点或者直接抛异常在能量里表现为阈值归零导致计算爆炸。规避方法很简单要么写0.5*cm要么写1.0/2*cm总之保证除法的某一侧是浮点。我自己写的时候习惯把所有系数都写成带小数点的形式1.0、2.0、0.5从源头上杜绝这个问题。3.2 输出环节G4BestUnit 与手动换算代码内部用内部单位没问题但输出给人看的时候就不能这样了。谁也不想在日志里看到能量 0.000661657然后自己去数小数点。Geant4 提供了G4BestUnit它会把数值自动换算到该类量最合适的量级上并附上单位符号。#include G4UnitsTable.hh G4double energy 661.657*keV; G4double length 5.0*cm; G4double mass 0.51099895000*MeV; G4cout G4BestUnit(energy, Energy) G4endl; // 输出 661.657 keV G4cout G4BestUnit(length, Length) G4endl; // 输出 5 cm G4cout G4BestUnit(mass, Mass) G4endl; // 输出 510.999 keVG4BestUnit的第二个参数是“类别名”是一个字符串必须和单位表里注册的类别完全一致。常见的类别名包括Length、Time、Energy、Mass、Angle、Electric Charge、Magnetic flux density、Density、Pressure、Velocity等等。这个字符串写错了不会崩溃但会退化成打印原始数值很容易被忽略。我一般的做法是在项目里定义几个常量字符串需要的地方直接引用避免手打出错。它还有一个针对三维向量的重载接受G4ThreeVector会把三个分量都换算一遍输出成类似(1,2,3) cm的形式。做几何调试的时候特别方便直接打印一个位置向量就能看出量级对不对。如果你不想要自动换算想要固定单位输出那就手动除G4cout thickness halfZ/mm mm G4endl; G4cout density density/(g/cm3) g/cm3 G4endl; G4cout field fieldValue/tesla T G4endl;手动除的好处是单位固定方便后续用脚本解析日志。做批量扫描的时候我倾向用这种方式因为自动换算是根据数值大小动态选单位的同一列数据可能这行是 keV 那行是 MeV脚本解析起来要额外处理。固定单位就没有这个问题代价是可能看到 0.000661657 这样的数字但那在脚本里完全无所谓。3.3 宏命令里的单位字符串是怎么被解析的Geant4 的命令交互体系对单位有一层专门的支持。如果你自己写过自定义命令应该知道G4UIcommand那一套。带单位的数值参数有专门的类G4UIcmdWithADoubleAndUnit用于单个带单位的浮点数G4UIcmdWith3VectorAndUnit用于三维向量。这两个类相比普通版本多做了两件事一是告诉用户在宏文件里可以带单位写二是在取值时自动把字符串解析并换算成内部单位。// 在命令类的构造函数里 fSetThicknessCmd new G4UIcmdWithADoubleAndUnit(/det/setThickness, this); fSetThicknessCmd-SetGuidance(Set the thickness of the sensitive layer); fSetThicknessCmd-SetParameterName(Thickness, false); fSetThicknessCmd-SetUnitCategory(Length); fSetThicknessCmd-SetDefaultUnit(mm); fSetThicknessCmd-SetRange(Thickness0.); fSetThicknessCmd-AvailableForStates(G4State_PreInit, G4State_Idle); // 在 SetNewValue 里 void MyDetectorMessenger::SetNewValue(G4UIcommand* command, G4String newValue) { if (command fSetThicknessCmd) { fDetector-SetThickness(fSetThicknessCmd-GetNewDoubleValue(newValue)); } }关键在GetNewDoubleValue。它接收宏文件里原始输入的那个字符串比如5 mm或者0.5 cm或者干脆2.5然后把它解析成一个内部单位的G4double。如果用户写了单位就按写的单位换算如果没写单位就按SetDefaultUnit里设的默认单位处理。这两个设置是配套的SetUnitCategory决定了候选单位列表用户在图形界面里能看到一个下拉框SetDefaultUnit决定了裸数字怎么解释。这里有个我踩过的坑值得提醒SetDefaultUnit如果忘了设裸数字会被当成内部基本单位处理。也就是你写/det/setThickness 5用户心里想的是 5 mm实际得到的也是 5 mm看起来没问题但只要有人改成/det/setThickness 0.5想表达 0.5 mm还是没问题。真正的麻烦在于这个命令和别的命令混用的时候有的命令默认 cm 有的默认 mm用户根本记不住。统一设置默认单位并且在SetGuidance里写清楚是唯一可靠的做法。SetRange那一行也值得提一句它做的是取值范围的检查但检查的是换算之后的值。所以如果用户写-1 mm会被拦下来写0.001 m也是 1 mm 能通过。这个机制把单位换算和参数校验串在了一起用起来很省心。3.4 自定义单位注册一次全局可用框架自带的单位覆盖了绝大多数场景但总有需要补充的时候。比如你在做磁场相关的模拟习惯用高斯而不是特斯拉或者你在处理剂量习惯用戈瑞。这时候可以用G4UnitDefinition注册新单位。#include G4UnitsTable.hh int main(int argc, char** argv) { // 注册必须在创建 RunManager 之前或者至少在第一次使用该单位之前 new G4UnitDefinition(kilogauss, kG, Magnetic flux density, 0.1*tesla); new G4UnitDefinition(milligauss, mG, Magnetic flux density, 1e-3*gauss); new G4UnitDefinition(mrad_dose, mrad, Dose, 1e-5*gray); // ... }构造函数有四个参数完整名称、符号、类别名、相对于基本单位的数值。注册之后这个单位就出现在单位表里了G4BestUnit在挑“最合适单位”的时候会把它纳入候选宏命令的SetUnitCategory也会把它列进下拉框。也就是说注册一次输入输出两边都自动生效。有一个细节需要注意符号名不要和已有单位冲突。比如上面例子里我用了mrad表示毫拉德剂量但mrad在角度类别里已经存在毫弧度。虽然类别不同、内部是两个独立条目但在某些打印场景下可能引起混淆。更稳妥的做法是用一个不冲突的符号比如mrd。另外G4UnitDefinition是 new 出来不 delete 的这在框架里是惯例因为它的生命周期和程序一致不算内存泄漏。调试单位表本身也有工具G4UnitDefinition::PrintUnitsTable()会在标准输出里打一张完整的表包含所有已注册的单位、符号、类别和数值。我在接手别人的项目时经常先加这么一行跑一次看看有没有人偷偷注册了什么奇怪的单位也顺便确认版本兼容性。4. 常见坑与排查实录4.1 症状速查从结果反推单位错误单位类的错误有个共同特征它不会让程序崩溃只会让结果“看起来不太对”。所以排查的第一步是建立一个症状到原因的映射。下面这张表是我这些年积累下来的命中率相当高。观察到的症状最可能的原因快速验证方法几何体积、重量与设计值差 10 或 1000 倍尺寸被当成 mm但心里想的是 cm 或 m打印对应的G4VSolid::GetCubicVolume()再除以cm3材料密度导致射程整体偏移密度用了kg/m3但没换算或漏了g/cm3的除法打印material-GetDensity()/(g/cm3)能量沉积谱整体平移或缩放能量阈值、入射能量、事例产生子里的能量单位不一致逐个打印所有能量相关的G4double除以MeV飞行时间结果差了千分之一量级手写了3e8当光速搜代码库里的3e8、2.998e8这类魔数次级粒子几乎不产生产生阈值设得过大比如把1*keV写成1*MeV打印G4ProductionCuts的当前值磁场偏转角度完全不对磁场值用了高斯但没乘gauss或者单位类别写错打印fieldValue/tesla看量级是否合理采样结果出现大量 0 或者异常某处整数除法把系数算成了 0搜形如a/b*unit且 a、b 都是整数的写法用这张表的时候有个技巧不要直接去读代码找错误而是先把所有关键量打印出来除以期望单位再和设计值对照。人眼对“这个数应该是 2.33 但我看到 2330”这种量级差异极其敏感比逐行读代码快得多。4.2 命名空间冲突与单字母常量的遮蔽这是 Geant4 单位系统里最有“特色”的一类坑。前面说过单位常量里有大量单字母和短名字m米、s秒、g克、cm、mm、ns、eV。它们全是全局作用域里的常量。这意味着你在任何地方声明一个同名的局部变量都会把这层全局常量遮蔽掉。举个最常见的例子写循环的时候习惯用i、j、k偶尔会用m、n表示维度。如果你在函数里写了G4int m 3;那这个函数体里所有的m指的都是这个整数你后面再写2.0*m想表示 2 米得到的其实是 6。编译器不会报错因为这个整数能隐式转换成G4double结果就是把一个长度算成了 6 mm。更阴险的是g。材料定义里经常要遍历有人会写for (G4int g 0; ...)然后循环体里恰好有一句density 2.7*g/cm3本意是 2.7 克每立方厘米实际算出来是 2.7 乘循环变量再除以 cm³。这种 bug 的表现是密度随循环次数变化但你要是不打印出来根本发现不了。规避这类问题只有一条路不要用单字母做变量名。这不是风格问题是安全问题。维度用nDim循环变量用idx或iLoop反正不要让任何局部变量和单位常量重名。同样不要写using namespace CLHEP;因为那会把几十个短名字一次性灌进全局作用域遮蔽风险成倍增加。还有一类冲突来自using namespace std;和单位常量的共存。比如s是秒std::string里有s吗没有。但std::chrono里有std::chrono::s如果你同时打开了using namespace std::chrono那s就有歧义了。这个组合不常见但确实存在尤其是有人在里面混用了计时工具。4.3 分析输出里的单位陷阱用G4AnalysisManager往 ntuple 或者 CSV 里填数据的时候单位问题会以一种更隐蔽的方式出现填充的值全是内部单位但列名里没写单位或者写了错误的单位。比如你填一个能量列代码是analysisManager-FillNtupleDColumn(0, energy);这里energy是内部单位下的G4double也就是以 MeV 为单位的数值。如果你的分析脚本默认按 keV 处理那所有结果都会差 1000 倍而且脚本作者看到数值范围是 0 到 2会以为这是 MeV 的量级根本不会怀疑。我的做法是在填充之前就显式换算把这层转换固定在代码里而不是留给分析脚本// 明确以 keV 为单位输出 analysisManager-FillNtupleDColumn(0, energy/keV); analysisManager-FillNtupleDColumn(1, length/cm); analysisManager-FillNtupleDColumn(2, time/ns);然后在列名或者 ntuple 的说明里把单位写死Edep [keV]、TrackLength [cm]、Time [ns]。这样即使换个脚本作者看一眼列名也知道该怎么处理。张量类型的列更容易出问题比如FillNtupleDColumn填充G4ThreeVector的三个分量时三个分量必须用同一个单位换算不能这一个除mm那一个除cm。还有一个细节G4BestUnit在正式的数据分析输出里最好别用。它的输出是给人看的格式化字符串解析起来麻烦而且同一列可能混杂不同单位的字符串脚本处理起来要写额外的正则。调试日志里用它很舒服正式输出还是老老实实固定单位。4.4 我自己踩过的几个具体坑说几个真实的例子都是我自己写错过的比抽象地讲原理更有说服力。第一个是双重换算。我写了一个工具函数签名是G4double ComputeRange(G4double energy)约定传进来的是内部单位。调用方传了10*MeV进去函数内部第一行写了energy energy*MeV;理由是“保险起见再乘一次”。结果能量被放大了一百倍射程算出来完全不对。这个 bug 找了两天才定位到因为函数本身逻辑没问题调用方也没问题错在函数内部那句“保险起见”。教训是内部单位只在边界上做一次转换函数内部绝不再乘。如果实在不放心就在函数注释里写清楚参数的单位约定而不是靠代码里再乘一次。第二个是宏文件里的默认单位。有个命令是设置探测器偏移量的我设了SetUnitCategory(Length)但忘了设SetDefaultUnit。宏文件里写的是/det/setOffset 0.5我以为这是 0.5 mm实际也是 0.5 mm因为默认就是内部单位看起来没问题。后来有人改成0.05想表达更小的偏移还是没问题。真正的坑在于另一个命令那个我设了SetDefaultUnit(cm)于是同一份宏文件里两个命令的裸数字含义不一样用户写了一串裸数字整个几何就歪了。现在的规矩是所有带单位的命令都必须显式设默认单位并且在 guidance 里写出来。第三个是G4BestUnit的类别字符串拼写。我写过G4BestUnit(field, Magnetic Field)少了flux density几个字程序不报错直接打印出原始数值1.5没有单位。我在日志里看到这个1.5还以为磁场单位丢了查了半天单位注册最后发现是类别名对不上。正确写法是Magnetic flux density。这类错误只能靠对照PrintUnitsTable()的输出确认没有别的好办法。第四个是版本迁移。把一个老项目从旧版本升上来代码里保留了using namespace CLHEP;同时又包含了新的G4SystemOfUnits.hh。在某些编译配置下能过某些配置下报歧义。当时的做法是逐个文件把using namespace CLHEP;删掉遇到报错的地方再看是哪个常量找不到从CLHEP::前缀改成裸名。整个过程花了大半天但一劳永逸。5. 一套可以直接抄的工程化做法5.1 在项目里把单位约定写成规矩单位系统能不能用对很大程度上不取决于你懂不懂原理而取决于项目里有没有形成硬性规矩。我在自己的项目里固定了几条写进了代码规范文档新加入的人第一天就要看。第一条所有物理量必须显式带单位常量任何裸的数字字面量都视为缺陷review 时直接打回。例外只有一种情况这个数字本身是无量纲的系数比如效率、分支比、角度正弦值。第二条禁止在任何代码里出现物理常量的魔数包括但不限于光速、普朗克常量、玻尔兹曼常量、电子质量。一律从G4PhysicalConstants.hh取。这条规矩执行起来最容易因为搜一下3e8、1.6e-19、6.6e-34这类数字就能查出来。第三条所有输出到日志或者数据文件的值必须在同一行或列名里写明单位。日志里写thickness 5 mm不要写thickness 5。数据文件里的列名统一用方括号标注单位。第四条几何构造函数、材料构造函数这类“单位敏感”的接口在调用处必须能一眼看出每个参数的单位。做法是参数名里带单位后缀比如halfX_mm或者把字面量写成5.0*cm这种自解释形式。我个人偏好后者因为参数名带后缀在重构时容易忘改。第五条禁止using namespace CLHEP;。这条是硬性的没有例外。这些规矩听起来琐碎但实际效果非常明显。我统计过自己维护的一个中等规模项目引入这套规矩之后与单位相关的调试时间从每月几个小时降到了接近零。因为大部分单位错误都是在 review 阶段被发现的根本没机会进入运行阶段。5.2 一个自检函数启动时打印关键量除了规矩还可以加一层运行时自检。我在每个项目的初始化阶段加了一个PrintUnitSanityCheck()函数把几个关键量打印出来配合设计值人工对一遍。这个函数不到三十行但救过我好几次。#include G4SystemOfUnits.hh #include G4PhysicalConstants.hh #include G4UnitsTable.hh void PrintUnitSanityCheck() { G4cout Unit sanity check G4endl; G4cout c_light c_light/(mm/ns) mm/ns G4endl; G4cout 1 cm cm/mm mm G4endl; G4cout 1 MeV MeV/keV keV G4endl; G4cout 1 GeV GeV/MeV MeV G4endl; G4cout electron m electron_mass_c2/MeV MeV G4endl; // 把探测器关键尺寸和材料密度也打出来和设计值对照 G4cout det halfX GetDetectorHalfX()/mm mm G4endl; G4cout det density GetDetectorDensity()/(g/cm3) g/cm3 G4endl; G4cout prod cut GetProductionCut()/keV keV G4endl; // 需要时打印完整单位表调试版本迁移时特别有用 // G4UnitDefinition::PrintUnitsTable(); G4cout G4endl; }这个函数的价值在于把“隐性假设”变成了“显性输出”。比如c_light那一行如果哪天有人不小心把物理常量头文件换成老版本的这一行会立刻暴露出数值差异。探测器尺寸和密度那两行更直接一旦有人改了构造函数里的单位启动日志里马上就能看出来。最后一行的PrintUnitsTable()平时注释掉在做版本升级或者怀疑单位表被污染的时候打开。它会打印出当前进程里所有已注册单位的完整列表包括你自定义的那些。我一般在接手新项目、升级 Geant4 版本、或者引入新的第三方模块之后会跑一次确认没有意外的单位注册。另外一个很实用的小习惯在几何和材料的构建代码里紧跟每一个new G4Box、new G4Material之后加一行G4cout打印关键参数用G4BestUnit格式化。初期调试阶段这行日志的价值极高因为它把“我写的”和“框架理解的”直接摆在一起对照任何单位误解都会在启动时暴露而不是等到跑完几万个事例之后才发现结果不对。真正让我把这套习惯坚持下来的是有一次跑一个需要十几个小时的批量模拟结果出来之后发现能量沉积谱整体偏了一个常数因子。回头查了两个小时最后定位到材料密度那里g/cm3写成了g/m3密度小了一百万倍射程长了几个数量级整个探测器的响应全变了。如果启动时打一行密度