Skip to content

Commit

Permalink
feat: pass additional surface variables into the atmospheric model (#93)
Browse files Browse the repository at this point in the history
  • Loading branch information
zhihong-tan authored Nov 17, 2023
1 parent 7761886 commit 740edc2
Show file tree
Hide file tree
Showing 3 changed files with 82 additions and 4 deletions.
32 changes: 32 additions & 0 deletions full/atm_land_ice_flux_exchange.F90
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,7 @@ module atm_land_ice_flux_exchange_mod
id_u_star, id_b_star, id_q_star, id_u_flux, id_v_flux, &
id_t_surf, id_t_flux, id_r_flux, id_q_flux, id_slp, &
id_t_atm, id_u_atm, id_v_atm, id_wind, &
id_thv_atm, id_thv_surf, &
id_t_ref, id_rh_ref, id_u_ref, id_v_ref, id_wind_ref, &
id_del_h, id_del_m, id_del_q, id_rough_scale, &
id_t_ca, id_q_surf, id_q_atm, id_z_atm, id_p_atm, id_gust, &
Expand Down Expand Up @@ -571,6 +572,9 @@ subroutine atm_land_ice_flux_exchange_init(Time, Atm, Land, Ice, atmos_ice_bound
allocate( land_ice_atmos_boundary%shflx(is:ie,js:je) )!miz
allocate( land_ice_atmos_boundary%lhflx(is:ie,js:je) )!miz
#endif
allocate( land_ice_atmos_boundary%wind(is:ie,js:je) )
allocate( land_ice_atmos_boundary%thv_atm(is:ie,js:je) )
allocate( land_ice_atmos_boundary%thv_surf(is:ie,js:je) )
allocate( land_ice_atmos_boundary%rough_mom(is:ie,js:je) )
allocate( land_ice_atmos_boundary%frac_open_sea(is:ie,js:je) )
! initialize boundary values for override experiments (mjh)
Expand Down Expand Up @@ -598,6 +602,9 @@ subroutine atm_land_ice_flux_exchange_init(Time, Atm, Land, Ice, atmos_ice_bound
land_ice_atmos_boundary%shflx=0.0
land_ice_atmos_boundary%lhflx=0.0
#endif
land_ice_atmos_boundary%wind=0.0
land_ice_atmos_boundary%thv_atm=0.0
land_ice_atmos_boundary%thv_surf=0.0
land_ice_atmos_boundary%rough_mom=0.01
land_ice_atmos_boundary%frac_open_sea=0.0

Expand Down Expand Up @@ -679,6 +686,7 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
ex_rough_mom, ex_rough_heat, ex_rough_moist, &
ex_rough_scale,&
ex_q_star, &
ex_thv_atm, ex_thv_surf, &
ex_cd_q, &
ex_ref, ex_ref_u, ex_ref_v, ex_u10, &
ex_ref2, &
Expand Down Expand Up @@ -1129,6 +1137,7 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
!$OMP ex_gust,ex_flux_t,ex_flux_tr,ex_flux_lw, &
!$OMP ex_flux_u,ex_flux_v,ex_cd_m,ex_cd_t,ex_cd_q, &
!$OMP ex_wind,ex_u_star,ex_b_star,ex_q_star, &
!$OMP ex_thv_atm,ex_thv_surf, &
!$OMP ex_dhdt_surf,ex_dedt_surf,ex_dfdtr_surf, &
!$OMP ex_drdt_surf,ex_dhdt_atm,ex_dfdtr_atm, &
!$OMP ex_dtaudu_atm, ex_dtaudv_atm,dt,ex_land, &
Expand All @@ -1146,6 +1155,7 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
ex_flux_t(is:ie), ex_flux_tr(is:ie,isphum), ex_flux_lw(is:ie), ex_flux_u(is:ie), ex_flux_v(is:ie), &
ex_cd_m(is:ie), ex_cd_t(is:ie), ex_cd_q(is:ie), &
ex_wind(is:ie), ex_u_star(is:ie), ex_b_star(is:ie), ex_q_star(is:ie), &
ex_thv_atm(is:ie), ex_thv_surf(is:ie), &
ex_dhdt_surf(is:ie), ex_dedt_surf(is:ie), ex_dfdtr_surf(is:ie,isphum), ex_drdt_surf(is:ie), &
ex_dhdt_atm(is:ie), ex_dfdtr_atm(is:ie,isphum), ex_dtaudu_atm(is:ie), ex_dtaudv_atm(is:ie), &
dt, &
Expand All @@ -1165,6 +1175,7 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
ex_flux_t, ex_flux_tr(:,isphum), ex_flux_lw, ex_flux_u, ex_flux_v, &
ex_cd_m, ex_cd_t, ex_cd_q, &
ex_wind, ex_u_star, ex_b_star, ex_q_star, &
ex_thv_atm, ex_thv_surf, &
ex_dhdt_surf, ex_dedt_surf, ex_dfdtr_surf(:,isphum), ex_drdt_surf, &
ex_dhdt_atm, ex_dfdtr_atm(:,isphum), ex_dtaudu_atm, ex_dtaudv_atm, &
dt, &
Expand All @@ -1174,6 +1185,7 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
endif
#endif


! call mpp_clock_end(fluxClock)
zrefm = 10.0
zrefh = z_ref_heat
Expand Down Expand Up @@ -1395,6 +1407,14 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%u_ref, 'ATM', ex_ref_u , xmap_sfc, complete=.false.) !bqx
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%v_ref, 'ATM', ex_ref_v , xmap_sfc, complete=.true.) !bqx

#ifndef use_AM3_physics
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%shflx, 'ATM', ex_flux_t , xmap_sfc)
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%lhflx, 'ATM', ex_flux_tr(:,isphum), xmap_sfc)
#endif
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%wind, 'ATM', ex_wind , xmap_sfc)
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%thv_atm, 'ATM', ex_thv_atm , xmap_sfc)
call fms_xgrid_get_from_xgrid (Land_Ice_Atmos_Boundary%thv_surf, 'ATM', ex_thv_surf , xmap_sfc)

#ifdef use_AM3_physics
if (do_forecast) then
call fms_xgrid_get_from_xgrid (Ice%t_surf, 'OCN', ex_t_surf, xmap_sfc)
Expand Down Expand Up @@ -1581,6 +1601,10 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Land_Ice_Atmos_Boundar
!------- moisture scale -----------
used = fms_diag_send_data ( id_q_star, Land_Ice_Atmos_Boundary%q_star, Time )

!------- surf and atm virtual potential temperature -----------
used = fms_diag_send_data ( id_thv_atm, Land_Ice_Atmos_Boundary%thv_atm, Time )
used = fms_diag_send_data ( id_thv_surf, Land_Ice_Atmos_Boundary%thv_surf, Time )

!-----------------------------------------------------------------------
!------ diagnostics for fields at bottom atmospheric level ------

Expand Down Expand Up @@ -3433,6 +3457,14 @@ subroutine diag_field_init ( Time, atmos_axes, land_axes, land_pe )
fms_diag_register_diag_field ( mod_name, 'q_star', atmos_axes, Time, &
'moisture scale', 'kg water/kg air' )

id_thv_atm = &
fms_diag_register_diag_field ( mod_name, 'thv_atm', atmos_axes, Time, &
'surface air virtual potential temperature', 'K')

id_thv_surf = &
fms_diag_register_diag_field ( mod_name, 'thv_surf', atmos_axes, Time, &
'surface virtual potential temperature', 'K')

id_u_flux = &
fms_diag_register_diag_field ( mod_name, 'tau_x', atmos_axes, Time, &
'zonal wind stress', 'pa' )
Expand Down
18 changes: 17 additions & 1 deletion shared/surface_flux.F90
Original file line number Diff line number Diff line change
Expand Up @@ -210,6 +210,7 @@ subroutine surface_flux_1d ( &
flux_t, flux_q, flux_r, flux_u, flux_v, &
cd_m, cd_t, cd_q, &
w_atm, u_star, b_star, q_star, &
thv_atm, thv_surf, &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm, dtaudv_atm, &
dt, land, seawater, avail )
Expand Down Expand Up @@ -250,6 +251,8 @@ subroutine surface_flux_1d ( &
u_star, & !< Turbulent velocity scale
b_star, & !< Turbulent buoyant scale
q_star, & !< Turbulent moisture scale
thv_atm, & ! Surface air theta_v
thv_surf, & ! Surface theta_v
cd_m, & !< Momentum exchange coefficient
cd_t, & ! Heat exchange coefficient
cd_q !< Moisture exchange coefficient
Expand All @@ -262,7 +265,7 @@ subroutine surface_flux_1d ( &

! ---- local vars ----------------------------------------------------------
real, dimension(size(t_atm(:))) :: &
thv_atm, th_atm, tv_atm, thv_surf, &
th_atm, tv_atm, & ! thv_atm and thv_surf are moved to output
e_sat, e_sat1, q_sat, q_sat1, p_ratio, &
t_surf0, t_surf1, u_dif, v_dif, &
rho_drag, drag_t, drag_m, drag_q, rho, &
Expand Down Expand Up @@ -448,6 +451,8 @@ subroutine surface_flux_1d ( &
q_star = 0.0
q_surf = 0.0
w_atm = 0.0
thv_atm = 0.0
thv_surf = 0.0
endwhere

! calculate d(stress component)/d(atmos wind component)
Expand Down Expand Up @@ -477,6 +482,7 @@ subroutine surface_flux_0d ( &
flux_t_0, flux_q_0, flux_r_0, flux_u_0, flux_v_0, &
cd_m_0, cd_t_0, cd_q_0, &
w_atm_0, u_star_0, b_star_0, q_star_0, &
thv_atm_0, thv_surf_0, &
dhdt_surf_0, dedt_surf_0, dedq_surf_0, drdt_surf_0, &
dhdt_atm_0, dedq_atm_0, dtaudu_atm_0, dtaudv_atm_0, &
dt, land_0, seawater_0, avail_0 )
Expand Down Expand Up @@ -518,6 +524,8 @@ subroutine surface_flux_0d ( &
u_star_0, & !< Turbulent velocity scale
b_star_0, & !< Turbulent buoyant scale
q_star_0, & !< Turbulent moisture scale
thv_atm_0, & ! Surface air theta_v
thv_surf_0, & ! Surface theta_v
cd_m_0, & !< Momentum exchange coefficient
cd_t_0, & ! Heat exchange coefficient
cd_q_0 !< Moisture exchange coefficient
Expand All @@ -536,6 +544,7 @@ subroutine surface_flux_0d ( &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm,dtaudv_atm, &
w_atm, u_star, b_star, q_star, &
thv_atm, thv_surf, &
cd_m, cd_t, cd_q
real, dimension(1) :: q_surf

Expand Down Expand Up @@ -571,6 +580,7 @@ subroutine surface_flux_0d ( &
flux_t, flux_q, flux_r, flux_u, flux_v, &
cd_m, cd_t, cd_q, &
w_atm, u_star, b_star, q_star, &
thv_atm, thv_surf, &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm, dtaudv_atm, &
dt, land, seawater, avail )
Expand All @@ -593,6 +603,8 @@ subroutine surface_flux_0d ( &
b_star_0 = b_star(1)
q_star_0 = q_star(1)
q_surf_0 = q_surf(1)
thv_atm_0 = thv_atm(1)
thv_surf_0 = thv_surf(1)
cd_m_0 = cd_m(1)
cd_t_0 = cd_t(1)
cd_q_0 = cd_q(1)
Expand All @@ -607,6 +619,7 @@ subroutine surface_flux_2d ( &
flux_t, flux_q, flux_r, flux_u, flux_v, &
cd_m, cd_t, cd_q, &
w_atm, u_star, b_star, q_star, &
thv_atm, thv_surf, &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm, dtaudv_atm, &
dt, land, seawater, avail )
Expand Down Expand Up @@ -648,6 +661,8 @@ subroutine surface_flux_2d ( &
u_star, & !< Turbulent velocity scale
b_star, & !< Turbulent buoyant scale
q_star, & !< Turbulent moisture scale
thv_atm, & ! Surface air theta_v
thv_surf, & ! Surface theta_v
cd_m, & !< Momentum exchange coefficient
cd_t, & ! Heat exchange coefficient
cd_q !< Moisture exchange coefficient
Expand All @@ -666,6 +681,7 @@ subroutine surface_flux_2d ( &
flux_t(:,j), flux_q(:,j), flux_r(:,j), flux_u(:,j), flux_v(:,j), &
cd_m(:,j), cd_t(:,j), cd_q(:,j), &
w_atm(:,j), u_star(:,j), b_star(:,j), q_star(:,j), &
thv_atm(:,j), thv_surf(:,j), &
dhdt_surf(:,j), dedt_surf(:,j), dedq_surf(:,j), drdt_surf(:,j), &
dhdt_atm(:,j), dedq_atm(:,j), dtaudu_atm(:,j), dtaudv_atm(:,j), &
dt, land(:,j), seawater(:,j), avail(:,j) )
Expand Down
36 changes: 33 additions & 3 deletions simple/flux_exchange.F90
Original file line number Diff line number Diff line change
Expand Up @@ -59,6 +59,7 @@ module flux_exchange_mod
id_u_star, id_b_star, id_q_star, id_u_flux, id_v_flux, &
id_t_surf, id_t_flux, id_q_flux, id_r_flux, &
id_t_atm, id_u_atm, id_v_atm, id_wind, &
id_thv_atm, id_thv_surf, &
id_t_ref, id_rh_ref, id_u_ref, id_v_ref, id_q_ref, &
id_del_h, id_del_m, id_del_q, id_albedo, id_gust, &
id_t_ca, id_q_surf, id_q_atm, id_z_atm, id_p_atm, &
Expand Down Expand Up @@ -152,7 +153,8 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Boundary )
albedo, albedo_vis_dir, albedo_nir_dir, &
albedo_vis_dif, albedo_nir_dif, &
del_m, del_h, del_q, land_frac, &
ref, ref2, t_ref, qs_ref, qs_ref_cmip
ref, ref2, t_ref, qs_ref, qs_ref_cmip, &
thv_atm, thv_surf

logical, dimension(is:ie,js:je) :: mask, seawater, avail
real :: zrefm, zrefh
Expand Down Expand Up @@ -256,6 +258,7 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Boundary )
flux_t, flux_q, flux_lw, flux_u, flux_v, &
cd_m, cd_t, cd_q, wind, &
u_star, b_star, q_star, &
thv_atm, thv_surf, &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm, dtaudv_atm, &
dt, Land%mask(:,:,1), seawater, avail )
Expand All @@ -282,6 +285,14 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Boundary )
Boundary%b_star = b_star
Boundary%q_star = q_star

! additional boundary variables for use in ncep-edmf
Boundary%shflx = flux_t
Boundary%lhflx = flux_q
Boundary%wind = wind
Boundary%thv_atm = thv_atm
Boundary%thv_surf = thv_surf


!=======================================================================
!-------------------- diagnostics section ------------------------------

Expand All @@ -303,6 +314,8 @@ subroutine sfc_boundary_layer ( dt, Time, Atm, Land, Ice, Boundary )
if ( id_u_star > 0 ) used = fms_diag_send_data ( id_u_star, u_star, Time )
if ( id_b_star > 0 ) used = fms_diag_send_data ( id_b_star, b_star, Time )
if ( id_q_star > 0 ) used = fms_diag_send_data ( id_q_star, q_star, Time )
if ( id_thv_atm > 0 ) used = fms_diag_send_data ( id_thv_atm, thv_atm, Time )
if ( id_thv_surf > 0 ) used = fms_diag_send_data ( id_thv_surf, thv_surf, Time )
if ( id_t_atm > 0 ) used = fms_diag_send_data ( id_t_atm, Atm%t_bot, Time )
if ( id_u_atm > 0 ) used = fms_diag_send_data ( id_u_atm, Atm%u_bot, Time )
if ( id_v_atm > 0 ) used = fms_diag_send_data ( id_v_atm, Atm%v_bot, Time )
Expand Down Expand Up @@ -791,6 +804,9 @@ subroutine flux_exchange_init ( Time, Atm, Land, Ice, &
allocate( land_ice_atmos_boundary%shflx(is:ie,js:je) )!miz
allocate( land_ice_atmos_boundary%lhflx(is:ie,js:je) )!miz
#endif
allocate( land_ice_atmos_boundary%wind(is:ie,js:je) )
allocate( land_ice_atmos_boundary%thv_atm(is:ie,js:je) )
allocate( land_ice_atmos_boundary%thv_surf(is:ie,js:je) )
allocate( land_ice_atmos_boundary%rough_mom(is:ie,js:je) )
allocate( land_ice_atmos_boundary%frac_open_sea(is:ie,js:je) )

Expand Down Expand Up @@ -818,6 +834,9 @@ subroutine flux_exchange_init ( Time, Atm, Land, Ice, &
land_ice_atmos_boundary%shflx=0.0
land_ice_atmos_boundary%lhflx=0.0
#endif
land_ice_atmos_boundary%wind=0.0
land_ice_atmos_boundary%thv_atm=0.0
land_ice_atmos_boundary%thv_surf=0.0
land_ice_atmos_boundary%rough_mom=0.01
land_ice_atmos_boundary%frac_open_sea=0.0

Expand Down Expand Up @@ -935,10 +954,18 @@ subroutine diag_field_init ( Time, atmos_axes )
fms_diag_register_diag_field ( mod_name, 'b_star', atmos_axes, Time, &
'buoyancy scale', 'm/s2' )

id_q_star = &
fms_diag_register_diag_field ( mod_name, 'q_star', atmos_axes, Time, &
id_q_star = &
fms_diag_register_diag_field ( mod_name, 'q_star', atmos_axes, Time, &
'moisture scale', 'kg water/kg air' )

id_thv_atm = &
fms_diag_register_diag_field ( mod_name, 'thv_atm', atmos_axes, Time, &
'surface air virtual potential temperature', 'K')

id_thv_surf = &
fms_diag_register_diag_field ( mod_name, 'thv_surf', atmos_axes, Time, &
'surface virtual potential temperature', 'K')

id_u_flux = &
fms_diag_register_diag_field ( mod_name, 'tau_x', atmos_axes, Time, &
'zonal wind stress', 'pa' )
Expand Down Expand Up @@ -1156,6 +1183,7 @@ subroutine surface_flux_2d ( &
flux_t, flux_q, flux_r, flux_u, flux_v, &
cd_m, cd_t, cd_q, &
w_atm, u_star, b_star, q_star, &
thv_atm, thv_surf, &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm, dtaudv_atm, &
dt, land, seawater, avail )
Expand All @@ -1173,6 +1201,7 @@ subroutine surface_flux_2d ( &
dhdt_surf, dedt_surf, dedq_surf, drdt_surf, &
dhdt_atm, dedq_atm, dtaudu_atm,dtaudv_atm, &
w_atm, u_star, b_star, q_star, &
thv_atm, thv_surf, &
cd_m, cd_t, cd_q
real, intent(inout), dimension(:,:) :: q_surf
real, intent(in) :: dt
Expand All @@ -1189,6 +1218,7 @@ subroutine surface_flux_2d ( &
flux_t(:,j), flux_q(:,j), flux_r(:,j), flux_u(:,j), flux_v(:,j), &
cd_m(:,j), cd_t(:,j), cd_q(:,j), &
w_atm(:,j), u_star(:,j), b_star(:,j), q_star(:,j), &
thv_atm(:,j), thv_surf(:,j), &
dhdt_surf(:,j), dedt_surf(:,j), dedq_surf(:,j), drdt_surf(:,j), &
dhdt_atm(:,j), dedq_atm(:,j), dtaudu_atm(:,j), dtaudv_atm(:,j), &
dt, land(:,j), seawater(:,j), avail(:,j) )
Expand Down

0 comments on commit 740edc2

Please sign in to comment.