aeroptproperties.F

路径

LMDZ.MARS\libf\phymars\aeroptproperties.F

所属目录/模块

libf\phymars

文件定位

该文件定义气溶胶光学属性计算模块 aeroptproperties_mod，导出唯一例程 aeroptproperties，由辐射调用入口 callradite（行 424）在 aeropacity 之前调用。它的职责是：在每个网格箱内，根据气溶胶颗粒的有效半径 reffrad 和有效方差 nueffrad，对预先按离散半径表（radiustab）给出的单粒子散射参数做对数正态尺寸分布的积分（Gauss-Legendre 求积），从而得到该网格箱实际尺寸分布下的有效散射参数——消光效率比 Q/Qref、单次散射反照率 omega、不对称因子 g，分可见光（nsun=2 带）与红外（nir=4 带）以及参考波长输出。

为避免每个网格箱都重做昂贵的卷积积分，模块在 firstcall 建立一张 (reff, nueff) 对数插值网格（refftabsize=100 × nuefftabsize=100），把每个网格节点的卷积结果缓存在 *grid 数组里（checkgrid 标记是否算过），随后对每个真实网格箱只做线性（常方差）或双线性（变方差）插值。

定义的符号

依赖的模块

调用的关键例程

输入

输出

共享状态与副作用

• 全部插值网格缓存、Gauss 中间量、pi、firstcall、out_qwg 均为 SAVE+THREADPRIVATE（行 41–221），在 firstcall 一次性 allocate 并初始化。
• checkgrid 跨调用累积：某 (reff,nueff) 网格节点一旦算过即标 .true.，后续调用复用缓存，因此本例程的输出与历史调用顺序无关但缓存状态会增长。
• 标准输出副作用：firstcall 打印插值网格各散射体/通道的最小/最大半径；除零风险时打印警告并 abort_physic。
• 测试副作用：out_qwg=.true.（源码默认 .false.）时对 out_iaer=2 气溶胶用 write_output 写 qextvis/omegvis/gvis/qextir/omegir/gir/omegvisref/omegirref 诊断（#ifndef MESOSCALE）。
• 只读 dimradmars_mod 的离散光学参数表；不修改它们。

核心逻辑

1. firstcall（行 233–345）：allocate 全部 SAVE 数组；checkgrid=.false.、varyingnueff=.false.；pi=2*asin(1)；对每个气溶胶 iaer、每个通道域 idomain（1=可见，2=红外）建立对数半径插值轴 refftab（从 refftabmin=radiusm-radiusr*radgaus(ngau) 到 refftabmax=最大离散半径，按体积比 logvratgrid 等比展开）和指数方差轴 nuefftab（exp 从 nuefftabmin=-4.6 到 nuefftabmax=0）。
2. 逐气溶胶主循环（行 347–1293）：
   • 单粒径捷径（行 349–370）：若 nsize(iaer,1)=nsize(iaer,2)=1，光学属性均匀，直接把离散表第 1 个半径的值铺到所有 ig,lg。
   • 否则进入双通道域循环（行 372）：
     • 1.1–1.2 Gauss 节点插值（行 375–504）：取半径中点 radiusm、半幅 radiusr；对 10 个 Gauss 节点分别在 radiusm±drad 处用 minloc 定位离散半径区间、线性内插（kint）得到 a 组（radiusm-drad）和 b 组（radiusm+drad）的单粒子参数 qsqref/omeg/g/qref（可见或红外）。
     • 常方差分支（.NOT.varyingnueff，行 507–871）：对每层每格点算 reff 网格索引 grid_i 与权重 kx；对 j=grid_i,grid_i+1 两节点若 checkgrid 未算则做对数正态分布卷积（dista/distb 为分布权重，normd 为归一化截面），分别累加 qext/qscat/g/qref/qscatref 网格量并归一化，得 qsqref/omeg；最后用 k1=1-kx,k2=kx 线性插值写出 3D 输出。
     • 变方差分支（varyingnueff，行 873–1287）：额外算方差网格索引 grid_j 与权重 ky；对 (j,k) 四节点做同样卷积缓存；用双线性权重 k1..k4 写出 3D 输出。
3. 测试输出（行 1295–1350）：out_qwg 时写诊断场。

伪代码

if firstcall:
    allocate 所有 SAVE 网格/缓存数组; checkgrid=false; varyingnueff=false
    pi = 2*asin(1)
    for iaer, idomain in (1=VIS, 2=IR):
        建立对数半径轴 refftab[1..100] 和指数方差轴 nuefftab[1..100]
    firstcall = false

for iaer in 1..naerkind:
    if nsize(iaer,1)==1 and nsize(iaer,2)==1:
        把离散表第 1 半径值直接铺到所有 (ig,lg)   # 单粒径均匀
    else:
        for idomain in (VIS, IR):
            # 在 10 个 Gauss 节点 radiusm±drad 处线性内插单粒子参数 -> a/b 组
            if not varyingnueff(iaer):     # 常方差: 一维网格 + 线性插值
                for (lg,ig):
                    定位 grid_i, kx
                    for j in grid_i..grid_i+1:
                        if not checkgrid(j,1,iaer,idomain):
                            对数正态卷积 -> qext/qscat/g/qref 网格量, 归一化
                            checkgrid=true
                    线性插值 (k1,k2) 写 QVISsQREF3d/.../QREFir3d 等
            else:                          # 变方差: 二维网格 + 双线性插值
                for (lg,ig):
                    定位 grid_i,kx 和 grid_j,ky
                    for (j,k) in 四节点:
                        if not checkgrid: 卷积缓存
                    双线性插值 (k1..k4) 写 3D 输出

if out_qwg:  # 默认 false
    write_output 诊断 qext/omeg/g (#ifndef MESOSCALE)

参与的主题流程

写法特点

• 固定格式 Fortran（.F，第 6 列续行）。
• 全部网格缓存为 SAVE+THREADPRIVATE+ALLOCATABLE，firstcall 分配；OpenMP 下每线程各持一份。
• checkgrid 惰性缓存：只在首次命中某 (reff[,nueff]) 网格节点时做卷积，显著省时；varyingnueff(iaer)=.false. 进一步退化为一维网格（注释提示设为 false 计算更快）。
• Gauss-Legendre 10 点求积，radgaus/weightgaus 为硬编码 DATA。
• 硬编码常数：nuefftabmin=-4.6、nuefftabmax=0、refftabsize=100、nuefftabsize=100、ngau=10、out_iaer=2、normd 初值 1e-30。
• 对数正态分布参数按 Hansen 1974：sizedistk1=r_g、sizedistk2=sigma_g^2。
• 数值保护：grid_i/grid_j 越界钳制到 [1,tabsize-1]，normd==1e-30 触发 abort_physic。

复现要点

• 维度常数：nsun=2（不可改，见 dimradmars 注释）、nir=4、naerkind 由 callphys.def 决定。
• 插值网格：对数半径轴 refftab 按体积比 logvratgrid 等比展开，节点数 100；方差轴 nuefftab 指数分布于 exp(-4.6)..exp(0)，节点数 100。
• 卷积归一化：normd 为对数正态分布的几何截面加权积分；所有 *grid 量除以 normd 得加权平均；omeg=qscat/qext、qsqref=qext/qref。
• 单粒径捷径要求 nsize(iaer,1)=nsize(iaer,2)=1；否则必须走 Gauss 积分路径。
• 输入 radiustab/QVISsQREF/... 必须先由 suaer 等填充，否则插值结果无意义。
• out_qwg 默认 .false.，复现常规运行时不产生诊断输出。

待确认

• 行 822、1213 红外 qsqrefIRgrid 用 qrefVISgrid（可见光参考）作分母，而非 qrefIRgrid；红外卷积（行 729、745 等）也用 qrefVISb/qrefVISa 加权。推断：这是 LMDZ.MARS 辐射方案“红外 Q 相对可见光参考波长归一”的约定（与 QIRsQREF 命名一致，QIRsQREF=Qext/Qext(longrefvis)，见 dimradmars 注释行 149），并非笔误；但跨通道用可见光参考量是易误解点，复现时需保留该约定。
• 常方差分支卷积循环（行 543）用 nueffrad(1,1,iaer)（固定取第 1 格点第 1 层方差）而非当前 nueffrad(ig,lg,iaer)；推断这正是“常方差”假设的体现——整层大气用同一方差以省时，与文件头注释 varyingnueff 说明一致。

相关页面

• aeropacity_mod — 紧随本例程被 callradite 调用，用本例程输出的 QREFvis3d/QREFir3d 等算光学厚度
• updaterad / updatereffrad_mod — 提供 reffrad/nueffrad
• radiation 主题、dust-cycle 主题
• write_output_mod — 测试诊断输出后端
• callradite_mod - 本例程唯一调用方和辐射入口。
• dimradmars_mod - 离散光学参数表与维度常数来源。
• suaer - 填充 radiustab/QVISsQREF/QIRsQREF 等离散光学参数表。