Added alpha_fd coefficient for Beljaars TOFD in GSL drag

physics/GWD/drag_suite.F90

@@ -219,8 +219,8 @@ subroutine drag_suite_run(                                           &
      &           dusfc_ms,dvsfc_ms,dusfc_bl,dvsfc_bl,                   &
      &           dusfc_ss,dvsfc_ss,dusfc_fd,dvsfc_fd,                   &
      &           slmsk,br1,hpbl,                                        &
-     &           g, cp, rd, rv, fv, pi, imx, cdmbgwd, me, master,       &
-     &           lprnt, ipr, rdxzb, dx, gwd_opt,                        &
+     &           g, cp, rd, rv, fv, pi, imx, cdmbgwd, alpha_fd,         &    
+     &           me, master, lprnt, ipr, rdxzb, dx, gwd_opt,            &  
      &           do_gsl_drag_ls_bl, do_gsl_drag_ss, do_gsl_drag_tofd,   &
      &           dtend, dtidx, index_of_process_orographic_gwd,         &
      &           index_of_temperature, index_of_x_wind,                 &
@@ -327,7 +327,8 @@ subroutine drag_suite_run(                                           &
    integer, intent(in) :: gwd_opt
    logical, intent(in) :: lprnt
    integer, intent(in) :: KPBL(:)
-   real(kind=kind_phys), intent(in) :: deltim, G, CP, RD, RV, cdmbgwd(:)
+   real(kind=kind_phys), intent(in) :: deltim, G, CP, RD, RV,    &
+                                       cdmbgwd(:), alpha_fd
    real(kind=kind_phys), intent(inout) :: dtend(:,:,:)
    logical, intent(in) :: ldiag3d
    integer, intent(in) :: dtidx(:,:), index_of_temperature,      &
@@ -443,6 +444,7 @@ subroutine drag_suite_run(                                           &
    real(kind=kind_phys), dimension(im,km) :: utendform,vtendform
    real(kind=kind_phys)                 :: a1,a2,wsp
    real(kind=kind_phys)                 :: H_efold
+   real(kind=kind_phys), parameter      :: coeff_fd = 6.325e-3_kind_phys
 
 ! critical richardson number for wave breaking : ! larger drag with larger value
    real(kind=kind_phys), parameter       ::  ric     = 0.25
@@ -1362,8 +1364,10 @@ subroutine drag_suite_run(                                           &
             H_efold = 1500.
             DO k=kts,km
                wsp=SQRT(uwnd1(i,k)**2 + vwnd1(i,k)**2)
-               ! alpha*beta*Cmd*Ccorr*2.109 = 12.*1.*0.005*0.6*2.109 = 0.0759
-               var_temp = 0.0759*EXP(-(zl(i,k)/H_efold)**1.5)*a2*          &
+               ! Note:  In Beljaars et al. (2004):
+               ! alpha_fd*beta*Cmd*Ccorr*2.109 = 12.*1.*0.005*0.6*2.109 = 0.0759
+               ! lump beta*Cmd*Ccorr*2.109 into 1.*0.005*0.6*2.109 = coeff_fd ~ 6.325e-3_kind_phys
+               var_temp = alpha_fd*coeff_fd*EXP(-(zl(i,k)/H_efold)**1.5)*a2*       &
                                  zl(i,k)**(-1.2)*ss_taper(i) ! this is greater than zero
                !  Note:  This is a semi-implicit treatment of the time differencing
                !  per Beljaars et al. (2004, QJRMS)

physics/GWD/drag_suite.meta

@@ -501,6 +501,14 @@
   type = real
   kind = kind_phys
   intent = in
+[alpha_fd]
+  standard_name = alpha_coefficient_for_turbulent_orographic_form_drag
+  long_name = alpha coefficient for Beljaars et al turbulent orographic form drag
+  units = none
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = in
 [me]
   standard_name = mpi_rank
   long_name = rank of the current MPI task

physics/GWD/ugwpv1_gsldrag.F90

@@ -307,7 +307,7 @@ subroutine ugwpv1_gsldrag_run(me, master, im, levs, ak, bk, ntrac, lonr, dtp,
           fhzero, kdt, ldiag3d, lssav, flag_for_gwd_generic_tend, do_gsl_drag_ls_bl,    &
           do_gsl_drag_ss, do_gsl_drag_tofd, do_ugwp_v1, do_ugwp_v1_orog_only,           &
           do_ugwp_v1_w_gsldrag, gwd_opt, do_tofd, ldiag_ugwp, ugwp_seq_update,          &
-          cdmbgwd, jdat, nmtvr, hprime, oc, theta, sigma, gamma,                        &
+          cdmbgwd, alpha_fd, jdat, nmtvr, hprime, oc, theta, sigma, gamma,              &
           elvmax, clx, oa4, varss,oc1ss,oa4ss,ol4ss, dx,  xlat, xlat_d, sinlat, coslat, &
           area, rain, br1, hpbl, kpbl, slmsk,                                           &
           ugrs, vgrs, tgrs, q1, prsi, prsl, prslk, phii, phil,  del, tau_amf,           &
@@ -371,7 +371,7 @@ subroutine ugwpv1_gsldrag_run(me, master, im, levs, ak, bk, ntrac, lonr, dtp,
 ! SSO parameters and variables
     integer,                 intent(in) :: gwd_opt                         !gwd_opt  and nmtvr are "redundant" controls
     integer,                 intent(in) :: nmtvr
-    real(kind=kind_phys),    intent(in) :: cdmbgwd(:)                      ! for gsl_drag
+    real(kind=kind_phys),    intent(in) :: cdmbgwd(:), alpha_fd            ! for gsl_drag
 
     real(kind=kind_phys),    intent(in), dimension(:)       :: hprime, oc, theta, sigma, gamma
 

physics/GWD/ugwpv1_gsldrag.meta

@@ -459,6 +459,14 @@
   type = real
   kind = kind_phys
   intent = in
+[alpha_fd]
+  standard_name = alpha_coefficient_for_turbulent_orographic_form_drag
+  long_name = alpha coefficient for Beljaars et al turbulent orographic form drag
+  units = none
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = in
 [jdat]
   standard_name = date_and_time_of_forecast_in_united_states_order
   long_name = current forecast date and time

physics/GWD/unified_ugwp.F90

@@ -250,7 +250,7 @@ subroutine unified_ugwp_run(me, master, im, levs, ak,bk, ntrac, dtp, fhzero, kdt
          dvsfc_ss,dusfc_fd,dvsfc_fd,dtaux2d_ms,dtauy2d_ms,dtaux2d_bl,dtauy2d_bl,       &
          dtaux2d_ss,dtauy2d_ss,dtaux2d_fd,dtauy2d_fd,dudt_ngw,dvdt_ngw,dtdt_ngw,       &
          br1,hpbl,slmsk, do_tofd, ldiag_ugwp, ugwp_seq_update,                         &
-         cdmbgwd, jdat, xlat, xlat_d, sinlat, coslat, area,                            &
+         cdmbgwd, alpha_fd, jdat, xlat, xlat_d, sinlat, coslat, area,                  &
          ugrs, vgrs, tgrs, q1, prsi, prsl, prslk, phii, phil,                          &
          del, kpbl, dusfcg, dvsfcg, gw_dudt, gw_dvdt, gw_dtdt, gw_kdis,                &
          tau_tofd, tau_mtb, tau_ogw, tau_ngw, zmtb, zlwb, zogw,                        &
@@ -287,7 +287,7 @@ subroutine unified_ugwp_run(me, master, im, levs, ak,bk, ntrac, dtp, fhzero, kdt
     real(kind=kind_phys),    intent(in),    dimension(:,:)  :: del, ugrs, vgrs, tgrs, prsl, prslk, phil
     real(kind=kind_phys),    intent(in),    dimension(:,:)  :: prsi, phii
     real(kind=kind_phys),    intent(in),    dimension(:,:)  :: q1
-    real(kind=kind_phys),    intent(in) :: dtp, fhzero, cdmbgwd(:)
+    real(kind=kind_phys),    intent(in) :: dtp, fhzero, cdmbgwd(:), alpha_fd
     integer, intent(in) :: jdat(:)
     logical, intent(in) :: do_tofd, ldiag_ugwp, ugwp_seq_update
 
@@ -497,8 +497,8 @@ subroutine unified_ugwp_run(me, master, im, levs, ak,bk, ntrac, dtp, fhzero, kdt
                  dvsfcg,dusfc_ms,dvsfc_ms,dusfc_bl,dvsfc_bl,         &
                  dusfc_ss,dvsfc_ss,dusfc_fd,dvsfc_fd,                &
                  slmsk,br1,hpbl,con_g,con_cp,con_rd,con_rv,          &
-                 con_fvirt,con_pi,lonr,                              &
-                 cdmbgwd,me,master,lprnt,ipr,rdxzb,dx,gwd_opt,       &
+                 con_fvirt,con_pi,lonr,cdmbgwd,alpha_fd,             &
+                 me,master,lprnt,ipr,rdxzb,dx,gwd_opt,               &
                  do_gsl_drag_ls_bl,do_gsl_drag_ss,do_gsl_drag_tofd,  &
                  dtend, dtidx, index_of_process_orographic_gwd,      &
                  index_of_temperature, index_of_x_wind,              &

physics/GWD/unified_ugwp.meta

@@ -685,6 +685,14 @@
   type = real
   kind = kind_phys
   intent = in
+[alpha_fd]
+  standard_name = alpha_coefficient_for_turbulent_orographic_form_drag
+  long_name = alpha coefficient for Beljaars et al turbulent orographic form drag
+  units = none
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = in
 [jdat]
   standard_name = date_and_time_of_forecast_in_united_states_order
   long_name = current forecast date and time