https://agupubs.onlinelibrary.wiley.com/share/YNSG74MV8B8BAAUMCHN3?target=10.1029/2022JD036708
那会没把提出的水汽稳定度修正函数加到CAS-ESM,当时对CAS-ESM模式还没这么熟悉,也想着师兄能帮我，但是师兄太忙了。现在我对模式已经比较熟悉啦，之前又把师兄加COARE的工作重复了一下，万事俱备只欠东风了。
试验设计：
为了探究温度稳定度修正函数psit和比湿稳定度修正函数对湍流潜热通量的模拟差异，我们采用CAS-ESM2.0全耦合模式，设置了两组实验，除了湍流通量方案中稳定度修正函数以外，两组实验的其他设置都是相同的。我们将采用Psit的COARE算法的CAS-ESM耦合模式作为控制实验（CTRL），我们将采用Psiq的COARE算法的CAS-ESM耦合模式定义为实验一（EXP1）。两组实验都是初始积分（startup）开始,为了防止非平衡初值或扰动的条件下模式调整到平衡态的过程对实验结果影响，我们的两组实验均积分10年，从1995年1月1日到2004年12月31日，并选取2000年到2004年的数据进行分析。
技术细节写在这里：
我已经编写了模式的运行脚本，我也将COARE算法放入到了CAS-ESM中，现在只要在module_coare_model.f90中增加我的psiq的计算子程序再调用就可以了。
real function psiq_26(zet)
!% computes specific humidity structure function (ma,2022)
real(r8) :: zet
real(r8) :: dzet, psi, x, psik, psic, f
dzet = min(50.,0.35*zet) !% stable
psiq_26 = -0.6*zet
if(zet < 0.) then !% unstable
x = (1.-15.*zet)**0.5
psik = 2.*log((1.+x)/2.)
x = (1.-34.15*zet)**0.3333
psic = 1.5*log((1.+x+x*x)/3.) - sqrt(3.)*atan((1.+2.*x)/sqrt(3.)) + &
4.*atan(1.)/sqrt(3.)
f = zet*zet/(1.+zet*zet)
psiq_26 = (1.-f)*psik + f*psic
endif
end function psiq_26在稳定的时候调用我们的比湿稳定度修正函数，计算潜热通量。下面是我的调试用的脚本，每次只要./tiaoshi.sh就可以自动进行模式运行试验了。我把这次试验叫HIST_coare_psiq_ri。
#case试验名称
CASE_NAME='HIST_coare_psiq_ri'
#试验开始年月日YYYYMMDD
start_ymd='20000101'
stop_option='nyears'
stop_n='1'
#判断和删除残留文件夹
if [ -d "/data/chengxl/CAS-ESM2.0-test1/scripts/${CASE_NAME}" ]; then
rm -r /data/chengxl/CAS-ESM2.0-test1/scripts/${CASE_NAME}
fi
if [ -d "/data/chengxl/CAS-ESM2.0-test1/run/${CASE_NAME}" ]; then
rm -r /data/chengxl/CAS-ESM2.0-test1/run/${CASE_NAME}
fi
#创建case试验和编译
cd /data/chengxl/CAS-ESM2.0-test1/scripts
./create_newcase -case ${CASE_NAME} -compset HIST_C6_C -res fd14_licom -mach huan_default
cd /data/chengxl/CAS-ESM2.0-test1/scripts/${CASE_NAME}
./configure -case
./${CASE_NAME}.huan_default.build
#运行试验前的设置
cd /data/chengxl/CAS-ESM2.0-test1/run/${CASE_NAME}/run
#cp /data/chengxl/CAS-ESM2.0-test1/run_demo/atm_in .
#cp /data/chengxl/CAS-ESM2.0-test1/run_demo/drv_in .
#输出变量管理
sed -i "12 a\\nhtfrq=0,0" atm_in
sed -i "13 a\\fincl1 = 'lhf_tbf','tau_tbf' ,'shf_tbf','n2','pbl_h','h_wave','tau_coare','hsb_coare','hlb_coare', 'LHFLX', 'SHFLX','qsss'" atm_in
sed -i "14 a\\fincl2 = 'lhf_tbf','tau_tbf','shf_tbf','tau_coare','hsb_coare','hlb_coare','h_wave' ,'h_wave_ln','pbl_h','ustar','obklen','n2','LHFLX','SHFLX','ri'" atm_in
#运行时间控制
sed -i "1,$ s/continue/startup/g" drv_in
sed -i "1,$ s/start_ymd = 00010101/start_ymd =${start_ymd}/g" drv_in
sed -i "1,$ s/ndays/${stop_option}/g" drv_in
sed -i "1,$ s/stop_n = 5/stop_n =${stop_n}/g" drv_in
sed -i "1,$ s/restart_n = 5/restart_n =${stop_n}/g" drv_in
cp /data/chengxl/CAS-ESM2.0-test1/run_demo/run.slurm .
sbatch run.slurm
sleep 2
squeue
sleep 2 经过二十四小时左右的运行，模式结果全部都输出好了，没有错误。我们对模式的输出进行分析。
我将我的边界层模块的变量都放在了h1文件中。我们首先查看12月的都在。
(base) [chengxl@login01 run]$ ls B20TR_C35_y100.iap.h1.2000-*
B20TR_C35_y100.iap.h1.2000-01.nc B20TR_C35_y100.iap.h1.2000-05.nc B20TR_C35_y100.iap.h1.2000-09.nc
B20TR_C35_y100.iap.h1.2000-02.nc B20TR_C35_y100.iap.h1.2000-06.nc B20TR_C35_y100.iap.h1.2000-10.nc
B20TR_C35_y100.iap.h1.2000-03.nc B20TR_C35_y100.iap.h1.2000-07.nc B20TR_C35_y100.iap.h1.2000-11.nc
B20TR_C35_y100.iap.h1.2000-04.nc B20TR_C35_y100.iap.h1.2000-08.nc B20TR_C35_y100.iap.h1.2000-12.nc
我们首先将其转化为一年平均的，因为是地球系统模式，所以分析年平均是有意义的，时间太短反而科学性有待考量。我们用一句语句来解决：
ls B20TR_C35_y100.iap.h1.2000-* | xargs -I{} cdo yearmean {} CASESM-COARE-psiq.2000ym.nccdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo).nc
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.27s 50MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.31s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.39s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.31s 50MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.30s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.44s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.37s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.19s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.32s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.49s 49MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.31s 50MB].
cdo: /public/software/apps/anaconda3/5.3.0/lib/libuuid.so.1: no version information available (required by cdo)
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable date_written!
Warning (cdfCheckVars): Unsupported data type (char/string), skipped variable time_written!
cdo yearmean: Processed 2720139 values from 29 variables over 1 timestep [0.23s 49MB].
成功得到了我们想要的数据CASESM-COARE-psiq.2000ym.nc
因为我比较习惯于使用本地主机上的Python还有panoply来分析我的结果，所以我先用ftp把我的资料从超算上下载到本地。
我们用panoply来简单看一下潜热的计算结果。
先保存一下，hlb_coare_in_CASESM-COARE-psiq.2000ym.png
然后我们再看下如果是原先用psit计算出来的潜热通量
我们看到计算出来的结果是合理的。比原先的潜热模拟值要大，这也是我们预料之内的。因为稳定度参数随着稳定度的增加减少量是变小的，这样的话分母变小，系数是相对变大的。
哦对了，因为极地区域没有经过处理，所以也是偏大的。这一点还存在一点问题。下次修改代码的时候把极地的通量采用原先的计算结果。就是在海冰的时候不该使用COARE算法。
除了这个问题，其他没有问题了。我们开始用Python来进行一些简单的数据处理。就是和OAflux潜热，以及Psiq,Psit,之间的对比
ctrl_data = xr.open_dataset(r'H:\fluxdeepl\CAS_ESM_psi_test\CASESM-COARE-psit.2000ym.nc', engine='netcdf4')
exp1_data = xr.open_dataset(r'H:\fluxdeepl\CAS_ESM_psi_test\CASESM-COARE-psiq.2000ym.nc', engine='netcdf4')
oaflux_data = xr.open_dataset(r'H:\fluxdeepl\CAS_ESM2_flux_test\latent256X128.nc', engine='netcdf4')
lhf_ctrl = ctrl_data.hlb_coare.data
lhf_ctrl = np.squeeze(lhf_ctrl, axis = 0)
lhf_exp1 = exp1_data.hlb_coare.data
lhf_exp1 = np.squeeze(lhf_exp1, axis = 0)
lhf_oa = oaflux_data.lhtfl.data[192:204,:,:]
lhf_oa = np.nanmean(lhf_oa,axis = 0)
exp1_oa = lhf_exp1 - lhf_oa
ctrl_oa = lhf_ctrl - lhf_oa
exp1_ctrl = lhf_exp1 - lhf_ctrl
write_to_nc_1(lhf_oa,'oa.nc')
write_to_nc_1(lhf_exp1,'psiq.nc')
write_to_nc_1(lhf_ctrl,'psit.nc')
write_to_nc_1(exp1_oa,'psiq_oa.nc')
write_to_nc_1(ctrl_oa,'psit_oa.nc')
write_to_nc_1(exp1_ctrl,'psiq_psit.nc')
然后我们就用panpoly来画图了，懒得用NCL或者python画了。
1.潜热通量对水汽稳定度修正函数响应较大。
2.热带海气边界层湍流潜热通量的模拟有所改进。
