nuclea.F

路径

LMDZ.MARS\libf\phymars\nuclea.F

所属目录 / 模块

libf\phymars

文件定位

nuclea.F 定义 nuclea_mod，提供水冰云完整微物理中的异质成核率计算。它不负责水汽、尘埃或 CCN tracer 的更新；唯一主调用方 improvedclouds_mod.F90 先把尘埃矩展开成水云半径 bin 上的 n_aer/m_aer，再调用 nuclea(ph2o, zt, satu, n_aer, rate) 得到每个 bin 的成核率，随后用 exp(-rate*microtimestep) 把被激活的尘埃数/质量转入 CCN tracer。

源码注释说明其物理基础为 Pruppacher & Klett (1978) 的固体基底上水冰成核，并由 Keese (1989) 修订；当前 Mars GCM 版本还加入了 temp_dependent_m 温度依赖接触参数分支。

定义的符号

依赖的模块

调用的关键例程

输入

输出

共享状态与副作用

• 本文件自身只保存 firstcall，并声明为 OpenMP THREADPRIVATE。
• firstcall 为真时会打印当前接触参数模式：常数 mteta、Maattaanen et al. 2014 线性拟合，或 Maattaanen et al. 2014 tanh 拟合。
• mteta、rad_cld、vo1 等关键输入来自 microphys_h 的线程私有共享状态；rad_cld 和 vo1 由 improvedclouds_mod.F90:148-182 初始化，mteta 由 conf_phys.F:863-872 或温度依赖分支决定。
• 本文件不读写 restart、NetCDF、诊断文件，也不直接更新 tracer；副作用限于写输出数组和首次打印说明。

核心逻辑

1. 从 microphys_h::mteta 复制本地双精度接触参数 mtetalocal。
2. 如果 temp_dependent_m=.true.，则按 cloud_adapt_ts 选择温度依赖公式：
   • cloud_adapt_ts=.false.：mtetalocal = min(0.0044*temp + 0.1831, 0.97)。
   • cloud_adapt_ts=.true.：mtetalocal = 0.469 + ((0.972-0.469)*tanh((temp/158.282)**4.244))。
3. 首次调用时打印接触参数模式提示。
4. 若 sat <= 1，所有 bin 的 nucrate 直接置零。
5. 若 sat > 1，先计算水分子数密度 nh2o = ph2o/(kbz*temp)、临界胚胎半径 rstar、临界胚胎分子数 gstar。
6. 计算平面基底近似形状因子 fshapesimple = (2+m)(1-m)^2/4。
7. 对每个水云半径 bin：
   • 若 n_ccn(i) < 1e-10，跳过昂贵计算并令成核率为 0。
   • 若 rad_cld(i) > 3000*rstar，使用平面近似形状因子。
   • 否则调用 fshape(mtetalocal, rad_cld(i)/rstar) 计算有限曲率形状因子。
   • 用 fistar、deltaf 和指数项计算 nucrate(i)；deltaf 被限制在 [-100,100] 以避免指数溢出。

伪代码

nuclea(ph2o, temp, sat, n_ccn):
    mtetalocal = mteta
    if temp_dependent_m:
        if cloud_adapt_ts:
            mtetalocal = tanh fit from Maattaanen et al. 2014
        else:
            mtetalocal = min(linear fit, 0.97)

    print contact-parameter mode on first call

    if sat <= 1:
        nucrate(:) = 0
        return

    nh2o = ph2o / (kbz * temp)
    rstar = 2 * sig_h2o(temp) * vo1 / (rgp * temp * log(sat))
    gstar = 4 * nav * pi * rstar^3 / (3 * vo1)
    fshapesimple = (2 + m) * (1 - m)^2 / 4

    for each cloud bin i:
        if n_ccn(i) is negligible:
            nucrate(i) = 0
            continue
        if rad_cld(i) is much larger than rstar:
            zefshape = fshapesimple
        else:
            zefshape = fshape(m, rad_cld(i)/rstar)
        fistar = 4/3*pi*sig_h2o(temp)*rstar^2*zefshape
        deltaf = clamp((2*desorp - surfdif - fistar)/(kbz*temp), -100, 100)
        nucrate(i) = prefactor * exp(deltaf)

参与的主题流程

写法特点

• 固定格式 Fortran 源码，模块名小写，主例程名小写。
• nucrate 输出为 REAL，但内部大多数物理量用 DOUBLE PRECISION 计算。
• n_ccn(i) < 1e-10 分支在无尘核 bin 上跳过成核率公式，避免无意义计算。
• deltaf 在指数前被限幅到 [-100,100]；当限到 -100 时源码直接给 nucrate=0。

复现要点

• 直接调用 nuclea 前必须保证 improvedclouds firstcall 已经初始化 microphys_h::rad_cld 和 vo1；否则半径 bin 和分子体积可能未定义。
• mteta 默认值在 conf_phys.F 中为 0.95；若 temp_dependent_m=.false.，还会读取 mteta 配置覆盖。
• cloud_adapt_ts=.true. 要求 temp_dependent_m=.true.，该一致性检查在 conf_phys.F:889-894。
• sat <= 1 时成核率全为零；调用方 improvedclouds 也只在 satu >= 1 分支中调用本例程。
• nucrate 只影响成核概率；实际 dust/CCN tracer 更新发生在 improvedclouds_mod.F90:296-308。

待确认

• fshape 的输入 rap=rad_cld/rstar 在公式中是否要求额外下限保护；当前源码没有对 rstar 或 rap 做保护，依赖 sat > 1 使 log(sat) 有效。
• rgp 与 vo1 混用在 rstar 公式中的量纲细节未在源码注释中展开，本页按源码公式记录。

相关页面

• phymars 目录：本模块所属目录。
• improvedclouds_mod：唯一主调用方，负责把成核率转成 dust/CCN tracer 更新。
• growthrate：完整水云微物理中的冰晶增长阻力例程。
• microphys_h：提供 nbin_cld、rad_cld、mteta、vo1 和分子常数。
• water-cycle-config：记录 temp_dependent_m、cloud_adapt_ts、mteta 等水云微物理配置。
• nucleaco2：CO2 冰成核例程，结构上与本水冰成核例程相近但用于 CO2 云路径。