LMDZ.MARS Wiki

dyn1d/testphys1d.F90

路径

LMDZ.MARS\libf\phymars\dyn1d\testphys1d.F90

文件定位

testphys1d.F90 定义 1D 单柱物理测试主程序 testphys1d。它不运行 3D 动力核心，而是在 ngrid=1 的单个物理列上：

  1. 检查 run.def 是否存在。
  2. 读取 startfiles_1D 并判断 start1D.txt、startfi.nc 是否可用。
  3. 调用 init_testphys1d_mod 初始化时间、压力、温度、风、tracer、地表、土壤、坡面和水 profile 控制状态。
  4. 必要时用 phyredem 写出新的 startfi.nc。
  5. 在 ndt 个物理步中调用 physiq_mod 的 physiq，并用 1D 专用的风场松弛、水 profile 强制/松弛和表压更新推进状态。
  6. 如果启用 startfiles_1D，调用 writerestart1D_mod 写出 restart1D.txt。

源码注释给出的编译示例是 makegcm -p mars -d 25 testphys1d；运行需要 testphys1d.def、callphys.def、带 INCLUDEDEF=callphys.def 的 run.def，以及 z2sig.def 等垂直层定义文件。

定义的符号

符号 类型 行号 作用
testphys1d program 1-268 1D 火星物理单柱主程序，初始化状态并循环调用 physiq。
gr_fi_dyn dummy subroutine 274-286 为编译 testphys1d 和使用 writediagfi 提供 1D 到动力网格的最小映射。

文件共 289 行。主程序内没有自定义 module；所有物理、tracer、并行和 I/O 功能都通过 use 引入。

关键状态和维度

变量 维度/单位 来源/用途
ngrid=1 单格点 固定 1D 物理列。
nlayer=llm 垂直层数 来自 dimensions.h。
odpref=610. Pa 尘埃光学厚度参考压力，传给初始化。
ndt, dttestphys 步数、秒 init_testphys1d 从配置计算；主循环 do idt=1,ndt 使用。
day0, day, time sol、sol 分数 物理时间和日内时间，传给 physiq 并逐步更新。
psurf, plev, play Pa 表压、界面压力、层中压力；主循环按 dpsurf 更新后用 ap/bp/aps/bps 重算。
u, v, temp (nlayer) 风和温度剖面；每步累加 du/dv/dtemp。
q, dq, dqdyn allocatable (1,nlayer,nq) tracer 混合比和 tendency；q 每步累加 dttestphys*dq。
ctrl_h2ovap, ctrl_h2oice logical 控制 1D 水汽/水冰参考廓线强制或松弛。
qref_h2ovap, qref_h2oice (nlayer) 水 profile 参考值，由初始化读取。

初始化流程

1. 命令行和并行占位

主程序先调用 parse_args()。若编译了 CPP_XIOS，随后调用 mod_const_mpi 的 init_const_mpi 和 parallel_lmdz 的 init_parallel，为 XIOS/并行输出路径准备 communicator 和 rank 状态。

2. run.def 和 start 文件检查

  1. inquire(file='run.def', exist=there) 检查运行目录下必须存在 run.def。缺失时打印说明并 error stop。
  2. startfiles_1D 默认 .false.，通过 getin('startfiles_1D',startfiles_1D) 改写。
  3. 如果 startfiles_1D=.true.，分别检查 start1D.txt 和 startfi.nc。
  4. 如果有 start1D.txt 但没有 startfi.nc，程序直接 error stop，因为 qsurf、tsurf 等物理 restart 状态可能无法一致初始化。

3. 调用 1D 初始化

主程序调用：

call init_testphys1d('start1D.txt','startfi.nc',therestart1D,therestartfi,ngrid,nlayer,odpref,nq,q, &
                     time,psurf,u,v,temp,ndt,ptif,pks,dttestphys,zqsat,dq,dqdyn,day0,day,gru,grv,w, &
                     play,plev,latitude,longitude,cell_area,                                        &
                     ctrl_h2ovap,relaxtime_h2ovap,qref_h2ovap,ctrl_h2oice,relaxtime_h2oice,qref_h2oice)

该调用负责分配 q/dq/dqdyn/zqsat，读取 traceur.def、profile、restart、垂直坐标和物理配置；详细分支见 init_testphys1d_mod。

4. 生成 startfi.nc

如果没有 startfi.nc，主程序先调用 physdem0("startfi.nc",...) 写头部/维度/控制信息，再调用 physdem1("startfi.nc",...) 写地表、土壤、坡面、tracer surface、云量、tauscaling 等物理状态。这个 startfi.nc 会在第一次 physiq 调用时被物理模块读取。

时间步主循环

主循环是 do idt = 1,ndt。每步执行：

  1. 若 idt==ndt，设 lastcall=.true.；首步前 firstcall=.true.。
  2. 用 aps/bps、psurf、rcp、cpp 和 temp 计算 Exner 相关量 s/h，再垂直积分得到 geopotential phi。
  3. 调用：
call physiq(1,llm,nq,firstcall,lastcall,day,time,dttestphys,plev,play,phi,u,v,temp,q,w,du,dv,dtemp,dq,dpsurf)

physiq 输出风、温度、tracer 和表压 tendency。

  1. 1D 专用风场更新不是当前启用的 Coriolis 方案，而是把风向背景风 gru/grv 以 1.e4 s 时间尺度松弛：
du = du + (gru - u)/1.e4
dv = dv + (grv - v)/1.e4
  1. 若 water 为真且 ctrl_h2ovap 或 ctrl_h2oice 为真，则改写或追加水汽/水冰 tendency：
    • relaxtime < 0：强制到参考 profile，即 (qref-q)/dttestphys。
    • relaxtime >= 0：在物理 tendency 基础上减去 (q-qref)/relaxtime。
  2. 更新时间：time = time + dttestphys/daysec；若 time > 1，减 1 并令 day=day+1。
  3. 更新 prognostic 状态：
    • u = u + dttestphys*du
    • v = v + dttestphys*dv
    • temp = temp + dttestphys*dtemp
    • psurf = psurf + dttestphys*dpsurf(1)
    • plev = ap + psurf*bp
    • play = aps + psurf*bps
    • q = q + dttestphys*dq

Restart 和输出

文件/输出 条件 写入者 内容
startfi.nc therestartfi=.false. physdem0、physdem1 physiq 第一次调用所需物理 restart 状态。
restart1D.txt startfiles_1D=.true. writerestart1D psurf/pa/preff、逐 tracer 垂直廓线、u/v、tsurf/temp。
标准输出 始终 主程序 缺文件提示、startfiles_1D 值和最终完成信息。

主程序 use write_output_mod, only: write_output，但当前源码中没有直接调用 write_output；具体诊断输出发生在 physiq 及其下游物理模块。

gr_fi_dyn dummy 例程

gr_fi_dyn(nfield,ngrid,im,jm,pfi,pdyn) 是为 1D 编译路径提供的最小动力-物理映射：

  1. 若 ngrid /= 1，直接 error stop 'gr_fi_dyn error: in 1D ngrid should be 1!!!'。
  2. 执行 pdyn(1,1,1:nfield)=pfi(1,1:nfield)。

它没有通用水平重映射能力，只是让需要 gr_fi_dyn 符号的诊断/输出路径在 1D 中可链接。

调用关系

方向 符号/页面 说明
被调 init_testphys1d_mod 时间循环前的主初始化入口。
被调 phyredem 无 startfi.nc 时写物理 restart。
被调 physiq_mod 每个物理步调用的 Mars 物理主入口。
被调 writerestart1D_mod 可选写 restart1D.txt。
被调 mod_const_mpi、parallel_lmdz CPP_XIOS 分支的 MPI/并行初始化。
依赖 comvert_mod 使用 ap/bp/aps/bps/pa/preff 构造压力和 restart 实参。
依赖 tracer_mod 使用水 tracer 索引和 noms 写 restart。
依赖 profile_temp_mod、read_profile_mod 通过初始化例程间接使用温度和 tracer profile。

复现要点

  1. 运行目录必须有 run.def；否则在任何初始化前停止。
  2. run.def 需要包含 INCLUDEDEF=callphys.def，否则 getin 无法读取物理配置。
  3. startfiles_1D=.true. 时允许两个 start 文件都存在、两个都缺失、或只存在 startfi.nc；但只存在 start1D.txt 而缺 startfi.nc 会停止。
  4. ndt 是初始化例程返回的物理步数，主程序按步数循环，不再把它解释为 sol 数。
  5. 1D 风场当前向 gru/grv 松弛，源码中的 Coriolis 增量块被注释掉；复现 3D 风场反馈时不能直接等同。
  6. 水汽/水冰 profile 控制发生在 physiq 之后、状态更新之前，因此会覆盖或修正物理过程给出的 dq。
  7. 表压由 dpsurf 逐步累加，并立即重算 plev/play；长时间 1D run 的质量守恒解释需要单独检查。
  8. restart1D.txt 只有在 startfiles_1D=.true. 时写出；默认 .false. 时不会写 1D 文本 restart。

待确认和风险

相关页面