improvedco2clouds_mod.F90

路径

LMDZ.MARS\libf\phymars\improvedco2clouds_mod.F90

所属目录/模块

libf/phymars

文件定位

CO2 云微物理核心模块。基于 improvedclouds_mod.F（水冰云）改编，包含 CO2 冰晶成核（尘埃、流星体、水冰三类 CCN）、质量增长/升华和 CCN 释放逻辑。被 co2cloud_mod 在每个微时间步调用。

定义的符号

依赖的模块

调用的关键例程

输入

输出

共享状态与副作用

• SAVE+THREADPRIVATE：sigma_ice（水冰分布方差）、dev3、meteor（流星体通量数据）、pression_meteor、firstcall。
• firstcall：初始化 sigma_ice/dev3；meteo_flux=true 时从 Meteo_flux_Plane.dat 读取流星体通量数据（130 层 × 100 bin）。
• 数据文件依赖：Meteo_flux_Plane.dat（meteo_flux=true 时必需，位于 datadir）。

核心逻辑

1. firstcall 初始化（行 246-289）：sigma_ice = sqrt(log(1+nuice_sed))；读取流星体通量文件。

2. 初始化（行 293-311）：zt = pteff + sum_subpdt*microtimestep，zq = pqeff + sum_subpdq*microtimestep（低于 threshold 归零），保存初始值 zq0。

3. 计算饱和（行 315-316）：co2sat → zqsat（饱和蒸汽压）；tcondco2 → tcond（凝结温度，输入含 co2 + co2_ice）。

4. 流星体 CCN 插值（行 320-332）：将流星体通量从文件气压层插值到模式层，meteor_ccn = meteor(idx_min,:)/masse。

5. 主循环（行 336-588）——对每个格点和层：

   4.1 成核（satu >= 1 时，行 351-495）：

   • 尘埃 CCN：展开尘埃 log-normal 分布为 binned n_aer/m_aer，调用 nucleaco2（接触因子 mtetaco2=0.95），概率 Proba = 1-exp(-microtimestep*rate)，更新 ccnco2/dust tracer。
   • 流星体 CCN（meteo_flux=true）：对 meteor_ccn 调用 nucleaco2，更新 ccnco2_meteor 和 ccnco2 总量。
   • 水冰 CCN（co2useh2o=true）：展开水冰 CCN 分布（用 sigma_ice 和 dev3），调用 nucleaco2（接触因子 mteta），按 ratioh2o_ccn 分配质量到 ccnco2_h2o_mass_ice/ccn，从水冰/CCN 中扣除。

   4.2 冰增长（ccnco2_number*tauscaling >= 1 时，行 503-555）：

   • 更新 riceco2/rhocloudco2。
   • 计算温度依赖潜热：lw = l0 + l1*T + l2*T² + l3*T³ + l4*T⁴（Azreg-Aïnou 2005）。
   • 能量限制：facteurmax = |tcond-zt| * cpp/lw。
   • 调用 massflowrateco2 → Ic_rice（单粒子质量增长率 kg/s）。
   • dMice = ccnco2_number * Ic_rice * microtimestep * tauscaling。
   • 裁剪：dMice ∈ [-facteurmax, -co2_ice] ∩ [-co2, facteurmax]。
   • 更新：co2_ice += dMice，co2 -= dMice，subpdtcloudco2 = dMice*lw/cpp/microtimestep。

   4.3 CCN 释放（co2_ice < threshold_2 或 ccnco2_number*tauscaling < 1 时，行 558-586）：

   • 尘埃 CCN 归还：dust += ccnco2，ccnco2 = 0。
   • 流星体 CCN：从 dust 中扣除流星体部分（避免重复计数），ccnco2_meteor = 0。
   • 水冰 CCN：ccnco2_h2o 归还到 ccn_mass/h2o_ice/ccn_number。
   • CO2 气相回收：co2 += co2_ice，co2_ice = 0。

6. 倾向计算（行 592-627）：subpdqcloudco2 = (zq - zq0) / microtimestep，对所有相关 tracer。

伪代码

improvedCO2clouds(ngrid, nlay, microtimestep, ...):

  if firstcall:
    sigma_ice = sqrt(log(1+nuice_sed))
    if meteo_flux: 读取 Meteo_flux_Plane.dat

  zt = pteff + sum_subpdt * microtimestep
  zq = pqeff + sum_subpdq * microtimestep (threshold 保护)
  zq0 = zq

  zqsat = co2sat(zt)
  tcond = tcondco2(pplay, co2 + co2_ice)

  if meteo_flux: 插值流星体通量到模式层

  for l, ig:
    pco2 = co2 * (mmean/mco2*1e3) * pplay
    satu = pco2 / zqsat
    rho_ice_co2T = density_co2_ice(zt)

    // --- 成核 (satu >= 1) ---
    尘埃 CCN: n_aer/m_aer = log-normal 展开, nucleaco2(mtetaco2), Proba=1-exp(-rate*dt)
    流星体 CCN (meteo_flux): nucleaco2, 更新 ccnco2_meteor
    水冰 CCN (co2useh2o): nucleaco2(mteta), 按 ratioh2o_ccn 分配

    // --- 冰增长 (ccnco2_number >= 1) ---
    lw = l0 + l1*T + l2*T² + l3*T³ + l4*T⁴
    facteurmax = |tcond-zt| * cpp/lw
    Ic_rice = massflowrateco2(pplay, zt, satu, riceco2, mmean)
    dMice = ccnco2_number * Ic_rice * dt * tauscaling
    dMice = clip(dMice, [-facteurmax, facteurmax])
    co2_ice += dMice, co2 -= dMice
    subpdtcloudco2 = dMice * lw / cpp / dt

    // --- CCN 释放 (co2_ice ≈ 0) ---
    dust += ccnco2, ccnco2 = 0
    流星体/水冰 CCN 归还
    co2 += co2_ice, co2_ice = 0

  subpdqcloudco2 = (zq - zq0) / microtimestep

参与的主题流程

写法特点

• 自由格式 Fortran 90（.F90 后缀）。
• SAVE + THREADPRIVATE：sigma_ice、dev3、meteor、pression_meteor、firstcall。
• firstcall 模式：读取流星体通量文件。
• 概率成核：Proba = 1-exp(-microtimestep*rate) —— 将连续成核率转换为离散时间步的概率。
• 三类 CCN：尘埃（mtetaco2=0.95）、流星体、水冰（mteta），每类独立调用 nucleaco2。
• 能量限制：facteurmax = |tcond-zt| * cpp/lw 确保凝结/升华不违反热力学约束。
• 温度依赖潜热：lw = l0 + l1*T + l2*T² + l3*T³ + l4*T⁴（Azreg-Aïnou 2005），4 阶多项式。
• CCN 释放：CO2 冰完全升华时，尘埃/流星体/水冰 CCN 归还原 tracer——保证质量守恒。
• 硬编码常数：threshold=1e-30、threshold_2=1e-13、l0-l4（潜热系数）。

复现要点

• tcondco2 输入为 co2 + co2_ice（行 316）——包含冰相，与 co2cloud_mod 中仅传 co2_vap 不同。
• co2sat 输入为 zt（有效温度 + 倾向），tcondco2 输入为 pplay 和 tracer 值。
• 潜热 lw 使用 Azreg-Aïnou (2005) 4 阶多项式——与 co2condens 中的常数 latcond=5.9e5 不同（这里是温度依赖的）。
• facteurmax 使用 cpp/lw（而非 latcond）——单位分析：|ΔT| * cpp/lw 为无量纲质量比。
• 流星体 CCN 释放时从 dust 中扣除（行 566-567）——因为成核时 ccnco2 += meteor_ccn 但 dust 未减，释放时需修正。
• 水冰 CCN 释放时分别归还 ccn_mass（CCN 部分）和 h2o_ice（冰部分）。

待确认

• time_phylmdz_mod::daysec 在 use 列表中但代码中未直接使用——可能为遗留。
• nbin_meteor=100 与 nbinco2_cld=100 相同——流星体 bin 与微物理 bin 一一对应。
• mteta vs mtetaco2：水冰 CCN 用 mteta（来自 microphys_h），尘埃用 mtetaco2=0.95——两个接触参数。

复现风险

• Meteo_flux_Plane.dat 缺失（meteo_flux=true 时）将导致 abort_physic。
• 概率成核在 microtimestep 很大时可能过高估计成核数。
• 能量限制 facteurmax 使用 cpp（干空气比热）而非湿空气——在 CO2 丰富大气中近似合理。
• 潜热多项式在极端温度下可能外推出非物理值。

相关页面

• co2cloud_mod — CO2 云主调度 - 调用本模块
• co2-saturation-helpers — CO2 饱和与凝结温度 - tcondco2/co2sat
• density_co2_ice — CO2 冰密度 - 温度依赖密度
• updaterad — 冰粒/尘埃半径更新 - updaterice_microco2/updaterccnCO2
• nucleaco2 — CO2 异质成核率 - 尘/流星/水冰三类 CCN 成核率
• massflowrateco2 — CO2 冰粒质量传输率 - 单粒子质量增长/升华通量
• CO2 循环主题页