calcstormfract_mod.F90

路径

LMDZ.MARS\libf\phymars\calcstormfract_mod.F90

所属目录/模块

libf\phymars

文件定位

该文件定义火箭式尘暴（rocket dust storm）亚格点面积分数计算模块。它只提供一个例程 calcstormfract：从每个水平格点的 stormdust_mass 垂直廓线中取第 2 层到顶层的最大绝对质量混合比，用参考质量混合比 5.e-4 归一化，再把结果限制在 [1.e-2, 0.6]。输出 totstormfract(ngrid) 表示每个网格中被风暴尘埃占据的面积分数。

physiq 在辐射计算前、且 rdstorm 开启时调用本例程。后续 callradite / aeropacity / rocketduststorm 使用 totstormfract 把格点平均 stormdust 转换为风暴内部浓度或风暴内部光学厚度。

定义的符号

依赖的模块

调用的关键例程

输入

输出

共享状态与副作用

• 只读 tracer_mod 中的 igcm_stormdust_mass。
• 不写 module 变量、common block、文件或日志。
• 无 SAVE、无 THREADPRIVATE、无 firstcall 初始化。
• 无显式错误处理；若 nlayer < 2，切片 pq(ig,2:nlayer,...) 不具备正常物理含义。复现 1 层测试时需额外确认。

核心逻辑

1. 设置三个硬编码参数：
   • mmr_ref = 5.e-4：参考质量混合比，源码注释说明对应 tau=10 的 OMEGA 本地尘暴。
   • fracmax = 0.6：最大风暴面积分数。
   • fracmin = 1.e-2：最小风暴面积分数。
2. 对每个水平格点 ig：
   • 读取 pq(ig,2:nlayer,igcm_stormdust_mass)，即跳过第 1 层，从第 2 层到顶层取风暴尘埃质量混合比。
   • 取该垂直切片的最大值，再取绝对值：abs(maxval(...))。
   • 除以 mmr_ref 得到未限制的面积分数。
   • 先用 max(..., fracmin) 施加下限 0.01。
   • 再用 min(..., fracmax) 施加上限 0.6。

源码行 56 的注释写有 why from 2nd layer?，说明跳过第 1 层的理由在源码内未解释。

伪代码

calcstormfract(ngrid, nlayer, nq, pq, totstormfract):
    mmr_ref = 5.e-4
    fracmin = 1.e-2
    fracmax = 0.6

    for ig in 1..ngrid:
        storm_mmr = abs(maxval(pq[ig, 2:nlayer, igcm_stormdust_mass]))
        fraction = storm_mmr / mmr_ref
        fraction = max(fraction, fracmin)
        fraction = min(fraction, fracmax)
        totstormfract[ig] = fraction

参与的主题流程

写法特点

• .F90 后缀但代码排版保留固定格式风格缩进。
• 纯函数式辅助例程，无分配、无 I/O、无诊断。
• nq 只用于 pq 形状声明，例程主体不直接读取。
• totstormfract 永远至少为 0.01，即没有 stormdust 的格点也不会返回 0；这是为了后续除以 totstormfract 时避免除零。
• totstormfract 最大为 0.6，对应最多 60% 网格面积被风暴尘埃占据。
• 第 1 层被排除在 maxval 之外，源码只留下疑问注释，未说明物理理由。

复现要点

• 必须使用 stormdust_mass，不是 stormdust_number 或背景 dust_mass。
• 垂直范围按源码为 2:nlayer；不要改成 1:nlayer，除非另有物理依据。
• 下限和上限顺序为先 max(fracmin) 再 min(fracmax)。
• 低 stormdust 或无 stormdust 时返回 0.01，会影响 aeropacity 中 aerosol(...stormdust...)/totstormfract 以及 rocketduststorm 中 stormdust/totstormfract。
• abs(maxval(...)) 只对最大值取绝对值；如果数组中全是负值，则结果是最接近 0 的负值取绝对值，而不是最大绝对值。正常 tracer 应非负，但复现异常输入时应按源码字面执行。

待确认

• 为什么 maxval 跳过第 1 层（源码注释直接写 why from 2nd layer?）。
• abs(maxval(...)) 是否有意而不是 maxval(abs(...))；正常非负 tracer 下二者等价。
• 1D 或测试配置若 nlayer=1 时该切片如何处理，当前源码没有保护。

相关页面

• rocketduststorm_mod：火箭式尘暴垂直输运主模块，使用 totstormfract 做风暴内部浓度转换。
• dust-cycle：尘埃循环主题页。
• dust_coagulation_mod：包含 stormdust 双矩 tracer 在凝并中的处理。
• callsedim_mod：包含 stormdust 沉降路径。
• aeropacity_mod：在 section 5 中使用 totstormfract 把 stormdust 光学厚度折算到风暴面积内部浓度。
• callradite_mod.md：辐射调度入口。