Skip to content

Commit

Permalink
Merge branch 'dev/master' of github.com:ESMG/MOM6 into dev/master
Browse files Browse the repository at this point in the history
  • Loading branch information
kshedstrom committed May 13, 2016
2 parents 4d5d490 + d22c2b8 commit 8cc50a1
Show file tree
Hide file tree
Showing 4 changed files with 100 additions and 55 deletions.
4 changes: 2 additions & 2 deletions src/ALE/MOM_ALE.F90
Original file line number Diff line number Diff line change
Expand Up @@ -702,8 +702,8 @@ subroutine ALE_remap_scalar(CS, G, nk_src, h_src, s_src, h_dst, s_dst, all_cells
if (present(old_remap)) use_remapping_core_w = old_remap
n_points = nk_src

!$OMP parallel default(none) shared(CS,G,h_src,s_src,h_dst,s_dst,use_remapping_core_w &
!$OMP ,ignore_vanished_layers, nk_src,dx ) &
!$OMP parallel default(none) shared(CS,G,h_src,s_src,h_dst,s_dst &
!$OMP ,ignore_vanished_layers, use_remapping_core_w, nk_src,dx ) &
!$OMP firstprivate(n_points)
!$OMP do
do j = G%jsc,G%jec
Expand Down
8 changes: 4 additions & 4 deletions src/core/MOM_isopycnal_slopes.F90
Original file line number Diff line number Diff line change
Expand Up @@ -133,11 +133,11 @@ subroutine calc_isoneutral_slopes(G, GV, h, e, tv, dt_kappa_smooth, &

!$OMP parallel do default(none) shared(nz,is,ie,js,je,use_EOS,G,GV,pres,T,S, &
!$OMP IsdB,tv,h,h_neglect,e,dz_neglect, &
!$OMP h_neglect2,present_N2_u,G_Rho0,N2_u) &
!$OMP h_neglect2,present_N2_u,G_Rho0,N2_u,slope_x) &
!$OMP private(drdiA,drdiB,drdkL,drdkR,pres_u,T_u,S_u, &
!$OMP drho_dT_u,drho_dS_u,hg2A,hg2B,hg2L,hg2R,haA, &
!$OMP haB,haL,haR,dzaL,dzaR,wtA,wtB,wtL,wtR,drdz, &
!$OMP drdx,mag_grad2,Slope,slope2_Ratio,slope_x)
!$OMP drdx,mag_grad2,Slope,slope2_Ratio)
do j=js,je ; do K=nz,2,-1
if (.not.(use_EOS)) then
drdiA = 0.0 ; drdiB = 0.0
Expand Down Expand Up @@ -220,11 +220,11 @@ subroutine calc_isoneutral_slopes(G, GV, h, e, tv, dt_kappa_smooth, &
! Calculate the meridional isopycnal slope.
!$OMP parallel do default(none) shared(nz,is,ie,js,je,use_EOS,G,GV,pres,T,S, &
!$OMP IsdB,tv,h,h_neglect,e,dz_neglect, &
!$OMP h_neglect2,present_N2_v,G_Rho0,N2_v) &
!$OMP h_neglect2,present_N2_v,G_Rho0,N2_v,slope_y) &
!$OMP private(drdjA,drdjB,drdkL,drdkR,pres_v,T_v,S_v, &
!$OMP drho_dT_v,drho_dS_v,hg2A,hg2B,hg2L,hg2R,haA, &
!$OMP haB,haL,haR,dzaL,dzaR,wtA,wtB,wtL,wtR,drdz, &
!$OMP drdy,mag_grad2,Slope,slope2_Ratio,slope_y)
!$OMP drdy,mag_grad2,Slope,slope2_Ratio)
do j=js-1,je ; do K=nz,2,-1
if (.not.(use_EOS)) then
drdjA = 0.0 ; drdjB = 0.0
Expand Down
60 changes: 30 additions & 30 deletions src/parameterizations/vertical/MOM_energetic_PBL.F90
Original file line number Diff line number Diff line change
Expand Up @@ -374,47 +374,47 @@ subroutine energetic_PBL(h_3d, u_3d, v_3d, tv, fluxes, dt, Kd_int, G, GV, CS, &
reset_diags = .false. ! This is the second call to mixedlayer.

if (reset_diags) then
!$OMP parallel default(none) shared(is,ie,js,je,CS)
!!OMP parallel default(none) shared(is,ie,js,je,CS)
if (CS%TKE_diagnostics) then
!$OMP do
!!OMP do
do j=js,je ; do i=is,ie
CS%diag_TKE_wind(i,j) = 0.0 ; CS%diag_TKE_MKE(i,j) = 0.0
CS%diag_TKE_conv(i,j) = 0.0 ; CS%diag_TKE_forcing(i,j) = 0.0
CS%diag_TKE_mixing(i,j) = 0.0 ; CS%diag_TKE_mech_decay(i,j) = 0.0
CS%diag_TKE_conv_decay(i,j) = 0.0
enddo ; enddo
endif
!$OMP end parallel
!!OMP end parallel
endif


!$OMP parallel do default(none) shared(js,je,nz,is,ie,h_3d,u_3d,v_3d,tv,dt, &
!$OMP CS,G,GV,sfc_connected,fluxes,IdtdR0, &
!$OMP TKE_forced,debug,H_neglect,dSV_dT, &
!$OMP dSV_dS,I_dtmrho,C1_3,h_tt_min,vonKar, &
!$OMP max_itt,Kd_int) &
!$OMP private(h,u,v,T,S,Kd,mech_TKE_k,conv_PErel_k, &
!$OMP U_Star,absf,mech_TKE,conv_PErel,nstar_k, &
!$OMP h_sum,I_hs,h_bot,hb_hs,T0,S0,num_itts, &
!$OMP pres,dMass,dPres,dT_to_dPE,dS_to_dPE, &
!$OMP dT_to_dColHt,dS_to_dColHt,Kddt_h, &
!$OMP b_den_1,dT_to_dPE_a,dT_to_dColHt_a, &
!$OMP dS_to_dColHt_a,htot,uhtot,vhtot, &
!$OMP Idecay_len_TKE,exp_kh,nstar_FC,tot_TKE, &
!$OMP TKE_reduc,dTe_t2,dSe_t2,dTe,dSe,dt_h, &
!$OMP Convectively_stable,sfc_disconnect,b1, &
!$OMP c1,dT_km1_t2,dS_km1_t2,dTe_term, &
!$OMP dSe_term,MKE2_Hharm,vstar,h_tt, &
!$OMP Kd_guess0,Kddt_h_g0,dPEc_dKd_Kd0, &
!$OMP PE_chg_max,dPEa_dKd_g0,PE_chg_g0, &
!$OMP MKE_src,dPE_conv,Kddt_h_max,Kddt_h_min, &
!$OMP TKE_left_max,TKE_left_min,Kddt_h_guess, &
!$OMP TKE_left_itt,dPEa_dKd_itt,PE_chg_itt, &
!$OMP MKE_src_itt,Kddt_h_itt,dPEc_dKd,PE_chg, &
!$OMP dMKE_src_dK,TKE_left,use_Newt, &
!$OMP dKddt_h_Newt,Kddt_h_Newt,Kddt_h_next, &
!$OMP dKddt_h,Te,Se,Hsfc_used,dS_to_dPE_a, &
!$OMP dMKE_max,TKE_here)
!!OMP parallel do default(none) shared(js,je,nz,is,ie,h_3d,u_3d,v_3d,tv,dt, &
!!OMP CS,G,GV,fluxes,IdtdR0, &
!!OMP TKE_forced,debug,H_neglect,dSV_dT, &
!!OMP dSV_dS,I_dtmrho,C1_3,h_tt_min,vonKar, &
!!OMP max_itt,Kd_int) &
!!OMP private(i,j,k,h,u,v,T,S,Kd,mech_TKE_k,conv_PErel_k, &
!!OMP U_Star,absf,mech_TKE,conv_PErel,nstar_k, &
!!OMP h_sum,I_hs,h_bot,hb_hs,T0,S0,num_itts, &
!!OMP pres,dMass,dPres,dT_to_dPE,dS_to_dPE, &
!!OMP dT_to_dColHt,dS_to_dColHt,Kddt_h, &
!!OMP b_den_1,dT_to_dPE_a,dT_to_dColHt_a, &
!!OMP dS_to_dColHt_a,htot,uhtot,vhtot, &
!!OMP Idecay_len_TKE,exp_kh,nstar_FC,tot_TKE, &
!!OMP TKE_reduc,dTe_t2,dSe_t2,dTe,dSe,dt_h, &
!!OMP Convectively_stable,sfc_disconnect,b1, &
!!OMP c1,dT_km1_t2,dS_km1_t2,dTe_term, &
!!OMP dSe_term,MKE2_Hharm,vstar,h_tt, &
!!OMP Kd_guess0,Kddt_h_g0,dPEc_dKd_Kd0, &
!!OMP PE_chg_max,dPEa_dKd_g0,PE_chg_g0, &
!!OMP MKE_src,dPE_conv,Kddt_h_max,Kddt_h_min, &
!!OMP TKE_left_max,TKE_left_min,Kddt_h_guess, &
!!OMP TKE_left_itt,dPEa_dKd_itt,PE_chg_itt, &
!!OMP MKE_src_itt,Kddt_h_itt,dPEc_dKd,PE_chg, &
!!OMP dMKE_src_dK,TKE_left,use_Newt, &
!!OMP dKddt_h_Newt,Kddt_h_Newt,Kddt_h_next, &
!!OMP dKddt_h,Te,Se,Hsfc_used,dS_to_dPE_a, &
!!OMP dMKE_max,sfc_connected,TKE_here)
do j=js,je
! Copy the thicknesses and other fields to 2-d arrays.
do k=1,nz ; do i=is,ie
Expand Down
83 changes: 64 additions & 19 deletions src/parameterizations/vertical/MOM_vert_friction.F90
Original file line number Diff line number Diff line change
Expand Up @@ -151,6 +151,10 @@ module MOM_vert_friction
! to calculate the viscous coupling between layers
! except near the bottom. Otherwise the arithmetic
! mean thickness is used except near the bottom.
real :: harm_BL_val ! A scale to determine when water is in the boundary
! layers based solely on harmonic mean thicknesses
! for the purpose of determining the extent to which
! the thicknesses used in the viscosities are upwinded.
logical :: direct_stress ! If true, the wind stress is distributed over the
! topmost Hmix_stress of fluid and KVML may be very small.
logical :: dynamic_viscous_ML ! If true, use the results from a dynamic
Expand Down Expand Up @@ -626,11 +630,15 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)
real :: topfn ! A function which goes from 1 at the top to 0 much more
! than Htbl into the interior.
real :: z2 ! The distance from the bottom, normalized by Hbbl, nondim.
real :: z2_wt ! A nondimensional (0-1) weight used when calculating z2.
real :: h_neglect ! A thickness that is so small it is usually lost
! in roundoff and can be neglected, in H.
real :: H_to_m, m_to_H ! Unit conversion factors.

real :: h_arith ! The arithmetic mean thickness, in m or kg m-2.
real :: I_valBL ! The inverse of a scaling factor determining when water is
! still within the boundary layer, as determined by the sum
! of the harmonic mean thicknesses.
logical, dimension(SZIB_(G)) :: do_i, do_i_shelf
logical :: do_any_shelf

Expand All @@ -644,6 +652,7 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)
h_neglect = GV%H_subroundoff
H_to_m = GV%H_to_m ; m_to_H = GV%m_to_H
I_Hbbl(:) = 1.0 / (CS%Hbbl * GV%m_to_H + h_neglect)
I_valBL = 0.0 ; if (CS%harm_BL_val > 0.0) I_valBL = 1.0 / CS%harm_BL_val

if (CS%debug .or. (CS%id_hML_u > 0)) then
allocate(hML_u(G%IsdB:G%IedB,G%jsd:G%jed)) ; hML_u(:,:) = 0.0
Expand All @@ -662,11 +671,11 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)
endif

!$OMP parallel do default(none) shared(G,GV,CS,visc,Isq,ieq,nz,u,h,fluxes,hML_u, &
!$OMP h_neglect,dt,m_to_H) &
!$OMP h_neglect,dt,m_to_H,I_valBL) &
!$OMP firstprivate(i_hbbl) &
!$OMP private(do_i,kv_bbl,bbl_thick,z_i,h_harm,h_arith,hvel,z2, &
!$OMP botfn,zh,Dmin,zcol,a,do_any_shelf,do_i_shelf, &
!$OMP a_shelf,Ztop_min,I_HTbl,hvel_shelf,topfn,h_ml)
!$OMP a_shelf,Ztop_min,I_HTbl,hvel_shelf,topfn,h_ml,z2_wt)
do j=G%Jsc,G%Jec
do I=Isq,Ieq ; do_i(I) = (G%mask2dCu(I,j) > 0) ; enddo

Expand Down Expand Up @@ -710,9 +719,15 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)

hvel(I,k) = h_arith
if (u(I,j,k) * (h(i+1,j,k)-h(i,j,k)) > 0) then
z2 = max(zh(I), max(zcol(i),zcol(i+1)) + Dmin(I)) * I_Hbbl(I)
botfn = 1.0 / (1.0 + 0.09*z2*z2*z2*z2*z2*z2)
hvel(I,k) = (1.0-botfn)*h_arith + botfn*h_harm(I,k)
if (zh(I) * I_Hbbl(I) < CS%harm_BL_val) then
hvel(I,k) = h_harm(I,k)
else
z2_wt = 1.0 ; if (zh(I) * I_Hbbl(I) < 2.0*CS%harm_BL_val) &
z2_wt = max(0.0, min(1.0, zh(I) * I_Hbbl(I) * I_valBL - 1.0))
z2 = z2_wt * (max(zh(I), max(zcol(i),zcol(i+1)) + Dmin(I)) * I_Hbbl(I))
botfn = 1.0 / (1.0 + 0.09*z2*z2*z2*z2*z2*z2)
hvel(I,k) = (1.0-botfn)*h_arith + botfn*h_harm(I,k)
endif
endif

endif ; enddo ! i loop
Expand Down Expand Up @@ -750,9 +765,15 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)

hvel_shelf(I,k) = hvel(I,k)
if (u(I,j,k) * (h(i+1,j,k)-h(i,j,k)) > 0) then
z2 = max(zh(I), Ztop_min(I) - min(zcol(i),zcol(i+1))) * I_HTbl(I)
topfn = 1.0 / (1.0 + 0.09*z2**6)
hvel_shelf(I,k) = min(hvel(I,k), (1.0-topfn)*h_arith + topfn*h_harm(I,k))
if (zh(I) * I_HTbl(I) < CS%harm_BL_val) then
hvel_shelf(I,k) = min(hvel(I,k), h_harm(I,k))
else
z2_wt = 1.0 ; if (zh(I) * I_HTbl(I) < 2.0*CS%harm_BL_val) &
z2_wt = max(0.0, min(1.0, zh(I) * I_HTbl(I) * I_valBL - 1.0))
z2 = z2_wt * (max(zh(I), Ztop_min(I) - min(zcol(i),zcol(i+1))) * I_HTbl(I))
topfn = 1.0 / (1.0 + 0.09*z2**6)
hvel_shelf(I,k) = min(hvel(I,k), (1.0-topfn)*h_arith + topfn*h_harm(I,k))
endif
endif
endif ; enddo
enddo
Expand All @@ -766,8 +787,11 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)

if (do_any_shelf) then
do K=1,nz+1 ; do I=Isq,Ieq ; if (do_i_shelf(I)) then
CS%a_u(I,j,K) = fluxes%frac_shelf_u(I,j) * a_shelf(I,k) + &
CS%a_u(I,j,K) = fluxes%frac_shelf_u(I,j) * a_shelf(I,K) + &
(1.0-fluxes%frac_shelf_u(I,j)) * a(I,K)
! This is Alistair's suggestion, but it destabilizes the model. I do not know why. RWH
! CS%a_u(I,j,K) = fluxes%frac_shelf_u(I,j) * max(a_shelf(I,K), a(I,K)) + &
! (1.0-fluxes%frac_shelf_u(I,j)) * a(I,K)
elseif (do_i(I)) then
CS%a_u(I,j,K) = a(I,K)
endif ; enddo ; enddo
Expand All @@ -788,11 +812,11 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)

! Now work on v-points.
!$OMP parallel do default(none) shared(G,GV,CS,visc,is,ie,Jsq,Jeq,nz,v,h,fluxes,hML_v, &
!$OMP h_neglect,dt,m_to_H) &
!$OMP h_neglect,dt,m_to_H,I_valBL) &
!$OMP firstprivate(i_hbbl) &
!$OMP private(do_i,kv_bbl,bbl_thick,z_i,h_harm,h_arith,hvel,z2, &
!$OMP botfn,zh,Dmin,zcol1,zcol2,a,do_any_shelf,do_i_shelf, &
!$OMP a_shelf,Ztop_min,I_HTbl,hvel_shelf,topfn,h_ml)
!$OMP a_shelf,Ztop_min,I_HTbl,hvel_shelf,topfn,h_ml,z2_wt)
do J=Jsq,Jeq
do i=is,ie ; do_i(i) = (G%mask2dCv(i,J) > 0) ; enddo

Expand Down Expand Up @@ -837,10 +861,16 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)

hvel(i,k) = h_arith
if (v(i,J,k) * (h(i,j+1,k)-h(i,j,k)) > 0) then
z2 = max(zh(i), max(zcol1(i),zcol2(i)) + Dmin(i)) * I_Hbbl(i)
botfn = 1.0 / (1.0 + 0.09*z2*z2*z2*z2*z2*z2)
hvel(i,k) = (1.0-botfn)*h_arith + botfn*h_harm(i,k)
endif
if (zh(i) * I_Hbbl(i) < CS%harm_BL_val) then
hvel(i,k) = h_harm(i,k)
else
z2_wt = 1.0 ; if (zh(i) * I_Hbbl(i) < 2.0*CS%harm_BL_val) &
z2_wt = max(0.0, min(1.0, zh(i) * I_Hbbl(i) * I_valBL - 1.0))
z2 = z2_wt * (max(zh(i), max(zcol1(i),zcol2(i)) + Dmin(i)) * I_Hbbl(i))
botfn = 1.0 / (1.0 + 0.09*z2*z2*z2*z2*z2*z2)
hvel(i,k) = (1.0-botfn)*h_arith + botfn*h_harm(i,k)
endif
endif

endif ; enddo ; enddo ! i & k loops
endif
Expand Down Expand Up @@ -876,10 +906,16 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)

hvel_shelf(i,k) = hvel(i,k)
if (v(i,j,k) * (h(i,j+1,k)-h(i,j,k)) > 0) then
z2 = max(zh(i), Ztop_min(i) - min(zcol1(i),zcol2(i))) * I_HTbl(i)
topfn = 1.0 / (1.0 + 0.09*z2**6)
hvel_shelf(i,k) = min(hvel(i,k), (1.0-topfn)*h_arith + topfn*h_harm(i,k))
endif
if (zh(i) * I_HTbl(i) < CS%harm_BL_val) then
hvel_shelf(i,k) = min(hvel(i,k), h_harm(i,k))
else
z2_wt = 1.0 ; if (zh(i) * I_HTbl(i) < 2.0*CS%harm_BL_val) &
z2_wt = max(0.0, min(1.0, zh(i) * I_HTbl(i) * I_valBL - 1.0))
z2 = z2_wt * (max(zh(i), Ztop_min(i) - min(zcol1(i),zcol2(i))) * I_HTbl(i))
topfn = 1.0 / (1.0 + 0.09*z2**6)
hvel_shelf(i,k) = min(hvel(i,k), (1.0-topfn)*h_arith + topfn*h_harm(i,k))
endif
endif
endif ; enddo
enddo
call find_coupling_coef(a_shelf, hvel_shelf, do_i_shelf, h_harm, &
Expand All @@ -894,6 +930,9 @@ subroutine vertvisc_coef(u, v, h, fluxes, visc, dt, G, GV, CS)
do K=1,nz+1 ; do i=is,ie ; if (do_i_shelf(i)) then
CS%a_v(i,J,K) = fluxes%frac_shelf_v(i,J) * a_shelf(i,k) + &
(1.0-fluxes%frac_shelf_v(i,J)) * a(i,K)
! This is Alistair's suggestion, but it destabilizes the model. I do not know why. RWH
! CS%a_v(i,J,K) = fluxes%frac_shelf_v(i,J) * max(a_shelf(i,K), a(i,K)) + &
! (1.0-fluxes%frac_shelf_v(i,J)) * a(i,K)
elseif (do_i(i)) then
CS%a_v(i,J,K) = a(i,K)
endif ; enddo ; enddo
Expand Down Expand Up @@ -1435,6 +1474,12 @@ subroutine vertvisc_init(MIS, Time, G, GV, param_file, diag, ADp, dirs, &
call get_param(param_file, mod, "HARMONIC_VISC", CS%harmonic_visc, &
"If true, use the harmonic mean thicknesses for \n"//&
"calculating the vertical viscosity.", default=.false.)
call get_param(param_file, mod, "HARMONIC_BL_SCALE", CS%harm_BL_val, &
"A scale to determine when water is in the boundary \n"//&
"layers based solely on harmonic mean thicknesses for \n"//&
"the purpose of determining the extent to which the \n"//&
"thicknesses used in the viscosities are upwinded.", &
default=0.0, units="nondim")
call get_param(param_file, mod, "DEBUG", CS%debug, default=.false.)

if (GV%nkml < 1) &
Expand Down

0 comments on commit 8cc50a1

Please sign in to comment.