From da007a6e8c5d16b5de8337fd88a6ebedd8f344e6 Mon Sep 17 00:00:00 2001 From: Michael Levy Date: Fri, 10 Apr 2015 16:30:28 -0600 Subject: [PATCH] Allow non-zero nonlocal term at surface Steve G requested the ability to match the non-local term to incoming flux values at the surface -- this commit is a preliminary step in that direction in that it adds a kpp_init flag lnonzero_surf_nonlocal that, when set to true, returns G(0) instead of Cs*G(0) as gamma / Q. Note that this still requires an option to allow G(1) = 1, which currently isn't available. --- src/shared/cvmix_kpp.F90 | 103 ++++++++++++++------------------------- 1 file changed, 37 insertions(+), 66 deletions(-) diff --git a/src/shared/cvmix_kpp.F90 b/src/shared/cvmix_kpp.F90 index 642f6caa6..7d381e81d 100644 --- a/src/shared/cvmix_kpp.F90 +++ b/src/shared/cvmix_kpp.F90 @@ -73,7 +73,6 @@ module cvmix_kpp public :: cvmix_kpp_compute_shape_function_coeffs public :: cvmix_kpp_compute_kOBL_depth public :: cvmix_kpp_compute_enhanced_diff - public :: cvmix_kpp_compute_nonlocal public :: cvmix_kpp_compute_nu_at_OBL_depth_LMD94 interface cvmix_coeffs_kpp @@ -116,6 +115,9 @@ module cvmix_kpp real(cvmix_r8) :: vonkarman ! von Karman constant real(cvmix_r8) :: Cstar ! coefficient for nonlinear transport + real(cvmix_r8) :: nonlocal_coeff ! Cs from Eq (20) in LMD94 + ! Default value comes from paper, but + ! some users may set it = 1. ! For velocity scale function, _m => momentum and _s => scalar (tracer) real(cvmix_r8) :: zeta_m ! parameter for computing vel scale func @@ -190,7 +192,7 @@ subroutine cvmix_init_kpp(ri_crit, minOBLdepth, maxOBLdepth, minVtsqr, & Cv, interp_type, interp_type2, MatchTechnique, & old_vals, lEkman, lMonOb, lnoDGat1, & lavg_N_or_Nsqr, lenhanced_diff, & - CVmix_kpp_params_user) + lnonzero_surf_nonlocal, CVmix_kpp_params_user) ! !DESCRIPTION: ! Initialization routine for KPP mixing. @@ -219,7 +221,8 @@ subroutine cvmix_init_kpp(ri_crit, minOBLdepth, maxOBLdepth, minVtsqr, & lMonOb, & lnoDGat1, & lavg_N_or_Nsqr, & - lenhanced_diff + lenhanced_diff, & + lnonzero_surf_nonlocal ! !OUTPUT PARAMETERS: type(cvmix_kpp_params_type), intent(inout), target, optional :: & @@ -229,6 +232,9 @@ subroutine cvmix_init_kpp(ri_crit, minOBLdepth, maxOBLdepth, minVtsqr, & !BOC real(cvmix_r8) :: zm, zs, a_m, a_s, c_m, c_s + real(cvmix_r8) :: Cstar_loc, vonkar_loc, surf_layer_ext_loc + real(cvmix_r8) :: nonlocal_coeff + if (present(ri_crit)) then if (ri_crit.lt.cvmix_zero) then @@ -275,16 +281,18 @@ subroutine cvmix_init_kpp(ri_crit, minOBLdepth, maxOBLdepth, minVtsqr, & print*, "ERROR: vonkarman can not be negative." stop 1 end if - call cvmix_put_kpp('vonkarman', vonkarman, CVmix_kpp_params_user) + vonkar_loc = vonkarman else - call cvmix_put_kpp('vonkarman', 0.4_cvmix_r8, CVmix_kpp_params_user) + vonkar_loc = 0.4_cvmix_r8 end if + call cvmix_put_kpp('vonkarman', vonkar_loc, CVmix_kpp_params_user) if (present(Cstar)) then - call cvmix_put_kpp('Cstar', Cstar, CVmix_kpp_params_user) + Cstar_loc = Cstar else - call cvmix_put_kpp('Cstar', 10, CVmix_kpp_params_user) + Cstar_loc = real(10,cvmix_r8) end if + call cvmix_put_kpp('Cstar', Cstar_loc, CVmix_kpp_params_user) if (present(zeta_m)) then if (zeta_m.ge.cvmix_zero) then @@ -332,11 +340,12 @@ subroutine cvmix_init_kpp(ri_crit, minOBLdepth, maxOBLdepth, minVtsqr, & print*, "surf_layer_ext must be between 0 and 1, inclusive." stop 1 end if - call cvmix_put_kpp('surf_layer_ext', surf_layer_ext, & - CVmix_kpp_params_user) + surf_layer_ext_loc = surf_layer_ext else - call cvmix_put_kpp('surf_layer_ext', 0.1_cvmix_r8, CVmix_kpp_params_user) + surf_layer_ext_loc = 0.1_cvmix_r8 end if + call cvmix_put_kpp('surf_layer_ext', surf_layer_ext_loc, & + CVmix_kpp_params_user) if (present(Cv)) then ! Use scalar Cv parameter @@ -469,6 +478,17 @@ subroutine cvmix_init_kpp(ri_crit, minOBLdepth, maxOBLdepth, minVtsqr, & call cvmix_put_kpp('lenhanced_diff', .true., CVmix_kpp_params_user) end if + ! By default, assume that G(0) = 0 for nonlocal term + nonlocal_coeff = (Cstar_loc*vonkar_loc* & + (vonkar_loc*surf_layer_ext_loc*c_s)** & + (cvmix_one/real(3,cvmix_r8))) + if (present(lnonzero_surf_nonlocal)) then + if (lnonzero_surf_nonlocal) then + nonlocal_coeff = real(1,cvmix_r8) + end if + end if + call cvmix_put_kpp('nonlocal_coeff',nonlocal_coeff,CVmix_kpp_params_user) + !EOC end subroutine cvmix_init_kpp @@ -618,7 +638,7 @@ subroutine cvmix_coeffs_kpp_low(Mdiff_out, Tdiff_out, Sdiff_out, zw, zt, & real(cvmix_r8) :: dMshapeAt1, dTshapeAt1, dSshapeAt1 ! [MTS]shapeAtS is value of shape function at sigma = S - real(cvmix_r8) :: MshapeAtS, TshapeAtS, SshapeAtS + real(cvmix_r8) :: MshapeAtS, TshapeAtS, SshapeAtS, GAtS ! [MTS]diff_OBL is value of diffusivity at OBL depth ! d[MTS]diff_OBL is value of derivative of diffusivity at OBL depth @@ -897,10 +917,10 @@ subroutine cvmix_coeffs_kpp_low(Mdiff_out, Tdiff_out, Sdiff_out, zw, zt, & ! (3c) Compute nonlocal term at each cell interface if ((.not.lstable).and.(kw.le.kwup)) then - call cvmix_kpp_compute_nonlocal(Tshape2, sigma(kw), Tnonlocal(kw), & - CVmix_kpp_params_user) - call cvmix_kpp_compute_nonlocal(Sshape2, sigma(kw), Snonlocal(kw), & - CVmix_kpp_params_user) + GAtS = cvmix_math_evaluate_cubic(Tshape2, sigma(kw)) + Tnonlocal(kw) = CVmix_kpp_params_user%nonlocal_coeff*GAtS + GAtS = cvmix_math_evaluate_cubic(Sshape2, sigma(kw)) + Snonlocal(kw) = CVmix_kpp_params_user%nonlocal_coeff*GAtS end if ! (3d) Diffusivity = OBL_depth * (turbulent scale) * G(sigma) @@ -1015,6 +1035,8 @@ subroutine cvmix_put_kpp_real(varname, val, CVmix_kpp_params_user) CVmix_kpp_params_out%surf_layer_ext = val case ('Cv') CVmix_kpp_params_out%Cv = val + case ('nonlocal_coeff') + CVmix_kpp_params_out%nonlocal_coeff = val case DEFAULT print*, "ERROR: ", trim(varname), " not a valid choice!" stop 1 @@ -1530,57 +1552,6 @@ end subroutine cvmix_kpp_compute_enhanced_diff !BOP -! !IROUTINE: cvmix_kpp_compute_nonlocal -! !INTERFACE: - - subroutine cvmix_kpp_compute_nonlocal(shape_fun, sigma, nonlocal, & - CVmix_kpp_params_user) - -! !DESCRIPTION: -! Compute the nonlocal transport contribution to vertical turbulent fluxes. -! Note that Large, et al., refer to $\gamma_x$ as the non-local term, while -! this routine computes $K_x\gamma_x/$[surface forcing] -!\\ -!\\ - -! !INPUT PARAMETERS: - type(cvmix_kpp_params_type), intent(in), optional, target :: & - CVmix_kpp_params_user - real(cvmix_r8), dimension(4), intent(in) :: shape_fun - real(cvmix_r8), intent(in) :: sigma - -! !OUTPUT PARAMETERS: - real(cvmix_r8), intent(out) :: nonlocal - - ! Local variables - type(cvmix_kpp_params_type), pointer :: CVmix_kpp_params_in - - ! Constants from params - real(cvmix_r8) :: Cstar, vonkar, c_s, surf_layer_ext - - real(cvmix_r8) :: GatS -!EOP -!BOC - - CVmix_kpp_params_in => CVmix_kpp_params_saved - if (present(CVmix_kpp_params_user)) then - CVmix_kpp_params_in => CVmix_kpp_params_user - end if - vonkar = cvmix_get_kpp_real('vonkarman', CVmix_kpp_params_in) - Cstar = cvmix_get_kpp_real('Cstar', CVmix_kpp_params_in) - surf_layer_ext = cvmix_get_kpp_real('surf_layer_ext', CVmix_kpp_params_in) - c_s = cvmix_get_kpp_real('c_s', CVmix_kpp_params_in) - - GatS = cvmix_math_evaluate_cubic(shape_fun, sigma) - - nonlocal = GatS*(Cstar*vonkar*(vonkar*surf_layer_ext*c_s)** & - (cvmix_one/real(3,cvmix_r8))) -! EOC - - end subroutine cvmix_kpp_compute_nonlocal - -!BOP - ! !IROUTINE: cvmix_kpp_compute_OBL_depth_wrap ! !INTERFACE: