lmdz_atke_exchange_coeff.F90

路径

LMDZ.MARS\libf\phymars\lmdz_atke_exchange_coeff.F90

所属目录 / 模块

libf\phymars

文件定位

lmdz_atke_exchange_coeff.F90 定义 lmdz_atke_exchange_coeff 模块，是 ATKE（Atelier TKE）1.5 阶湍流闭合方案的核心计算模块。它被 vdifc_mod.F 调用，用于计算湍流交换系数 Km（动量）和 Kh（热量），以及湍流动能（TKE）的演化。

本模块实现了多种 TKE 数值格式（iflag_atke = 0-4）和混合长度参数化（iflag_atke_lmix = 1-3），支持稳定和不稳定条件下的湍流计算。

定义的符号

依赖的模块

调用的关键例程

输入

输出

共享状态与副作用

• 本文件没有 save 变量、COMMON、THREADPRIVATE、文件 I/O 或配置读取。
• 副作用是原地修改 tke 数组。
• 所有配置参数来自 lmdz_atke_turbulence_ini 模块。

核心逻辑

atke_compute_km_kh

1. 位温计算：

   • 使用参考压力 preff=610 Pa 计算位温
   • 若 atke_ok_virtual 为真，考虑水汽对浮力的影响

2. 层间距计算：

   • 计算层中点间距 dz_lay 和界面间距 dz_interf

3. 梯度 Richardson 数和稳定性函数：

   • 计算 Brunt-Väisälä 频率平方 N²
   • 计算风切变平方 shear2
   • 计算 Richardson 数 Ri = N² / shear2
   • 不稳定条件（Ri < 0）：Sm 使用 atan 函数，Prandtl 使用 atan 函数
   • 稳定条件（Ri ≥ 0）：Sm 线性递减，Prandtl 使用指数+线性公式

4. 混合长度计算（4 种方案）：

   • iflag_atke_lmix = 0：Blackadar 公式，忽略稳定度影响
   • iflag_atke_lmix = 1：Blackadar + Deardoff 浮力长度尺度
   • iflag_atke_lmix = 2：添加风切变影响（Grisogono & Belusic 2008）
   • iflag_atke_lmix = 3：添加风切变影响（Grisogono 2010）

5. TKE 计算（5 种数值格式）：

   • iflag_atke = 0：稳态解（dtke/dt = 0）
   • iflag_atke = 1：全隐式格式（二阶多项式求解）
   • iflag_atke = 2：半隐式格式（正定性保证）
   • iflag_atke = 3：MAR 格式（Deleersnijder 1992）
   • iflag_atke = 4：Arpege 格式（V. Masson 方法）

6. 边界条件：

   • 顶层 tke(:,nlay+1) = 0
   • 地表 tke(:,1) = ctkes * ustar²

7. TKE 垂直扩散（可选）：

   • 若 atke_ok_vdiff 为真，调用 atke_vdiff_tke

8. 交换系数计算：

   • Km = viscom + l_exchange * Sm * sqrt(tke)
   • Kh = viscoh + l_exchange * Sh * sqrt(tke)

atke_vdiff_tke

1. 计算 TKE 的涡扩散系数 Ke
2. 构造三对角矩阵系数 ak, bk, ck
3. 使用隐式格式求解 TKE 的垂直扩散

伪代码

atke_compute_km_kh(ngrid, nlay, dtime, wind_u, wind_v, temp, qvap,
                   play, pinterf, z_lay, z_interf, cdrag_uv, tke,
                   Km_out, Kh_out):

    ! 1. 位温计算
    theta = temp * (preff/play)^(rd/rcpd)
    if atke_ok_virtual:
        theta = theta * (1 + rvap/(RD/RV)) / (1 + rvap)

    ! 2. 层间距计算
    dz_lay = z_lay(ilay) - z_lay(ilay-1)
    dz_interf = z_interf(ilay) - z_interf(ilay-1)

    ! 3. Richardson 数和稳定性函数
    N2 = g * (theta(ilay) - theta(ilay-1)) / thetam / dz_lay
    shear2 = (du/dz)² + (dv/dz)²
    Ri = N2 / max(shear2, 1e-10)

    if Ri < 0:  ! 不稳定
        Sm = (2/π) * (cinf-cn) * atan(-Ri/Ri0) + cn
        Prandtl = -(2/π) * (pr_asym-pr_neut) * atan(Ri/Ri1) + pr_neut
    else:  ! 稳定
        Sm = max(smmin, cn * (1 - Ri/Ric))
        Prandtl = pr_neut * exp(-pr_slope/pr_neut * Ri + Ri/pr_neut) + Ri * pr_slope

    Sh = Sm / Prandtl

    ! 4. 混合长度
    l_exchange = kappa*l0*z_interf / (kappa*z_interf + l0)
    if N2 > 0:
        lstrat = clmix * sqrt(tke) / sqrt(N2)
        l_exchange = 1 / (1/l_exchange + 1/lstrat)

    ! 5. TKE 计算
    if iflag_atke == 0:  ! 稳态
        tke = cepsilon * l² * Sm * shear2 * (1 - Ri/Prandtl)
    elif iflag_atke == 1:  ! 全隐式
        solve quadratic equation for qq
        tke = 0.5 * qq²
    elif iflag_atke == 2:  ! 半隐式
        qq = (qq + source) / (1 + loss)
        tke = 0.5 * qq²
    elif iflag_atke == 3:  ! MAR
        qq = (qq²/dtime + qq*netsource) / (qq/dtime + netloss)
        tke = 0.5 * qq²
    elif iflag_atke == 4:  ! Arpege
        qq = (source + qq*(1+dtime*qq/...)) / (1 + 2*qq*dtime/...)
        tke = 0.5 * qq²

    ! 6. 边界条件
    tke(:,nlay+1) = 0
    tke(:,1) = ctkes * ustar²

    ! 7. TKE 垂直扩散（可选）
    if atke_ok_vdiff:
        atke_vdiff_tke(...)

    ! 8. 交换系数
    Km = viscom + l_exchange * Sm * sqrt(tke)
    Kh = viscoh + l_exchange * Sh * sqrt(tke)

参与的主题流程

写法特点

• .F90 自由格式模块，使用 Fortran 90 的数组语法。
• 实现了 5 种 TKE 数值格式（iflag_atke = 0-4），通过配置参数选择。
• 实现了 4 种混合长度参数化（iflag_atke_lmix = 0-3），支持不同稳定度条件。
• 使用 abort_physic 进行错误处理。
• atke_vdiff_tke 使用隐式三对角格式求解 TKE 垂直扩散。

复现要点

• 必须确保 lmdz_atke_turbulence_ini 模块中的配置参数已正确初始化。
• iflag_atke 和 iflag_atke_lmix 的选择会影响 TKE 和交换系数的计算结果。
• 分子粘性 viscom 和 viscoh 会被加到湍流交换系数上。
• 地表 TKE 使用 ctkes * ustar² 参数化。
• 顶层 TKE 强制为 0。

待确认

• preff=610 Pa 的选择是否对所有火星条件都适用。
• 稳定条件下 Ri >= Prandtl 的检查（行 147-150）是否过于严格。
• atke_vdiff_tke 中 TKE 扩散系数 Ke 的计算是否正确。

相关页面

• lmdz_atke_turbulence_ini：ATKE 配置参数模块。
• vdifc_mod.md：垂直扩散主模块，调用本模块。
• pbl_parameters_mod：PBL 诊断参数计算。