Skip to content

+Add fluxes%tau_mag and forces%tau_mag #384

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
marshallward merged 2 commits into from
Jun 28, 2023
Merged

+Add fluxes%tau_mag and forces%tau_mag #384

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
51de820
+Add fluxes%tau_mag and forces%tau_mag
Hallberg-NOAA Jun 21, 2023
f9764f5
Merge branch 'dev/gfdl' into add_tau_mag
marshallward Jun 28, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
44 changes: 29 additions & 15 deletions config_src/drivers/FMS_cap/MOM_surface_forcing_gfdl.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -624,13 +624,16 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
endif

! Set the wind stresses and ustar.
if (associated(fluxes%ustar) .and. associated(fluxes%ustar_gustless)) then
if (associated(fluxes%ustar) .and. associated(fluxes%ustar_gustless) .and. associated(fluxes%tau_mag)) then
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, ustar=fluxes%ustar, &
gustless_ustar=fluxes%ustar_gustless)
elseif (associated(fluxes%ustar)) then
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, ustar=fluxes%ustar)
elseif (associated(fluxes%ustar_gustless)) then
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, gustless_ustar=fluxes%ustar_gustless)
mag_tau=fluxes%tau_mag, gustless_ustar=fluxes%ustar_gustless)
else
if (associated(fluxes%ustar)) &
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, ustar=fluxes%ustar)
if (associated(fluxes%ustar_gustless)) &
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, gustless_ustar=fluxes%ustar_gustless)
if (associated(fluxes%tau_mag)) &
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, mag_tau=fluxes%tau_mag)
endif

if (coupler_type_initialized(fluxes%tr_fluxes) .and. &
Expand Down Expand Up @@ -674,7 +677,8 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS, dt_
real, dimension(SZI_(G),SZJ_(G)) :: &
rigidity_at_h, & ! Ice rigidity at tracer points [L4 Z-1 T-1 ~> m3 s-1]
net_mass_src, & ! A temporary of net mass sources [R Z T-1 ~> kg m-2 s-1].
ustar_tmp ! A temporary array of ustar values [Z T-1 ~> m s-1].
ustar_tmp, & ! A temporary array of ustar values [Z T-1 ~> m s-1].
tau_mag_tmp ! A temporary array of surface stress magnitudes [R Z L T-2 ~> Pa]

real :: I_GEarth ! The inverse of the gravitational acceleration [T2 Z L-2 ~> s2 m-1]
real :: Kv_rho_ice ! (CS%Kv_sea_ice / CS%density_sea_ice) [L4 Z-2 T-1 R-1 ~> m5 s-1 kg-1]
Expand Down Expand Up @@ -778,12 +782,13 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS, dt_
! Set the wind stresses and ustar.
if (wt1 <= 0.0) then
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux=forces%taux, tauy=forces%tauy, &
ustar=forces%ustar, tau_halo=1)
ustar=forces%ustar, mag_tau=forces%tau_mag, tau_halo=1)
else
call extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux=forces%taux, tauy=forces%tauy, &
ustar=ustar_tmp, tau_halo=1)
ustar=ustar_tmp, mag_tau=tau_mag_tmp, tau_halo=1)
do j=js,je ; do i=is,ie
forces%ustar(i,j) = wt1*forces%ustar(i,j) + wt2*ustar_tmp(i,j)
forces%tau_mag(i,j) = wt1*forces%tau_mag(i,j) + wt2*tau_mag_tmp(i,j)
enddo ; enddo
endif

Expand Down Expand Up @@ -880,7 +885,7 @@ end subroutine convert_IOB_to_forces
!! Ice_ocean_boundary_type into optional argument arrays, including changes of units, sign
!! conventions, and putting the fields into arrays with MOM-standard sized halos.
subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy, ustar, &
gustless_ustar, tau_halo)
gustless_ustar, mag_tau, tau_halo)
type(ice_ocean_boundary_type), &
target, intent(in) :: IOB !< An ice-ocean boundary type with fluxes to drive
!! the ocean in a coupled model
Expand All @@ -900,6 +905,9 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
real, dimension(SZI_(G),SZJ_(G)), &
optional, intent(out) :: gustless_ustar !< The surface friction velocity without
!! any contributions from gustiness [Z T-1 ~> m s-1].
real, dimension(SZI_(G),SZJ_(G)), &
optional, intent(inout) :: mag_tau !< The magintude of the wind stress at tracer points
!! including subgridscale variability and gustiness [R Z L T-2 ~> Pa]
integer, optional, intent(in) :: tau_halo !< The halo size of wind stresses to set, 0 by default.

! Local variables
Expand All @@ -916,7 +924,7 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
real :: tau_mag ! magnitude of the wind stress [R Z L T-2 ~> Pa]
real :: stress_conversion ! A unit conversion factor from Pa times any stress multiplier [R Z L T-2 Pa-1 ~> 1]

logical :: do_ustar, do_gustless
logical :: do_ustar, do_gustless, do_tau_mag
integer :: wind_stagger ! AGRID, BGRID_NE, or CGRID_NE (integers from MOM_domains)
integer :: i, j, is, ie, js, je, ish, ieh, jsh, jeh, Isqh, Ieqh, Jsqh, Jeqh, i0, j0, halo

Expand All @@ -929,7 +937,7 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
IRho0 = US%L_to_Z / CS%Rho0
stress_conversion = US%Pa_to_RLZ_T2 * CS%wind_stress_multiplier

do_ustar = present(ustar) ; do_gustless = present(gustless_ustar)
do_ustar = present(ustar) ; do_gustless = present(gustless_ustar) ; do_tau_mag = present(mag_tau)

wind_stagger = CS%wind_stagger
if ((IOB%wind_stagger == AGRID) .or. (IOB%wind_stagger == BGRID_NE) .or. &
Expand Down Expand Up @@ -1022,13 +1030,13 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
endif ! endif for extracting wind stress fields with various staggerings
endif

if (do_ustar .or. do_gustless) then
if (do_ustar .or. do_tau_mag .or. do_gustless) then
! Set surface friction velocity directly or as a function of staggering.
! ustar is required for the bulk mixed layer formulation and other turbulent mixing
! parametizations. The background gustiness (for example with a relatively small value
! of 0.02 Pa) is intended to give reasonable behavior in regions of very weak winds.
if (associated(IOB%stress_mag)) then
if (do_ustar) then ; do j=js,je ; do i=is,ie
if (do_ustar .or. do_tau_mag) then ; do j=js,je ; do i=is,ie
gustiness = CS%gust_const
if (CS%read_gust_2d) then
if ((wind_stagger == CGRID_NE) .or. &
Expand All @@ -1038,7 +1046,10 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
(G%mask2dBu(I,J-1) + G%mask2dBu(I-1,J))) > 0.0)) ) &
gustiness = CS%gust(i,j)
endif
ustar(i,j) = sqrt(gustiness*IRho0 + IRho0*US%Pa_to_RLZ_T2*IOB%stress_mag(i-i0,j-j0))
if (do_tau_mag) &
mag_tau(i,j) = gustiness + US%Pa_to_RLZ_T2*IOB%stress_mag(i-i0,j-j0)
if (do_ustar) &
ustar(i,j) = sqrt(gustiness*IRho0 + IRho0*US%Pa_to_RLZ_T2*IOB%stress_mag(i-i0,j-j0))
enddo ; enddo ; endif
if (CS%answer_date < 20190101) then
if (do_gustless) then ; do j=js,je ; do i=is,ie
Expand All @@ -1062,6 +1073,7 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
if (CS%read_gust_2d) gustiness = CS%gust(i,j)
endif
if (do_ustar) ustar(i,j) = sqrt(gustiness*IRho0 + IRho0 * tau_mag)
if (do_tau_mag) mag_tau(i,j) = gustiness + tau_mag
if (CS%answer_date < 20190101) then
if (do_gustless) gustless_ustar(i,j) = sqrt(US%L_to_Z*tau_mag / CS%Rho0)
else
Expand All @@ -1074,6 +1086,7 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
gustiness = CS%gust_const
if (CS%read_gust_2d .and. (G%mask2dT(i,j) > 0.0)) gustiness = CS%gust(i,j)
if (do_ustar) ustar(i,j) = sqrt(gustiness*IRho0 + IRho0 * tau_mag)
if (do_tau_mag) mag_tau(i,j) = gustiness + tau_mag
if (CS%answer_date < 20190101) then
if (do_gustless) gustless_ustar(i,j) = sqrt(US%L_to_Z*tau_mag / CS%Rho0)
else
Expand All @@ -1095,6 +1108,7 @@ subroutine extract_IOB_stresses(IOB, index_bounds, Time, G, US, CS, taux, tauy,
if (CS%read_gust_2d) gustiness = CS%gust(i,j)

if (do_ustar) ustar(i,j) = sqrt(gustiness*IRho0 + IRho0 * tau_mag)
if (do_tau_mag) mag_tau(i,j) = gustiness + tau_mag
if (CS%answer_date < 20190101) then
if (do_gustless) gustless_ustar(i,j) = sqrt(US%L_to_Z*tau_mag / CS%Rho0)
else
Expand Down
4 changes: 4 additions & 0 deletions config_src/drivers/mct_cap/mom_surface_forcing_mct.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -774,6 +774,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
((G%mask2dBu(I,J) + G%mask2dBu(I-1,J-1)) + (G%mask2dBu(I,J-1) + G%mask2dBu(I-1,J))) )
if (CS%read_gust_2d) gustiness = CS%gust(i,j)
endif
forces%tau_mag(i,j) = gustiness + tau_mag
forces%ustar(i,j) = sqrt(gustiness*Irho0 + Irho0*tau_mag)
enddo ; enddo

Expand All @@ -799,6 +800,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
do j=js,je ; do i=is,ie
gustiness = CS%gust_const
if (CS%read_gust_2d .and. (G%mask2dT(i,j) > 0.0)) gustiness = CS%gust(i,j)
forces%tau_mag(i,j) = gustiness + G%mask2dT(i,j) * sqrt(taux_at_h(i,j)**2 + tauy_at_h(i,j)**2)
forces%ustar(i,j) = sqrt(gustiness*Irho0 + Irho0 * G%mask2dT(i,j) * &
sqrt(taux_at_h(i,j)**2 + tauy_at_h(i,j)**2))
enddo ; enddo
Expand All @@ -820,8 +822,10 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
G%mask2dCv(i,J)*forces%tauy(i,J)**2) / (G%mask2dCv(i,J-1) + G%mask2dCv(i,J))

if (CS%read_gust_2d) then
forces%tau_mag(i,j) = CS%gust(i,j) + sqrt(taux2 + tauy2)
forces%ustar(i,j) = sqrt(CS%gust(i,j)*Irho0 + Irho0*sqrt(taux2 + tauy2))
else
forces%tau_mag(i,j) = CS%gust_const + sqrt(taux2 + tauy2)
forces%ustar(i,j) = sqrt(CS%gust_const*Irho0 + Irho0*sqrt(taux2 + tauy2))
endif
enddo ; enddo
Expand Down
4 changes: 4 additions & 0 deletions config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -851,6 +851,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
((G%mask2dBu(I,J) + G%mask2dBu(I-1,J-1)) + (G%mask2dBu(I,J-1) + G%mask2dBu(I-1,J))) )
if (CS%read_gust_2d) gustiness = CS%gust(i,j)
endif
forces%tau_mag(i,j) = gustiness + tau_mag
forces%ustar(i,j) = sqrt(gustiness*Irho0 + Irho0*tau_mag)
enddo ; enddo
call pass_vector(forces%taux, forces%tauy, G%Domain, halo=1)
Expand All @@ -876,6 +877,7 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
do j=js,je ; do i=is,ie
gustiness = CS%gust_const
if (CS%read_gust_2d .and. (G%mask2dT(i,j) > 0.0)) gustiness = CS%gust(i,j)
forces%tau_mag(i,j) = gustiness + G%mask2dT(i,j) * sqrt(taux_at_h(i,j)**2 + tauy_at_h(i,j)**2)
forces%ustar(i,j) = sqrt(gustiness*Irho0 + Irho0 * G%mask2dT(i,j) * &
sqrt(taux_at_h(i,j)**2 + tauy_at_h(i,j)**2))
enddo ; enddo
Expand All @@ -897,8 +899,10 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS)
G%mask2dCv(i,J)*forces%tauy(i,J)**2) / (G%mask2dCv(i,J-1) + G%mask2dCv(i,J))

if (CS%read_gust_2d) then
forces%tau_mag(i,j) = CS%gust(i,j) + sqrt(taux2 + tauy2)
forces%ustar(i,j) = sqrt(CS%gust(i,j)*Irho0 + Irho0*sqrt(taux2 + tauy2))
else
forces%tau_mag(i,j) = CS%gust(i,j) + sqrt(taux2 + tauy2)
forces%ustar(i,j) = sqrt(CS%gust_const*Irho0 + Irho0*sqrt(taux2 + tauy2))
endif
enddo ; enddo
Expand Down
61 changes: 44 additions & 17 deletions config_src/drivers/solo_driver/MOM_surface_forcing.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -407,10 +407,16 @@ subroutine wind_forcing_const(sfc_state, forces, tau_x0, tau_y0, day, G, US, CS)
if (associated(forces%ustar)) then ; do j=js,je ; do i=is,ie
forces%ustar(i,j) = sqrt( US%L_to_Z * ( mag_tau + CS%gust(i,j) ) / CS%Rho0 )
enddo ; enddo ; endif
if (associated(forces%tau_mag)) then ; do j=js,je ; do i=is,ie
forces%tau_mag(i,j) = mag_tau + CS%gust(i,j)
enddo ; enddo ; endif
else
if (associated(forces%ustar)) then ; do j=js,je ; do i=is,ie
forces%ustar(i,j) = sqrt( US%L_to_Z * ( mag_tau + CS%gust_const ) / CS%Rho0 )
enddo ; enddo ; endif
if (associated(forces%tau_mag)) then ; do j=js,je ; do i=is,ie
forces%tau_mag(i,j) = mag_tau + CS%gust_const
enddo ; enddo ; endif
endif

call callTree_leave("wind_forcing_const")
Expand Down Expand Up @@ -522,9 +528,11 @@ subroutine wind_forcing_gyres(sfc_state, forces, day, G, US, CS)
! set the friction velocity
if (CS%answer_date < 20190101) then
do j=js,je ; do i=is,ie
forces%tau_mag(i,j) = CS%gust_const + sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2)))
forces%ustar(i,j) = sqrt(US%L_to_Z * ((CS%gust_const/CS%Rho0) + &
sqrt(0.5*(forces%tauy(i,j-1)*forces%tauy(i,j-1) + forces%tauy(i,j)*forces%tauy(i,j) + &
forces%taux(i-1,j)*forces%taux(i-1,j) + forces%taux(i,j)*forces%taux(i,j)))/CS%Rho0) )
sqrt(0.5*(forces%tauy(i,J-1)*forces%tauy(i,J-1) + forces%tauy(i,J)*forces%tauy(i,J) + &
forces%taux(I-1,j)*forces%taux(I-1,j) + forces%taux(I,j)*forces%taux(I,j)))/CS%Rho0) )
enddo ; enddo
else
call stresses_to_ustar(forces, G, US, CS)
Expand Down Expand Up @@ -725,11 +733,12 @@ subroutine wind_forcing_from_file(sfc_state, forces, day, G, US, CS)
if (.not.read_Ustar) then
if (CS%read_gust_2d) then
do j=js,je ; do i=is,ie
forces%ustar(i,j) = sqrt((CS%gust(i,j) + &
sqrt(temp_x(i,j)*temp_x(i,j) + temp_y(i,j)*temp_y(i,j))) * US%L_to_Z / CS%Rho0)
forces%tau_mag(i,j) = CS%gust(i,j) + sqrt(temp_x(i,j)**2 + temp_y(i,j)**2)
forces%ustar(i,j) = sqrt(forces%tau_mag(i,j) * US%L_to_Z / CS%Rho0)
enddo ; enddo
else
do j=js,je ; do i=is,ie
forces%tau_mag(i,j) = CS%gust_const + sqrt(temp_x(i,j)**2 + temp_y(i,j)**2)
forces%ustar(i,j) = sqrt(US%L_to_Z * (CS%gust_const/CS%Rho0 + &
sqrt(temp_x(i,j)*temp_x(i,j) + temp_y(i,j)*temp_y(i,j)) / CS%Rho0) )
enddo ; enddo
Expand Down Expand Up @@ -772,15 +781,19 @@ subroutine wind_forcing_from_file(sfc_state, forces, day, G, US, CS)
if (.not.read_Ustar) then
if (CS%read_gust_2d) then
do j=js,je ; do i=is,ie
forces%ustar(i,j) = sqrt((CS%gust(i,j) + &
sqrt(0.5*((forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2) + &
(forces%taux(i-1,j)**2 + forces%taux(i,j)**2))) ) * US%L_to_Z / CS%Rho0 )
forces%tau_mag(i,j) = CS%gust(i,j) + &
sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2)))
forces%ustar(i,j) = sqrt( forces%tau_mag(i,j) * US%L_to_Z / CS%Rho0 )
enddo ; enddo
else
do j=js,je ; do i=is,ie
forces%tau_mag(i,j) = CS%gust_const + &
sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2)))
forces%ustar(i,j) = sqrt(US%L_to_Z * ( (CS%gust_const/CS%Rho0) + &
sqrt(0.5*((forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2) + &
(forces%taux(i-1,j)**2 + forces%taux(i,j)**2)))/CS%Rho0))
sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2)))/CS%Rho0))
enddo ; enddo
endif
endif
Expand All @@ -792,6 +805,9 @@ subroutine wind_forcing_from_file(sfc_state, forces, day, G, US, CS)
if (read_Ustar) then
call MOM_read_data(filename, CS%Ustar_var, forces%ustar(:,:), &
G%Domain, timelevel=time_lev, scale=US%m_to_Z*US%T_to_s)
do j=js,je ; do i=is,ie
forces%tau_mag(i,j) = US%Z_to_L * CS%Rho0 * forces%ustar(i,j)**2
enddo ; enddo
endif

CS%wind_last_lev = time_lev
Expand All @@ -815,6 +831,7 @@ subroutine wind_forcing_by_data_override(sfc_state, forces, day, G, US, CS)
! Local variables
real :: temp_x(SZI_(G),SZJ_(G)) ! Pseudo-zonal wind stresses at h-points [R Z L T-2 ~> Pa].
real :: temp_y(SZI_(G),SZJ_(G)) ! Pseudo-meridional wind stresses at h-points [R Z L T-2 ~> Pa].
real :: ustar_tmp(SZI_(G),SZJ_(G)) ! The pre-override value of ustar [Z T-1 ~> m s-1]
integer :: i, j

call callTree_enter("wind_forcing_by_data_override, MOM_surface_forcing.F90")
Expand All @@ -840,17 +857,25 @@ subroutine wind_forcing_by_data_override(sfc_state, forces, day, G, US, CS)
if (CS%read_gust_2d) then
call data_override(G%Domain, 'gust', CS%gust, day, scale=US%Pa_to_RLZ_T2)
do j=G%jsc,G%jec ; do i=G%isc,G%iec
forces%ustar(i,j) = sqrt((sqrt(temp_x(i,j)**2 + temp_y(i,j)**2) + &
CS%gust(i,j)) * US%L_to_Z / CS%Rho0)
forces%tau_mag(i,j) = sqrt(temp_x(i,j)**2 + temp_y(i,j)**2) + CS%gust(i,j)
forces%ustar(i,j) = sqrt( forces%tau_mag(i,j) * US%L_to_Z / CS%Rho0 )
enddo ; enddo
else
do j=G%jsc,G%jec ; do i=G%isc,G%iec
forces%tau_mag(i,j) = sqrt(temp_x(i,j)**2 + temp_y(i,j)**2) + CS%gust_const
! forces%ustar(i,j) = sqrt( forces%tau_mag(i,j) * US%L_to_Z / CS%Rho0 )
forces%ustar(i,j) = sqrt(US%L_to_Z * (sqrt(temp_x(i,j)**2 + temp_y(i,j)**2)/CS%Rho0 + &
CS%gust_const/CS%Rho0))
enddo ; enddo
endif

! Give the data override the option to modify the newly calculated forces%ustar.
ustar_tmp(:,:) = forces%ustar(:,:)
call data_override(G%Domain, 'ustar', forces%ustar, day, scale=US%m_to_Z*US%T_to_s)
! Only reset values where data override of ustar has occurred
do j=G%jsc,G%jec ; do i=G%isc,G%iec ; if (ustar_tmp(i,j) /= forces%ustar(i,j)) then
forces%tau_mag(i,j) = US%Z_to_L * CS%Rho0 * forces%ustar(i,j)**2
endif ; enddo ; enddo

call pass_vector(forces%taux, forces%tauy, G%Domain, To_All)

Expand All @@ -875,15 +900,17 @@ subroutine stresses_to_ustar(forces, G, US, CS)

if (CS%read_gust_2d) then
do j=js,je ; do i=is,ie
forces%ustar(i,j) = sqrt( (CS%gust(i,j) + &
sqrt(0.5*((forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2) + &
(forces%taux(i-1,j)**2 + forces%taux(i,j)**2))) ) * I_rho )
forces%tau_mag(i,j) = CS%gust(i,j) + &
sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2)))
forces%ustar(i,j) = sqrt( forces%tau_mag(i,j) * I_rho )
enddo ; enddo
else
do j=js,je ; do i=is,ie
forces%ustar(i,j) = sqrt( (CS%gust_const + &
sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2))) ) * I_rho )
forces%tau_mag(i,j) = CS%gust_const + &
sqrt(0.5*((forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2) + &
(forces%taux(I-1,j)**2 + forces%taux(I,j)**2)))
forces%ustar(i,j) = sqrt( forces%tau_mag(i,j) * I_rho )
enddo ; enddo
endif

Expand Down
5 changes: 3 additions & 2 deletions config_src/drivers/solo_driver/user_surface_forcing.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -88,9 +88,10 @@ subroutine USER_wind_forcing(sfc_state, forces, day, G, US, CS)
! is always positive.
if (associated(forces%ustar)) then ; do j=js,je ; do i=is,ie
! This expression can be changed if desired, but need not be.
forces%ustar(i,j) = G%mask2dT(i,j) * sqrt((CS%gust_const + &
forces%tau_mag(i,j) = G%mask2dT(i,j) * (CS%gust_const + &
sqrt(0.5*(forces%taux(I-1,j)**2 + forces%taux(I,j)**2) + &
0.5*(forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2))) * (US%L_to_Z/CS%Rho0))
0.5*(forces%tauy(i,J-1)**2 + forces%tauy(i,J)**2)))
forces%ustar(i,j) = G%mask2dT(i,j) * sqrt(forces%tau_mag(i,j) * (US%L_to_Z/CS%Rho0))
enddo ; enddo ; endif

end subroutine USER_wind_forcing
Expand Down
Loading