calchim_mod.F90

路径

LMDZ.MARS\libf\aeronomars\calchim_mod.F90

所属目录/模块

libf\aeronomars

文件定位

calchim_mod.F90 定义 calchim_mod，是 Mars 物理包中把 advected tracers、温度、风场、太阳几何、尘/冰表面积和化学开关整理后送入光化学/热层化学例程的调度模块。核心例程 calchim 在 firstcall 阶段识别可参与化学的 tracer 索引，选择在线或离线光解初始化路径，并在每个物理网格列上把质量混合比转换为体积混合比、调用 photochemistry、perosat、可选 chemthermos/deposition，最后把更新后的 mole fraction 转回 dqchim 倾向。

本文件还保存 ichemistry、zdqchim、zdqschim 三个模块状态；ini_calchim_mod/end_calchim_mod 负责分配和释放化学倾向缓存。

定义的符号

依赖的模块

calchim 内部：

调用的关键例程

输入

输出

共享状态与副作用

• ichemistry、zdqchim、zdqschim 为模块级 SAVE 且 THREADPRIVATE；ini_calchim_mod 分配 zdqchim(ngrid,nlayer,nq) 和 zdqschim(ngrid,nq) 并清零，end_calchim_mod 释放。
• calchim 内部的 firstcall、depos、ionchem、deutchem、jonline、unichim、output、nbq、niq、chemthermod 以及 tracer 索引缓存均为 SAVE/THREADPRIVATE。firstcall 后不会重新扫描 tracer 列表。
• firstcall 默认开关为：unichim=.true.、jonline=.true.、ionchem=.false.、deutchem=.false.、depos=.false.、output=.true.。
• calchim 会在 firstcall 打印找到的 tracer 数量和 niq 索引；缺少必需 tracer 时多处调用 abort_physic。
• 若 photochem .and. jonline，调用 init_photolysis；若 photochem .and. .not.jonline，调用 read_phototable 并产生文件读取/广播副作用。
• 若 .not.unichim，firstcall 调 chemthermos_readini 读取运行目录下的热层反应速率文件；逐列还会调用 chemthermos。
• 若 photochem .and. output .and. ngrid>1，通过 write_output 和 wstats 写诊断/统计字段。

核心逻辑

1. firstcall 初始化开关与 tracer 列表：设定 unichim/jonline/ionchem/deutchem/depos/output 默认值，分配 niq(nq)，按 igcm_* 检查化学需要的 tracer。CO2、CO、O、O1D、O2、O3、H、H2、OH、HO2、H2O2、N2、H2O 等缺失会中止；CH4 缺失只打印并忽略。
2. 按 tracer 存在性开启扩展化学：若 HDO 存在则 deutchem=.true.，随后要求 OD、D、HD、DO2、HDO2 等氘相关物种；若 O2+ 存在则 ionchem=.true.，随后要求 CO2+、O+、NO+、CO+、C+、N2+、N+、H+、HCO2+、HCO+、H2O+、H3O+、OH+、e- 等离子化学物种。
3. 光解/速率初始化：若 photochem 为真，jonline=true 时调用 init_photolysis，且 deutchem 会先 nphot=nphot+2；jonline=false 时调用 read_phototable。若 .not.unichim，读取 chemthermos_reactionrates.def。
4. 逐列准备化学输入：对每个 ig,l，计算 zq=pq+pdq*ptimestep，把化学 tracer 从质量混合比换算为 mole fraction：zycol=zq*mmean/mmol；同时更新温度/风场、把气压转 mbar、数密度转 cm-3、高度转 km、表面积转 cm2/cm3，并通过 temp_elect(...,1) 得到电子温度。
5. 选择热层切换层：unichim=true 时 lswitch=nlayer+1；否则首次遇到 pplay(ig,l)<10. 的层作为 lswitch，用于分隔低层光化学和高层热层化学。
6. 执行化学分支：photochem 为真时设置反应数：基础二体反应 31、三体反应 6、quench/heterogeneous 9；ionchem 额外加 64 个二体反应并设 nphotion=18；deutchem 额外加 35 个二体反应；nb_phot_max=nphot+nphotion+nquench 后调用 photochemistry。
7. 后处理与可选过程：光化学后保存 jo3/jh2o/iter，调用 perosat 处理 H2O2 凝结；.not.unichim 时以 chemthermod=3 调 chemthermos；depos 为真时调用 deposition。
8. 回写倾向：每层用 maxloc(zycol(l,:)) 找到主量物种 iqmax；其它化学物种按 (new_mass_mixing_ratio-zq)/ptimestep 形成 dqchim，主量物种的 dqchim 设为这些倾向的负和。
9. 诊断输出：若 photochem .and. output .and. ngrid>1，输出光解率、迭代数、NO/O2 夜辉发射率，并写统计。

伪代码

calchim(...):
  if firstcall:
    set unichim=true, jonline=true, ionchem=false, deutchem=false, depos=false, output=true
    allocate niq(nq)
    scan igcm_* tracer indices; append required species into niq or abort if mandatory species missing
    if HDO exists: deutchem=true; require OD/D/HD/DO2/HDO2
    if O2+ exists: ionchem=true; require ion species and electron tracer
    if photochem:
      if jonline:
        if deutchem: nphot = nphot + 2
        call init_photolysis()
      else:
        call read_phototable()
    if not unichim: call chemthermos_readini()
    firstcall=false

  zero dqchim and dqschim
  for each grid column ig:
    for each layer l:
      for each chemistry tracer iq in niq:
        zq = pq + pdq*ptimestep
        zycol(l,iq) = zq * mmean / mmol(iq)
      ztemp = pt + pdt*ptimestep
      zpress = pplay/100
      zdens = zpress/(kb*1e4*ztemp)
      zlocal = zzlay/1000
      ztelec = temp_elect(zlocal,ztemp,1)
      convert dust/ice surface area from m2/m3 to cm2/cm3
      choose lswitch = nlayer+1 if unichim else first pplay < 10 Pa

    if photochem:
      set reaction counters from ionchem/deutchem/nphot/nquench
      call photochemistry(..., zycol, ...)
      save jo3/jh2o/iter
      call perosat(..., zycol, dqcloud, dqscloud)
      if not unichim: call chemthermos(..., zycol, ...)
      save emission fields

    if depos: call deposition(..., zycol, ...)

    for each layer:
      iqmax = index of max zycol(l,:)
      for each chemistry tracer iq except iqmax:
        dqchim(ig,l,iq) = (zycol(l,iq)*mmol(iq)/mmean(ig,l) - zq(ig,l,iq))/ptimestep
        dqchim(ig,l,iqmax) -= dqchim(ig,l,iq)

  if photochem and output and ngrid > 1:
    write_output jo3, jh2o, iter_chemistry, emission_no, emission_o2
    wstats jo3, emission_no, emission_o2, mmean

参与的主题流程

写法特点

• .F90 自由格式模块，但包含大量 SAVE/THREADPRIVATE 缓存；firstcall 后 tracer 索引和开关不会随运行中配置变化重置。
• 开关默认值在源码中硬编码，且 firstcall 前有一段较早的光解/热层初始化代码被注释；活跃初始化位于 tracer 检查之后（765-780 行）。
• 质量守恒采用“每层最大 mole fraction 物种补偿”的写法：其它化学物种显式计算倾向，主量物种取负和。
• temp_elect 第三个参数固定为 1，即使用 Viking 路径；源码注释说明 2 对应 MAVEN。
• .not.unichim 的切换压力阈值为 pplay < 10.（Pa），而传给化学例程的 zpress 已转为 mbar。
• depos 默认 .false.，且下游 deposition.F 当前也不实际修改 zycol。
• ig_vl1 的 Viking-1 诊断网格计算代码被注释；只有 depos=true 且实际调用 deposition 时才可能用到该变量，当前默认路径不会使用。
• 疑似源码细节：氘相关缺失检查中若对应 tracer 不存在，若干 else 分支仍检查 i_oplus/i_d/i_hd 等非本物种条件并报 HDO/HD/DO2/HDO2 缺失（449-520 行）；这里按源码事实记录，是否为复制粘贴错误待开发者确认。

复现要点

• 运行中需提供基础化学 tracer：CO2、CO、O、O1D、O2、O3、H、H2、OH、HO2、H2O2、N2、H2O；若 photochem=true，N、N2D、NO、NO2 缺失也会中止。
• CH4 缺失不会中止，只打印 “no ch4, ignored”。
• HDO tracer 存在会开启 deutchem，因此必须同步提供 OD、D、HD、DO2、HDO2 等氘相关 tracer。
• O2+ tracer 存在会开启 ionchem，因此必须同步提供所有离子化学 tracer 和电子 tracer。
• 默认 jonline=true，光解走 init_photolysis；若源码改为 jonline=false，运行需要可读取离线光解查找表。
• 默认 unichim=true，不会读取 chemthermos_reactionrates.def，也不会调用 chemthermos；若改为 .false.，运行目录需包含该速率文件。
• 默认 output=true，但诊断输出只在 photochem 为真且 ngrid>1 时发生。

待确认

• abort_physic 的宿主模块未在本例程 use 列表中显式声明；是否依赖旧式外部过程或间接可见性，待编译链路确认。
• ichemistry 在本文件声明但未赋值；其配置来源和调用方如何使用该频率字段需要结合 conf_phys/physiq 继续核验。
• 氘 tracer 缺失检查中若干 else 分支的条件与报错文字不完全匹配，疑似复制粘贴问题，待源码维护者确认。
• .not.unichim、jonline=false、depos=true 这三个路径默认都关闭；本页按源码路径记录，实际生产配置是否启用待运行配置确认。

相关页面

• aeronomars/index：aeronomars 目录总览，本文件位于“光化学”分类。
• chemistrydata：离线光解查找表读取方，供 calchim 在 jonline=false 时调用。
• chemthermos_readini：热层反应速率读取方，供 .not.unichim firstcall 调用。
• deposition：可选干沉降路径，受 depos 开关控制。
• iono_h：提供 temp_elect 电子温度函数。
• conc_mod：提供 mmean，参与质量混合比/体积混合比换算。
• phys_state_var_init_mod：初始化/释放阶段调用 ini_calchim_mod 和 end_calchim_mod。
• physiq_mod：物理主调度中与化学/热层相关的上游入口。
• photochemistry_mod：主光化学积分文件页；说明 photochemistry 的在线/离线光解、反应速率、隐式积分、jo3/jh2o/iter 和夜辉诊断生成。
• photolysis_mod：在线光解初始化和共享谱表状态；calchim firstcall 在 jonline=true 时调用 init_photolysis，并可能先对 nphot 加 2。
• photolysis_online.F：在线光解率计算文件页。
• chemthermos：热层化学文件页，说明 .not.unichim 路径中的数密度矩阵、jthermcalc_e107/paramfoto_compact 调用和 lswitch:nlayer 写回。
• perosat：H2O2 饱和/凝结文件页，说明 dqcloud/dqscloud 的写入逻辑。