diff --git a/README b/README index bc5df3e3c4..507b11775b 100644 --- a/README +++ b/README @@ -1,4 +1,4 @@ -WRF Model Version 4.7.1 +WRF Model Version 4.8.0 https://www2.mmm.ucar.edu/wrf/users/ diff --git a/Registry/Registry.EM_COMMON b/Registry/Registry.EM_COMMON index 025600ae50..d0065b6063 100644 --- a/Registry/Registry.EM_COMMON +++ b/Registry/Registry.EM_COMMON @@ -1155,6 +1155,11 @@ state real maxwidth ij misc 1 - h "m state real ztop_plume ij misc 1 - h "ztop_plume" "Height of tallest plume" "m" state real excess_h ij misc 1 - h "excess_h" "Max heat excess of plume initialization" "K" state real excess_q ij misc 1 - h "excess_q" "Max moisture excess of plume initialization" "kg/kg" +state real maxMF_dd ij misc 1 - h "maxMF_dd" "Maximum mass-flux of downdrafts" "m/s * area" +state real maxwidth_dd ij misc 1 - h "maxwidth_dd" "Maximum downdraft width" "m" +state real maxtkeprod ij misc 1 - h "maxtkeprod" "Maximum tke production by cloud top cooling" "m2/s2/hr" +state real cldtop_cooling ij misc 1 - h "cldtop_cooling" "Radiative cooling at cloud top" "K/hr" +state real ent_eff ij misc 1 - h "ent_eff" "Entrainment efficiency for cloud tops" "unitless" #FogDES variables state real fgdp ij misc 1 - - "fgdp" "Accumulated fog deposition" "mm" @@ -2438,8 +2443,8 @@ rconfig integer nssl_3moment namelist,physics 1 0 rconfig integer nssl_density_on namelist,physics 1 -1 rh "NSSL graupel/hail density flag" "" "" rconfig integer nssl_ssat_output namelist,physics 1 0 rh "NSSL ssat output flag" "" "" -rconfig integer tempo_aerosolaware namelist,physics 1 1 rh "TEMPO aerosol-aware flag" "" "" -rconfig integer tempo_hailaware namelist,physics 1 1 rh "TEMPO hail-aware flag" "" "" +rconfig integer tempo_aerosolaware namelist,physics 1 0 rh "TEMPO aerosol-aware flag" "" "" +rconfig integer tempo_hailaware namelist,physics 1 0 rh "TEMPO hail-aware flag" "" "" rconfig integer CCNTY namelist,physics 1 2 rh "Aerosol background type for NTU microphysics" "" "" @@ -2649,6 +2654,8 @@ rconfig real oml_gamma namelist,physics 1 0.1 rconfig real oml_relaxation_time namelist,physics 1 0. h "oml_relaxation_time" "Relaxation time of mixed layer ocean model back to original values" "s" rconfig integer isftcflx namelist,physics 1 0 h "isftcflx" "switch to control sfc fluxes" "" rconfig integer iz0tlnd namelist,physics 1 0 h "iz0tlnd" "switch to control land thermal roughness length" "" +rconfig logical kim_tofd namelist,physics max_domains .false. rh "kim_tofd" "turbulent form drag (tofd) option for kim gravity wave drag" "" +rconfig real tofd_factor namelist,physics max_domains .003 rh "tofd_factor" "factor in kim tofd scheme" "" rconfig real shadlen namelist,physics 1 25000. - "shadow_length" "maximum length of orographic shadow" "m" rconfig integer slope_rad namelist,physics max_domains 0 - "slope_rad" "1: use slope-dependent radiation, 0:not" "" rconfig integer topo_shading namelist,physics max_domains 0 - "topo_shading" "1: apply topographic shading to radiation, 0:not" "" @@ -2892,13 +2899,8 @@ rconfig integer damp_opt namelist,dynamics 1 3 rconfig integer rad_nudge namelist,dynamics 1 0 irh "rad_nudge" "" "" rconfig integer gwd_opt namelist,dynamics max_domains 0 irh "gwd_opt" "" "" rconfig integer gwd_diags namelist,dynamics 1 0 irh "gwd_diags" "switch to turn on extra gwd diagnostics if available for given gwd_opt" "" -<<<<<<< HEAD -rconfig real gwd_dx_factor namelist,dynamics max_domains 2. rh "gwd_dx_factor" "effective grid factor for kim gravity wave drag" "" -rconfig logical gwd_if_nonhyd namelist,dynamics max_domains .true. rh "gwd_if_nonhyd" "non_hydrostatic option for kim gravity wave drag" "" -======= -rconfig logical kim_tofd namelist,dynamics max_domains .true. rh "kim_tofd" "turbulent form drag (tofd) option for kim gravity wave drag" "" -rconfig real tofd_factor namelist,dynamics max_domains .003 rh "tofd_factor" "factor in kim tofd scheme" "" ->>>>>>> cc524de3ed762a43176cf5318b47656d4a284e7c +rconfig real gwd_dx_factor namelist,dynamics max_domains 2. rh "gwd_dx_factor" "effective grid factor for kim gravity wave drag" "" +rconfig logical gwd_if_nonhyd namelist,dynamics max_domains .true. rh "gwd_if_nonhyd" "non_hydrostatic option for kim gravity wave drag" "" rconfig real zdamp namelist,dynamics max_domains 5000. h "zdamp" "" "" rconfig real dampcoef namelist,dynamics max_domains 0.2 h "dampcoef" "" "" rconfig real khdif namelist,dynamics max_domains 0 h "khdif" "" "" @@ -3227,7 +3229,7 @@ package ysuscheme bl_pbl_physics==1 - - package myjpblscheme bl_pbl_physics==2 - state:tke_pbl,el_pbl package gfsscheme bl_pbl_physics==3 - - package qnsepblscheme bl_pbl_physics==4 - state:tke_pbl,el_pbl,massflux_EDKF,entr_EDKF,detr_EDKF,thl_up,thv_up,rv_up,rt_up,rc_up,u_up,v_up,frac_up,rc_mf -package mynnpblscheme bl_pbl_physics==5 - scalar:qke_adv;state:qke,sh3d,sm3d,tsq,qsq,cov,el_pbl,tke_pbl +package mynnpblscheme bl_pbl_physics==5 - scalar:qke_adv;state:qke,sh3d,sm3d,tsq,qsq,cov,cldfra_bl,qc_bl,qi_bl,el_pbl,tke_pbl package acmpblscheme bl_pbl_physics==7 - - package boulacscheme bl_pbl_physics==8 - state:el_pbl,tke_pbl,wu_tur,wv_tur,wt_tur,wq_tur package camuwpblscheme bl_pbl_physics==9 - state:tauresx2d,tauresy2d,tpert2d,qpert2d,wpert2d,tke_pbl,smaw3d,wsedl3d,turbtype3d @@ -3239,11 +3241,10 @@ package kepsscheme bl_pbl_physics==17 - scalar:tke_ad package mrfscheme bl_pbl_physics==99 - - package tkebudget tke_budget==1 - state:qSHEAR,qBUOY,qDISS,qWT,dqke -package mynn_dmp_edmf1 bl_mynn_edmf==1 - state:ztop_plume,maxmf,maxwidth,excess_h,excess_q -package mynn_dmp_edmf2 bl_mynn_edmf==2 - state:ztop_plume,maxmf,maxwidth,excess_h,excess_q +package mynn_dmp_edmf1 bl_mynn_edmf==1 - state:ztop_plume,maxmf,maxwidth,excess_h,excess_q,maxmf_dd,maxwidth_dd,maxtkeprod,cldtop_cooling,ent_eff +package mynn_dmp_edmf2 bl_mynn_edmf==2 - state:ztop_plume,maxmf,maxwidth,excess_h,excess_q,maxmf_dd,maxwidth_dd,maxtkeprod,cldtop_cooling,ent_eff package mynn_3Doutput1 bl_mynn_output==1 - state:edmf_a,edmf_w,edmf_thl,edmf_qt,edmf_ent,edmf_qc,sub_thl3D,sub_sqv3D,det_thl3D,det_sqv3D package mynn_3Doutput2 bl_mynn_output==2 - state:edmf_a,edmf_w,edmf_thl,edmf_qt,edmf_ent,edmf_qc,sub_thl3D,sub_sqv3D,det_thl3D,det_sqv3D -package pbl_cloud icloud_bl==1 - state:cldfra_bl,qc_bl,qi_bl package sms_3dtke km_opt==5 - state:gamu,gamv,nlflux,dlk,l_diss,elmin,xkmv_meso @@ -3552,6 +3553,8 @@ halo HALO_EM_COUPLE_B dyn_em 48:ph_1,ph_2,w_1,w_2,t_1,t_2,u_1,u_2,v_1,v_2 period PERIOD_EM_COUPLE_B dyn_em 3:ph_1,ph_2,w_1,w_2,t_1,t_2,u_1,u_2,v_1,v_2,\ moist,chem,tracer,scalar +halo HALO_EM_CFBM dyn_em 8:znt,t2,q2,psfc,rainc,rainnc,u10,v10 + #cyl for trajectory halo HALO_EM_F dyn_em 24:muu,muv,mut halo HALO_EM_F_1 dyn_em 24:muus,muvs diff --git a/Registry/registry.cordex b/Registry/registry.cordex index fced987075..cda0dafd45 100644 --- a/Registry/registry.cordex +++ b/Registry/registry.cordex @@ -88,10 +88,10 @@ state real TOTWSGSMAX ij misc 1 - rh{19} "t state real TOTUGSMAX ij misc 1 - rh{19} "totugsmax" "Total Eastward maximum near-surface wind speed of gust" "ms-1" state real TOTVGSMAX ij misc 1 - rh{19} "totvgsmax" "Total Northward maximum near-surface wind speed of gust" "ms-1" state real TOTWSGSPERCEN ij misc 1 - rh{19} "totwsgspercen" "Percentage of time steps where grid point got total wind gust" "%" -state real WSZ100 ij misc 1 - h{19} "wsz100" "100m wind speed" "ms-1" +state real WSZ100 ij misc 1 - h{19} "wsz100" "100 m wind speed" "ms-1" state real UZ100 ij misc 1 - h{19} "uz100" "Eastward 100 m wind speed" "ms-1" state real VZ100 ij misc 1 - h{19} "vz100" "Northward 100 m wind speed" "ms-1" -state real WSZ50 ij misc 1 - h{19} "wsz50" "50m wind speed" "ms-1" +state real WSZ50 ij misc 1 - h{19} "wsz50" "50 m wind speed" "ms-1" state real UZ50 ij misc 1 - h{19} "uz50" "Eastward 50 m wind speed" "ms-1" state real VZ50 ij misc 1 - h{19} "vz50" "Northward 50 m wind speed" "ms-1" state real TAZ50 ij misc 1 - h{19} "taz50" "50 m air temperature" "K" @@ -180,10 +180,6 @@ state real IVTMEAN ij misc 1 - rh{19} "i state real ZEROISOTHERM ij misc 1 - h{19} "zeroisotherm" "height above ground of the 0-isotherm" "m" state real prhf ij misc 1 - rh{18} "prhf" "precipitation flux" "kgm-2s-1" state real COLMAX ij misc 1 - h{18} "colmax" "instantaneous maximum radar reflectivity in the column" "dBz" - -# Instantaneous values. If one require them into the output, modify the second '-' by 'h{19}' -state real IUT ij misc 1 - - "iut" "instantaneous vertically integrated eastward transport of water vapour" "kgm-1s-1" -state real IVT ij misc 1 - - "ivt" "instantaneous vertically integrated northward transport of water vapour" "kgm-1s-1" endif ifdef CDXWRF_2 @@ -314,7 +310,32 @@ state real WBACZLH ij misc 1 - rh{19} "w state real WBACZMH ij misc 1 - rh{19} "wbaczmh" "W.B. mid level acc. ver. convergence of QH" "mm" state real WBACZHH ij misc 1 - rh{19} "wbaczhh" "W.B. high level acc. ver. convergence of QH" "mm" -# Instantaneous values. If one require them into the output, modify the second '-' by 'h{19}' +# instantaneous values. If one require them into the output, modify the second '-' by 'h{19}' +ifdef INSTVALS +state real QVTTEND ikj misc 1 - h{19} "qvttend" "inter time-step water vapor tendency" "kgkg-1s-1" +state real QCTTEND ikj misc 1 - h{19} "qcttend" "inter time-step cloud tendency" "kgkg-1s-1" +state real QRTTEND ikj misc 1 - h{19} "qrttend" "inter time-step rain tendency" "kgkg-1s-1" +state real QSTTEND ikj misc 1 - h{19} "qsttend" "inter time-step snow tendency" "kgkg-1s-1" +state real QITTEND ikj misc 1 - h{19} "qittend" "inter time-step ice tendency" "kgkg-1s-1" +state real QHTTEND ikj misc 1 - h{19} "qhttend" "inter time-step hail tendency" "kgkg-1s-1" +state real QGTTEND ikj misc 1 - h{19} "qgttend" "inter time-step graupel tendency" "kgkg-1s-1" +state real QV_HADV ikj misc 1 - h{19} "qv_hadv" "instantaneous QV Horizontal advection" "kgkg-1" +state real QC_HADV ikj misc 1 - h{19} "qc_hadv" "instantaneous QC Horizontal advection" "kgkg-1" +state real QR_HADV ikj misc 1 - h{19} "qr_hadv" "instantaneous QR Horizontal advection" "kgkg-1" +state real QS_HADV ikj misc 1 - h{19} "qs_hadv" "instantaneous QS Horizontal advection" "kgkg-1" +state real QI_HADV ikj misc 1 - h{19} "qi_hadv" "instantaneous QI Horizontal advection" "kgkg-1" +state real QH_HADV ikj misc 1 - h{19} "qh_hadv" "instantaneous QH Horizontal advection" "kgkg-1" +state real QG_HADV ikj misc 1 - h{19} "qg_hadv" "instantaneous QG Horizontal advection" "kgkg-1" +state real QV_ZADV ikj misc 1 - h{19} "qv_zadv" "instantaneous QV Vertical advection" "kgkg-1" +state real QC_ZADV ikj misc 1 - h{19} "qc_zadv" "instantaneous QC Vertical advection" "kgkg-1" +state real QR_ZADV ikj misc 1 - h{19} "qr_zadv" "instantaneous QR Vertical advection" "kgkg-1" +state real QS_ZADV ikj misc 1 - h{19} "qs_zadv" "instantaneous QS Vertical advection" "kgkg-1" +state real QI_ZADV ikj misc 1 - h{19} "qi_zadv" "instantaneous QI Vertical advection" "kgkg-1" +state real QH_ZADV ikj misc 1 - h{19} "qh_zadv" "instantaneous QH Vertical advection" "kgkg-1" +state real QG_ZADV ikj misc 1 - h{19} "qg_zadv" "instantaneous QG Vertical advection" "kgkg-1" +state real CDXDIABH ikj misc 1 - h{19} "cdxdiabh" "diabatic heating from Micro-Physics" "K" +endif +ifndef INSTVALS state real QVTTEND ikj misc 1 - - "qvttend" "inter time-step water vapor tendency" "kgkg-1s-1" state real QCTTEND ikj misc 1 - - "qcttend" "inter time-step cloud tendency" "kgkg-1s-1" state real QRTTEND ikj misc 1 - - "qrttend" "inter time-step rain tendency" "kgkg-1s-1" @@ -322,53 +343,56 @@ state real QSTTEND ikj misc 1 - - "q state real QITTEND ikj misc 1 - - "qittend" "inter time-step ice tendency" "kgkg-1s-1" state real QHTTEND ikj misc 1 - - "qhttend" "inter time-step hail tendency" "kgkg-1s-1" state real QGTTEND ikj misc 1 - - "qgttend" "inter time-step graupel tendency" "kgkg-1s-1" -state real QV_HADV ikj misc 1 - - "qv_hadv" "Instantaneous QV Horizontal advection" "kgkg-1" -state real QC_HADV ikj misc 1 - - "qc_hadv" "Instantaneous QC Horizontal advection" "kgkg-1" -state real QR_HADV ikj misc 1 - - "qr_hadv" "Instantaneous QR Horizontal advection" "kgkg-1" -state real QS_HADV ikj misc 1 - - "qs_hadv" "Instantaneous QS Horizontal advection" "kgkg-1" -state real QI_HADV ikj misc 1 - - "qi_hadv" "Instantaneous QI Horizontal advection" "kgkg-1" -state real QH_HADV ikj misc 1 - - "qh_hadv" "Instantaneous QH Horizontal advection" "kgkg-1" -state real QG_HADV ikj misc 1 - - "qg_hadv" "Instantaneous QG Horizontal advection" "kgkg-1" -state real QV_ZADV ikj misc 1 - - "qv_zadv" "Instantaneous QV Vertical advection" "kgkg-1" -state real QC_ZADV ikj misc 1 - - "qc_zadv" "Instantaneous QC Vertical advection" "kgkg-1" -state real QR_ZADV ikj misc 1 - - "qr_zadv" "Instantaneous QR Vertical advection" "kgkg-1" -state real QS_ZADV ikj misc 1 - - "qs_zadv" "Instantaneous QS Vertical advection" "kgkg-1" -state real QI_ZADV ikj misc 1 - - "qi_zadv" "Instantaneous QI Vertical advection" "kgkg-1" -state real QH_ZADV ikj misc 1 - - "qh_zadv" "Instantaneous QH Vertical advection" "kgkg-1" -state real QG_ZADV ikj misc 1 - - "qg_zadv" "Instantaneous QG Vertical advection" "kgkg-1" +state real QV_HADV ikj misc 1 - - "qv_hadv" "instantaneous QV Horizontal advection" "kgkg-1" +state real QC_HADV ikj misc 1 - - "qc_hadv" "instantaneous QC Horizontal advection" "kgkg-1" +state real QR_HADV ikj misc 1 - - "qr_hadv" "instantaneous QR Horizontal advection" "kgkg-1" +state real QS_HADV ikj misc 1 - - "qs_hadv" "instantaneous QS Horizontal advection" "kgkg-1" +state real QI_HADV ikj misc 1 - - "qi_hadv" "instantaneous QI Horizontal advection" "kgkg-1" +state real QH_HADV ikj misc 1 - - "qh_hadv" "instantaneous QH Horizontal advection" "kgkg-1" +state real QG_HADV ikj misc 1 - - "qg_hadv" "instantaneous QG Horizontal advection" "kgkg-1" +state real QV_ZADV ikj misc 1 - - "qv_zadv" "instantaneous QV Vertical advection" "kgkg-1" +state real QC_ZADV ikj misc 1 - - "qc_zadv" "instantaneous QC Vertical advection" "kgkg-1" +state real QR_ZADV ikj misc 1 - - "qr_zadv" "instantaneous QR Vertical advection" "kgkg-1" +state real QS_ZADV ikj misc 1 - - "qs_zadv" "instantaneous QS Vertical advection" "kgkg-1" +state real QI_ZADV ikj misc 1 - - "qi_zadv" "instantaneous QI Vertical advection" "kgkg-1" +state real QH_ZADV ikj misc 1 - - "qh_zadv" "instantaneous QH Vertical advection" "kgkg-1" +state real QG_ZADV ikj misc 1 - - "qg_zadv" "instantaneous QG Vertical advection" "kgkg-1" endif ##state real H_TENDENCY ikj misc 1 - - "h_tendency" "Horizontal tendency" "" ##state real Z_TENDENCY ikj misc 1 - - "z_tendency" "Vertical tendency" "" +endif ### L. Fita, January 2018. CIMA ### These variables would make too much output for a climate run. Include them if you need instantaneaous values (change '##INST##state' by 'state') ifdef INSTVALS -state real CLT ij misc 1 - h{19} "clt" "TOTAL CLOUDINESS" "1" -state real CLL ij misc 1 - h{19} "cll" "LOW-LEVEL CLOUDINESS (p >= 68000 Pa)" "1" -state real CLM ij misc 1 - h{19} "clm" "MID-LEVEL CLOUDINESS (44000 <= p < 68000 Pa)" "1" -state real CLH ij misc 1 - h{19} "clh" "HIGH-LEVEL CLOUDINESS (p < 44000 Pa)" "1" -state real WSGS ij misc 1 - h{19} "wsgs" "near-surface wind speed of gust" "ms-1" -state real USGS ij misc 1 - h{19} "usgs" "Eastward near-surface wind speed of gust" "ms-1" -state real VSGS ij misc 1 - h{19} "vsgs" "Northward near-surface wind speed of gust" "ms-1" -state integer GUSTPOINT ij misc 1 - h{19} "gustpoint" "whether grid point got wind gust (1) or not (0)" "1" -state real TOTWSGS ij misc 1 - h{19} "totwsgs" "total (TKE + H. pr) near-surface wind speed of gust" "ms-1" -state real TOTUSGS ij misc 1 - h{19} "totusgs" "total Eastward near-surface wind speed of gust" "ms-1" -state real TOTVSGS ij misc 1 - h{19} "totvsgs" "total Northward near-surface wind speed of gust" "ms-1" -state integer TOTGUSTPOINT ij misc 1 - h{19} "totgustpoint" "whether grid point got total wind gust (1) or not (0)" "1" -state real POTEVAPO ij misc 1 - h{19} "potevapo" "potential evapotranspiration" "kgm-2s-1" -state real POTEVAPOGEN ij misc 1 - h{19} "potevapogen" "generic potential evapotranspiration" "kgm-2s-1" -state real CDXCAPE ij misc 1 - h{19} "cdxcape" "convective available potential energy" "Jkg-1" -state real CIN ij misc 1 - h{19} "cin" "convective inhibition" "Jkg-1" -state real LFCP ij misc 1 - h{19} "lfcp" "PRESSURE LEVEL FREE CONVECTION" "Pa" -state real LFCZ ij misc 1 - h{19} "lfcz" "HEIGHT LEVEL FREE CONVECTION" "m" -state real LI ij misc 1 - h{19} "li" "LIFTED INDEX" "1" -state integer fog ij misc 1 - h{19} "fog" "Whether there is fog (1: yes [vis < 1km]; 0: not)" "-" -state real fogvisblty ij misc 1 - h{19} "fogvisblty" "visibility inside fog" "km" -state real tds ij misc 1 - h{19} "tds" "surface dew point temperature" "K" -state real tws ij misc 1 - h{19} "tws" "surface wet-bulb temperature" "K" -state real CDXDIABH ij misc 1 - h{19} "cdxdiabh" "diabatic heating from Micro-Physics" "K" +state real CLT ij misc 1 - h{19} "clt" "instantaneous TOTAL CLOUDINESS" "1" +state real CLL ij misc 1 - h{19} "cll" "instantaneous LOW-LEVEL CLOUDINESS (p >= 68000 Pa)" "1" +state real CLM ij misc 1 - h{19} "clm" "instantaneous MID-LEVEL CLOUDINESS (44000 <= p < 68000 Pa)" "1" +state real CLH ij misc 1 - h{19} "clh" "instantaneous HIGH-LEVEL CLOUDINESS (p < 44000 Pa)" "1" +state real WSGS ij misc 1 - h{19} "wsgs" "instantaneous near-surface wind speed of gust" "ms-1" +state real USGS ij misc 1 - h{19} "usgs" "instantaneous Eastward near-surface wind speed of gust" "ms-1" +state real VSGS ij misc 1 - h{19} "vsgs" "instantaneous Northward near-surface wind speed of gust" "ms-1" +state integer GUSTPOINT ij misc 1 - h{19} "gustpoint" "instantaneous whether grid point got wind gust (1) or not (0)" "1" +state real TOTWSGS ij misc 1 - h{19} "totwsgs" "instantaneous total (TKE + H. pr) near-surface wind speed of gust" "ms-1" +state real TOTUSGS ij misc 1 - h{19} "totusgs" "instantaneous total Eastward near-surface wind speed of gust" "ms-1" +state real TOTVSGS ij misc 1 - h{19} "totvsgs" "instantaneous total Northward near-surface wind speed of gust" "ms-1" +state integer TOTGUSTPOINT ij misc 1 - h{19} "totgustpoint" "instantaneous whether grid point got total wind gust (1) or not (0)" "1" +state real POTEVAPO ij misc 1 - h{19} "potevapo" "instantaneous potential evapotranspiration" "kgm-2s-1" +state real POTEVAPOGEN ij misc 1 - h{19} "potevapogen" "instantaneous generic potential evapotranspiration" "kgm-2s-1" +ifdef CDXWRF_1 +state real CDXCAPE ij misc 1 - h{19} "cdxcape" "instantaneous convective available potential energy" "Jkg-1" +state real CIN ij misc 1 - h{19} "cin" "instantaneous convective inhibition" "Jkg-1" +state real LFCP ij misc 1 - h{19} "lfcp" "instantaneous PRESSURE LEVEL FREE CONVECTION" "Pa" +state real LFCZ ij misc 1 - h{19} "lfcz" "instantaneous HEIGHT LEVEL FREE CONVECTION" "m" +state real LI ij misc 1 - h{19} "li" "instantaneous LIFTED INDEX" "1" +state real IUT ij misc 1 - h{19} "iut" "instantaneous vertically integrated eastward transport of water vapour" "kgm-1s-1" +state real IVT ij misc 1 - h{19} "ivt" "instantaneous vertically integrated northward transport of water vapour" "kgm-1s-1" +endif +ifdef CDXWRF_2 +state integer fog ij misc 1 - h{19} "fog" "instantaneous Whether there is fog (1: yes [vis < 1km]; 0: not)" "-" +state real fogvisblty ij misc 1 - h{19} "fogvisblty" "instantaneous visibility inside fog" "km" +state real tds ij misc 1 - h{19} "tds" "instantaneous surface dew point temperature" "K" +state real tws ij misc 1 - h{19} "tws" "instantaneous surface wet-bulb temperature" "K" +endif endif - -# Water budget related. Only to work if CDXWRF = 3. - diff --git a/arch/configure.defaults b/arch/configure.defaults index b0dc9ad502..a11b427927 100644 --- a/arch/configure.defaults +++ b/arch/configure.defaults @@ -2092,9 +2092,8 @@ RWORDSIZE = CONFIGURE_RWORDSIZE AMDARCH = -march=znver3 AMDMATHLIB = -fveclib=AMDLIBM -AMDLDFLAGS = -Wl,-mllvm -Wl,-enable-loop-reversal -Wl,-mllvm -Wl,-enable-gather -Wl,-mllvm -Wl,-vectorize-noncontigous-memory-aggressively +AMDLDFLAGS = -Wl,-mllvm -Wl,-enable-loop-reversal -Wl,-mllvm -Wl,-enable-gather -Wl,-mllvm -Wl,-vectorize-non-contiguous-memory-aggressively AMDFCFLAGS = - PROMOTION = #-fdefault-real-8 ARCH_LOCAL = -DNONSTANDARD_SYSTEM_SUBR -DWRF_USE_CLM -DRPC_TYPES=2 CFLAGS_LOCAL = -w -c -m64 -Ofast -ffast-math $(AMDARCH) diff --git a/chem/module_optical_averaging.F b/chem/module_optical_averaging.F index 0d0b3752c1..dcbf5c5f4d 100644 --- a/chem/module_optical_averaging.F +++ b/chem/module_optical_averaging.F @@ -5084,19 +5084,23 @@ subroutine mieaer( & ! check for refr and refi outside of lookup table range to prevent unstable extrapolation 'binterp' below if (refr < refrtabsw(1,ns)) then refr = refrtabsw(1,ns) - write(*,*) 'Warning: refr is smaller than lookup table range and reset to minimum bound at SW band ', ns + if (number_bin_col(m,klevel).ge.1.e-10) & + write(*,*) 'Warning: refr is smaller than lookup table range and reset to minimum bound at SW band ', ns endif if (refr > refrtabsw(prefr,ns)) then refr = refrtabsw(prefr,ns) - write(*,*) 'Warning: refr is larger than lookup table range and reset to maximum bound at SW band ', ns + if (number_bin_col(m,klevel).ge.1.e-10) & + write(*,*) 'Warning: refr is larger than lookup table range and reset to maximum bound at SW band ', ns endif if (refi < refitabsw(1,ns)) then refi = refitabsw(1,ns) - write(*,*) 'Warning: refi is smaller than lookup table range and reset to minimum bound at SW band ', ns + if (number_bin_col(m,klevel).ge.1.e-10) & + write(*,*) 'Warning: refi is smaller than lookup table range and reset to minimum bound at SW band ', ns endif if (refi > refitabsw(prefi,ns)) then + if (number_bin_col(m,klevel).ge.1.e-10 .and. abs(refi).gt.1.e-12) & + write(*,*) 'Warning: refi is larger than lookup table range and reset to maximum bound at SW band ', ns refi = refitabsw(prefi,ns) - write(*,*) 'Warning: refi is larger than lookup table range and reset to maximum bound at SW band ', ns endif ! interpolate coefficients linear in refractive index @@ -5359,19 +5363,23 @@ subroutine mieaer( & ! check for refr and refi outside of lookup table range to prevent unstable extrapolation 'binterp' below if (refr < refrtablw(1,ns)) then refr = refrtablw(1,ns) - write(*,*) 'Warning: refr is smaller than lookup table range and reset to minimum bound at LW band ', ns + if (number_bin_col(m,klevel).ge.1.e-10) & + write(*,*) 'Warning: refr is smaller than lookup table range and reset to minimum bound at LW band ', ns endif if (refr > refrtablw(prefr,ns)) then refr = refrtablw(prefr,ns) - write(*,*) 'Warning: refr is larger than lookup table range and reset to maximum bound at LW band ', ns + if (number_bin_col(m,klevel).ge.1.e-10) & + write(*,*) 'Warning: refr is larger than lookup table range and reset to maximum bound at LW band ', ns endif - if (refi < refitablw(1,ns)) then + if (refi < refitablw(1,ns)) then refi = refitablw(1,ns) - write(*,*) 'Warning: refi is smaller than lookup table range and reset to minimum bound at LW band ', ns + if (number_bin_col(m,klevel).ge.1.e-10) & + write(*,*) 'Warning: refi is smaller than lookup table range and reset to minimum bound at LW band ', ns endif if (refi > refitablw(prefi,ns)) then + if (number_bin_col(m,klevel).ge.1.e-10 .and. abs(refi).gt.1.e-12) & + write(*,*) 'Warning: refi is larger than lookup table range and reset to maximum bound at LW band ', ns refi = refitablw(prefi,ns) - write(*,*) 'Warning: refi is larger than lookup table range and reset to maximum bound at LW band ', ns endif ! interpolate coefficients linear in refractive index diff --git a/doc/README.cfbm b/doc/README.cfbm new file mode 100644 index 0000000000..244efd4ce6 --- /dev/null +++ b/doc/README.cfbm @@ -0,0 +1,35 @@ + +WRF should be compiled with CMake in order to use the Community Fire Behavior model: +git clone --recurse-submodules https://github.com/wrf-model/WRF +./configure_new -- -DENABLE_CFBM=ON +./compile_new + +To activate the model you need to set ifire = 1 in the fire block of the WRF namelist. +This instructs WRF to read the CFBM namelist: namelist.cfbm + +In a given WRF simulation, CFBM can only be activated in one domain. + +It only works with Lambert-Conformal projections at this time. + +This is an example of the namelist.cfbm: + + &time + dt = 2 ! WRF time step [s] for the nest CFBM is active + interval_output = 2 ! Frequency [s] saving fire output +/ + + &atm + / + + &fire + fire_num_ignitions = 1, + fire_ignition_ros1 = 1.2, + fire_ignition_start_lat1 = 39.68, + fire_ignition_start_lon1 = -103.58, + fire_ignition_end_lat1 = 39.68, + fire_ignition_end_lon1 = -103.58, + fire_ignition_radius1 = 2000.0, + fire_ignition_start_time1 = 6.0, + fire_ignition_end_time1 = 60.0, +/ + diff --git a/doc/README.cordex b/doc/README.cordex index 55e4a5cff0..dfd8a99aee 100644 --- a/doc/README.cordex +++ b/doc/README.cordex @@ -1,5 +1,5 @@ -******* WRF module to compute atmospheric model diagnostics required by CORDEX on WRFV3.9.1.1 -* L. Fita, CIMA. December 2017 +******* WRF module to compute atmospheric model diagnostics required by CORDEX during model integration +* L. Fita, CIMA. May 2026 See for more detail wiki page: http://wiki.cima.fcen.uba.ar/mediawiki/index.php/CDXWRF @@ -10,7 +10,7 @@ This module computes that variables required by CORDEX which usually are post-pr Efficiency constrains make necessary to introduce a system of compilation in order to do not overload WRF execution. Long simulations are usually used in CORDEX experiments, thus model temporal efficiency is even more important. In order to achieve this, a new pre-compilation variable has been added. According to the value of this giving variable and accordingly modifications in the 'Registry/registry.cordex' file certain variables are computed and provided at the output. This is not very 'WRF'-friendly, since usually all options are directly available with just a single compilation, but it has been shown that the large amount of variables introduced by this module, it slows down too much (up to 40% of additional time) the running time of the model. See below for more details (INSTALLATION section) -* Subroutines/Functions +* Main subroutines/Functions - module_diag_cordex: Main subroutine to compute CORDEX required variables - module_diagvar_cordex: subroutine with the specific computation of CORDEX required variables - registry.cordex: required definition of the new CORDEX variables @@ -55,7 +55,7 @@ Efficiency constrains make necessary to introduce a system of compilation in ord - zmla: pbl height following a generic method [m] - colmax: high-frequency maximum radar reflectivity in the column [dBz] (on auxhist8, wrfhfcdx) -** Only if 'INSTVALS' modifications are made +** Only if 'instvals' is activated when executing configure - cape: Convective Available Potential Energy [Jkg-1] - cin: Convective inhibition [Jkg-1] - zlfc: Height at the Level of free convection [m] @@ -127,7 +127,7 @@ Efficiency constrains make necessary to introduce a system of compilation in ord - mrfsomean: mean total frozen water content of soil layer [kgm-2] - mrfsosmean: mean frozen water content in Upper Portion of soil column (0-10cm) [kgm-2] -* Only if 'INSTVALS' modifications are made +* Only if 'instvals' are activated - clt: total cloud cover [%] - cll: low-level cloud cover [%] - clm: mid-level cloud cover [%] @@ -178,7 +178,7 @@ Efficiency constrains make necessary to introduce a system of compilation in ord - wbacf{l/m/h}, wbacf[c/r/s/i/g/h]{l/m/h}: Water-budget vertically integrated accumulated horizontal advection for water vapour, cloud, rain, snow, ice, graupel, hail at low, medium and high levels (same as cloudiness) [mm] - wbacz{l/m/h}, wbacz[c/r/s/i/g/h]{l/m/h}: Water-budget vertically integrated accumulated vertical advection for water vapour, cloud, rain, snow, ice, graupel, hail at low, medium and high levels (same as cloudiness) [mm] -** Only if 'INSTVALS' modifications are made +** Only if 'instvals' is activated - fog: whether point is inside fog (vis < 1km) [1] - vis: visibility inside fog [km] - tds: 2m dew point temperature [K] @@ -192,56 +192,23 @@ Efficiency constrains make necessary to introduce a system of compilation in ord * INSTALLATION These steps must be followed prior the re-compilation of the WRF model (it requires a clean -a) - 1.- Getting the code from the fork in a $INSTDIR - $ cd $INSTDIR - $ git clone git@github.com:LluisFB/WRF.git - $ cd WRF + 1.- In order to compile the module user needs to choose a level of output CDXWRF=[CDXWRF_N] that meets the requirements. Amuont of output is incremental, so a given number for [CDXWRF_N] will provide the variables associated to the previous values. So for example [CDXWRF_N] = 3, will provide all the variables for values 1, 2 and 3 + + 2.- When executing the configuration of the compilation, user needs to provide a given value to the variable cdxwrf as follows: + $ ./configure -cordex [CDXWRF_N] - 2.- Getting submodules - $ git submodule update --init --recursive + Accordingly to the value [CXWRF_N] given to the pre-compilation variable CDXWRF one obtains: + - Without adding the variable: all CORDEX 'Core' variables + - CDXWRF=1 CORDEX 'Tier' variables: clivg, clivh, zmla, [cape/cin/zlfc/plfc/lidx]{min/max/mean} + - CDXWRF=2 The same as with CDXWRF=1 and additional variables: ua, va, ws, ta, press, zg, hur, hus, tfog, fogvisbltymin, fogvisbltymax, fogvisbltymean, tdsmin, tdsmax, tdsmean + - CDXWRF=3 The same as with CDXWRF=1,2 and additional residence-time variables + - CDXWRF=4 The same as with CDXWRF=1,2 and 3 and the Water-Budget relarted ones: wbacdiabh, wbacpw, wbacpw[c/r/s/i/g/h], wbacf, wbacf[c/r/s/i/g/h], wbacz, wbacz[c/r/s/i/g/h], wbacdiabh{l/m/h}, wbacpw{l/m/h}, wbacpw[c/r/s/i/g/h]{l/m/h}, wbacf{l/m/h}, wbacf[c/r/s/i/g/h]{l/m/h}, wbacz{l/m/h}, wbacz[c/r/s/i/g/h]{l/m/h} - 3.- Check to the right branch - $ git checkout cdxwrf - - ** NOTE **: assuming that you got the right submodules [e.g. noahMP] (when downloading the fork) !! - - 4.- Clean the code (in order to avoid to run again configure one can make a copy of the 'configure.wrf' and recover it after the clean, otherwise it is erased) - $ cp configure.wrf configure.cordex.wrf - $ ./clean -a - $ cp configure.cordex.wrf configure.wrf - - 5.- edit the `configure.wrf' and add the line (after the line -DNETCDF and/or -DCLWRFGHG) - -DCORDEXDIAG \ - 6.- Set up (or not) the pre-compilation variable CDXWRF (after the line -DCORDEXDIAG) - -DCDXWRF=[value] \ - - Accordingly to the value given to the pre-compilation variable CDXWRF one obtains: - - Without adding the variable: all CORDEX 'Core' variables - - CDXWRF=1 CORDEX 'Tier' variables: clivg, clivh, zmla, [cape/cin/zlfc/plfc/lidx]{min/max/mean} - - CDXWRF=2 The same as with CDXWRF=1 and additional variables: ua, va, ws, ta, press, zg, hur, hus, tfog, fogvisbltymin, fogvisbltymax, fogvisbltymean, tdsmin, tdsmax, tdsmean - - CDXWRF=3 The same as with CDXWRF=1,2 and additional residence-time variables - - CDXWRF=4 The same as with CDXWRF=1,2 and 3 and the Water-Budget relarted ones: wbacdiabh, wbacpw, wbacpw[c/r/s/i/g/h], wbacf, wbacf[c/r/s/i/g/h], wbacz, wbacz[c/r/s/i/g/h], wbacdiabh{l/m/h}, wbacpw{l/m/h}, wbacpw[c/r/s/i/g/h]{l/m/h}, wbacf{l/m/h}, wbacf[c/r/s/i/g/h]{l/m/h}, wbacz{l/m/h}, wbacz[c/r/s/i/g/h]{l/m/h} - Simultanesouly, one needs to: - a.- make a copy of Registry/registry.cordex - $ cp Registry/registry.cordex Registry/registry.cordex_comp - b.- modify the Registry/registry.cordex_comp accordingly to the value of CDXWRF: - - User must provide a level for CDXWRF. - > Adding CDXWRF=1, one needs to remove the comment ##CDXWRF1## at the beginning of the line of the definition of certain variables - > Adding CDXWRF=2, one needs to remove the comment ##CDXWRF1## and ##CDXWRF2## at the beginning of the line of the definition of certain variables - > Adding CDXWRF=3, one needs to remove the comment ##CDXWRF1##, ##CDXWRF2## and ##CDXWRF3## at the beginning of the line of the definition of certain variables - > Adding CDXWRF=4, one needs to remove the comment ##CDXWRF1##, ##CDXWRF2##, ##CDXWRF3## and ##CDXWRF4## at the beginning of the line of the definition of certain variables - - Example for CDXWRF=3 - $ sed -i 's/##CDXWRF1##//g' Registry/registry.cordex_comp - $ sed -i 's/##CDXWRF2##//g' Registry/registry.cordex_comp - $ sed -i 's/##CDXWRF3##//g' Registry/registry.cordex_comp - - - Additionally, one can also get the instantaneous values for the variables which only certain statistics (accumulation, minimum, mean, ...) are provided. In order to get them, one need to: - a.- Search in 'phys/module_diagnostics_driver.F' and 'phys/module_diag_cordex.F' the lines of code marked with 'INSTVALS' and change accordingly. - b.- Modify Registry/registry.cordex accordingly (removing ##INST## at the beginning of the line of the definition of certain variables, and adding 'h9' to certain others) - c.- re-compile - - 7.- compile as always - · ./compile em_real >& compile.log + * NOTE: Additionally, one can also get the instantaneous values for the variables which only certain statistics (accumulation, minimum, mean, ...) are provided. In order to get them, one need to add -instvals argument when configuring the compilation. + $ ./configure -cordex [CDXWRF_N] -instvals + + 3.- compile as always + ./compile em_real >& compile.log NOTE: after any change into a `Registry' related file, one needs to before the compilation refresh entirely the code throughout $ ./clean -a diff --git a/dyn_em/module_first_rk_step_part1.F b/dyn_em/module_first_rk_step_part1.F index d7a1557124..2f1809c815 100644 --- a/dyn_em/module_first_rk_step_part1.F +++ b/dyn_em/module_first_rk_step_part1.F @@ -56,7 +56,7 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & #ifdef DM_PARALLEL USE module_dm, ONLY : local_communicator, mytask, ntasks, ntasks_x, ntasks_y, local_communicator_periodic, wrf_dm_maxval USE module_comm_dm, ONLY : halo_em_phys_a_sub,halo_em_fdda_sfc_sub,halo_pwp_sub,halo_em_chem_e_3_sub, & - halo_em_chem_e_5_sub, halo_em_hydro_noahmp_sub + halo_em_chem_e_5_sub, halo_em_hydro_noahmp_sub, halo_em_tracer_e_5_sub, halo_em_cfbm_sub #if ( WRFPLUS == 1 ) USE module_comm_dm, ONLY : halo_em_phys_a_bl_surf_sub #endif @@ -64,10 +64,11 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & USE module_utility #ifdef WRF_CFBM - ! CFBM + ! CFBM modules use advance_mod, only : Advance_state - use wrf_mod, only : Interp_wrf2dvar_to_cfbm, Interp_wrfwinds_to_cfbm - use wrf_mod, only : Provide_atm_feedback + use wrf_mod, only : Interp_wrfwinds_to_cfbm, Provide_atm_feedback + use interp_mod, only : VINTERP_WINDS_FROM_3D_WINDS, VINTERP_WINDS_FROM_10M_WINDS + use coupling_mod, only : Interp_horizontal #endif IMPLICIT NONE @@ -1270,7 +1271,6 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & & ,tke_budget=config_flags%tke_budget & & ,bl_mynn_cloudpdf=config_flags%bl_mynn_cloudpdf & & ,bl_mynn_mixlength=config_flags%bl_mynn_mixlength & - & ,icloud_bl=config_flags%icloud_bl & & ,qc_bl=grid%qc_bl,qi_bl=grid%qi_bl,cldfra_bl=grid%cldfra_bl & & ,bl_mynn_edmf=config_flags%bl_mynn_edmf & ,bl_mynn_edmf_dd=config_flags%bl_mynn_edmf_dd & @@ -1293,7 +1293,9 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & & ,qcg=grid%qcg, grav_settling=config_flags%grav_settling & ! & ,K_m=grid%K_m, K_h=grid%K_h, K_q=grid%K_q & & ,vdfg=grid%vdfg,maxwidth=grid%maxwidth,maxMF=grid%maxmf & - & ,ztop_plume=grid%ztop_plume & + & ,ztop_plume=grid%ztop_plume,maxmf_dd=grid%maxmf_dd & + & ,maxtkeprod=grid%maxtkeprod,maxwidth_dd=grid%maxwidth_dd & + & ,cldtop_cooling=grid%cldtop_cooling,ent_eff=grid%ent_eff & & ,excess_h=grid%excess_h,excess_q=grid%excess_q & & ,spp_pbl=config_flags%spp_pbl & & ,pattern_spp_pbl=grid%pattern_spp_pbl & @@ -1408,38 +1410,68 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & #ifdef WRF_CFBM case (1) ! Interpolate from Atm to fire grid - call Interp_wrf2dvar_to_cfbm (grid%znt, ims, ime, jms, jme, & - grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & - grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & - grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fz0) - - call Interp_wrf2dvar_to_cfbm (grid%t2, ims, ime, jms, jme, & +#ifdef DM_PARALLEL +# include "HALO_EM_CFBM.inc" +#endif + call Interp_horizontal (grid%znt, grid%fire_state%proj, ims, ime, jms, jme, & grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & - grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & - grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fire_t2) + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%fz0) - call Interp_wrf2dvar_to_cfbm (grid%q2, ims, ime, jms, jme, & + call Interp_horizontal (grid%t2, grid%fire_state%proj, ims, ime, jms, jme, & grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & - grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & - grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fire_q2) + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%fire_t2) - call Interp_wrf2dvar_to_cfbm (grid%psfc, ims, ime, jms, jme, & + call Interp_horizontal (grid%q2, grid%fire_state%proj, ims, ime, jms, jme, & grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & - grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & - grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fire_psfc) + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%fire_q2) - call Interp_wrf2dvar_to_cfbm (grid%rainc + grid%rainnc, ims, ime, jms, jme, & + call Interp_horizontal (grid%psfc, grid%fire_state%proj, ims, ime, jms, jme, & grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & - grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & - grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fire_rain) + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%fire_psfc) - call Interp_wrfwinds_to_cfbm (grid%u_phy, grid%v_phy, grid%z_at_w, ims, ime, kms, kme, jms, jme, & + call Interp_horizontal (grid%rainc + grid%rainnc, grid%fire_state%proj, ims, ime, jms, jme, & grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & - grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & - grid%fire_state%kfds, grid%fire_state%kfde, & - grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fz0, & - grid%cfbm_namelist%fire_lsm_zcoupling, grid%cfbm_namelist%fire_lsm_zcoupling_ref, & - grid%cfbm_namelist%fire_wind_height, grid%fire_state%uf, grid%fire_state%vf) + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%fire_rain) + + Select_vinterp_opt: select case (grid%cfbm_namelist%wind_vinterp_opt) + case (VINTERP_WINDS_FROM_3D_WINDS) + call Interp_wrfwinds_to_cfbm (grid%u_phy, grid%v_phy, grid%z_at_w, ims, ime, kms, kme, jms, jme, & + grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & + grid%fire_state%ifps, grid%fire_state%ifpe, grid%fire_state%jfps, grid%fire_state%jfpe, & + grid%fire_state%kfds, grid%fire_state%kfde, & + grid%fire_state%lats, grid%fire_state%lons, grid%fire_state%proj, grid%fire_state%fz0, & + grid%cfbm_namelist%fire_lsm_zcoupling, grid%cfbm_namelist%fire_lsm_zcoupling_ref, & + grid%cfbm_namelist%fire_wind_height, grid%fire_state%uf, grid%fire_state%vf) + + case (VINTERP_WINDS_FROM_10M_WINDS) + call Interp_horizontal (grid%u10, grid%fire_state%proj, ims, ime, jms, jme, & + grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%uf) + + call Interp_horizontal (grid%v10, grid%fire_state%proj, ims, ime, jms, jme, & + grid%fire_state%ifms, grid%fire_state%ifme, grid%fire_state%jfms, grid%fire_state%jfme, & + grid%fire_state%num_tiles, grid%fire_state%i_start, grid%fire_state%i_end, grid%fire_state%j_start, & + grid%fire_state%j_end, grid%cfbm_namelist%hinterp_opt, grid%fire_state%lats, grid%fire_state%lons, & + grid%fire_state%vf) + + call grid%fire_state%Apply_wafs () + + case default + call wrf_error_fatal ('The vertical wind interpolation option does not exist in CFBM') + + end select Select_vinterp_opt call Advance_state (grid = grid%fire_state, config_flags = grid%cfbm_namelist) @@ -1475,6 +1507,12 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & grid%rthfrten, grid%rqvfrten) ! theta and Qv tendencies end do !$OMP END PARALLEL DO + + if (config_flags%tracer_opt == 3) then +#ifdef DM_PARALLEL +# include "HALO_EM_TRACER_E_5.inc" +#endif + end if else call wrf_error_fatal ('num_tiles in WRF and CFBM should match') endif diff --git a/dyn_em/module_initialize_real.F b/dyn_em/module_initialize_real.F index 275157b4e4..9fd521a74b 100644 --- a/dyn_em/module_initialize_real.F +++ b/dyn_em/module_initialize_real.F @@ -2744,7 +2744,8 @@ SUBROUTINE init_domain_rk ( grid & .OR. (config_flags%mp_physics .EQ. TEMPO .AND. config_flags%tempo_aerosolaware .EQ. 1) & .OR. config_flags%mp_physics .EQ. RCON_MP_SCHEME) .and. & config_flags%wif_input_opt .EQ. 0 ) THEN - CALL wrf_error_fatal ('wif_input_opt=0 but mp_physics=28 or mp_physics=29' ) + CALL wrf_error_fatal ( & + 'wif_input_opt=0 but mp_physics=28, 29, or 88' ) END IF if_thompsonaero_3d !========================================================================================= diff --git a/dyn_em/start_em.F b/dyn_em/start_em.F index b5c886714d..524f338dd3 100644 --- a/dyn_em/start_em.F +++ b/dyn_em/start_em.F @@ -20,11 +20,12 @@ SUBROUTINE start_domain_em ( grid, allowed_to_read & USE module_dm, ONLY : wrf_dm_min_real, wrf_dm_max_real, wrf_dm_maxval, & ntasks_x, ntasks_y, & local_communicator_periodic, local_communicator, mytask, ntasks + #else USE module_dm, ONLY : wrf_dm_min_real, wrf_dm_max_real #endif USE module_comm_dm - USE module_llxy, ONLY : proj_cassini + USE module_llxy, ONLY : proj_cassini, PROJ_LC USE module_physics_init USE NoahmpGroundwaterInitMod, ONLY: NoahmpGroundwaterInitMain USE module_lightning_driver, ONLY : lightning_init @@ -145,6 +146,7 @@ SUBROUTINE start_domain_em ( grid, allowed_to_read & ifms, ifme, jfms, jfme, kfms, kfme, & ifps, ifpe, jfps, jfpe, kfps, kfpe real :: cen_lat, cen_lon, truelat1, truelat2, stand_lon + logical :: is_cfbm_initialized = .false. #endif @@ -2246,6 +2248,20 @@ SUBROUTINE start_domain_em ( grid, allowed_to_read & ! Community Fire Behavior model case (1) #ifdef WRF_CFBM + Init_cfbm: if (.not. is_cfbm_initialized) then +#ifdef DM_PARALLEL + ! Broadcast the namelist + call grid%cfbm_namelist%Broadcast_nml (grid%communicator) +#endif + + if (abs (grid%dt - grid%cfbm_namelist%dt) > 1.0e-5) then + call wrf_error_fatal ('CFBM dt should match WRF dt for the fire domain') + end if + + if (config_flags%map_proj /= PROJ_LC) then + call wrf_error_fatal ('CFBM only supports Lambert Conformal projections at this moment') + end if + ! Get fire dims call get_ijk_from_subgrid (grid, ifds, ifde, jfds, jfde, kfds, kfde, & ifms, ifme, jfms, jfme, kfms, kfme, & @@ -2280,7 +2296,11 @@ SUBROUTINE start_domain_em ( grid, allowed_to_read & sr_x = config_flags%sr_x, sr_y = config_flags%sr_y, & map_proj = config_flags%map_proj, cen_lat = cen_lat, cen_lon = cen_lon, truelat1 = truelat1, & truelat2 = truelat2, stand_lon = stand_lon, nfuel_cat = grid%nfuel_cat, zsf = grid%zsf, & - dzdxf = grid%dzdxf, dzdyf = grid%dzdyf) + dzdxf = grid%dzdxf, dzdyf = grid%dzdyf, communicator = grid%communicator) + + is_cfbm_initialized = .true. + + end if Init_cfbm #else call wrf_error_fatal ('CFBM requires compilation with CMake') #endif diff --git a/inc/version_decl b/inc/version_decl index a9eda5589f..72dcf84562 100644 --- a/inc/version_decl +++ b/inc/version_decl @@ -1 +1 @@ - CHARACTER (LEN=*), PARAMETER :: release_version = 'V4.7.1' + CHARACTER (LEN=*), PARAMETER :: release_version = 'V4.8.0' diff --git a/main/depend.common b/main/depend.common index f78d1d008d..3b3dc21149 100644 --- a/main/depend.common +++ b/main/depend.common @@ -2531,6 +2531,9 @@ module_microphysics_driver.o: \ module_mp_cammgmp_driver.o \ module_irrigation.o \ module_mp_fast_sbm.o \ + module_calc_lpi_new.o \ + module_ltng_pe.o \ + module_ltng_strokes.o \ ../frame/module_driver_constants.o \ ../frame/module_state_description.o \ ../frame/module_wrf_error.o \ diff --git a/phys/CMakeLists.txt b/phys/CMakeLists.txt index 75e7798beb..5fe11c1b2e 100644 --- a/phys/CMakeLists.txt +++ b/phys/CMakeLists.txt @@ -412,6 +412,7 @@ target_sources( # Shared physics physics_mmm/bl_gwdo.F90 physics_mmm/bl_ysu.F90 + physics_mmm/bl_shinhong.F90 physics_mmm/cu_ntiedtke.F90 physics_mmm/module_libmassv.F90 physics_mmm/mp_radar.F90 diff --git a/phys/GFL b/phys/GFL index ec9b3ba1ba..873b414670 160000 --- a/phys/GFL +++ b/phys/GFL @@ -1 +1 @@ -Subproject commit ec9b3ba1ba21a7c6f677a7a896a8676240ddaa31 +Subproject commit 873b414670323ae766ec4fe1f4b6b5ab891358ed diff --git a/phys/MYNN-EDMF b/phys/MYNN-EDMF index 43e0992eca..8f8e8b5902 160000 --- a/phys/MYNN-EDMF +++ b/phys/MYNN-EDMF @@ -1 +1 @@ -Subproject commit 43e0992ecab47892c6ec70c46d6b27af77ff2555 +Subproject commit 8f8e8b59021249949c04d1686edc5a02ca325627 diff --git a/phys/MYNN-SFC b/phys/MYNN-SFC index 9dd8e52ead..b66187f54a 160000 --- a/phys/MYNN-SFC +++ b/phys/MYNN-SFC @@ -1 +1 @@ -Subproject commit 9dd8e52ead0bc67e4da0cd41fea31956a3af9d54 +Subproject commit b66187f54a0c169e36592552634867474e603ac0 diff --git a/phys/Makefile b/phys/Makefile index 6988b5cddf..93cc5fcc37 100644 --- a/phys/Makefile +++ b/phys/Makefile @@ -493,7 +493,7 @@ submodules : fi @if [ \( ! -f module_mp_tempo_driver.F90 \) -o \( ! -f module_mp_tempo_main.F90 \) -o \ \( ! -f module_mp_tempo_cfgs.F90 \) -o \( ! -f module_mp_tempo_aerosols.F90 \) -o \ - \( ! -f module_mp_tempo_ml.F90 \) -o \( ! -f modulemp_tempo_diags.F90 \) -o \ + \( ! -f module_mp_tempo_ml.F90 \) -o \( ! -f module_mp_tempo_diags.F90 \) -o \ \( ! -f module_mp_tempo_utils.F90 \) -o \( ! -f module_mp_tempo_params.F90 \) ] ; then \ echo Pulling in TEMPO submodule ; \ ( cd .. ; git submodule update --init --recursive ) ; \ diff --git a/phys/TEMPO b/phys/TEMPO index 2f68b2ce13..881614f68a 160000 --- a/phys/TEMPO +++ b/phys/TEMPO @@ -1 +1 @@ -Subproject commit 2f68b2ce13a3f76fd126b787c0b97c1f4046f611 +Subproject commit 881614f68aa65de8f383d56139162049baf600fa diff --git a/phys/fire_behavior b/phys/fire_behavior index e3b13645c7..eb7758055a 160000 --- a/phys/fire_behavior +++ b/phys/fire_behavior @@ -1 +1 @@ -Subproject commit e3b13645c7e6638b6831179c068a3984726246b0 +Subproject commit eb7758055a5d436f022a634cae780f5512c42bb5 diff --git a/phys/module_diag_cordex.F b/phys/module_diag_cordex.F index d9e3a74089..3b36a7111b 100644 --- a/phys/module_diag_cordex.F +++ b/phys/module_diag_cordex.F @@ -337,15 +337,17 @@ SUBROUTINE cordex_output_calc( , clt, cll, clm, clh, & iut, ivt, & wsgs, ugustwind, vgustwind, gustpoint, & - totwsgs, totugustwind, totvgustwind, gustpoint, & - potevap, potevapogen, & + totwsgs, totugustwind, totvgustwind, totgustpoint, & + potevap, potevapgen, & #if CDXWRF>=1 cape, cin, zlfc, plfc, lidx, & #endif #if CDXWRF>=2 - fog, fogvisblty, & + fog, vis, & tds, tws & #endif +#if CDXWRF>=4 +#endif #endif !! End of INSTVALS #endif @@ -758,6 +760,19 @@ SUBROUTINE cordex_output_calc( #else REAL, DIMENSION(ims:ime,kms:kme,jms:jme) :: ta, press, zg, hur #endif +#if CDXWRF>=3 +! if WRF could output rank 5 or higher variables... +! REAL, DIMENSION(ims:ime,jms:jme,Nhtasrng,Nhhursrng), & +! INTENT(inout) :: tashurstreshighres +! REAL, DIMENSION(ims:ime,jms:jme,Nltasrng,Nlhursrng), & +! INTENT(inout) :: tashurstreslowres + REAL, DIMENSION(ims:ime,Nhtashursrng,jms:jme), & + INTENT(inout) :: tashurstreshighres + REAL, DIMENSION(ims:ime,Nltashursrng,jms:jme), & + INTENT(inout) :: tashurstreslowres + REAL, DIMENSION(ims:ime,Nwbdswssrng,jms:jme), & + INTENT(inout) :: wbdswsstres +#endif #if CDXWRF>=4 REAL, DIMENSION(ims:ime,jms:jme), INTENT(inout) :: wbacdh REAL, DIMENSION(ims:ime,jms:jme), INTENT(inout) :: wbacpw, wbacpwc, wbacpwr, wbacpws, & @@ -777,18 +792,6 @@ SUBROUTINE cordex_output_calc( wbaczlc, wbaczmc, wbaczhc, wbaczlr, wbaczmr, wbaczhr, wbaczls, wbaczms, wbaczhs, & wbaczli, wbaczmi, wbaczhi, wbaczlg, wbaczmg, wbaczhg, wbaczlh, wbaczmh, wbaczhh #endif -#if CDXWRF>=3 -! REAL, DIMENSION(ims:ime,jms:jme,Nhtasrng,Nhhursrng), & -! INTENT(inout) :: tashurstreshighres -! REAL, DIMENSION(ims:ime,jms:jme,Nltasrng,Nlhursrng), & -! INTENT(inout) :: tashurstreslowres - REAL, DIMENSION(ims:ime,Nhtashursrng,jms:jme), & - INTENT(inout) :: tashurstreshighres - REAL, DIMENSION(ims:ime,Nltashursrng,jms:jme), & - INTENT(inout) :: tashurstreslowres - REAL, DIMENSION(ims:ime,Nwbdswssrng,jms:jme), & - INTENT(inout) :: wbdswsstres -#endif #endif REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: uas, vas, wss, wbds, hurs, huss REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: psl @@ -819,34 +822,35 @@ SUBROUTINE cordex_output_calc( REAL, DIMENSION(ims:ime,jms:jme), INTENT(inout) :: rsdt, rsut, rlut REAL, DIMENSION(ims:ime,jms:jme), INTENT(inout) :: tas_hm, qvs_hm -!! INSTVALS: Instantaneous values (uncomment and recompile) -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: clt, cll, clm, clh -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: wsgs, ugustwind, vgustwind -! INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(out) :: gustpoint -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: totwsgs, totugustwind, totvgustwind -! INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(out) :: totgustpoint -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: potevap, potevapgen -#ifdef CDXWRF +#ifdef INSTVALS +!! INSTVALS: Instantaneous values + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: clt, cll, clm, clh + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: wsgs, ugustwind, vgustwind + INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(out) :: gustpoint + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: totwsgs, totugustwind, totvgustwind + INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(out) :: totgustpoint + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: potevap, potevapgen #if CDXWRF>=1 -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: cape, cin, zlfc, plfc, lidx -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: iut, ivt + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: cape, cin, zlfc, plfc, lidx + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: iut, ivt #endif #if CDXWRF>=2 -! INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(out) :: fog -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: vis -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: tds -! REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: tws + INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(out) :: fog + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: vis + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: tds + REAL, DIMENSION(ims:ime,jms:jme), INTENT(out) :: tws #endif +#if CDXWRF>=3 #endif -!! INSTVALS: comment these lines +#if CDXWRF>=4 +#endif +#else +!! INSTVALS: without output REAL, DIMENSION(ims:ime,jms:jme) :: clt, cll, clm, clh REAL, DIMENSION(ims:ime,jms:jme) :: wsgs, ugustwind, vgustwind INTEGER, DIMENSION(ims:ime,jms:jme) :: gustpoint REAL, DIMENSION(ims:ime,jms:jme) :: totwsgs, totugustwind, totvgustwind INTEGER, DIMENSION(ims:ime,jms:jme) :: totgustpoint - REAL, DIMENSION(ims:ime,jms:jme) :: wsgsb01, ugustwindb01, vgustwindb01 - REAL, DIMENSION(ims:ime,jms:jme) :: wsgshp, ugustwindhp, vgustwindhp - INTEGER, DIMENSION(ims:ime,jms:jme) :: gustpointb01, gustpointhp REAL, DIMENSION(ims:ime,jms:jme) :: potevap, potevapgen #if CDXWRF>=1 REAL, DIMENSION(ims:ime,jms:jme) :: cape, cin, zlfc, plfc, lidx @@ -859,10 +863,10 @@ SUBROUTINE cordex_output_calc( REAL, DIMENSION(ims:ime,jms:jme) :: tws #endif #if CDXWRF>=3 - INTEGER :: itas, ihurs, iwbds, iwss #endif #if CDXWRF>=4 #endif +#endif !! End of INSTVALS ! LOCAL VAR @@ -882,11 +886,16 @@ SUBROUTINE cordex_output_calc( REAL, DIMENSION(kms:kme) :: qvar, zagl, ri, qsolid REAL, DIMENSION(ims:ime,kms:kme+1,jms:jme) :: stzg REAL, DIMENSION(ims:ime,jms:jme) :: smoothp + REAL, DIMENSION(ims:ime,jms:jme) :: wsgsb01, ugustwindb01, vgustwindb01 + REAL, DIMENSION(ims:ime,jms:jme) :: wsgshp, ugustwindhp, vgustwindhp + INTEGER, DIMENSION(ims:ime,jms:jme) :: gustpointb01, gustpointhp + #ifdef CDXWRF #if CDXWRF>=1 #endif #if CDXWRF>=3 + INTEGER :: itas, ihurs, iwbds, iwss REAL :: uuas, vvas REAL, DIMENSION(ims:ime,jms:jme) :: wbdsv, wssv #endif @@ -981,7 +990,7 @@ SUBROUTINE cordex_output_calc( ! air-temperature CALL var_ta(thp(i,:,j), press(i,:,j), dimz, ta(i,:,j)) -#if CXWRF>= 2 +#if CDXWRF>= 2 ! 2m relative humidity CALL var_hurs(t2(i,j), ps(i,j), q2(i,j), hurs(i,j)) @@ -2313,6 +2322,14 @@ SUBROUTINE cordex_output_calc( !WRITE(msg, *)'i:',i,' i+1', i+1, ' j:', j, ' j+1:', j+1, ' ilims:', i_start(ij), i_end(ij), ' jlims:', j_start(ij), j_end(ij) !CALL wrf_debug(10,' ' // TRIM(msg) ) +#if CDXWRF< 2 + ! 2m relative humidity + CALL var_hurs(t2(i,j), ps(i,j), q2(i,j), hurs(i,j)) + + ! 10m Earth rotated winds + CALL var_uasvas(u10(i,j), v10(i,j), sina(i,j), cosa(i,j), uas(i,j), vas(i,j)) +#endif + ! Height above surface at half-mass levels zsfc = zg(i,:,j)-hgt(i,j) diff --git a/phys/module_diagnostics_driver.F b/phys/module_diagnostics_driver.F index faeb284357..f168889250 100644 --- a/phys/module_diagnostics_driver.F +++ b/phys/module_diagnostics_driver.F @@ -1307,7 +1307,7 @@ SUBROUTINE diagnostics_driver ( grid, config_flags, & #if CDXWRF>=4 #ifdef INSTVALS !! INSTVALS: Instantaneous values (uncomment and re-compile) - grid%cdxdiabh = grid%diab_h + grid%cdxdiabh = grid%h_diabatic #endif #endif #endif @@ -1587,13 +1587,14 @@ SUBROUTINE diagnostics_driver ( grid, config_flags, & !! INSTVALS: Instantaneous values ,CLT=grid%clt & ,CLL=grid%cll, CLM=grid%clm, CLH=grid%clh & + ,IUT=grid%iut, IVT=grid%ivt & ,WSGS=grid%wsgs, UGUSTWIND=grid%usgs & ,VGUSTWIND=grid%vsgs, GUSTPOINT=grid%gustpoint & ,TOTWSGS=grid%totwsgs, TOTUGUSTWIND=grid%totusgs & ,TOTVGUSTWIND=grid%totvsgs & ,TOTGUSTPOINT=grid%totgustpoint & - ,POTEVAPO=grid%potevapo & - ,POTEVAPOGEN=grid%potevapogen & + ,POTEVAP=grid%potevapo & + ,POTEVAPGEN=grid%potevapogen & #ifdef CDXWRF #if CDXWRF>=1 ,CAPE=grid%cdxcape & @@ -1601,9 +1602,11 @@ SUBROUTINE diagnostics_driver ( grid, config_flags, & ,PLFC=grid%lfcp, LIDX=grid%li & #endif #if CDXWRF>=2 - ,FOG=grid%fog, FOGVISBLTY=grid%fogvisblty & + ,FOG=grid%fog, VIS=grid%fogvisblty & ,TDS=grid%tds, TWS=grid%tws & #endif +#if CDXWRF>=4 +#endif #endif !! End of INSTVALS #endif diff --git a/phys/module_microphysics_driver.F b/phys/module_microphysics_driver.F index 34337f6238..df2e400644 100644 --- a/phys/module_microphysics_driver.F +++ b/phys/module_microphysics_driver.F @@ -71,8 +71,8 @@ SUBROUTINE microphysics_driver( & ,f_qns,f_qnr,f_qng,f_qnc,f_qnn,f_qh,f_qnh & , f_qzr,f_qzi,f_qzs,f_qzg,f_qzh & ,f_qvolg,f_qvolh & - ,f_qic,f_qip,f_qid & - ,f_qnic,f_qnip,f_qnid & + ,f_qic,f_qip,f_qid & + ,f_qnic,f_qnip,f_qnid & ,f_qir,f_qib & ! for P3 ,f_qi2,f_qni2,f_qir2,f_qib2 & ! for P3 ,f_qvoli,f_qaoli & ! for Jensen ISHMAEL @@ -87,9 +87,8 @@ SUBROUTINE microphysics_driver( & ,hail,ice2 & ! for mp_gsfcgce !NUWRF JJS 20110525 vvvvv ,phys_tot, physc, physe, physd, physs, physm, physf& ! for gsfcgce - ,acphys_tot, acphysc, acphyse, acphysd & ! for gsfcgce - ,acphyss, acphysm, acphysf & ! for gsfcgce - + ,acphys_tot, acphysc, acphyse, acphysd & ! for gsfcgce + ,acphyss, acphysm, acphysf & ! for gsfcgce ,re_cloud_gsfc, re_rain_gsfc, re_ice_gsfc & ,re_snow_gsfc, re_graupel_gsfc, re_hail_gsfc & ! cloud effective radius ,precr3d, preci3d, precs3d, precg3d, prech3d & @@ -1480,12 +1479,12 @@ SUBROUTINE microphysics_driver( & if (config_flags%nwp_diagnostics == 1) then if (present(hail_maxk1)) then if (allocated(tempo_driver_diags%max_hail_diameter_sfc)) then - hail_maxk1(i,j) = tempo_driver_diags%max_hail_diameter_sfc(i,j) + hail_maxk1(i,j) = 1.e-3 * tempo_driver_diags%max_hail_diameter_sfc(i,j) endif endif if (present(hail_max2d)) then if (allocated(tempo_driver_diags%max_hail_diameter_column)) then - hail_max2d(i,j) = tempo_driver_diags%max_hail_diameter_column(i,j) + hail_max2d(i,j) = 1.e-3 * tempo_driver_diags%max_hail_diameter_column(i,j) endif endif endif diff --git a/phys/module_mp_ntu.F b/phys/module_mp_ntu.F index 5629d63415..e647a5da98 100644 --- a/phys/module_mp_ntu.F +++ b/phys/module_mp_ntu.F @@ -5692,7 +5692,6 @@ SUBROUTINE LARGE_DT(DT,TK1D,QV1D,P1D,RHO,QC1D,QR1D,QI1D,QS1D, & VTNRG = SQRT(VTN0*VTN0+0.04*VTNR*VTNG) VTARG = SQRT(VTAX*VTAX+0.04*VTAR*VTAG) RHOGW = (QG1D+QR1D)/(QG1D/RHOG+QR1D/RHOW+ISMALL) - RATIO = (RHOGW/DNGRM)**THRD STOKE = (RHOW*ABS(VTQ0)*MVDR**2.)/(9.*MUA*MVDG) STOKE = MAX(1.5,MIN(10.,STOKE)) ERG = 5.5E-1*LOG10(2.51*STOKE) ! parameterization based on Cober and List, 1993 [JAS] @@ -5705,6 +5704,7 @@ SUBROUTINE LARGE_DT(DT,TK1D,QV1D,P1D,RHO,QC1D,QR1D,QI1D,QS1D, & DNGRM = 3.E2*(-5.E5*MVDR*0.6*VTQRG/TC1D)**0.44 ENDIF DNGRM = MIN(MAX(DNGRM,RHOIMIN),RHOG0) + RATIO = (RHOGW/DNGRM)**THRD IF (MVDG.GE.2.E-4.AND.MVDG.GE.MVDR) THEN QRMrg = QRMR1*ERG*VTQRG*NG1D*(GG3*GR4+2.*GG2*GR5+GR6) NRMrg = NRMR1*ERG*VTNRG*NG1D*(GG3+2.*GG2*GR2+GR3) diff --git a/phys/module_mp_udm.F b/phys/module_mp_udm.F index 2a98ac582c..40d235faef 100644 --- a/phys/module_mp_udm.F +++ b/phys/module_mp_udm.F @@ -944,6 +944,7 @@ subroutine udm2d(t1, q1 & psmlt(:) = 0. ; pgmlt(:) = 0. ; pseml(:) = 0. ; pgeml(:) = 0. prevp(:) = 0. ; psevp(:) = 0. ; pgevp(:) = 0. denqr(:) = 0. ; denqs(:) = 0. ; denqi(:) = 0. + denqh(:) = 0. ; denqg(:) = 0. vtr(:) = 0. ; vts(:) = 0. ; vtg(:) = 0. ; vti(:) = 0. qrpath = 0. ; qspath = 0. ; qgpath = 0. ; qipath = 0. precip_r = 0. ; precip_s = 0. ; precip_g = 0. ; precip_i = 0. @@ -1955,7 +1956,7 @@ subroutine udm2d(t1, q1 & do k = ktop, kts, -1 if(ifice(k)) then if(supsat(k)>0 .and. (.not.ifsat(k))) then - supice(k) = satrdt(k)-prevp(k)-pidep(k)-psdep(k)-pgdep(k) + supice(k) = satrdt(k)-prevp(k)-pidep(k)-psdep(k)-pgdep(k)-phdep(k) if(slimsk(i)==0) then ni0(k) = 1000.*exp(-0.2*tcelci(k)-5.) else if(slimsk(i)==1) then diff --git a/phys/module_pbl_driver.F b/phys/module_pbl_driver.F index 0455148925..07712a3308 100644 --- a/phys/module_pbl_driver.F +++ b/phys/module_pbl_driver.F @@ -51,7 +51,7 @@ SUBROUTINE pbl_driver( & ,dqke,qWT,qSHEAR,qBUOY,qDISS,tke_budget & ,bl_mynn_closure,bl_mynn_cloudpdf & ,bl_mynn_mixlength & - ,icloud_bl,qc_bl,qi_bl,cldfra_bl & + ,qc_bl,qi_bl,cldfra_bl & ,bl_mynn_edmf,bl_mynn_edmf_dd & ,bl_mynn_edmf_mom,bl_mynn_edmf_tke & ,bl_mynn_mixscalars,bl_mynn_mixaerosols & @@ -64,6 +64,8 @@ SUBROUTINE pbl_driver( & ,vdfg & ,maxwidth,maxMF,ztop_plume & ,excess_h,excess_q & + ,maxmf_dd,maxwidth_dd,maxtkeprod & + ,cldtop_cooling,ent_eff & ,spp_pbl,pattern_spp_pbl & ! EEPS ,pek,pep,pek_adv,pep_adv & @@ -603,7 +605,6 @@ SUBROUTINE pbl_driver( & bl_mynn_cloudmix, & bl_mynn_mixqt, & bl_mynn_ess, & - icloud_bl, & tke_budget REAL, OPTIONAL, INTENT(IN) :: bl_mynn_closure @@ -615,8 +616,9 @@ SUBROUTINE pbl_driver( & REAL, DIMENSION( ims:ime , jms:jme ), OPTIONAL, & & INTENT(INOUT):: vdfg REAL, OPTIONAL, DIMENSION( ims:ime , jms:jme ), & - & INTENT(OUT) :: maxwidth,maxMF,ztop_plume,& - excess_h,excess_q + & INTENT(OUT) :: maxwidth,maxMF,ztop_plume, & + maxmf_dd,maxtkeprod,cldtop_cooling,ent_eff, & + maxwidth_dd,excess_h,excess_q REAL, OPTIONAL, DIMENSION( ims:ime, kms:kme, jms:jme ), & INTENT(INOUT) :: qnwfa_curr,qnifa_curr,qnbca_curr @@ -1728,14 +1730,15 @@ SUBROUTINE pbl_driver( & &sub_thl3D=sub_thl3D,sub_sqv3D=sub_sqv3D, & &det_thl3D=det_thl3D,det_sqv3D=det_sqv3D, & &maxwidth=maxwidth,maxMF=maxMF, & - &ztop_plume=ztop_plume, & + &ztop_plume=ztop_plume,maxmf_dd=maxmf_dd, & + &maxtkeprod=maxtkeprod,cldtop_cooling=cldtop_cooling, & + &ent_eff=ent_eff,maxwidth_dd=maxwidth_dd, & &excess_h=excess_h,excess_q=excess_q, & &RTHRATEN=RTHRATEN, & &bl_mynn_tkeadvect=bl_mynn_tkeadvect, & &tke_budget=tke_budget, & &bl_mynn_cloudpdf=bl_mynn_cloudpdf, & &bl_mynn_mixlength=bl_mynn_mixlength, & - &icloud_bl=icloud_bl, & &bl_mynn_edmf=bl_mynn_edmf, & &bl_mynn_edmf_dd=bl_mynn_edmf_dd, & &bl_mynn_edmf_mom=bl_mynn_edmf_mom, & diff --git a/phys/module_ra_rrtmg_swf.F b/phys/module_ra_rrtmg_swf.F index 9ba21d34ec..1260707330 100644 --- a/phys/module_ra_rrtmg_swf.F +++ b/phys/module_ra_rrtmg_swf.F @@ -11689,7 +11689,7 @@ SUBROUTINE RRTMG_SWRAD_FAST( & else ! da=6.2831853071795862*(julian-1)/365. ! eot=(0.000075+0.001868*cos(da)-0.032077*sin(da) & -! -0.014615*cos(2*da)-0.04089*sin(2*da))*(229.18) +! -0.014615*cos(2*da)-0.040849*sin(2*da))*(229.18) xt24 = mod(xtime+radt*0.5,1440.)+eot tloctm = gmt + xt24/60. + xlong(i,j)/15. hrang = 15. * (tloctm-12.) * degrad diff --git a/phys/module_radiation_driver.F b/phys/module_radiation_driver.F index 1421cbd34f..00cfd82f25 100644 --- a/phys/module_radiation_driver.F +++ b/phys/module_radiation_driver.F @@ -3524,9 +3524,9 @@ SUBROUTINE calc_coszen(ims,ime,jms,jme,its,ite,jts,jte, & integer :: i,j real :: da,eot,xt24,tloctm,xxlat - da=6.2831853071795862*(julian-1)/365. + da=6.2831853071795862*(julian-0.5)/365. eot=(0.000075+0.001868*cos(da)-0.032077*sin(da) & - -0.014615*cos(2*da)-0.04089*sin(2*da))*(229.18) + -0.014615*cos(2*da)-0.040849*sin(2*da))*(229.18) xt24=mod(xtime,1440.)+eot do j=jts,jte do i=its,ite diff --git a/phys/module_sf_urban.F b/phys/module_sf_urban.F index c49b9c31fe..b63208ef55 100644 --- a/phys/module_sf_urban.F +++ b/phys/module_sf_urban.F @@ -103,8 +103,9 @@ MODULE module_sf_urban ! !===Huang and Wang, 2025/12/28, urban nature-based solutions for single layer UCM (including trees)=== - INTEGER :: TREEOPTION ! Urban Tree&Vegetation Option ! add by yuqi - REAL :: fvg_data = 0.0 ! vegetated ground fraction from URBPARM + INTEGER :: TREEOPTION ! Urban Tree&Vegetation Option ! urban-nbs + INTEGER :: TTscheme ! Tree temperature solution [0: diagnostic, 1: prognostic] + REAL, ALLOCATABLE, DIMENSION(:) :: FVG_TBL ! vegetated ground fraction from URBPARM REAL, ALLOCATABLE, DIMENSION(:) :: RTREE_TBL ! Tree crown radius [normalized] REAL, ALLOCATABLE, DIMENSION(:) :: RTREE_M_TBL ! Tree crown radius in meter REAL, ALLOCATABLE, DIMENSION(:) :: RW_M_TBL ! Road width in meter @@ -117,6 +118,8 @@ MODULE module_sf_urban REAL, ALLOCATABLE, DIMENSION(:) :: CAPVG_TBL ! Heat capacity of vegetated ground REAL, ALLOCATABLE, DIMENSION(:) :: EPSVG_TBL ! Emissivity of vegetated ground REAL, ALLOCATABLE, DIMENSION(:) :: EPST_TBL ! Emissivity of tree canopy + REAL, ALLOCATABLE, DIMENSION(:) :: ALBVG_TBL ! Emissivity of vegetated ground + REAL, ALLOCATABLE, DIMENSION(:) :: ALBT_TBL ! Emissivity of tree canopy !===end urban nature-based solutions processes=== CHARACTER (LEN=256) , PRIVATE :: mesg @@ -138,7 +141,7 @@ MODULE module_sf_urban ! NCAR/RAP mukul@ucar.edu ! ! Purpose: -! Calculate surface temeprature, fluxes, canopy air temperature, and canopy wind +! Calculate surface temperature, fluxes, canopy air temperature, and canopy wind ! ! Subroutines: ! module_sf_urban @@ -151,8 +154,8 @@ MODULE module_sf_urban ! ! Input Data from WRF [MKS unit]: ! -! UTYPE [-] : Urban type. 1=Commercial/Industrial; 2=High-intensity residential; -! : 3=low-intensity residential +! UTYPE [-] : Urban type. 1=low-intensity residential; 2=High-intensity residential; +! : 3= Commercial/Industrial ! TA [K] : Potential temperature at 1st wrf level (absolute temp) ! QA [kg/kg] : Mixing ratio at 1st atmospheric level ! UA [m/s] : Wind speed at 1st atmospheric level @@ -180,7 +183,7 @@ MODULE module_sf_urban ! OMG [rad] : solar hour angle ! XLAT [deg] : latitude ! DELT [sec] : Time step -! ZNT [m] : Roughnes length +! ZNT [m] : Roughness length ! ! Output Data to WRF [MKS unit]: ! @@ -302,7 +305,8 @@ MODULE module_sf_urban ! Kusaka et al. (2001) Bound.-Layer Meteor., vol.101, p329-358 ! ! History: -! 2025/11, modified by Yuqi Huang, Chenghao Wang (Univ. Oklahoma) [urban nature-based solutions, including trees] +! 2025/11, modified by Yuqi Huang, Chenghao Wang (Univ. Oklahoma) [urban nature-based solutions] +! 2024/08, modified by Yuqi Huang, Chenghao Wang (Univ. Oklahoma) [detailed code annotation] ! 2014/10, modified by Jiachuan Yang (ASU) ! 2006/06 modified by H. Kusaka (Univ. Tsukuba), M. Tewari ! 2005/10/26, modified by Fei Chen, Mukul Tewari @@ -333,16 +337,15 @@ SUBROUTINE urban(LSOLAR, & ! L lp_urb,hgt_urb,frc_urb,lb_urb,zo_check, & ! O CMCR,TGR,TGRL,SMR,CMGR_URB,CHGR_URB,jmonth, & ! H DRELR,DRELB,DRELG,FLXHUMR,FLXHUMB,FLXHUMG, & - ! added by Chenghao Wang - TVG,TT,XXXVG,TVGL,SMG,CMCG,FLXHUMVG,FLXHUMT, & ! H (add by huang) + TVG,TT,XXXVG,TVGL,SMG,CMCG,FLXHUMVG,FLXHUMT, & ! H (add in urban-nbs) lf_urb_s, z0_urb, vegfrac_in & ! I (distributed aerodynamics) ) IMPLICIT NONE - REAL, PARAMETER :: CP=0.24 ! heat capacity of dry air [cgs unit] ! cal/g/c (yuqi) - REAL, PARAMETER :: EL=583. ! latent heat of vaporization [cgs unit] ! cal/g (yuqi) + REAL, PARAMETER :: CP=0.24 ! heat capacity of dry air [cgs unit] ! cal/g/c + REAL, PARAMETER :: EL=583. ! latent heat of vaporization [cgs unit] ! cal/g REAL, PARAMETER :: SIG=8.17E-11 ! stefun bolzman constant [cgs unit] ! cal/cm^2/min/k^4 - REAL, PARAMETER :: SIG_SI=5.67E-8 ! [MKS unit] ! w/m2/k^4 (yuqi) + REAL, PARAMETER :: SIG_SI=5.67E-8 ! [MKS unit] ! w/m2/k^4 REAL, PARAMETER :: AK=0.4 ! kalman const. [-] REAL, PARAMETER :: PI=3.14159 ! pi [-] REAL, PARAMETER :: TETENA=7.5 ! const. of Tetens Equation [-] @@ -369,7 +372,6 @@ SUBROUTINE urban(LSOLAR, & ! L !------------------------------------------------------------------------------- REAL :: TVGP, TTP ! Temperature at Pervious Time Step REAL :: ALBGE ! Effective Surface Albedo of Ground - REAL :: ALBVG ! Surface Albedo of Vegetated Ground REAL :: XI ! Tree Shadow Related Variable REAL :: ANGLE, TAN1, TAN2 ! Reference Angles Between Tree and Roof REAL :: SDT ! Direct Short Wave Radiation incident on Tree @@ -391,7 +393,7 @@ SUBROUTINE urban(LSOLAR, & ! L REAL, DIMENSION(1:num_road_layers) :: ZSOILG ! Total Depth of Vegetated Ground (Negative) [m], (Accumulated From Surface) REAL, DIMENSION(1:num_road_layers) :: SMGP ! Soil Moisture at Each Layer on Vegetated Ground at Pervious Time Step REAL, DIMENSION(1:num_road_layers) :: ETG ! Vegetated Canopy Evapotranspiration - REAL, DIMENSION(1:num_road_layers) :: DZVG ! Depth of Each Single Layer (Positive)[m] ! This Could be an Input Variable (yuqi) + REAL, DIMENSION(1:num_road_layers) :: DZVG ! Depth of Each Single Layer (Positive)[m] REAL :: QS0VG, DQS0VGDTVG, QS0T, DQS0TDTT ! Saturation QS and derivative for vegetated ground/tree REAL :: EPVG ! Potential Evaporation For Vegetated Ground [kg/m2/s] REAL :: EDIRG ! Direct Evaporation From Vegetated Ground @@ -415,8 +417,8 @@ SUBROUTINE urban(LSOLAR, & ! L REAL :: ALPHAT ! sensible heat transfer coefficient for urban trees REAL :: W2 ! Deep Soil Moisture REAL :: RS_LEAF ! Leaf Stomatal Resistance - REAL :: ELET ! Latent Heat of Tree Leaf ! w/m/m (yuqi) - REAL :: AVAILABLE_WATER_SUM, POTENTIAL_ELET ! newly added by yuqi + REAL :: ELET ! Latent Heat of Tree Leaf ! w/m/m + REAL :: AVAILABLE_WATER_SUM, POTENTIAL_ELET REAL :: QS0C ! Saturation specific humidity at canyon air temperature CHARACTER(LEN=256) :: WARN_MSG ! A character string for the warning message REAL, DIMENSION(1:num_road_layers) :: SROOT, SROOTP ! Root Water Uptake (m/s) @@ -424,8 +426,8 @@ SUBROUTINE urban(LSOLAR, & ! L REAL(Kind = 8) :: DRS_LEAF, M,D,DD,FE,DFE,S,DM,KK,L,N REAL(Kind = 8) :: DELEVGDTVG,DELEVGDTT,DELETDTB,DELETDTG,DELETDTVG,DELETDTT,DTCDTVG,DTCDTT ! derivatives wrt temperature (latent heat) REAL :: DQCDTVG,DQCDTT,ELEVG,KVG,G0VG,DG0BDTVG,DG0BDTT,DG0GDTVG,DG0GDTT,DG0VGDTB,DG0VGDTG,DG0VGDTVG,DG0VGDTT ! derivatives wrt temperature - REAL(Kind = 8) :: WF,A,B,C,E,FF,GG,H,VGF,I,J,O,P,CHARA,det,det_tol,det_scale !,RVGT - REAL(Kind = 8) :: A_use,FF_use,KK_use,lambda,relax + REAL(Kind = 8) :: WF,A,B,C,E,FF,GG,H,VGF,I,J,O,P,CHARA,det,det_tol,det_scale + REAL(Kind = 8) :: A_use,FF_use,KK_use,lambda,relax,W2_num,W2_den,den_qc REAL :: max_step REAL :: SSOILG ! Soil Heat Flux output from SHFLX W/m2 REAL :: TGE, TGEP, QS0GE, ALPHAGE, FLXTHVG, FLXTHT ! effective ground values @@ -433,12 +435,11 @@ SUBROUTINE urban(LSOLAR, & ! L !----------- fixed universe parameters for urban vegetation --------------- INTEGER, PARAMETER :: NVG = 4 ! Total Layer of Vegetated Ground REAL, PARAMETER :: CPP=1004.5 ! heat capacity of dry air [J/K/kg] - REAL, PARAMETER :: ALBT = 0.2 ! Leaf Surface Albedo REAL, PARAMETER :: DLEAF = 0.1 ! Unit = m REAL, PARAMETER :: TREF = 298.0 ! Unit = K REAL, PARAMETER :: EREF = 3167.0 ! Unit = K - !----------- local variables from read_para (Yuqi) --------------- + !----------- local variables from read_para (urban-nbs) --------------- REAL :: FVG ! Vegetated Fraction for Ground REAL :: Z0VG ! Roughness Length for Momentum of Vegetated Ground REAL :: Z0HVG ! Roughness Length for Heat of Vegetated Ground @@ -450,8 +451,10 @@ SUBROUTINE urban(LSOLAR, & ! L REAL :: RW_M ! road width in meters REAL :: DTREE ! Distance of Tree Crown Center From Wall (RW/2) REAL :: HTREE ! Height of Tree Crown Center - REAL :: LAI_VEG ! Vegetated Ground Leaf Areal Index - REAL :: LAI_TREE ! Leaf Area Index of Trees (Could be an input) + REAL :: ALBT ! Leaf Surface Albedo + REAL :: ALBVG ! Vegetated ground albedo + REAL :: LAI_TREE ! Leaf Area Index of Trees + REAL :: LAI_VEG ! Vegetated Ground Leaf Areal Index ! -------------------------------------End Urban Vegetation Variable Define----------------------------------------------------------------------------------------------------- @@ -478,8 +481,8 @@ SUBROUTINE urban(LSOLAR, & ! L ! I: input variables from LSM to Urban !------------------------------------------------------------------------------- - INTEGER, INTENT(IN) :: UTYPE ! urban type [1=Commercial/Industrial, 2=High-intensity residential, - ! 3=low-intensity residential] + INTEGER, INTENT(IN) :: UTYPE ! urban type [1=low-intensity residential, 2=High-intensity residential, + ! 3=Commercial/Industrial] INTEGER, INTENT(IN) :: jmonth! current month REAL, INTENT(IN) :: TA ! potential temp at 1st atmospheric level [K] REAL, INTENT(IN) :: QA ! mixing ratio at 1st atmospheric level [kg/kg] @@ -511,7 +514,7 @@ SUBROUTINE urban(LSOLAR, & ! L !------------------------------------------------------------------------------- REAL, INTENT(INOUT) :: mh_urb ! mean building height [m] REAL, INTENT(INOUT) :: stdh_urb ! standard deviation of building height [m] - REAL, INTENT(INOUT) :: hgt_urb ! area weighted mean building height [m] ! if hgt_urb>0 use NUDAPT else use default urban parameter table (yuqi) + REAL, INTENT(INOUT) :: hgt_urb ! area weighted mean building height [m] ! if hgt_urb>0 use NUDAPT else use default urban parameter table REAL, INTENT(INOUT) :: lp_urb ! plan area fraction [-] REAL, INTENT(INOUT) :: frc_urb ! urban fraction [-] REAL, INTENT(INOUT) :: lb_urb ! building surface to plan area ratio [-] @@ -607,7 +610,7 @@ SUBROUTINE urban(LSOLAR, & ! L REAL :: RLMO_URB REAL :: AKANDA_URBAN - REAL :: TH2X !m + REAL :: TH2X !m INTEGER :: BOUNDR, BOUNDB, BOUNDG INTEGER :: CH_SCHEME, TS_SCHEME @@ -671,7 +674,7 @@ SUBROUTINE urban(LSOLAR, & ! L REAL :: XXX2, PSIM2, PSIH2, XXX10, PSIM10, PSIH10 REAL :: PSIX, PSIT, PSIX2, PSIT2, PSIX10, PSIT10 - REAL :: TRP, TBP, TGP, TCP, QCP, TST, QST ! TCP is local variable ? (yuqi) + REAL :: TRP, TBP, TGP, TCP, QCP, TST, QST REAL :: TSP, CHS_LOCAL, CHS2_LOCAL REAL :: WDR,HGT2,BW,DHGT @@ -696,14 +699,14 @@ SUBROUTINE urban(LSOLAR, & ! L ! REAL :: DQS0GRDTGR, ETR, ECR,RAIN1, RAINDR, DEW, ETAR, BETGR REAL :: DQS0GRDTGR, ECR,RAIN1, RAINDR, DEW, ETAR, BETGR ! REAL :: DF1, RGR, RGRR, RCH, RR1, RR2, YY, ZZ1, SSOILR - REAL :: DF1, RGR, RGRR, RCH, YY, ZZ1, SSOILR, RVGT !(RVGT add by yuqi) + REAL :: DF1, RGR, RGRR, RCH, YY, ZZ1, SSOILR, RVGT !(RVGT add by Yuqi) REAL :: DRRDTGR, DHRDTGR, DELERDTGR, DG0RDTGR, DFDVT real,parameter :: SHDFAC = 0.80 ! Vegetated area fraction of green roof vegetation real,parameter :: ALBV = 0.20 ! green roof albedo real,parameter :: EPSV = 0.93 ! green roof emissivity real,parameter :: LAI = 1.50 ! leaf area index on green roof real,parameter :: CMCMAX = 0.5E-3 ! Maximum canopy interception capacity - real,parameter :: SMCREF = 0.329 ! Reference soil moisture + real,parameter :: SMCREF = 0.329 ! Reference soil moisture ! below variables should be soil type dependent (Yuqi) real,parameter :: SMCDRY = 0.066 ! Residual soil moisture real,parameter :: SMCWLT = 0.084 ! Wilting point real,parameter :: SMCMAX = 0.439 ! Saturated soil moisture @@ -748,9 +751,9 @@ SUBROUTINE urban(LSOLAR, & ! L CALL read_param(UTYPE,ZR,SIGMA_ZED,Z0C,Z0HC,ZDC,SVF,R,RW,HGT, & AH,CAPR,CAPB,CAPG,AKSR,AKSB,AKSG,ALBR,ALBB, & ALBG,EPSR,EPSB,EPSG,Z0R,Z0B,Z0G,Z0HB,Z0HG, & - BETR,BETB,BETG,TRLEND,TBLEND,TGLEND, & - RTREE,RTREE_M,RW_M,HTREE,DTREE,LAI_TREE,LAI_VEG, & ! added by yuqi - Z0VG,Z0HVG,CAPVG,EPSVG,EPST, & ! added by yuqi + BETR,BETB,BETG,TRLEND,TBLEND,TGLEND,ALBVG,ALBT,FVG, & + RTREE,RTREE_M,RW_M,HTREE,DTREE,LAI_TREE,LAI_VEG, & ! urban-nbs + Z0VG,Z0HVG,CAPVG,EPSVG,EPST, & ! urban-nbs !for BEP NUMDIR, STREET_DIRECTION, STREET_WIDTH, & BUILDING_WIDTH, NUMHGT, HEIGHT_BIN, & @@ -758,7 +761,7 @@ SUBROUTINE urban(LSOLAR, & ! L !end BEP BOUNDR,BOUNDB,BOUNDG,CH_SCHEME,TS_SCHEME, & AKANDA_URBAN,ALH) - FVG = FVG_DATA + ! Glotfelty, 2012/07/05, NUDAPT Modification @@ -799,7 +802,7 @@ SUBROUTINE urban(LSOLAR, & ! L !WRITE_MESSAGE(mesg) !Calculate Building Width and Street Width Based on BEP formulation - if(lb_urb.gt.lp_urb)THEN ! building surface to plan area ratio > plain area fraction (yuqi) + if(lb_urb.gt.lp_urb)THEN ! building surface to plan area ratio > plain area fraction BW=2.*hgt_urb*lp_urb/(lb_urb-lp_urb) SW=2.*hgt_urb*lp_urb*((frc_urb/lp_urb)-1.)/(lb_urb-lp_urb) !write(mesg,*) 'Building Width',BW @@ -854,10 +857,10 @@ SUBROUTINE urban(LSOLAR, & ! L beta_macd = 1.0 - ZDC = ZR * ( 1.0 + ( alpha_macd ** ( -R ) ) * ( R - 1.0 ) ) ! e.q 23 in Macdonald's (1998) (yuqi) + ZDC = ZR * ( 1.0 + ( alpha_macd ** ( -R ) ) * ( R - 1.0 ) ) ! e.q 23 in Macdonald's (1998) Z0C = ZR * ( 1.0 - ZDC/ZR ) * & - exp (-(0.5 * beta_macd * Cd / (VonK**2) * ( 1.0-ZDC/ZR) * lambda_f )**(-0.5)) ! e.q 26 in Macdonald's (1998) (yuqi) + exp (-(0.5 * beta_macd * Cd / (VonK**2) * ( 1.0-ZDC/ZR) * lambda_f )**(-0.5)) ! e.q 26 in Macdonald's (1998) if(zo_check.eq.1)THEN write(mesg,*) 'Roughness Length NUDAPT',Z0C @@ -1017,7 +1020,7 @@ SUBROUTINE urban(LSOLAR, & ! L ZC=0.7*ZR XLB=0.4*(ZR-ZDC) ! BB formulation from Inoue (1963) -! #ifdef DOUBLE_PRECISION ! dlog is newly added (yuqi) - code follows HRLDAS (2025 summer) +! #ifdef DOUBLE_PRECISION ! dlog is newly added (Yuqi) - code follows HRLDAS (2025 summer) ! BB = 0.4 * ZR / ( XLB * dlog( ( ZR - ZDC ) / Z0C ) ) ! #else BB = 0.4 * ZR / ( XLB * alog( ( ZR - ZDC ) / Z0C ) ) @@ -1037,6 +1040,11 @@ SUBROUTINE urban(LSOLAR, & ! L ! SHADOW = .false. ! SHADOW has not been activated in previous versions of SLUCM (chenghao) - it's not turned on + TAU = EXP(-0.61 * LAI_TREE) + IF (TREEOPTION .eq. 1 ) THEN + CALL VF_TREE_ANALYTICAL(RW, HGT, HTREE, RTREE, TAU, VFGS_TREE, VFGW_TREE, VFGT_TREE, VFWS_TREE, VFWG_TREE, VFWW_TREE, VFWT_TREE, VFTS_TREE, VFTG_TREE, VFTW_TREE) + END IF + IF (SSG > 0.0) THEN IF(SHADOW) THEN ! shadow effects, calculate SLX @@ -1084,7 +1092,7 @@ SUBROUTINE urban(LSOLAR, & ! L SG1=SX*VFGS*(1.-ALBG) SB1=SX*VFWS*(1.-ALBB) SG2=SB1*ALBB/(1.-ALBB)*VFGW*(1.-ALBG) - SB2=SG1*ALBG/(1.-ALBG)*VFWG*(1.-ALBB) + SB1*ALBB*VFWW ! add one term in the new version (yuqi) + SB2=SG1*ALBG/(1.-ALBG)*VFWG*(1.-ALBB) + SB1*ALBB*VFWW ! add one term in the new version (urban-nbs) ELSEIF (TREEOPTION .ne. 1) THEN ! shadow effects model without urban vegetations @@ -1093,17 +1101,17 @@ SUBROUTINE urban(LSOLAR, & ! L SG1=SD*(RW-SLX)/RW*(1.-ALBG)+SQ*VFGS*(1.-ALBG) SB1=SD*SLX/W*(1.-ALBB)+SQ*VFWS*(1.-ALBB) SG2=SB1*ALBB/(1.-ALBB)*VFGW*(1.-ALBG) - SB2=SG1*ALBG/(1.-ALBG)*VFWG*(1.-ALBB) + SB1*ALBB*VFWW ! newly added last term (Huang) - ELSE ! shadow effects model with urban vegetations (Huang) + SB2=SG1*ALBG/(1.-ALBG)*VFWG*(1.-ALBB) + SB1*ALBB*VFWW ! newly added last term (urban-nbs) + ELSE ! shadow effects model with urban vegetations (urban-nbs) ! ------------------------------------------------------------------------------------------------------------------ ! tree module only activated under shadow case ! ------------------------------------------------------------------------------------------------------------------ - IF ( HTREE+RTREE .gt. ZR) THEN ! add by Huang + IF ( HTREE+RTREE .gt. ZR) THEN ! add in urban-nbs FATAL_ERROR("tree crown protrude roof height, adjust tree size") END IF ALBGE = (1. - FVG) * ALBG + FVG * ALBVG - XI = SLX/HGT ! check if XI matches with our ucm, (yuqi) + XI = SLX/HGT ! check if XI matches with offline UCM (urban-nbs) ANGLE = DTREE**2. + (HGT - HTREE)**2. - RTREE**2. TAN1 = (RTREE * (HGT-HTREE) + DTREE * SQRT(ANGLE))/((HGT - HTREE) * SQRT(ANGLE) - RTREE * DTREE) @@ -1120,13 +1128,11 @@ SUBROUTINE urban(LSOLAR, & ! L ELSE ! XI <= tan2 in the sun - completely illuminated SDT = SD * (2. * RTREE * SQRT(1. + XI**2.))/(2. * PI * RTREE) END IF - print*, 'THEATAZ is:', THEATAZ ! add by yuqi for testing purpose + print*, 'THEATAZ is:', THEATAZ ! urban-nbs for testing purpose - TAU = EXP(-0.61 * LAI_TREE) + CALL SHADOW_TREE(XI, HGT, RW, HTREE, RTREE, DTREE, CHI_SHADED, ETA_SHADED, CHI_SHADED_TREE, ETA_SHADED_TREE) - ! print*, 'road weidth, building height, tree height, tree radius are:', RW, HGT, HTREE, RTREE ! add by yuqi for testing purpose - CALL VF_TREE_ANALYTICAL(RW, HGT, HTREE, RTREE, TAU, VFGS_TREE, VFGW_TREE, VFGT_TREE, VFWS_TREE, VFWG_TREE, VFWW_TREE, VFWT_TREE, VFTS_TREE, VFTG_TREE, VFTW_TREE) - + ! print*, 'road weidth, building height, tree height, tree radius are:', RW, HGT, HTREE, RTREE ! add for urban-nbs for testing purpose SDG = SD * (RW - CHI_SHADED + CHI_SHADED_TREE * TAU)/RW SDB = SD * XI * (HGT - ETA_SHADED +ETA_SHADED_TREE * TAU)/W @@ -1151,7 +1157,7 @@ SUBROUTINE urban(LSOLAR, & ! L SR=SR1 ! units: [cal/cm2/s] SGR=SGR1 SG=SG1+SG2 - SVG=SVG1+SVG2 ! Huang + SVG=SVG1+SVG2 ! urban-nbs SB=SB1+SB2 IF (TREEOPTION == 1) THEN @@ -1347,14 +1353,12 @@ SUBROUTINE urban(LSOLAR, & ! L QS0GR=0.622*ES/(PS-0.378*ES) DQS0GRDTGR = DESDT*0.622*PS/((PS-0.378*ES)**2.) EPGR=RHOO*CHGR*UA*(QS0GR-QA) ! Potential evaporation [kg/m2/s] - ! print*, 'vegetated roof temperature is:', TGRP ! add by yuqi for testing purpose - ! print*, 'potential et from vegetated roof is:', EPGR ! add by yuqi for testing purpose IF (EPGR > 0.0) THEN ! Direct evaporation from soil on green roof CALL DIREVAP (EDIR,EPGR,SMRP(KZ),SHDFAC,SMCMAX,SMCDRY,FXEXP) ! CONTRIBUTED FROM BARE SOIL ! Evapotranspiration and canopy intercepted evaporation - CALL TRANSP (ETTR,ETR,ECR,SHDFAC,EPGR,CMCRP,CFACTR,CMCMAX,LAI,RSMIN,RSMAX,RGL,SX, & ! CONTRIBUTED FROM VEGETED all layers (yuqi) + CALL TRANSP (ETTR,ETR,ECR,SHDFAC,EPGR,CMCRP,CFACTR,CMCMAX,LAI,RSMIN,RSMAX,RGL,SX, & ! CONTRIBUTED FROM VEGETED all layers TGRP,TA,QA,SMRP,SMCWLT,SMCREF,CPP,PS,CHGR,EPSV,DELT,NROOT,NGR,DZGR, & ZSOILR,HS) ! Update moisture in soil layers @@ -1410,7 +1414,7 @@ SUBROUTINE urban(LSOLAR, & ! L DHRDTGR = RHO*CP*CHGR*UA*100. DELERDTGR = RHO*EL*CHGR*UA*BETGR*DQS0GRDTGR*100. DG0RDTGR = 2.*DF1/ DZGR(KZ) * ( 1.0 / 4.1868 ) * 1.E-4 - DFDVT = DRRDTGR - DHRDTGR - DELERDTGR - DG0RDTGR ! HERE USE ETP INSTEAD OF ACTUAL ET TO DERIVE DELEDTGR (YUQI) + DFDVT = DRRDTGR - DHRDTGR - DELERDTGR - DG0RDTGR DTGR = FV/DFDVT/ 6 TGR = TGRP - DTGR TGRP = TGR @@ -1448,7 +1452,7 @@ SUBROUTINE urban(LSOLAR, & ! L T1VC = TCP* (1.0+ 0.61 * QA) RLMO_URB=0.0 - CALL SFCDIF_URB(ZA,Z0C,T1VC,TH2V,UA,AKANDA_URBAN,CMC_URB,CHC_URB,RLMO_URB,CDC) ! CDC: momentum exchange coefficient for canopy (yuqi) + CALL SFCDIF_URB(ZA,Z0C,T1VC,TH2V,UA,AKANDA_URBAN,CMC_URB,CHC_URB,RLMO_URB,CDC) ! CDC: momentum exchange coefficient for canopy ALPHAC = RHO*CP*CHC_URB IF (CH_SCHEME == 1) THEN @@ -1730,8 +1734,7 @@ SUBROUTINE urban(LSOLAR, & ! L ELSE !------------------------------------------------------------------------------- - ! ENERGY BUDGET WITH URBAN TREE (Huang) - ! VFAC: VEGE FRACTION IN VEGETATED GROUND, SIMILAR TO SHDFAC, BUT SHOULD BE PARAMETERIZED (YUQI) + ! ENERGY BUDGET WITH URBAN TREE (urban-nbs) ! NEW VARIABLES: ! BHVG,Z0VG,Z0HVG,RIBVG,TVGP,XXXVG,ALPHAVG,CDVG,CHVG,SRUN1,SRUN2,SRUN3,ZSOILG(total depth),DZVG(depth of each single layer), ! QS0VG,DQS0VGDTVG,QS0T,DQS0TDTT,EPVG,EDIRG,ETTG,ETG,ECG,LAI_VEG,NVG,DRIPG,ETAG,BETVG, @@ -1751,7 +1754,7 @@ SUBROUTINE urban(LSOLAR, & ! L IF (TS_SCHEME == 1) THEN ! ------------ FIRST CALCULATE CHVG, CDVG, BETVG---------- - IF (CH_SCHEME == 1) THEN ! calculate heat exchange coefficient for vegetated ground (yuqi) + IF (CH_SCHEME == 1) THEN ! calculate heat exchange coefficient for vegetated ground Z=ZDC BHVG=LOG(Z0VG/Z0HVG)/0.4 RIBVG=(9.8*2/(TCP+TVGP))*(TCP-TVGP)*(Z+Z0VG)/(UC*UC) @@ -1762,39 +1765,64 @@ SUBROUTINE urban(LSOLAR, & ! L END IF CHVG=ALPHAVG/RHO/CP/UC - SRUN1 = 0.0 ! mainly following green roof's scheme below (yuqi) + SRUN1 = 0.0 ! mainly following green roof's scheme below (urban-nbs) SRUN2 = 0.0 SRUN3 = 0.0 DZVG = DZGR ! set identical green roof depth and vegetated ground layer depth (Yuqi) [m] KZ = 1 - ZSOILG (KZ) = - DZVG (KZ) ! dz: (positive) depth of each single layer; zsoil: (negative) accumulated depth from ground (yuqi) - DO KZ = 2,NVG ! nvg: total layer of vegetated ground, = 4 in here (yuqi) + ZSOILG (KZ) = - DZVG (KZ) ! dz: (positive) depth of each single layer; zsoil: (negative) accumulated depth from ground + DO KZ = 2,NVG ! nvg: total layer of vegetated ground, = 4 in here (urban-nbs) ZSOILG (KZ) = - DZVG(KZ) + ZSOILG (KZ -1) END DO - TTP = TCP + IF (TTscheme .eq. 0 ) THEN ! Option1: residual approach / diagnostic approach + TTP = TCP + ELSEIF (TTscheme .eq. 1 ) THEN ! Option2: resistance approach / prognostic approach + TTP = TT + ENDIF + + ES=6.11*EXP( (2.5*10**6/461.51)*(TTP-273.15)/(273.15*TTP) ) + QS0T=0.622*ES/(PS-0.378*ES) + RA_LEAF = 58 * (LAI_TREE**0.56) * SQRT(DLEAF/UC) + + ! Option1: mid-layer approach + !W2 = 0.5 * (SMG(NVG/2-1) + SMG(NVG/2)) + + ! Option2: weighted averaged approach + W2_num = 0.0 + W2_den = 0.0 + DO K = 1, 3 + W2_num = W2_num + SMG(K) * DZVG(K) + W2_den = W2_den + DZVG(K) + END DO + W2 = W2_num / W2_den + + CALL STOMATAL(RS_LEAF, RSMIN, SLEAF, LAI_TREE, W2, SMCWLT, SMCMAX, TCP, TTP, QCP, QS0T) + DO ITERATION=1,100 - ! print* , 'Current iteration:', ITERATION + + TBP = MAX(223.15, MIN(353.15, TBP)) + TGP = MAX(223.15, MIN(353.15, TGP)) + TVGP = MAX(223.15, MIN(353.15, TVGP)) + KZ=1 ES=6.11*EXP( (2.5*10**6/461.51)*(TVGP-273.15)/(273.15*TVGP) ) DESDT=(2.5*10**6/461.51)*ES/(TVGP**2) QS0VG=0.622*ES/(PS-0.378*ES) DQS0VGDTVG = DESDT*0.622*PS/((PS-0.378*ES)**2) EPVG=RHOO*CHVG*UC*(QS0VG-QCP) + + SMGP = MAX(SMGP, SMCWLT) !urban-nbs + SMGP = MIN(SMGP, SMCMAX) - - !print*, 'vegetated ground temperature is:', TVGP ! add by yuqi for testing purpose - !print*, 'potential et from vegetated ground is:', EPVG ! add by yuqi for testing purpose - !print*, 'soil moisture for vegetated ground is:', SMGP ! add by yuqi for testing purpose - - VFAC = 1.0 - EXP(-0.52 * LAI_VEG) ! following PhenologyMainMod.F90 in Noah_MP (yuqi) + VFAC = 1.0 - EXP(-0.52 * LAI_VEG) ! following PhenologyMainMod.F90 in Noah_MP (urban-nbs) IF (EPVG > 0.0) THEN - CALL DIREVAP (EDIRG,EPVG,SMGP(KZ),VFAC,SMCMAX,SMCDRY,FXEXP) ! CONTRIBUTED FROM BARE SOIL (YUQI) - CALL TRANSP (ETTG,ETG,ECG,VFAC,EPVG,CMCGP,CFACTR,CMCMAX,LAI_VEG,RSMIN,RSMAX,RGL,SX, & ! CONTRIBUTED FROM VEGETED, USE TC AND QC INSTEAD OF TA AND QA (YUQI) - TVGP,TCP,QCP,SMGP,SMCWLT,SMCREF,CPP,PS,CHVG,EPSVG,DELT,NROOT,NVG,DZVG, & ! NROOT IS NUMBER OF VEGE ROOT, SAME AS IN GREEN ROOF (YUQI) + CALL DIREVAP (EDIRG,EPVG,SMGP(KZ),VFAC,SMCMAX,SMCDRY,FXEXP) ! CONTRIBUTED FROM BARE SOIL (urban-nbs) + CALL TRANSP (ETTG,ETG,ECG,VFAC,EPVG,CMCGP,CFACTR,CMCMAX,LAI_VEG,RSMIN,RSMAX,RGL,SX, & ! CONTRIBUTED FROM VEGETED, USE TC AND QC INSTEAD OF TA AND QA (urban-nbs) + TVGP,TCP,QCP,SMGP,SMCWLT,SMCREF,CPP,PS,CHVG,EPSVG,DELT,NROOT,NVG,DZVG, & ! NROOT IS NUMBER OF VEGE ROOT, SAME AS IN GREEN ROOF (urban-nbs) ZSOILG,HS) - CALL SMFLX (SMGP,SMG,NVG,CMCGP,CMCG,DELT,RAIN,ZSOILG,SMCMAX,BEXP,SMCWLT,SROOTP,DKSAT,& ! NVG: LAYER OF VEGETATED GROUND (YUQI) + CALL SMFLX (SMGP,SMG,NVG,CMCGP,CMCG,DELT,RAIN,ZSOILG,SMCMAX,BEXP,SMCWLT,SROOTP,DKSAT,& ! NVG: LAYER OF VEGETATED GROUND (urban-nbs) DWSAT,VFAC,CMCMAX,SRUN1,SRUN2,SRUN3,EDIRG,ECG,ETG,DRIPG) ELSE DEW = - EPVG @@ -1805,8 +1833,10 @@ SUBROUTINE urban(LSOLAR, & ! L CALL SMFLX (SMGP,SMG,NVG,CMCGP,CMCG,DELT,RAINDR,ZSOILG,SMCMAX,BEXP,SMCWLT,SROOTP,DKSAT, & DWSAT,VFAC,CMCMAX,SRUN1,SRUN2,SRUN3,EDIRG,ECG,ETG,DRIPG) END IF - ! print*, 'Vegetated soil moisture=', SMG ! add by yuqi for testing purpose + SMG = MAX(SMG, SMCWLT) + SMG = MIN(SMG, SMCMAX) + EDIRG = EDIRG * 1000 ETTG = ETTG * 1000 ECG = ECG * 1000 @@ -1827,15 +1857,11 @@ SUBROUTINE urban(LSOLAR, & ! L ES=6.11*EXP( (2.5*10**6/461.51)*(TGP-273.15)/(273.15*TGP) ) DESDT=(2.5*10**6/461.51)*ES/(TGP**2) QS0G=0.622*ES/(PS-0.378*ES) - DQS0GDTG=DESDT*0.622*PS/((PS-0.378*ES)**2) + DQS0GDTG=DESDT*0.622*PS/((PS-0.378*ES)**2) ES=6.11*EXP( (2.5*10**6/461.51)*(TTP-273.15)/(273.15*TTP) ) - DESDT=(2.5*10**6/461.51)*ES/(TTP**2) - QS0T=0.622*ES/(PS-0.378*ES) - - + QS0T=0.622*ES/(PS-0.378*ES) EPSGE = (1 - FVG) * EPSG + FVG * EPSVG - TGEP = (1 - FVG) * TGP + FVG * TVGP EMITW = EPSB * SIG * TBP**4./60 @@ -1884,7 +1910,6 @@ SUBROUTINE urban(LSOLAR, & ! L RVG= RVG1+RVG2 ! Note: unit as cal/cm2/s RT = RT1 + RT2 ! Note: unit as cal/cm2/s RLEAF = RT * (2 * PI * RTREE * (1-TAU))/(2 * RTREE * LAI_TREE) *697.7*60. ! unit: w/m/m - ! print*,'RB, RVG, RT, RLEAF=', RB, RVG, RT, RLEAF ! add by yuqi for testing purpose DRBDTB1 = EPSB * (4 * EPSB * SIG * TB**3 * VFWW_TREE - 4 * SIG * TB**3)/60 @@ -1939,31 +1964,34 @@ SUBROUTINE urban(LSOLAR, & ! L ELEB=RHO*EL*CHB*UC*BETB*(QS0B-QCP)*100 ! Note: unit [cal/cm2/s] ELEG=RHO*EL*CHG*UC*BETG*(QS0G-QCP)*100 ! Note: unit [cal/cm2/s] - RA_LEAF = 58 * (LAI_TREE**0.56) * SQRT(DLEAF/UC) ! DLEAF = 0.1 (YUQI) - W2 = 0.5 * (SMGP(NVG/2-1) + SMGP(NVG/2)) - CALL STOMATAL(RS_LEAF, RSMIN, SLEAF, LAI_TREE, W2, SMCREF, SMCMAX, TCP, TT, QCP, QS0T) + CALL LELEAF(ELET, SLEAF, RLEAF, RA_LEAF, RS_LEAF, RHOO, TCP, QCP, QS0T) AVAILABLE_WATER_SUM = 0.0 DO K = 1, 3 - AVAILABLE_WATER_SUM = AVAILABLE_WATER_SUM + ((SMG(K)-SMCWLT) * DZVG(K)) ! unit: m (Yuqi) + AVAILABLE_WATER_SUM = AVAILABLE_WATER_SUM + ((SMG(K)-SMCWLT) * DZVG(K)) ! unit: m (urban-nbs) END DO IF ((RTREE > 0.) .AND. (LAI_TREE > 0.)) THEN POTENTIAL_ELET = (ELL * 1000 * AVAILABLE_WATER_SUM * FVG * RW / DELT) / & - (2.0 * RTREE * LAI_TREE) ! [J/kg]*[1000 kg/m3 - density]*[m]/[s] = [W/m2] + (2.0 * RTREE * LAI_TREE) ! [J/kg]*[1000 kg/m3 - density]*[m]/[s] = [W/m2] ELSE POTENTIAL_ELET = 0.0 END IF ELET = MIN(POTENTIAL_ELET, ELET) ! unit: W/m2 - CALL ROOT_UPTAKE(SROOT, NVG, SMG, SMCMAX, SMCDRY, RTREE_M, LAI_TREE, DZVG, ZSOILG, ELET, FVG, RW_M) ! SROOT unit: m/s (Yuqi) - + CALL ROOT_UPTAKE(SROOT, NVG, SMG, SMCMAX, SMCDRY, RTREE_M, LAI_TREE, DZVG, ZSOILG, ELET, FVG, RW_M) ! SROOT unit: m/s (urban-nbs) + SROOTP = SROOT ! update root water uptake (urban-nbs) - ! residual approach (Yuqi) - HT = (SLEAF + RLEAF - ELET - (TT - TTP)/DELT*640) * (1/697.7/60.) ! following e.q 28 in (Ryu, 2016), unit same as HB, HG in [cal/cm2/s] - ! ALPHAT = (RHO*CP*100.)/(1.27*RA_LEAF) ! unit [cal/cm2/s/K] ([g/cm3] * [cal/g/K] / [s/m] * [m/cm]) % commented out, no longer used for residual approach (Chenghao Wang) + IF (TTscheme .eq. 0 ) THEN ! Option1: residual approach / diagnostic approach + HT = (SLEAF + RLEAF - ELET - (TT - TTP)/DELT*640) * (1/697.7/60.) ! following e.q 28 in (Ryu, 2016), unit same as HB, HG in [cal/cm2/s] + ELSEIF (TTscheme .eq. 1 ) THEN ! Option2: resistance approach / prognostic approach + ALPHAT = (RHO * CP * 100) / (1.274 * RA_LEAF) ! unit [cal/cm2/s/K] ([g/cm3] * [cal/g/K] / [s/m] * [m/cm]) + HT = ALPHAT * (TTP - TCP) + ENDIF + + DTCDTB = (W * ALPHAB)/(RW * ALPHAC + RW * FVG * ALPHAVG + RW * (1 - FVG) * ALPHAG + W * ALPHAB) DTCDTG = (RW * ALPHAG * (1 - FVG))/(RW * ALPHAC + RW * FVG * ALPHAVG + RW * (1 - FVG) * ALPHAG + W * ALPHAB) DTCDTVG= (RW * ALPHAVG * FVG)/(RW * ALPHAC + RW * FVG * ALPHAVG + RW * (1 - FVG) * ALPHAG + W * ALPHAB) @@ -1980,11 +2008,11 @@ SUBROUTINE urban(LSOLAR, & ! L DHVGDTG=RHO*CP*CHVG*UC*(0-DTCDTG)*100 DHVGDTVG=RHO*CP*CHVG*UC*(1-DTCDTVG)*100 - - DQCDTB=W*ALPHAB*BETB*DQS0BDTB/(W * ALPHAB * BETB + RW * (1 - FVG) * ALPHAG * BETG + RW * ALPHAC) - DQCDTG=RW * (1 - FVG) * ALPHAG * BETG * DQS0GDTG/(W * ALPHAB * BETB + RW * (1 - FVG) * ALPHAG * BETG + RW * ALPHAC) + den_qc = W * ALPHAB * BETB + RW * (1 - FVG) * ALPHAG * BETG + RW * FVG * ALPHAVG * BETVG + RW * ALPHAC + DQCDTB=W*ALPHAB*BETB*DQS0BDTB/den_qc + DQCDTG=RW * (1 - FVG) * ALPHAG * BETG * DQS0GDTG/den_qc ! DQCDTVG=0 ! no QS0VG in QC2 - DQCDTVG=RW * FVG * ALPHAVG * BETVG * DQS0VGDTVG/(W * ALPHAB * BETB + RW * (1 - FVG) * ALPHAG * BETG + RW * ALPHAC) + DQCDTVG=RW * FVG * ALPHAVG * BETVG * DQS0VGDTVG/den_qc DELEBDTB=RHO*EL*CHB*UC*BETB*(DQS0BDTB-DQCDTB)*100 DELEBDTG=RHO*EL*CHB*UC*BETB*(0-DQCDTG)*100 @@ -1998,15 +2026,14 @@ SUBROUTINE urban(LSOLAR, & ! L ELEVG = ETAG * EL * 0.1 ! unit: cal/cm2/s - ! print*, 'actual et from vegetated ground=', ELEVG ! add by yuqi for testing purpose - G0B=AKSB*(TBP-TBL(1))/(DZB(1)/2) ! unit: cal/cm2/s (yuqi) + G0B=AKSB*(TBP-TBL(1))/(DZB(1)/2) ! unit: cal/cm2/s G0G=AKSG*(TGP-TGL(1))/(DZG(1)/2) CALL TDFCND (KVG,SMG(KZ), QUARTZ, SMCMAX ) - KVG = KVG * EXP(-2.0 * VFAC) - G0VG= KVG * (TVGP-TVGL(1))/(DZVG(1)/2.)/697.7/60 ! unit: cal/cm2/s (Yuqi) - ! print*, 'G0B, G0G, G0VG=', G0B, G0G, G0VG + !KVG = KVG * EXP(-2.0 * VFAC) ! comment this to allow heat to be transferred efficiently to deeper soil (urban-nbs, Yuqi) + G0VG= KVG * (TVGP-TVGL(1))/(DZVG(1)/2.)/697.7/60 ! unit: cal/cm2/s + DG0BDTB = 2 * AKSB/DZB(1) DG0BDTG = 0 @@ -2074,17 +2101,14 @@ SUBROUTINE urban(LSOLAR, & ! L TB = TBP + DTB TG = TGP + DTG TVG= TVGP+DTVG - - - ! print*, 'DTB, DTG, DTVG=', DTB, DTG, DTVG ! add by yuqi for testing purpose + TGE = (1 - FVG) * TG + FVG * TVG TBP = TB TGP = TG TVGP= TVG - - ! print*, 'TBP, TGP, TVGP=',TBP, TGP, TVGP ! add by yuqi for testing purpose - - ! After updating the temperatures (e.g., TBP = TB) (yuqi) + TGEP=TGE + + ! After updating the temperatures (e.g., TBP = TB) (urban-nbs) TBP = MAX(223.15, MIN(353.15, TBP)) !-50C to +80C TGP = MAX(223.15, MIN(353.15, TGP)) TVGP = MAX(223.15, MIN(353.15, TVGP)) @@ -2097,13 +2121,15 @@ SUBROUTINE urban(LSOLAR, & ! L TC1 = RW * ALPHAC + RW * FVG * ALPHAVG + RW * (1 - FVG) * ALPHAG + W * ALPHAB TC2 = RW * ALPHAC * TA + RW * FVG * ALPHAVG * TVGP + RW * (1 - FVG) * ALPHAG * TGP + W * ALPHAB * TBP + 2 * RTREE * LAI_TREE * (HT/100) TC = TC2 / TC1 - ! Chenghao: HT unit is cgs, keep it as is [RHO - g/cm3, CP - cal/g/K] - ! print*, 'ALPHAT, ALPHAC, ALPHAVG, TC1, TC2, tmp=', ALPHAT, ALPHAC, ALPHAVG,TC1, TC2, 2 * RTREE * LAI_TREE * HT * ALPHAT ! add by yuqi - + TC = MAX(223.15, MIN(353.15, TC)) + TCP = TC ! update TCP (urban-nbs) + + ! Chenghao: HT unit is cgs, keep it as is [RHO - g/cm3, CP - cal/g/K] + QC1 = W * ALPHAB * BETB + RW * (1 - FVG) * ALPHAG * BETG + RW * ALPHAC QC2 = W * ALPHAB * BETB * QS0B + RW * (1 - FVG) * ALPHAG * BETG * QS0G + RW * FVG * ETAG *CP / 1000 + RW * ALPHAC * QA + (2 * RTREE * LAI_TREE * CP * (ELET/697.7/60.)/100/EL) QC=QC2/QC1 - ! Chenghao: ELET - [W/m2], convert it to cgs to be consistent with all other terms (e.g., ALPHAB) + ! Chenghao: ELET - [W/m2], convert it to cgs to be consistent with all other terms (e.g., ALPHAB) ES=6.11*EXP( (2.5*10**6/461.51)*(TC-273.15)/(273.15*TC) ) @@ -2119,26 +2145,26 @@ SUBROUTINE urban(LSOLAR, & ! L WRITE_MESSAGE(WARN_MSG) END IF - ! print*, 'ELET, QS0B, QS0G, QS0VG,QC1,QC2=', ELET, QS0B, QS0G, QS0VG, QC1,QC2 ! add by yuqi - TT = TC ! update TT (yuqi) - !TTP = TT ! update TTP (yuqi) + IF (TTscheme .eq. 0 ) THEN ! Option1: residual approach / diagnostic approach + TT = TC ! update TT (urban-nbs) + ELSEIF (TTscheme .eq. 1 ) THEN ! Option2: resistance approach / prognostic approach + TT = TTP + 0.1 * (TC - TTP) + ENDIF DTC = TCP - TC - TCP = TC ! update TCP (yuqi) QCP = QC - !END IF - ! print*, 'updated TT ,TC and QC=', TT, TC, QC ! add by yuqi for testing purpose + RVGT= (SVG+RVG) * 697.7 * 60. RCH = RHOO*CPP*CHVG - RR1 = EPSVG*(TVGP**4) * 6.48E-8 / (PS* CHVG) + 1.0 ! FOLLOWING JIACHUAN'S SCHEME, BUT DIFFERENT FROM OSU CODE (YUQI) - IF (RAIN > 0.0) THEN ! May cause some problem sicne we can't use TA to linearize TC, we can do that for TR (yuqi) + RR1 = EPSVG*(TVGP**4) * 6.48E-8 / (PS* CHVG) + 1.0 ! FOLLOWING YANG'S SCHEME, BUT DIFFERENT FROM OSU CODE (urban-nbs) + IF (RAIN > 0.0) THEN ! May cause some problem sicne we can't use TA to linearize TC, we can do that for TR (urban-nbs) RR2 = RR1 + RAIN / 3600 * 4.218E+3 / RCH ELSE RR2 = RR1 END IF - YY = TVG + (RVGT / RCH - BETVG * EPVG * ELL/ RCH) / RR2 ! THIS TVG IS WRONG (yuqi) + YY = TVG + (RVGT / RCH - BETVG * EPVG * ELL/ RCH) / RR2 ZZ1 = KVG / (-0.5 * ZSOILG (KZ) * RCH * RR2 ) + 1.0 SRUN3 = SRUN3/ DELT SRUN2 = SRUN2 + SRUN3 @@ -2151,21 +2177,18 @@ SUBROUTINE urban(LSOLAR, & ! L .AND. ABS(DTC) < 0.000001) EXIT END DO - TT = TC ! add by yuqi - TTP = TT ! update TTP (yuqi) CALL multi_layer(num_wall_layers,BOUNDB,G0B,CAPB,AKSB,TBL,DZB,DELT,TBLEND) CALL multi_layer(num_road_layers,BOUNDG,G0G,CAPG,AKSG,TGL,DZG,DELT,TGLEND) - ! CALL multi_layer(num_road_layers,BOUNDG,G0VG,CAPVG,KVG,TVGL,DZVG*0.01,DELT,TGLEND) ! DZVG unit ([m] as input, needed to convert to[cm]) (Chenghao 2026/01/20) - CALL SHFLX (SSOILG,TVGL,SMG,SMCMAX,NVG,TVGP,DELT,YY,ZZ1,ZSOILG,TGLEND,ZBOT,SMCWLT,KVG,QUARTZ,CSOIL,CAPVG) ! Update temperature in vegetated ground soil layer (yuqi) - ! print*, 'TVGL=', TVGL ! add by yuqi for testing purpose - !print*, 'Soil heat flux=', SSOILG ! add by yuqi for testing purpose + ! CALL multi_layer(num_road_layers,BOUNDG,G0VG,CAPVG,KVG/697.7/60,TVGL,DZVG*0.01,DELT,TGLEND) ! DZVG unit ([m] as input, needed to convert to[cm]) (Chenghao 2026/01/20) + CALL SHFLX (SSOILG,TVGL,SMG,SMCMAX,NVG,TVGP,DELT,YY,ZZ1,ZSOILG,TGLEND,ZBOT,SMCWLT,KVG,QUARTZ,CSOIL,CAPVG) ! Update temperature in vegetated ground soil layer (urban-nbs) + ELSE !------------------------------------------------------------------------------- ! TB, TG by Force-Restore Method - ! Note, Conflict With Vegetated Ground Scheme (Huang) + ! Note, Conflict With Vegetated Ground Scheme (urban-nbs) !------------------------------------------------------------------------------- ES=6.11*EXP( (2.5*10**6/461.51)*(TBP-273.15)/(273.15*TBP) ) @@ -2177,10 +2200,10 @@ SUBROUTINE urban(LSOLAR, & ! L ES=6.11*EXP( (2.5*10**6/461.51)*(TTP-273.15)/(273.15*TTP) ) QS0T=0.622*ES/(PS-0.378*ES) - ES=6.11*EXP( (2.5*10**6/461.51)*(TVGP-273.15)/(273.15*TVGP) ) ! following green roof scheme (yuqi) + ES=6.11*EXP( (2.5*10**6/461.51)*(TVGP-273.15)/(273.15*TVGP) ) ! following green roof scheme (urban-nbs) QS0VG=0.622*ES/(PS-0.378*ES) - EPSGE = (1 - FVG) * EPSG + FVG * EPSVG ! DEFINE NEW EPSVG HERE (YUQI) + EPSGE = (1 - FVG) * EPSG + FVG * EPSVG ! DEFINE NEW EPSVG HERE (urban-nbs) TGEP = (1 - FVG) * TGP + FVG * TVGP QS0GE = (1 - FVG) * QS0G + FVG * QS0VG ALPHAGE = (1 - FVG) * ALPHAG + FVG * ALPHAVG @@ -2241,7 +2264,7 @@ SUBROUTINE urban(LSOLAR, & ! L RA_LEAF = 58 * (LAI_TREE**0.56) * SQRT(DLEAF/UC) HT = (SLEAF + RLEAF - ELET - (TT - TTP)/DELT*640)/697.7/60. ! HT in unit of cal/cm2/s (unit consistent with HB and HG) - W2 = 0.5 * (SMGP(NVG/2-1) + SMGP(NVG/2)) ! POTENTIAL ERROR (YUQI) + W2 = 0.5 * (SMGP(NVG/2-1) + SMGP(NVG/2)) ! POTENTIAL ERROR (Yuqi urban-nbs) CALL STOMATAL(RS_LEAF, RSMIN, SLEAF, LAI_TREE, W2, SMCREF, SMCMAX, TC, TT, QC, QS0T) CALL LELEAF(ELET, SLEAF, RLEAF, RA_LEAF, RS_LEAF, RHOO, TC, QC, QS0T) @@ -2256,12 +2279,12 @@ SUBROUTINE urban(LSOLAR, & ! L TC1 = RW * ALPHAC + RW * FVG * ALPHAVG + RW * (1 - FVG) * ALPHAG + W * ALPHAB TC2 = RW * ALPHAC * TA + RW * FVG * ALPHAVG * TVGP + RW * (1 - FVG) * ALPHAG * TGP + W * ALPHAB * TBP + 2 * RTREE * LAI_TREE * (HT/100) TC = TC2 / TC1 - ! Chenghao: HT unit is cgs, keep it as is [RHO - g/cm3, CP - cal/g/K] + ! Chenghao: HT unit is cgs, keep it as is [RHO - g/cm3, CP - cal/g/K] QC1 = W * ALPHAB * BETB + RW * (1 - FVG) * ALPHAG * BETG + RW * ALPHAC QC2 = W * ALPHAB * BETB * QS0B + RW * (1 - FVG) * ALPHAG * BETG * QS0G + RW * FVG * ETAG *CP / 1000 + RW * ALPHAC * QA + (2 * RTREE * LAI_TREE * CP * (ELET/697.7/60.)/100/EL) QC=QC2/QC1 - ! Chenghao: ELET - [W/m2], convert it to cgs to be consistent with all other terms (e.g., ALPHAB) + ! Chenghao: ELET - [W/m2], convert it to cgs to be consistent with all other terms (e.g., ALPHAB) TCP=TC QCP=QC @@ -2275,7 +2298,7 @@ SUBROUTINE urban(LSOLAR, & ! L FLXHUMB=ELEB/RHO/EL/100. ! ELEB is cgs unit, final unit m/s FLXTHG=HG/RHO/CP/100. FLXHUMG=ELEG/RHO/EL/100. ! ELEG is cgs unit, final unit m/s - FLXTHVG=HVG/RHO/CP/100. ! ADD BY (Huang) + FLXTHVG=HVG/RHO/CP/100. ! ADD IN urban-nbs FLXHUMVG=ELEVG/RHO/EL/100. ! ELEVG is cgs unit, final unit m/s FLXTHT=HT/RHO/CP/100. FLXHUMT=ELET*(1/697.7/60.)/RHO/EL/100. ! ELET is W/m2, final unit m/s @@ -2386,7 +2409,7 @@ SUBROUTINE urban(LSOLAR, & ! L PSIH = -5. * XXX ELSE X = (1.-16.*XXX)**0.25 -! #ifdef DOUBLE_PRECISION ! newly added in thie version (yuqi) +! #ifdef DOUBLE_PRECISION ! newly added in this version (nbs comment) ! Code follows HRLDAS here (2025 summer) (Yuqi) ! PSIM = 2.*DLOG((1.+X)/2.) + DLOG((1.+X*X)/2.) - 2.*ATAN(X) + PI/2. ! PSIH = 2.*DLOG((1.+X*X)/2.) @@ -2395,7 +2418,7 @@ SUBROUTINE urban(LSOLAR, & ! L PSIH = 2.*ALOG((1.+X*X)/2.) !#endif END IF -!#ifdef DOUBLE_PRECISION !newly added in this version (yuqi) +!#ifdef DOUBLE_PRECISION !newly added in this version (nbs comment) ! GZ1OZ0 = DLOG(Z/Z0) ! CD = 0.4**2./(DLOG(Z/Z0)-PSIM)**2. @@ -2416,7 +2439,7 @@ SUBROUTINE urban(LSOLAR, & ! L !m TS = TA + FLXTH/CH/UA ! surface potential temp (flux temp) !m QS = QA + FLXHUM/CH/UA ! surface humidity ! - !IF (distributed_aerodynamics_option) THEN ! comment this in the new version (yuqi) + !IF (distributed_aerodynamics_option) THEN ! comment this in the new version (urban-nbs comment) ! TS = TA + FLXTH / (ALPHAC / (RHO * CP)) ! surface potential temp (flux temp) ! QS = QA + FLXHUM / (ALPHAC / (RHO * CP)) ! surface humidity !ELSE @@ -2439,7 +2462,7 @@ SUBROUTINE urban(LSOLAR, & ! L X = (1.-16.*XXX2)**0.25 !#ifdef DOUBLE_PRECISION -! PSIM2 = 2.*DLOG((1.+X)/2.) + DLOG((1.+X*X)/2.) - 2.*ATAN(X) + 2.*ATAN(1.) ! newly added in this version (yuqi) +! PSIM2 = 2.*DLOG((1.+X)/2.) + DLOG((1.+X*X)/2.) - 2.*ATAN(X) + 2.*ATAN(1.) ! newly added in this version (nbs comment) ! PSIH2 = 2.*DLOG((1.+X*X)/2.) !#else PSIM2 = 2.*ALOG((1.+X)/2.) + ALOG((1.+X*X)/2.) - 2.*ATAN(X) + 2.*ATAN(1.) @@ -2450,7 +2473,7 @@ SUBROUTINE urban(LSOLAR, & ! L ! !#ifdef DOUBLE_PRECISION -! CHS2 = 0.4*UST/(DLOG(2./Z0H)-PSIH2) ! cenlin 03/09/2023: uncomment to allow CHS2 calc for offline hrldas-urban ! newly commented in this version following WRF v4.7.1 (yuqi) +! CHS2 = 0.4*UST/(DLOG(2./Z0H)-PSIH2) ! cenlin 03/09/2023: uncomment to allow CHS2 calc for offline hrldas-urban ! newly commented in this version following WRF v4.7.1 (urban-nbs comment) ! !#else ! CHS2 = 0.4*UST/(ALOG(2./Z0H)-PSIH2) ! cenlin 03/09/2023: uncomment to allow CHS2 calc for offline hrldas-urban @@ -2467,15 +2490,15 @@ SUBROUTINE urban(LSOLAR, & ! L ELSE X = (1.-16.*XXX10)**0.25 !#ifdef DOUBLE_PRECISION -! PSIM10 = 2.*DLOG((1.+X)/2.) + DLOG((1.+X*X)/2.) - 2.*ATAN(X) + 2.*ATAN(1.) ! newly added in this version (yuqi) -! PSIH10 = 2.*DLOG((1.+X*X)/2.) ! newly added in this version (yuqi) +! PSIM10 = 2.*DLOG((1.+X)/2.) + DLOG((1.+X*X)/2.) - 2.*ATAN(X) + 2.*ATAN(1.) ! newly added in this version (nbs comment) +! PSIH10 = 2.*DLOG((1.+X*X)/2.) ! newly added in this version (nbs comment) !#else PSIM10 = 2.*ALOG((1.+X)/2.) + ALOG((1.+X*X)/2.) - 2.*ATAN(X) + 2.*ATAN(1.) PSIH10 = 2.*ALOG((1.+X*X)/2.) !#endif END IF -!#ifdef DOUBLE_PRECISION ! newly added below in this version (yuqi) +!#ifdef DOUBLE_PRECISION ! newly added below in this version (urban-nbs comment) ! PSIX = DLOG(Z/Z0) - PSIM ! PSIT = DLOG(Z/Z0H) - PSIH @@ -2502,7 +2525,7 @@ SUBROUTINE urban(LSOLAR, & ! L ! TH2 = TS + (TA-TS)*(PSIT2/PSIT) ! potential temp at 2 m [K] ! TH2 = TS + (TA-TS)*(PSIT2/PSIT) ! Fei: this seems to be temp (not potential) at 2 m [K] !Fei: consistant with M-O theory - IF (distributed_aerodynamics_option) THEN ! comment this in the new version (yuqi) + IF (distributed_aerodynamics_option) THEN ! comment this in the new version (urban-nbs comment) CHS2_LOCAL = 0.4 * UST / (ALOG(2. / Z0H) - PSIH2) TH2 = TS + (TA - TS) * (CHS_LOCAL / CHS2_LOCAL) Q2 = QS + (QA - QS) * (CHS_LOCAL / CHS2_LOCAL) @@ -2549,7 +2572,7 @@ SUBROUTINE mos(XXX,ALPHA,CD,B1,RIB,Z,Z0,UA,TA,TSF,RHO) X=(1.-16.*XXX)**0.25 X0=(1.-16.*XXX0)**0.25 -!#ifdef DOUBLE_PRECISION ! newly added in the new version (yuqi) +!#ifdef DOUBLE_PRECISION ! newly added in the new version (nbs comment) ! PSIM=DLOG((Z+Z0)/Z0) & ! -DLOG((X+1.)**2.*(X**2.+1.)) & @@ -2564,7 +2587,7 @@ SUBROUTINE mos(XXX,ALPHA,CD,B1,RIB,Z,Z0,UA,TA,TSF,RHO) PSIM=ALOG((Z+Z0)/Z0) & -ALOG((X+1.)**2.*(X**2.+1.)) & +2.*ATAN(X) & - +ALOG((X0+1.)**2.*(X0**2.+1.)) & ! they replace first X as X0 in here (yuqi) + +ALOG((X0+1.)**2.*(X0**2.+1.)) & ! replace first X with X0 in the new version (urban-nbs comment) -2.*ATAN(X0) FAIH=1./SQRT(1.-16.*XXX) PSIH=ALOG((Z+Z0)/Z0)+0.4*B1 & @@ -2592,7 +2615,7 @@ SUBROUTINE mos(XXX,ALPHA,CD,B1,RIB,Z,Z0,UA,TA,TSF,RHO) XXX=0.714 !#ifdef DOUBLE_PRECISION -! PSIM=DLOG((Z+Z0)/Z0)+7.*XXX ! newly added in the new version (yuqi) +! PSIM=DLOG((Z+Z0)/Z0)+7.*XXX ! newly added in the new version (nbs comment) !#else PSIM=ALOG((Z+Z0)/Z0)+7.*XXX !#endif @@ -2600,7 +2623,7 @@ SUBROUTINE mos(XXX,ALPHA,CD,B1,RIB,Z,Z0,UA,TA,TSF,RHO) ELSE !#ifdef DOUBLE_PRECISION -! AL=DLOG((Z+Z0)/Z0) ! newly added in the new version (yuqi) +! AL=DLOG((Z+Z0)/Z0) ! newly added in the new version (nbs comment) !#else AL=ALOG((Z+Z0)/Z0) !#endif @@ -2610,7 +2633,7 @@ SUBROUTINE mos(XXX,ALPHA,CD,B1,RIB,Z,Z0,UA,TA,TSF,RHO) XXX=(AL+XKB-2.*RIB*7.*AL-SQRT(DD))/(2.*(RIB*7.**2-7.)) !#ifdef DOUBLE_PRECISION -! PSIM=DLOG((Z+Z0)/Z0)+7.*MIN(XXX,0.714) ! newly added in the new version (yuqi) +! PSIM=DLOG((Z+Z0)/Z0)+7.*MIN(XXX,0.714) ! newly added in the new version (nbs comment) !#else PSIM=ALOG((Z+Z0)/Z0)+7.*MIN(XXX,0.714) @@ -2645,7 +2668,7 @@ SUBROUTINE louis79(ALPHA,CD,RIB,Z,Z0,UA,RHO) !#ifdef DOUBLE_PRECISION -! A2=(0.4/DLOG(Z/Z0))**2. ! newly added in the new version (yuqi) +! A2=(0.4/DLOG(Z/Z0))**2. ! newly added in the new version (nbs comment) !#else A2=(0.4/ALOG(Z/Z0))**2. @@ -2687,7 +2710,7 @@ SUBROUTINE louis82(ALPHA,CD,RIB,Z,Z0,UA,RHO) REAL :: A2, FM, FH, CH, CHH !#ifdef DOUBLE_PRECISION -! A2=(0.4/DLOG(Z/Z0))**2. ! newly added in the new version (yuqi) +! A2=(0.4/DLOG(Z/Z0))**2. ! newly added in the new version (nbs comment) !#else A2=(0.4/ALOG(Z/Z0))**2. !#endif @@ -2800,9 +2823,9 @@ SUBROUTINE read_param(UTYPE, & ! in ZR,SIGMA_ZED,Z0C,Z0HC,ZDC,SVF,R,RW,HGT,AH, & ! out CAPR,CAPB,CAPG,AKSR,AKSB,AKSG,ALBR,ALBB,ALBG, & ! out EPSR,EPSB,EPSG,Z0R,Z0B,Z0G,Z0HB,Z0HG, & ! out - BETR,BETB,BETG,TRLEND,TBLEND,TGLEND, & ! out - RTREE,RTREE_M,RW_M,HTREE,DTREE,LAI_TREE,LAI_VEG, & ! out ! added by yuqi - Z0VG,Z0HVG,CAPVG,EPSVG,EPST, & ! out ! added by yuqi + BETR,BETB,BETG,TRLEND,TBLEND,TGLEND,ALBVG,ALBT,FVG, & ! out + RTREE,RTREE_M,RW_M,HTREE,DTREE,LAI_TREE,LAI_VEG, & ! out ! urban-nbs + Z0VG,Z0HVG,CAPVG,EPSVG,EPST, & ! out ! urban-nbs !for BEP NUMDIR, STREET_DIRECTION, STREET_WIDTH, & ! out BUILDING_WIDTH, NUMHGT, HEIGHT_BIN, & ! out @@ -2818,7 +2841,7 @@ SUBROUTINE read_param(UTYPE, & ! in SIGMA_ZED, & EPSR,EPSB,EPSG,Z0R,Z0B,Z0G,Z0HB,Z0HG, & BETR,BETB,BETG,TRLEND,TBLEND,TGLEND,RTREE,RTREE_M,RW_M,HTREE,& - DTREE,LAI_TREE,LAI_VEG,Z0VG,Z0HVG,CAPVG,EPSVG,EPST ! added by yuqi + DTREE,LAI_TREE,LAI_VEG,Z0VG,Z0HVG,CAPVG,EPSVG,EPST,ALBVG,ALBT,FVG ! added for urban-nbs REAL, INTENT(OUT) :: AKANDA_URBAN !for BEP INTEGER, INTENT(OUT) :: NUMDIR @@ -2870,7 +2893,7 @@ SUBROUTINE read_param(UTYPE, & ! in TRLEND= TRLEND_TBL(UTYPE) TBLEND= TBLEND_TBL(UTYPE) TGLEND= TGLEND_TBL(UTYPE) - RTREE= RTREE_TBL(UTYPE) ! added by yuqi + RTREE= RTREE_TBL(UTYPE) ! urban-nbs RTREE_M = RTREE_M_TBL(UTYPE) RW_M = RW_M_TBL(UTYPE) HTREE= HTREE_TBL(UTYPE) @@ -2882,15 +2905,15 @@ SUBROUTINE read_param(UTYPE, & ! in CAPVG= CAPVG_TBL(UTYPE) EPSVG= EPSVG_TBL(UTYPE) EPST = EPST_TBL(UTYPE) - - - + ALBVG= ALBVG_TBL(UTYPE) + ALBT = ALBT_TBL(UTYPE) BOUNDR= BOUNDR_DATA BOUNDB= BOUNDB_DATA BOUNDG= BOUNDG_DATA CH_SCHEME = CH_SCHEME_DATA TS_SCHEME = TS_SCHEME_DATA AKANDA_URBAN = AKANDA_URBAN_TBL(UTYPE) + FVG = FVG_TBL(UTYPE) !for BEP @@ -2940,7 +2963,7 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & INTEGER :: num_road_layers INTEGER :: dummy REAL :: DHGT, HGT, VFWS, VFGS - REAL :: fvg ! Yuqi + REAL :: fvg ! urban-nbs REAL, allocatable, dimension(:) :: ROOF_WIDTH REAL, allocatable, dimension(:) :: ROAD_WIDTH @@ -3085,11 +3108,13 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & if(allocate_status /= 0) FATAL_ERROR('Error allocating TGLEND_TBL in urban_param_init') ALLOCATE( FRC_URB_TBL(ICATE), stat=allocate_status ) if(allocate_status /= 0) FATAL_ERROR('Error allocating FRC_URB_TBL in urban_param_init') - ALLOCATE( RTREE_TBL(ICATE), stat=allocate_status ) ! added by yuqi + ALLOCATE( FVG_TBL(ICATE), stat=allocate_status ) + if(allocate_status /= 0) FATAL_ERROR('Error allocating FVG_TBL in urban_param_init') + ALLOCATE( RTREE_TBL(ICATE), stat=allocate_status ) if(allocate_status /= 0) FATAL_ERROR('Error allocating RTREE_TBL in urban_param_init') - ALLOCATE( RTREE_M_TBL(ICATE), stat=allocate_status ) ! added by yuqi + ALLOCATE( RTREE_M_TBL(ICATE), stat=allocate_status ) if(allocate_status /= 0) FATAL_ERROR('Error allocating RTREE_M_TBL in urban_param_init') - ALLOCATE( RW_M_TBL(ICATE), stat=allocate_status ) ! added by yuqi + ALLOCATE( RW_M_TBL(ICATE), stat=allocate_status ) if(allocate_status /= 0) FATAL_ERROR('Error allocating RW_M_TBL in urban_param_init') ALLOCATE( HTREE_TBL(ICATE), stat=allocate_status ) if(allocate_status /= 0) FATAL_ERROR('Error allocating HTREE_TBL in urban_param_init') @@ -3109,6 +3134,10 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & if(allocate_status /= 0) FATAL_ERROR('Error allocating EPSVG_TBL in urban_param_init') ALLOCATE( EPST_TBL(ICATE), stat=allocate_status ) if(allocate_status /= 0) FATAL_ERROR('Error allocating EPST_TBL in urban_param_init') + ALLOCATE( ALBVG_TBL(ICATE), stat=allocate_status ) + if(allocate_status /= 0) FATAL_ERROR('Error allocating ALBVG_TBL in urban_param_init') + ALLOCATE( ALBT_TBL(ICATE), stat=allocate_status ) + if(allocate_status /= 0) FATAL_ERROR('Error allocating ALBT_TBL in urban_param_init') ! ALLOCATE( ROOF_WIDTH(ICATE), stat=allocate_status ) ! if(allocate_status /= 0) FATAL_ERROR('Error allocating ROOF_WIDTH in urban_param_init') ! ALLOCATE( ROAD_WIDTH(ICATE), stat=allocate_status ) @@ -3287,8 +3316,7 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & else if (name == "DZGR") then read(string(indx+1:),*) dzgr(1:4) else if (name == "FVG") then - read(string(indx+1:),*) fvg - fvg_data = fvg + read(string(indx+1:),*) fvg_tbl(1:icate) else if (name == "RTREE") then read(string(indx+1:),*) rtree_tbl(1:icate) else if (name == "HTREE") then @@ -3305,7 +3333,7 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & read(string(indx+1:),*) z0hvg_tbl(1:icate) else if (name == "CAPVG") then read(string(indx+1:),*) capvg_tbl(1:icate) - ! Convert CAPVG_TBL from J m{-3} K{-1} to cal cm{-3} deg{-1} + ! Convert CAPVG_TBL from J m{-3} K{-1} to cal cm{-3} deg{-1} capvg_tbl = capvg_tbl * ( 1.0 / 4.1868 ) * 1.E-6 else if (name == "EPSVG") then read(string(indx+1:),*) epsvg_tbl(1:icate) @@ -3314,6 +3342,12 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & read(string(indx+1:),*) epst_tbl(1:icate) else if (name == "TREEOPTION") then read(string(indx+1:),*) treeoption + else if (name == "TTscheme") then + read(string(indx+1:),*) ttscheme + else if (name == "ALBVG") then + read(string(indx+1:),*) albvg_tbl(1:icate) + else if (name == "ALBT") then + read(string(indx+1:),*) albt_tbl(1:icate) !for BEP else if (name == "STREET PARAMETERS") then @@ -3413,7 +3447,7 @@ SUBROUTINE urban_param_init(DZR,DZB,DZG,num_soil_layers, & ! R: Normalized Roof Width (a.k.a. "building coverage ratio") R_TBL(LC) = ROOF_WIDTH(LC) / ( ROAD_WIDTH(LC) + ROOF_WIDTH(LC) ) - ! normalized rtree_tbl, dtree_tbl and htree_tbl ! added by yuqi to normalize tree dimension + ! normalized rtree_tbl, dtree_tbl and htree_tbl ! urban-nbs to normalize tree dimension RTREE_M_TBL(LC) = RTREE_TBL(LC) RTREE_TBL(LC) = RTREE_TBL(LC) / (ROAD_WIDTH(LC) + ROOF_WIDTH(LC)) DTREE_TBL(LC) = DTREE_TBL(LC) / (ROAD_WIDTH(LC) + ROOF_WIDTH(LC)) @@ -3529,7 +3563,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, CMCR_URB2D,TGR_URB2D,TGRL_URB3D,SMR_URB3D, & ! inout DRELR_URB2D,DRELB_URB2D,DRELG_URB2D, & ! inout FLXHUMR_URB2D, FLXHUMB_URB2D, FLXHUMG_URB2D, & ! inout - TVG_URB2D,XXXVG_URB2D,TVGL_URB3D,SMG_URB3D,TT_URB2D,CMCG_URB2D,FLXHUMVG_URB2D,FLXHUMT_URB2D, & ! inout (add by Yuqi Huang) + TVG_URB2D,XXXVG_URB2D,TVGL_URB3D,SMG_URB3D,TT_URB2D,CMCG_URB2D,FLXHUMVG_URB2D,FLXHUMT_URB2D, & ! inout (urban-nbs) A_U_BEP,A_V_BEP,A_T_BEP,A_Q_BEP, & ! inout multi-layer urban A_E_BEP,B_U_BEP,B_V_BEP, & ! inout multi-layer urban B_T_BEP,B_Q_BEP,B_E_BEP,DLG_BEP, & ! inout multi-layer urban @@ -3594,7 +3628,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: RN_URB2D REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TS_URB2D -!-----------add by Yuqi Huang ---------------- +!-----------add for urban-nbs ---------------- REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TVG_URB2D REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TT_URB2D REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: XXXVG_URB2D @@ -3689,7 +3723,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, !m !FS FRC_URB2D(I,J)=0. UTYPE_URB2D(I,J)=0 - distributed_aerodynamics_check: IF (distributed_aerodynamics_option) THEN ! comment this in the new version (yuqi) + distributed_aerodynamics_check: IF (distributed_aerodynamics_option) THEN ! comment this in the new version (urban-nbs comment) IF (IVGTYP(I, J) == ISURBAN) THEN UTYPE_URB2D(I, J) = 2 ELSE @@ -3702,7 +3736,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, IF(use_wudapt_lcz==0) THEN UTYPE_URB2D(I,J) = 2 ! for default. high-intensity ELSE - UTYPE_URB2D(I,J) = 5 ! for default. open mid-rise (Yuqi) + UTYPE_URB2D(I,J) = 5 ! for default. open mid-rise ENDIF ELSE IF( IVGTYP(I,J) == LCZ_1) THEN UTYPE_URB2D(I,J) = 1 @@ -3788,8 +3822,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, XXXR_URB2D(I,J)=0. XXXB_URB2D(I,J)=0. XXXG_URB2D(I,J)=0. - XXXC_URB2D(I,J)=0. - XXXVG_URB2D(I,J)=0. ! add by huang + XXXC_URB2D(I,J)=0. IF ( sf_urban_physics == 1 ) THEN DRELR_URB2D(I,J) = 0. @@ -3798,13 +3831,14 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, FLXHUMR_URB2D(I,J) = 0. FLXHUMB_URB2D(I,J) = 0. FLXHUMG_URB2D(I,J) = 0. - FLXHUMVG_URB2D(I,J) = 0. ! add by huang - FLXHUMT_URB2D(I,J) = 0. ! add by huang + FLXHUMVG_URB2D(I,J) = 0. ! add in urban-nbs + FLXHUMT_URB2D(I,J) = 0. ! add in urban-nbs CMCR_URB2D(I,J) = 0. - CMCG_URB2D(I,J) = 0. ! add by huang + CMCG_URB2D(I,J) = 0. ! add in urban-nbs TGR_URB2D(I,J)=TSURFACE0_URB(I,J)+0. - TVG_URB2D(I,J)=TSURFACE0_URB(I,J)+0. ! add by huang - TT_URB2D(I,J) =TSURFACE0_URB(I,J)+0. ! add by huang + TVG_URB2D(I,J)=TSURFACE0_URB(I,J)+0. ! add in urban-nbs + TT_URB2D(I,J) =TSURFACE0_URB(I,J)+0. ! add in urban-nbs + XXXVG_URB2D(I,J)=0. ! add in urban-nbs ENDIF TC_URB2D(I,J)=TSURFACE0_URB(I,J)+0. @@ -3832,7 +3866,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, TGRL_URB3D(I,3,J)=TLAYER0_URB(I,2,J)+0. TGRL_URB3D(I,4,J)=TLAYER0_URB(I,2,J)+(TLAYER0_URB(I,3,J)-TLAYER0_URB(I,2,J))*0.29 - TVGL_URB3D(I,1,J)=TLAYER0_URB(I,1,J)+0. ! the 4 lines below added by yuqi + TVGL_URB3D(I,1,J)=TLAYER0_URB(I,1,J)+0. ! 4 lines below - urban-nbs TVGL_URB3D(I,2,J)=0.5*(TLAYER0_URB(I,1,J)+TLAYER0_URB(I,2,J)) TVGL_URB3D(I,3,J)=TLAYER0_URB(I,2,J)+0. TVGL_URB3D(I,4,J)=TLAYER0_URB(I,2,J)+(TLAYER0_URB(I,3,J)-TLAYER0_URB(I,2,J))*0.29 @@ -3840,7 +3874,7 @@ SUBROUTINE urban_var_init(ISURBAN, TSURFACE0_URB,TLAYER0_URB,TDEEP0_URB,IVGTYP, SMR_URB3D(I,2,J)=0.2 SMR_URB3D(I,3,J)=0.2 SMR_URB3D(I,4,J)=0. - SMG_URB3D(I,1,J)=0.2 ! add by huang + SMG_URB3D(I,1,J)=0.2 ! add in urban-nbs SMG_URB3D(I,2,J)=0.2 SMG_URB3D(I,3,J)=0.2 SMG_URB3D(I,4,J)=0.0 @@ -4119,14 +4153,14 @@ SUBROUTINE SFCDIF_URB (ZLM,Z0,THZ0,THLM,SFCSPD,AKANDA,AKMS,AKHS,RLMO,CD,ZT_OUT,V ! OVER SOLID SURFACE (LAND, SEA-ICE). ! ---------------------------------------------------------------------- ! ---------------------------------------------------------------------- -! LECH'S SURFACE FUNCTIONS ! E.q A7&A8 (YUQI) +! LECH'S SURFACE FUNCTIONS ! E.q A7&A8 ! ---------------------------------------------------------------------- PSLMU (ZZ)= -0.96* log (1.0-4.5* ZZ) PSLMS (ZZ)= (ZZ/RFC) -2.076* (1. -1./ (ZZ +1.)) PSLHU (ZZ)= -0.96* log (1.0-4.5* ZZ) - PSLHS (ZZ)= ZZ * RFAC -2.076* (1. - exp(-1.2 * ZZ)) ! bug fixed in 4.7.1 (Yuqi) + PSLHS (ZZ)= ZZ * RFAC -2.076* (1. - exp(-1.2 * ZZ)) ! bug fixed in 4.7.1 ! ---------------------------------------------------------------------- -! PAULSON'S SURFACE FUNCTIONS ! E.q A5&A6 (YUQI) +! PAULSON'S SURFACE FUNCTIONS ! E.q A5&A6 ! ---------------------------------------------------------------------- PSPMU (XX)= -2.* log ( (XX +1.)*0.5) - log ( (XX * XX +1.)*0.5) +2.* ATAN (XX) - PIHF PSPMS (YY)= 5.* YY @@ -4143,13 +4177,13 @@ SUBROUTINE SFCDIF_URB (ZLM,Z0,THZ0,THLM,SFCSPD,AKANDA,AKMS,AKHS,RLMO,CD,ZT_OUT,V ! ---------------------------------------------------------------------- ! ZILFC = - CZIL * VKRM * SQVISC ! C.......ZT=Z0*ZTFC - ZU = Z0 ! momentum roughness length (yuqi), here set ZU = Z0 to initialize momentum roughness length - RDZ = 1./ ZLM ! ZLM is initial mixing length (yuqi) - CXCH = EXCM * RDZ ! used to limit the AKHS (kinematic heat transfer coefficient) and AKMS (kinematic momentum transfer coefficient) (yuqi) - DTHV = THLM - THZ0 ! Theta_V term in Obukhov Length, which is the potential temperature different of first layer atmosphere and Urban level (yuqi) + ZU = Z0 ! momentum roughness length, here set ZU = Z0 to initialize momentum roughness length + RDZ = 1./ ZLM ! ZLM is initial mixing length + CXCH = EXCM * RDZ ! used to limit the AKHS (kinematic heat transfer coefficient) and AKMS (kinematic momentum transfer coefficient) + DTHV = THLM - THZ0 ! Theta_V term in Obukhov Length, which is the potential temperature different of first layer atmosphere and Urban level ! ---------------------------------------------------------------------- -! BELJARS CORRECTION OF USTAR ! used to correct friction velocity (yuqi) +! BELJARS CORRECTION OF USTAR ! used to correct friction velocity ! ---------------------------------------------------------------------- DU2 = MAX (SFCSPD * SFCSPD,EPSU2) !cc If statements to avoid TANGENT LINEAR problems near zero @@ -4163,31 +4197,31 @@ SUBROUTINE SFCDIF_URB (ZLM,Z0,THZ0,THLM,SFCSPD,AKANDA,AKMS,AKHS,RLMO,CD,ZT_OUT,V ! ---------------------------------------------------------------------- ! ZILITINKEVITCH APPROACH FOR ZT ! ---------------------------------------------------------------------- - USTAR = MAX (SQRT (AKMS * SQRT (DU2+ WSTAR2)),EPSUST) ! give an initial value of friction velocity (here AKMS might from last iteration) (yuqi) + USTAR = MAX (SQRT (AKMS * SQRT (DU2+ WSTAR2)),EPSUST) ! give an initial value of friction velocity (here AKMS might from last iteration) ! ---------------------------------------------------------------------- ! KCL/TL Try Kanda approach instead (Kanda et al. 2007, JAMC) ! ZT = EXP (ZILFC * SQRT (USTAR * Z0))* Z0 - IF (PRESENT(VEGFRAC)) THEN ! commented out in HRLDAS (yuqi) + IF (PRESENT(VEGFRAC)) THEN ! commented out in HRLDAS ! Kawai et al. (2009) JAMC ZT = EXP (2.0-(AKANDA-0.9*VEGFRAC**0.29)*(SQVISC**2 * USTAR * Z0)**0.25)* Z0 ELSE - ZT = EXP (2.0-AKANDA*(SQVISC**2 * USTAR * Z0)**0.25)* Z0 ! thermal roughness length (yuqi) + ZT = EXP (2.0-AKANDA*(SQVISC**2 * USTAR * Z0)**0.25)* Z0 ! thermal roughness length END IF - ZSLU = ZLM + ZU ! effective momentum surface layer depth (yuqi) + ZSLU = ZLM + ZU ! effective momentum surface layer depth - ZSLT = ZLM + ZT ! effective thermal surface layer depth (yuqi) + ZSLT = ZLM + ZT ! effective thermal surface layer depth RLOGU = log (ZSLU / ZU) RLOGT = log (ZSLT / ZT) - RLMO = ELFC * AKHS * DTHV / USTAR **3 ! RLMO: 1/Obukhov length, AKHS is w'theta_v' term which is to be parameterized (yuqi) + RLMO = ELFC * AKHS * DTHV / USTAR **3 ! RLMO: 1/Obukhov length, AKHS is w'theta_v' term which is to be parameterized ! ---------------------------------------------------------------------- ! 1./MONIN-OBUKKHOV LENGTH-SCALE ! ---------------------------------------------------------------------- DO ITR = 1,ITRMX - ZETALT = MAX (ZSLT * RLMO,ZTMIN) ! before each iternation, first calculate ZETALT (which is dimensionless stability parameter, which is ζ=z/L) (yuqi) + ZETALT = MAX (ZSLT * RLMO,ZTMIN) ! before each iternation, first calculate ZETALT (which is dimensionless stability parameter, which is ζ=z/L) RLMO = ZETALT / ZSLT ZETALU = ZSLU * RLMO ZETAU = ZU * RLMO ! dimensionless stability parameter for momentum at surface layer, which is ζm=zm/L @@ -4195,7 +4229,7 @@ SUBROUTINE SFCDIF_URB (ZLM,Z0,THZ0,THLM,SFCSPD,AKANDA,AKMS,AKHS,RLMO,CD,ZT_OUT,V IF (ILECH .eq. 0) THEN - IF (RLMO .lt. 0.0)THEN ! unstable condition where ζ<0 (yuqi) + IF (RLMO .lt. 0.0)THEN ! unstable condition where ζ<0 XLU4 = 1. -16.* ZETALU XLT4 = 1. -16.* ZETALT XU4 = 1. -16.* ZETAU @@ -4211,7 +4245,7 @@ SUBROUTINE SFCDIF_URB (ZLM,Z0,THZ0,THLM,SFCSPD,AKANDA,AKMS,AKHS,RLMO,CD,ZT_OUT,V SIMM = PSPMU (XLU) - PSMZ + RLOGU PSHZ = PSPHU (XT) SIMH = PSPHU (XLT) - PSHZ + RLOGT - ELSE ! stable condition where ζ>0 (yuqi) + ELSE ! stable condition where ζ>0 ZETALU = MIN (ZETALU,ZTMAX) ZETALT = MIN (ZETALT,ZTMAX) PSMZ = PSPMS (ZETAU) @@ -4238,42 +4272,42 @@ SUBROUTINE SFCDIF_URB (ZLM,Z0,THZ0,THLM,SFCSPD,AKANDA,AKMS,AKHS,RLMO,CD,ZT_OUT,V END IF END IF ! ---------------------------------------------------------------------- -! BELJAARS CORRECTION FOR USTAR ! BELJAARS conflict with BelJARS (yuqi) +! BELJAARS CORRECTION FOR USTAR ! BELJAARS conflict with BelJARS ! ---------------------------------------------------------------------- - USTAR = MAX (SQRT (AKMS * SQRT (DU2+ WSTAR2)),EPSUST) ! update Ustar since AKMS has been updated by SIMM (yuqi) + USTAR = MAX (SQRT (AKMS * SQRT (DU2+ WSTAR2)),EPSUST) ! update Ustar since AKMS has been updated by SIMM !KCL/TL !ZT = EXP (ZILFC * SQRT (USTAR * Z0))* Z0 - IF (PRESENT(VEGFRAC)) THEN ! they comment this in the new version (yuqi) + IF (PRESENT(VEGFRAC)) THEN ! they comment this in the new version (urban-nbs comment) ! Kawai et al. (2009) JAMC ZT = EXP (2.0-(AKANDA-0.9*VEGFRAC**0.29)*(SQVISC**2 * USTAR * Z0)**0.25)* Z0 ELSE - ZT = EXP (2.0-AKANDA*(SQVISC**2 * USTAR * Z0)**0.25)* Z0 ! update ZT (yuqi) + ZT = EXP (2.0-AKANDA*(SQVISC**2 * USTAR * Z0)**0.25)* Z0 ! update ZT END IF ZSLT = ZLM + ZT RLOGT = log (ZSLT / ZT) USTARK = USTAR * VKRM - AKMS = MAX (USTARK / SIMM,CXCH) ! update AKMS and AKHS (yuqi) + AKMS = MAX (USTARK / SIMM,CXCH) ! update AKMS and AKHS AKHS = MAX (USTARK / SIMH,CXCH) ! IF (BTGH * AKHS * DTHV .ne. 0.0) THEN - WSTAR2 = WWST2* ABS (BTGH * AKHS * DTHV)** (2./3.) ! update WSTAR2 (yuqi) + WSTAR2 = WWST2* ABS (BTGH * AKHS * DTHV)** (2./3.) ! update WSTAR2 ELSE WSTAR2 = 0.0 END IF !----------------------------------------------------------------------- - RLMN = ELFC * AKHS * DTHV / USTAR **3 ! updata the 1/L (yuqi) + RLMN = ELFC * AKHS * DTHV / USTAR **3 ! updata the 1/L !----------------------------------------------------------------------- ! IF(ABS((RLMN-RLMO)/RLMA).LT.EPSIT) GO TO 110 !----------------------------------------------------------------------- - RLMA = RLMO * WOLD+ RLMN * WNEW ! performe a weighted average of the previous RLMO and the newly calculated RLMN (yuqi) + RLMA = RLMO * WOLD+ RLMN * WNEW ! performe a weighted average of the previous RLMO and the newly calculated RLMN !----------------------------------------------------------------------- - RLMO = RLMA ! assgin to the RLMO and finally get inversed MO length which is 1/L (yuqi) + RLMO = RLMA ! assgin to the RLMO and finally get inversed MO length which is 1/L END DO - CD = USTAR*USTAR/SFCSPD**2 ! calculate the drag coefficient (yuqi) - IF (PRESENT(ZT_OUT)) ZT_OUT = ZT ! they comment this in the new version (yuqi) + CD = USTAR*USTAR/SFCSPD**2 ! calculate the drag coefficient + IF (PRESENT(ZT_OUT)) ZT_OUT = ZT ! they comment this in the new version (urban-nbs comment) ! ---------------------------------------------------------------------- END SUBROUTINE SFCDIF_URB ! ---------------------------------------------------------------------- @@ -4374,7 +4408,7 @@ SUBROUTINE TRANSP (ETT,ET,EC,SHDFAC,ETP1,CMC,CFACTR,CMCMAX,LAI,RSMIN,RSMAX,RGL,S RCSOIL = MAX (RCSOIL,0.0001) - RC = RSMIN / (LAI * RCS * RCT * RCQ * RCSOIL) ! canopy resistance (yuqi) + RC = RSMIN / (LAI * RCS * RCT * RCQ * RCSOIL) ! canopy resistance DESDT = 0.622*SLV*EA/461.51/TA/TA/1013 DELTA = (SLV / CPP)* DESDT RR = (4.* EPSV *SIGMA * 287.04 / CPP)* (TA **4.)/ (TS * CH) + 1.0 @@ -4396,12 +4430,12 @@ SUBROUTINE TRANSP (ETT,ET,EC,SHDFAC,ETP1,CMC,CFACTR,CMCMAX,LAI,RSMIN,RSMAX,RGL,S DO K = 1, NROOT ET(K) = ETT1 * GX (K) / DENOM - ETT = ETT + ET (K) ! Evapotranspiration from all the layers (yuqi) + ETT = ETT + ET (K) ! Evapotranspiration from all the layers END DO IF (CMC > 0.0) THEN - EC = SHDFAC * ( ( CMC / CMCMAX ) ** CFACTR ) * ETP1 * 0.001 ! evaporation of precipitation intercepted by the canopy (yuqi) + EC = SHDFAC * ( ( CMC / CMCMAX ) ** CFACTR ) * ETP1 * 0.001 ! evaporation of precipitation intercepted by the canopy urban-nbs ELSE EC = 0.0 END IF @@ -4414,7 +4448,7 @@ END SUBROUTINE TRANSP ! ---------------------------------------------------------------------- SUBROUTINE SMFLX (SMCP,SMC,NSOIL,CMCP,CMC,DT,PRCP1,ZSOIL, & - & SMCMAX,BEXP,SMCWLT,SROOT,DKSAT,DWSAT, & ! added sroot variable (yuqi) + & SMCMAX,BEXP,SMCWLT,SROOT,DKSAT,DWSAT, & ! added sroot variable (urban-nbs) & SHDFAC,CMCMAX,RUNOFF1,RUNOFF2,RUNOFF3, & EDIR,EC,ET,DRIP) @@ -4441,17 +4475,17 @@ SUBROUTINE SMFLX (SMCP,SMC,NSOIL,CMCP,CMC,DT,PRCP1,ZSOIL, & ! ADD PRECIPITATION TO EXISTING CMC.IF RESULTING AMT EXCEEDS MAX CAPACITY, ! IT BECOMES DRIP AND WILL FALL TO THE GRND. ! ---------------------------------------------------------------------- - RHSCT = SHDFAC * PRCP1 * 0.001 /3600. - EC ! how many water currently in the canopy after rain (gain) and evaporation (loss) (yuqi) + RHSCT = SHDFAC * PRCP1 * 0.001 /3600. - EC ! how many water currently in the canopy after rain (gain) and evaporation (loss) (urban-nbs) DRIP = 0. TRHSCT = DT * RHSCT - EXCESS = CMCP + TRHSCT ! update the current canopy water (yuqi) + EXCESS = CMCP + TRHSCT ! update the current canopy water (urban-nbs) ! ---------------------------------------------------------------------- ! PCPDRP IS THE COMBINED PRCP1 AND DRIP (FROM CMCP) THAT GOES INTO THE ! SOIL ! ---------------------------------------------------------------------- IF (EXCESS > CMCMAX) DRIP = EXCESS - CMCMAX - PCPDRP = (1. - SHDFAC) * PRCP1 * 0.001 /3600. + DRIP / DT ! total surfce water influx (yuqi) + PCPDRP = (1. - SHDFAC) * PRCP1 * 0.001 /3600. + DRIP / DT ! total surfce water influx (urban-nbs) ! ---------------------------------------------------------------------- ! CALL SUBROUTINES SRT AND SSTEP TO SOLVE THE SOIL MOISTURE @@ -4511,7 +4545,7 @@ SUBROUTINE SRT (RHSTT,EDIR,ET,SMCP,NSOIL,PCPDRP,ZSOIL,DWSAT, & DDMAX (3) = (ZSOIL (2) - ZSOIL (3))* SMCAV DDMAX (3) = DDMAX (3)* (1.0- (SMCP (3) - SMCWLT)/ SMCAV) - DD = DDMAX(1)+DDMAX(2)+DDMAX(3) ! accumulated maximum holdable soil water for all three layers (yuqi) + DD = DDMAX(1)+DDMAX(2)+DDMAX(3) ! accumulated maximum holdable soil water for all three layers (urban-nbs) DT1 = DT/86400 KDT = 3.0 * DKSAT / PAR VAL = (1. - EXP ( - KDT * DT1)) @@ -4519,7 +4553,7 @@ SUBROUTINE SRT (RHSTT,EDIR,ET,SMCP,NSOIL,PCPDRP,ZSOIL,DWSAT, & PX = PCPDRP * DT IF (PX < 0.0) PX = 0.0 - INFMAX = (PX * (DDT / (PX + DDT)))/ DT ! identical to 'InfilRateMax' in Noah_MP schaake scheme (yuqi) + INFMAX = (PX * (DDT / (PX + DDT)))/ DT ! identical to 'InfilRateMax' in Noah_MP schaake scheme (urban-nbs) MXSMC = SMCP (1) CALL WDFCND (WDF,WCND,MXSMC,SMCMAX,BEXP,DKSAT,DWSAT) INFMAX = MAX (INFMAX,WCND) @@ -4527,7 +4561,7 @@ SUBROUTINE SRT (RHSTT,EDIR,ET,SMCP,NSOIL,PCPDRP,ZSOIL,DWSAT, & IF (PCPDRP > INFMAX) THEN - RUNOFF1 = PCPDRP - INFMAX ! runoff1 is excess infiltration runoff (yuqi) + RUNOFF1 = PCPDRP - INFMAX ! runoff1 is excess infiltration runoff (urban-nbs) PDDUM = INFMAX END IF END IF @@ -4539,8 +4573,8 @@ SUBROUTINE SRT (RHSTT,EDIR,ET,SMCP,NSOIL,PCPDRP,ZSOIL,DWSAT, & AI (1) = 0.0 BI (1) = WDF * DDZ / ( - ZSOIL (1) ) CI (1) = - BI (1) - DSMDZ = (SMCP (1) - SMCP (2) )/( - 0.5 * ZSOIL(2)) ! Eq.3.7.253 (i=1) in noah_mp decument (yuqi) - RHSTT (1) = (WDF * DSMDZ + WCND- PDDUM + EDIR + ET(1)+ SROOT(1))/ ZSOIL (1) ! should minus conductivity between two layers (yuqi) + DSMDZ = (SMCP (1) - SMCP (2) )/( - 0.5 * ZSOIL(2)) ! Eq.3.7.253 (i=1) in noah_mp decument (urban-nbs) + RHSTT (1) = (WDF * DSMDZ + WCND- PDDUM + EDIR + ET(1)+ SROOT(1))/ ZSOIL (1) ! should minus conductivity between two layers (urban-nbs) SSTT = WDF * DSMDZ + WCND+ EDIR + ET(1) ! DON'T KNOW WHAT IS THIS, maybe for testing purpose (yuqi) ! ---------------------------------------------------------------------- @@ -4567,7 +4601,7 @@ SUBROUTINE SRT (RHSTT,EDIR,ET,SMCP,NSOIL,PCPDRP,ZSOIL,DWSAT, & AI (K) = - WDF * DDZ / DENOM2 BI (K) = - ( AI (K) + CI (K) ) IF (K .eq. NSOIL-1) THEN - RUNOFF2 = 0.0 ! subsurface runoff for layer 2 (yuqi) + RUNOFF2 = 0.0 ! subsurface runoff for layer 2 END IF IF (K .ne. NSOIL-1) THEN WDF = WDF2 @@ -4610,8 +4644,8 @@ SUBROUTINE SSTEP (SMCP,SMC,CMCP,CMC,RHSTT,RHSCT,DT, & ! CREATE 'AMOUNT' VALUES OF VARIABLES TO BE INPUT TO THE ! TRI-DIAGONAL MATRIX ROUTINE. ! ---------------------------------------------------------------------- - DO K = 1,NSOIL-1 ! E.q. 3.7.265 in noah_mp document (yuqi) - RHSTT (K) = RHSTT (K) * DT ! update matric coefficient to solve the saturation excess water for each soil layer (yuqi) + DO K = 1,NSOIL-1 ! E.q. 3.7.265 in noah_mp document + RHSTT (K) = RHSTT (K) * DT ! update matric coefficient to solve the saturation excess water for each soil layer AI (K) = AI (K) * DT BI (K) = 1. + BI (K) * DT CI (K) = CI (K) * DT @@ -4640,7 +4674,7 @@ SUBROUTINE SSTEP (SMCP,SMC,CMCP,CMC,RHSTT,RHSCT,DT, & DDZ = - ZSOIL (1) DO K = 1,NSOIL-1 IF (K /= 1) DDZ = ZSOIL (K - 1) - ZSOIL (K) - SMCOUT (K) = SMCP (K) + CI (K) + WPLUS / DDZ ! previous soil moisture + updated soil water + soil water surplus from above layer (yuqi) + SMCOUT (K) = SMCP (K) + CI (K) + WPLUS / DDZ ! previous soil moisture + updated soil water + soil water surplus from above layer STOT = SMCOUT (K) IF (STOT > SMCMAX) THEN IF (K .eq. 1) THEN @@ -4653,14 +4687,14 @@ SUBROUTINE SSTEP (SMCP,SMC,CMCP,CMC,RHSTT,RHSCT,DT, & ELSE WPLUS = 0. END IF - SMC (K) = MAX ( MIN (STOT,SMCMAX),0.07 ) ! change 0.066 to 0.07 in order to avoid floating invalid error in root_uptake step (Yuqi) + SMC (K) = MAX ( MIN (STOT,SMCMAX),0.07 ) ! change 0.066 to 0.07 in order to avoid floating invalid error in root_uptake step (urban-nbs) END DO ! ---------------------------------------------------------------------- ! UPDATE CANOPY WATER CONTENT/INTERCEPTION (CMC). CONVERT RHSCT TO ! AN 'AMOUNT' VALUE AND ADD TO PREVIOUS CMC VALUE TO GET NEW CMC. ! ---------------------------------------------------------------------- - RUNOFF3 = WPLUS ! subsurface runoff for layer 3 (yuqi) + RUNOFF3 = WPLUS ! subsurface runoff for layer 3) CMC = CMCP + DT * RHSCT IF (CMC < 1.E-20) CMC = 0.0 CMC = MIN (CMC,CMCMAX) @@ -4736,7 +4770,7 @@ SUBROUTINE ROSR12 (P,A,B,C,D,DELTA,NSOIL) ! ---------------------------------------------------------------------- ! INITIALIZE EQN COEF C FOR THE LOWEST SOIL LAYER ! ---------------------------------------------------------------------- - C (NSOIL) = 0.0 ! follow E.q 3.7.294-3.7.301 in noah_mp document (yuqi) + C (NSOIL) = 0.0 ! follow E.q 3.7.294-3.7.301 in noah_mp document P (1) = - C (1) / B (1) DELTA (1) = D (1) / B (1) DO K = 2,NSOIL @@ -4754,7 +4788,7 @@ SUBROUTINE ROSR12 (P,A,B,C,D,DELTA,NSOIL) ! ---------------------------------------------------------------------- DO K = 2,NSOIL KK = NSOIL - K + 1 - P (KK) = P (KK) * P (KK +1) + DELTA (KK) ! E.q 3.7.301 in noah_mp document (yuqi) + P (KK) = P (KK) * P (KK +1) + DELTA (KK) ! E.q 3.7.301 in noah_mp document END DO ! ---------------------------------------------------------------------- END SUBROUTINE ROSR12 @@ -4778,7 +4812,7 @@ SUBROUTINE SHFLX (SSOIL,STC,SMC,SMCMAX,NSOIL,T1,DT,YY,ZZ1,ZSOIL, & REAL, INTENT(IN) :: DF1,DT,SMCMAX, SMCWLT, TBOT,YY, ZBOT,ZZ1, QUARTZ REAL, INTENT(IN) :: CSOIL, CAPR REAL, INTENT(INOUT) :: T1 - REAL, INTENT(OUT) :: SSOIL ! output : soil heat flux, w/m-2 (yuqi) + REAL, INTENT(OUT) :: SSOIL ! output : soil heat flux, w/m-2 REAL, DIMENSION(1:NSOIL), INTENT(IN) :: SMC,ZSOIL REAL, DIMENSION(1:NSOIL), INTENT(INOUT) :: STC REAL, DIMENSION(1:NSOIL) :: AI, BI, CI, STCF,RHSTS @@ -4912,12 +4946,12 @@ SUBROUTINE HRT (RHSTS,STC,SMC,SMCMAX,NSOIL,ZSOIL,YY,ZZ1, & IF (ITAVG) THEN CALL TBND (STC (K),TBOT,ZSOIL,ZBOT,K,NSOIL,TBK1) END IF - ELSE !(k=nsoil, for last layer) (yuqi) + ELSE !(k=nsoil, for last layer) ! ---------------------------------------------------------------------- ! SPECIAL CASE OF BOTTOM LAYER (CONCRETE ROOF) ! ---------------------------------------------------------------------- HCPCT = CAPR * 4.1868 * 1.E6 - DF1K = 3.24 ! default in NOAH_MP technical notebook page 47 (yuqi) + DF1K = 3.24 ! default in NOAH_MP technical notebook page 47 ! ---------------------------------------------------------------------- ! CALC THE VERTICAL TEMP GRADIENT THRU BOTTOM LAYER. ! ---------------------------------------------------------------------- @@ -5143,7 +5177,7 @@ END FUNCTION kanda_kawai_svf ! ---------------------------------------------------------------------- -! SHADOW_TREE : CALCULATE SHADOW CAST BY TREES AND WALLS (Huang) +! SHADOW_TREE : CALCULATE SHADOW CAST BY TREES AND WALLS (urban-nbs) ! ---------------------------------------------------------------------- SUBROUTINE SHADOW_TREE (XI, H, W, HTREE, RTREE, DTREE, CHI_SHADED, ETA_SHADED, & CHI_SHADED_TREE, ETA_SHADED_TREE) @@ -5203,7 +5237,7 @@ END SUBROUTINE SHADOW_TREE ! ---------------------------------------------------------------------- ! ---------------------------------------------------------------------- -! STOMATAL : CALCULATE LEAF STOMATAL RESISTANCE (Huang) +! STOMATAL : CALCULATE LEAF STOMATAL RESISTANCE (urban-nbs) ! ---------------------------------------------------------------------- SUBROUTINE STOMATAL(RS, RSMIN, SD_FACET, LAI, W2, WR, WS, TA, TS, QA, QSAT) @@ -5251,7 +5285,7 @@ END SUBROUTINE STOMATAL ! ---------------------------------------------------------------------- -! LELEAF : CALCULATE LATENT HEAT OF TREE LEAF (Huang) +! LELEAF : CALCULATE LATENT HEAT OF TREE LEAF (urban-nbs) ! ---------------------------------------------------------------------- SUBROUTINE LELEAF (LE, SW, LW, RA_LEAF, RS_LEAF, RA, TA, QA, QSAT) @@ -5281,7 +5315,7 @@ END SUBROUTINE LELEAF ! ---------------------------------------------------------------------- -! force_restore_tree (Huang) +! force_restore_tree (urban-nbs) ! Refer to Sang and Soon, BLM, 2007 ! ---------------------------------------------------------------------- @@ -5303,7 +5337,7 @@ END SUBROUTINE force_restore_tree ! ---------------------------------------------------------------------- -! ROOT WATER UPTAKE (Huang) +! ROOT WATER UPTAKE (urban-nbs) ! FOLLOWING NOAH_MP, C. He, P. Valayamkunnath, & refactor team (He et al. 2023) ! ---------------------------------------------------------------------- diff --git a/phys/module_surface_driver.F b/phys/module_surface_driver.F index 50cbc820be..22f96580fb 100644 --- a/phys/module_surface_driver.F +++ b/phys/module_surface_driver.F @@ -112,9 +112,8 @@ SUBROUTINE surface_driver( & & ,zi3d, watsat3d, csol3d, tkmg3d & !lake & ,tkdry3d, tksatu3d, LakeModel, lake_min_elev & !lake #if (EM_CORE==1) - ! & ,lakemask, lakeflag & !lake - & ,lakemask & !lake - , restart_flag & ! restart_flag + & ,lakemask & !lake + , restart_flag & ! restart_flag #endif ! KIM TOFD ,kim_tofd & ! kim tofd @@ -2477,6 +2476,7 @@ SUBROUTINE surface_driver( & its=i_start(ij),ite=i_end(ij),jts=j_start(ij), & jte=j_end(ij), kts=kts,kte=kte, & ustm=ustm,ck=ck,cka=cka,cd=cd,cda=cda, & + xice=xice,xice_threshold=xice_threshold, & sf_mynn_sfcflux_land=sf_mynn_sfcflux_land, & sf_mynn_sfcflux_water=sf_mynn_sfcflux_water, & isfflx=isfflx,shalwater_z0=shalwater_z0, & @@ -5374,8 +5374,9 @@ SUBROUTINE mynn_seaice_wrapper(U3D,V3D,T3D,QV3D, & ids,ide, jds,jde, kds,kde, & ims,ime, jms,jme, kms,kme, & its,ite, jts,jte, kts,kte, & - ustm,ck,cka,cd,cda,sf_mynn_sfcflux_land, & - sf_mynn_sfcflux_water,flag_lsm,restart, & + ustm,ck,cka,cd,cda, & + sf_mynn_sfcflux_land,sf_mynn_sfcflux_water, & + flag_lsm,restart, & shalwater_z0,water_depth,lakemask, & errmsg,errflg ) @@ -5565,6 +5566,7 @@ SUBROUTINE mynn_seaice_wrapper(U3D,V3D,T3D,QV3D, & ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme,& its=its,ite=ite,jts=jts,jte=jte,kts=kts,kte=kte,& ustm=ustm,ck=ck,cka=cka,cd=cd,cda=cda, & + xice=xice,xice_threshold=xice_threshold, & sf_mynn_sfcflux_land=sf_mynn_sfcflux_land, & sf_mynn_sfcflux_water=sf_mynn_sfcflux_water, & isfflx=isfflx,shalwater_z0=shalwater_z0, & @@ -5632,6 +5634,7 @@ SUBROUTINE mynn_seaice_wrapper(U3D,V3D,T3D,QV3D, & ims=ims,ime=ime, jms=jms,jme=jme, kms=kms,kme=kme, & its=its,ite=ite, jts=jts,jte=jte, kts=kts,kte=kte, & ustm=ustm_sea,ck=ck_sea,cka=cka_sea,cd=cd_sea,cda=cda_sea, & + xice=xice,xice_threshold=xice_threshold, & sf_mynn_sfcflux_land=sf_mynn_sfcflux_land, & sf_mynn_sfcflux_water=sf_mynn_sfcflux_water,isfflx=isfflx, & flag_lsm=flag_lsm,restart=restart, & diff --git a/phys/physics_mmm b/phys/physics_mmm index 63c806deef..550b5b4947 160000 --- a/phys/physics_mmm +++ b/phys/physics_mmm @@ -1 +1 @@ -Subproject commit 63c806deef6a5fe75f808866fcdc6a22596076ce +Subproject commit 550b5b494758ad2756394e93cb08b5709b8feb83 diff --git a/run/README.namelist b/run/README.namelist index 96207eeb61..dfe5dc727e 100644 --- a/run/README.namelist +++ b/run/README.namelist @@ -255,10 +255,10 @@ Namelist variables specifically for the WPS input for real: = 2: it uses an alternative way (less biased when compared against input data) to compute height in program real and pressure in model. - wif_input_opt = 0 ! = 1: option to process the Water Ice Friendly Aerosol input from metgrid for use with mp_physics=28,29 + wif_input_opt = 0 ! = 1: option to process the Water Ice Friendly Aerosol input from metgrid for use with mp_physics=28,29, or 88 with tempo_aerosolaware=1 = 2: since V4.4, option to use black carbon aerosol category with mp_physics=28,29, as well as its radiative effect. Must include file QNWFA_QNIFA_QNBCA_SIGMA_MONTHLY.dat during WPS - num_wif_levels = 30 ! number of levels in the Thompson Water Ice Friendly aerosols (mp_physic=28,29) + num_wif_levels = 30 ! number of levels in the Thompson Water Ice Friendly aerosols (mp_physic=28,29, or 88 with tempo_aerosolaware=1) p_top_requested = 5000 ! p_top (Pa) to use in the model vert_refine_fact = 1 ! vertical refinement factor for ndown, not used for concurrent vertical grid refinement vert_refine_method (max_dom) = 0 ! vertical refinement method @@ -530,7 +530,7 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t 2 ! continental clean (default), for ntu3m 3 ! continental average, for ntu3m 4 ! continental urban, for ntu3m - + = 88, TEMPO (Thompson-Eidhammer Microphysics Parameterization for Operations) = 95, Ferrier (old Eta) microphysics = 96, Madwrf = 97, Goddard GCE scheme (also uses gsfcgce_hail, gsfcgce_2ice) @@ -775,6 +775,7 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t (replacing bl_mynn_tkebudget in prior versions) bl_mynn_closure = 2.5: level 2.5 2.6: level 2.6 (default) + 2.7: level 2.7 (TTE) 3.0: level 3.0 bl_mynn_tkeadvect (max_dom) = .false., ! default off; = .true. do MYNN tke advection icloud_bl option to couple the subgrid-scale clouds from the PBL scheme (MYNN only) @@ -793,7 +794,7 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t 0: original (Sommeria and Deardorf 1977); 1: Kuwano et al 2010, similar to option 0, but uses resolved scale gradients as opposed to higher order moments ; - 2: from Chaboureau and Bechtold (2002, JAS, with mods, default) + 2: heavily modified Chaboureau and Bechtold (2002, JAS, default) bl_mynn_edmf (max_dom) = 1 ! activate mass-flux option in MYNN, 0: mass-flux option off Related (hidden) options: bl_mynn_edmf_mom (max_dom) = 1 ! - activates momentum transport in MYNN mass-flux scheme @@ -802,14 +803,28 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t bl_mynn_edmf_tke (max_dom) = 0, ! default; = 1 activates TKE transport in MYNN mass-flux scheme (assuming bl_mynn_edmf > 0) 0 no TKE transport (default);1: activate TKE transport - bl_mynn_output (max_dom) = 0, ! do not output extra arrays - 1 allocate and output extra 3D arrays from MYNN PBL + bl_mynn_output (max_dom) = 0, ! = 0: do not output extra arrays (default) + 1 ! = 1: allocate and output extra 3D arrays from MYNN-EDMF mass-flux component + bl_mynn_mixaerosols (max_dom) = 1, ! = 0: do not mix qnwfa, qnifa, or qnbca + ! = 1: mix aerosols (default) + bl_mynn_mixnumcon (max_dom) = 0, ! = 0: do not mix number concentrations (qnc, qni) (default) + ! = 1: mix number concentrations + bl_mynn_mixscalars (max_dom) = 0 ! = 0: off, 1: activate mixing of scalars (must set up 4D scalar array manually in WRF) + bl_mynn_ess (max_dom) ! = 0: use Chaboureau and Bechtold (2002) buoyancy-flux relationships to estimate effect static stability + ! = 1: use modified OGorman (2011, JAS) effective static stability + bl_mynn_mixqt (max_dom) = 0 ! mixing moisture species separately, 1: mixing total water (do not use) + bl_mynn_edmf_dd (max_dom) = 0, ! default = 0; 1 activates downdraft component in MYNN + bl_mynn_output (max_dom) = 0, ! default = 0; do not output extra arrays + ! =1:allocate and output extra 3D arrays from MYNN PBL for mass flux updrafts + ! =2:allocate and output extra 3D arrays from MYNN PBL for mass flux downdrafts bl_mynn_mixscalars (max_dom) = 0 ! off, 1: activate mixing of scalars bl_mynn_mixqt (max_dom) = 0 ! mixing moisture species separately, 1: mixing total water scalar_pblmix (max_dom) = 1 ! mix scalar fields consistent with PBL option (exch_h) tracer_pblmix (max_dom) = 1 ! mix tracer fields consistent with PBL option (exch_h) - shinhong_tke_diag (max_dom) = 0 ! diagnostic TKE and mixing length from Shin-Hong PBL + shinhong_nonlocal_flux (max_dom) = .true., ! default; =.false. to use YSU (gamma-type) option + shinhong_scu_mixing (max_dom) = .false., ! strato-cumulus mixing (topdown mxing in YSU pbl) + shinhong_ke_dissipation (max_dom) = .false., ! heating due to kinetic energy dissipation opt_thcnd option to treat thermal conductivity in Noah LSM = 1, original (default) = 2, McCumber and Pielke for silt loam and sandy loam @@ -826,7 +841,7 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t = 0, no cumulus = 1, Kain-Fritsch (new Eta) scheme = 2, Betts-Miller-Janjic scheme - = 3, Grell-Freitas ensemble scheme + = 3, Grell-Freitas-Li scheme = 4, Scale-aware GFS Simplified Arakawa-Schubert (SAS) scheme = 5, Grell 3D ensemble scheme = 6, Modifed Tiedtke scheme @@ -853,7 +868,7 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t This scheme only works with MYJ PBL scheme, and should not be combined with icloud=3 This scheme can also work with MYNN PBL scheme, but one should turn off EDMF (bl_mynn_edmf=0) - ishallow = 0, ! = 1 turns on shallow convection, used with Grell 3D ensemble schemes (cu_physics = 3 or 5) + ishallow = 0, ! = 1 turns on shallow convection, used with GFL and Grell 3D ensemble schemes (cu_physics = 3 or 5) clos_choice = 0, ! closure choice (place holder only) cu_diag (max_dom) = 0, ! additional t-averaged stuff for cu physics (cu_phys = 3, 5, 10, and 93 only) kf_edrates (max_dom) = 0, ! Add entrainment/detrainment rates and convective timescale output variables for KF-based @@ -1051,8 +1066,12 @@ use_rap_aero_icbc = .false. ! Set to .true. to ingest real-t seaice_thickness_default = 3.0 ! default value of seaice thickness for seaice_thickness_opt=0 tice2tsk_if2cold = .false. ! set Tice to Tsk to avoid unrealistically low sea ice temperatures iz0tlnd = 0 ! thermal roughness length for sfclay (0 = old, 1 = veg dependent Chen-Zhang Czil, - 2 = Zilitinkevitch (czil=0.1)) - for mynn sfc (0=Zilitinkevitch (def),1=Chen-Zhang,2=mod Yang,3=const zt) + ! 2 = Zilitinkevitch (czil=0.1)). + ! Note for mynn sfc, this iz0tlnd option has been replaced by sf_mynn_sfcflux_water and sf_mynn_sfcflux_land + sf_mynn_sfcflux_water = 1 ! MYNN-SFC option for the calculation of z0, zt and zq over water + ! 0:COARE3.0, 1:COARE3.5 (default), 2:Davis/COARE3.5 ,3:Davis/Garratt, 4:Taylor-Yelland + sf_mynn_sfcflux_land = 0 ! MYNN-SFC option for the calculation of zt and zq over land + ! 0:constant Czil (default), 1:variable Czil, 2:Yang (2008) mp_tend_lim = 10., ! limit on temp tendency from mp latent heating from radar data assimilation prec_acc_dt (max_dom) = 0., ! number of minutes in precipitation bucket - will add three new 2d output fields: prec_acc_c, prec_acc_nc and snow_acc_nc @@ -1720,8 +1739,12 @@ The following are for observation nudging: swap_pole_with_next_j = .false. ! T/F replace the entire j=1 (jds-1) with the values from j=2 (jds-2) actual_distance_average = .false. ! T/F average the field at each i location in the j-loop with a number of grid points based on a map-factor ratio gwd_opt (max_dom) = 0 ! for running without gravity wave drag - = 1 ; for running with gravity wave drag (Choi and Hong 2015) + = 1 ; for running with KIM gravity wave drag = 3 ; NOAA/GSL gravity wave drag and turbulent orographic form drag + kim_tofd = .false. ! .true., KIM turbulent orographic form drag for sfclay_physics=1 (TOFD) + tofd_factor = .003 ! a factor to compute TOFD in surface layer physics + gwd_if_nonhyd = .true. ! whether to incorporate the nonhydrostatic effect in KIM GWD + gwd_dx_factor = 2. ! effective grid spacing ratio in KIM GWD gwd_diags = 0 ; set to 1 to output diagnostics for gwd_opt = 3 sfs_opt (max_dom) = 0 ! nonlinear backscatter and anisotropy (NBA) off = 1 ! NBA1 using diagnostic stress terms (km_opt=2,3 for scalars) diff --git a/run/URBPARM.TBL b/run/URBPARM.TBL index 74894e0791..a4bb660c53 100644 --- a/run/URBPARM.TBL +++ b/run/URBPARM.TBL @@ -232,7 +232,7 @@ DZGR: 0.05 0.10 0.15 0.20 # (sf_urban_physics=1) # -FVG: 0.1 +FVG: 0.2,0.2,0.2 # # TREEOPTION [ 0: No urban trees - default, 1: Enable urban trees ] @@ -241,6 +241,13 @@ FVG: 0.1 TREEOPTION: 0 +# +# TTscheme [ 0: diagnostic approach - default, 1: prognostic approach ] +# (sf_urban_physics=1) +# + +TTscheme: 0 + # # RTREE: Tree crown radius [ m ] # HTREE: Tree crown center height [ m ] @@ -252,6 +259,8 @@ TREEOPTION: 0 # CAPVG: Heat capacity of vegetated ground [ J m{-3} K{-1} ] # EPSVG: Emissivity of vegetated ground [ - ] # EPST: Emissivity of tree canopy [ - ] +# ALBVG: Albedo of vegetated ground [ - ] +# ALBT: Albedo of tree canopy [ - ] # (sf_urban_physics=1, TREEOPTION=1) # @@ -265,6 +274,8 @@ Z0HVG: 0.001, 0.001, 0.001 CAPVG: 1.3E6, 1.3E6, 1.3E6 EPSVG: 0.93, 0.93, 0.93 EPST: 0.9, 0.9, 0.9 +ALBVG: 0.2, 0.2, 0.2 +ALBT: 0.15, 0.15, 0.15 # # FRC_URB: Fraction of the urban landscape which does not have natural diff --git a/run/URBPARM_LCZ.TBL b/run/URBPARM_LCZ.TBL index 629f97bcbb..b738aa3f83 100644 --- a/run/URBPARM_LCZ.TBL +++ b/run/URBPARM_LCZ.TBL @@ -231,7 +231,7 @@ DZGR: 0.05 0.10 0.15 0.20 # (sf_urban_physics=1) # -FVG: 0.1 +FVG: 0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2 # # TREEOPTION [ 0: No urban trees - default, 1: Enable urban trees ] @@ -240,6 +240,13 @@ FVG: 0.1 TREEOPTION: 0 +# +# TTscheme [ 0: diagnostic approach - default, 1: prognostic approach ] +# (sf_urban_physics=1) +# + +TTscheme: 0 + # # RTREE: Tree crown radius [ m ] # HTREE: Tree crown center height [ m ] @@ -251,6 +258,8 @@ TREEOPTION: 0 # CAPVG: Heat capacity of vegetated ground [ J m{-3} K{-1} ] # EPSVG: Emissivity of vegetated ground [ - ] # EPST: Emissivity of tree canopy [ - ] +# ALBVG: Albedo of vegetated ground [ - ] +# ALBT: Albedo of tree canopy [ - ] # (sf_urban_physics=1, TREEOPTION=1) # @@ -266,6 +275,8 @@ Z0HVG: 0.001, 0.001, 0.001,0.001, 0.001, 0.001,0.001, 0.001, 0.001,0.001, 0.001 CAPVG: 1.3E6, 1.3E6, 1.3E6,1.3E6, 1.3E6, 1.3E6,1.3E6, 1.3E6, 1.3E6,1.3E6, 1.3E6 EPSVG: 0.93, 0.93, 0.93,0.93, 0.93, 0.93,0.93, 0.93, 0.93,0.93, 0.93 EPST: 0.9, 0.9, 0.9,0.9, 0.9, 0.9,0.9, 0.9, 0.9,0.9, 0.9 +ALBVG: 0.2, 0.2, 0.2,0.2, 0.2, 0.2,0.2, 0.2, 0.2,0.2, 0.2 +ALBT: 0.15, 0.15, 0.15,0.15, 0.15, 0.15,0.15, 0.15, 0.15,0.15, 0.15 # # FRC_URB: Fraction of the urban landscape which does not have natural diff --git a/run/URBPARM_UZE.TBL b/run/URBPARM_UZE.TBL index 27534fb34f..38688ad1df 100644 --- a/run/URBPARM_UZE.TBL +++ b/run/URBPARM_UZE.TBL @@ -62,7 +62,7 @@ FRC_URB: 0.5, 0.6, 0.75 # (sf_urban_physics=1) # -FVG: 0.1 +FVG: 0.2,0.2,0.2 # # TREEOPTION [ 0: No urban trees - default, 1: Enable urban trees ] @@ -71,6 +71,13 @@ FVG: 0.1 TREEOPTION: 0 +# +# TTscheme [ 0: diagnostic approach - default, 1: prognostic approach ] +# (sf_urban_physics=1) +# + +TTscheme: 0 + # # RTREE: Tree crown radius [ m ] # HTREE: Tree crown center height [ m ] @@ -82,6 +89,8 @@ TREEOPTION: 0 # CAPVG: Heat capacity of vegetated ground [ J m{-3} K{-1} ] # EPSVG: Emissivity of vegetated ground [ - ] # EPST: Emissivity of tree canopy [ - ] +# ALBVG: Albedo of vegetated ground [ - ] +# ALBT: Albedo of tree canopy [ - ] # (sf_urban_physics=1, TREEOPTION=1) # @@ -95,6 +104,8 @@ Z0HVG: 0.001, 0.001, 0.001 CAPVG: 1.3E6, 1.3E6, 1.3E6 EPSVG: 0.93, 0.93, 0.93 EPST: 0.9, 0.9, 0.9 +ALBVG: 0.2, 0.2, 0.2 +ALBT: 0.15, 0.15, 0.15 # # CAPR: Heat capacity of roof [ J m{-3} K{-1} ] diff --git a/share/module_check_a_mundo.F b/share/module_check_a_mundo.F index c32a9fc5a2..c35d74d885 100644 --- a/share/module_check_a_mundo.F +++ b/share/module_check_a_mundo.F @@ -2774,6 +2774,24 @@ END FUNCTION bep_bem_ngr_u count_fatal_error = count_fatal_error + 1 END IF +!----------------------------------------------------------------------- +! CFBM checks +!----------------------------------------------------------------------- + + ! CFBM can be active in only one domain + oops = 0 + DO i = 1, model_config_rec % max_dom + IF (model_config_rec%ifire(i) == 1) THEN + oops = oops + 1 + END IF + END DO + + IF ( oops .GT. 1 ) THEN + wrf_err_message = '--- ERROR: CFBM can be active in only one domain' + CALL wrf_debug ( 0, TRIM( wrf_err_message ) ) + count_fatal_error = count_fatal_error + 1 + END IF + !----------------------------------------------------------------------- ! This WRF version does not support trajectories on a global domain !-----------------------------------------------------------------------