LMDZ.MARS Wiki

co2snow.F

路径

LMDZ.MARS\libf\phymars\co2snow.F

所属目录/模块

libf\phymars

文件定位

该文件定义 co2snow_mod 模块及其唯一例程 co2snow。它的职责是：在某列发生 CO2 凝结/升华时，根据落到地表的 CO2 冰通量更新地表（红外）发射率 pemisurf。新沉降的 CO2 雪由散射性小颗粒组成，会降低地表发射率；同时存在一个使发射率向参考值 emisref 弛豫（雪变质/metamorphism）的过程。这是一个纯辅助过程，不改变质量或温度，只改 pemisurf 一个输出。

文件头注释说明算法分两步：1 初始化，2 计算地表发射率，并引用 Forget & Pollack (1996, JGR Vol.101 No.E7) 作为物理依据。

注意：本例程被 co2condens（文件 co2condens_mod.F）的 CO2 snow/clouds 段按地表坡度 nslope 循环逐个坡面调用，每次只处理一个坡面的发射率（见调用点）。

定义的符号

符号 类型 行号 作用
co2snow_mod module 1 容器模块
co2snow subroutine 26 按 CO2 雪沉降更新地表发射率

依赖的模块

use 模块 only 列表 用途 待确认
surfdat_h iceradius, dtemisice iceradius(2)：CO2 雪平均散射半径（南/北半球，m）；dtemisice(2)：雪变质时间尺度（南/北半球，单位为火星日 sol） 否
geometry_mod latitude 格点纬度（rad），用于判断南/北半球选择 icap 否
time_phylmdz_mod daysec 一个火星日的秒数（s），把 dtemisice 的 sol 单位换算为秒 否

iceradius 与 dtemisice 的取值不在本文件内设定。默认值来源（推断：取决于运行方式）：

待确认：3D 正式运行时 dtemisice 的实际生效值（0.4 还是其他），取决于 start 文件中的 tab_cntrl。

调用的关键例程

被调用例程 所在模块/文件 调用位置 作用
无 — — 本例程不调用其他例程，仅含算术与 write(*,*) 诊断输出

输入

输入 来源 类型/维度 单位 含义
ngrid 调用方 integer, intent(in) — 大气列数
nlayer 调用方 integer, intent(in) — 大气层数
ptimestep 调用方 real, intent(in) s 物理时间步长
emisref(ngrid) 调用方 real, intent(in) — 无雪时的地表/冰参考发射率
condsub(ngrid) 调用方 logical, intent(in) — 该列是否存在 CO2 凝结或升华
pplev(ngrid,nlayer+1) 调用方 real, intent(in) Pa 层间压力
pcondicea(ngrid,nlayer) 调用方 real, intent(in) kg/m2/s 各层 CO2 凝结率
pcondices(ngrid) 调用方 real, intent(in) kg/m2/s 地表 CO2 凝结率
pfallice(ngrid,nlayer+1) 调用方 real, intent(in) kg/m2/s 下落的 CO2 冰通量

复现风险：pplev、pcondicea、pcondices 在当前实现中被声明为 intent(in) 但例程主体未使用（核心公式只用到 pfallice(ig,1)）。它们应是为接口兼容/历史保留。详见“写法特点”。

输出

输出 去向 类型/维度 单位 含义
pemisurf(ngrid) 调用方 real, intent(out) — 更新后的地表发射率

注意：pemisurf 声明为 intent(out)，但例程在 condsub(ig) 为真的分支里读取其入口值（pemisurf(ig)/emisref(ig)）参与计算。调用方 co2condens_mod.F 在调用前已用 pemisurf_tmp(:) = pemisurf(:,islope) 传入上一步的发射率，因此实际作为 in/out 使用。详见“写法特点”。

共享状态与副作用

核心逻辑

  1. firstcall（每线程一次）：用 Kscat(icap) = (0.001/3.) * alpha / iceradius(icap) 计算南/北半球散射系数，alpha = 0.45（硬编码 parameter），随后置 firstcall=.false.。
  2. 对每个格点 ig 循环：
    • 若 condsub(ig) 为假（无凝结/升华）：直接令 pemisurf(ig) = emisref(ig)（回到裸地/冰参考发射率），不做雪修正。
    • 若 condsub(ig) 为真：
      • 按纬度选半球：latitude(ig) < 0 → icap=2（南），否则 icap=1（北）。
      • 计算发射率变化率 zdemisurf，由两项相加：
        • 弛豫项 (emisref - pemisurf) / (dtemisice(icap) * daysec)：把发射率以时间尺度 dtemisice（sol，乘 daysec 转秒）拉回参考值，代表雪变质恢复。
        • 新雪沉降项：用一个积分形式更新——emisref * ((pemisurf/emisref)^(-3) + 3*Kscat(icap)*pfallice(ig,1)*ptimestep)^(-1/3) - pemisurf，再除以 ptimestep。源码注释标明这是为“数值安全”采用的积分形式。该项随地表落冰通量 pfallice(ig,1) 增大而压低发射率。
      • 更新 pemisurf(ig) = pemisurf(ig) + zdemisurf * ptimestep。
      • 若结果 < 0.1，打印告警（不做钳制，发射率不被强制抬回）。

复现风险：指数 (-1/3.) 中分子 -1 为整数、分母 3. 为实数，按 Fortran 运算 -1/3. 等于 -0.3333...（实数除法，因为有一个实数操作数）。(...)**(-3) 的 -3 为整数指数。复现实现时需保持该混合整型/实型写法对应的数值，不要误写成整数除法 -1/3 = 0。

伪代码

subroutine co2snow(ngrid, nlayer, ptimestep, emisref, condsub,
                   pplev, pcondicea, pcondices, pfallice, pemisurf):

    if firstcall:                       # 每 OpenMP 线程一次
        for icap in {1,2}:
            Kscat(icap) = (0.001/3.) * 0.45 / iceradius(icap)
        firstcall = false

    for ig in 1..ngrid:
        if not condsub(ig):
            pemisurf(ig) = emisref(ig)          # 无凝结：直接取参考发射率
        else:
            icap = (latitude(ig) < 0) ? 2 : 1   # 南/北半球

            relax = (emisref(ig) - pemisurf(ig)) / (dtemisice(icap) * daysec)

            snow  = ( emisref(ig) *
                      ( (pemisurf(ig)/emisref(ig))**(-3)
                        + 3*Kscat(icap)*pfallice(ig,1)*ptimestep )**(-1/3.)
                      - pemisurf(ig) ) / ptimestep

            zdemisurf = relax + snow
            pemisurf(ig) = pemisurf(ig) + zdemisurf * ptimestep

            if pemisurf(ig) < 0.1:
                print warning(ig, pemisurf(ig), zdemisurf*ptimestep)

参与的主题流程

主题 参与方式
CO2 循环 CO2 凝结后处理的一环：把地表 CO2 雪沉积对地表红外发射率的影响反馈给地表能量收支（发射率影响地表辐射冷却）

写法特点

复现要点

待确认

相关页面