Skip to content

Commit

Permalink
Merge branch 'dev/master' into dev/gfdl
Browse files Browse the repository at this point in the history
  • Loading branch information
Hallberg-NOAA committed Apr 24, 2018
2 parents 586ba3a + 77a97ba commit 3c70582
Show file tree
Hide file tree
Showing 17 changed files with 2,467 additions and 1,518 deletions.
2 changes: 1 addition & 1 deletion config_src/mct_driver/ocn_comp_mct.F90
Original file line number Diff line number Diff line change
Expand Up @@ -580,7 +580,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename )

if (debug .and. root_pe().eq.pe_here()) print *, "calling ocn_domain_mct"
call ocn_domain_mct(lsize, MOM_MCT_gsmap, MOM_MCT_dom)
call ocn_domain_mct(lsize*km, MOM_MCT_gsmap3d, MOM_MCT_dom3d) !TODO: this is not used
!call ocn_domain_mct(lsize*km, MOM_MCT_gsmap3d, MOM_MCT_dom3d) !TODO: this is not used

! Inialize mct attribute vectors

Expand Down
10 changes: 7 additions & 3 deletions src/core/MOM.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1466,7 +1466,8 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
real, dimension(:,:), pointer :: shelf_area
type(MOM_restart_CS), pointer :: restart_CSp_tmp => NULL()
type(group_pass_type) :: tmp_pass_uv_T_S_h, pass_uv_T_S_h
type(group_pass_type) :: tmp_pass_Kv_turb
! GMM, the following *is not* used. Should we delete it?
type(group_pass_type) :: tmp_pass_Kv_shear

real :: default_val ! default value for a parameter
logical :: write_geom_files ! If true, write out the grid geometry files.
Expand Down Expand Up @@ -2308,8 +2309,8 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &

call do_group_pass(pass_uv_T_S_h, G%Domain)

if (associated(CS%visc%Kv_turb)) &
call pass_var(CS%visc%Kv_turb, G%Domain, To_All+Omit_Corners, halo=1)
if (associated(CS%visc%Kv_shear)) &
call pass_var(CS%visc%Kv_shear, G%Domain, To_All+Omit_Corners, halo=1)

call cpu_clock_end(id_clock_pass_init)

Expand Down Expand Up @@ -2920,6 +2921,9 @@ subroutine MOM_end(CS)
call tracer_registry_end(CS%tracer_Reg)
call tracer_flow_control_end(CS%tracer_flow_CSp)

! GMM, the following is commented because it fails on Travis.
!if (associated(CS%diabatic_CSp)) call diabatic_driver_end(CS%diabatic_CSp)

if (CS%offline_tracer_mode) call offline_transport_end(CS%offline_CSp)

DEALLOC_(CS%uhtr) ; DEALLOC_(CS%vhtr)
Expand Down
13 changes: 8 additions & 5 deletions src/core/MOM_variables.F90
Original file line number Diff line number Diff line change
Expand Up @@ -173,7 +173,7 @@ module MOM_variables
!! coefficients, and related fields.
type, public :: vertvisc_type
real :: Prandtl_turb !< The Prandtl number for the turbulent diffusion
!! that is captured in Kd_turb.
!! that is captured in Kd_shear.
real, pointer, dimension(:,:) :: &
bbl_thick_u => NULL(), & !< The bottom boundary layer thickness at the
!! u-points, in m.
Expand Down Expand Up @@ -224,10 +224,13 @@ module MOM_variables
! Kd_extra_S is positive for salt fingering; Kd_extra_T
! is positive for double diffusive convection. These
! are only allocated if DOUBLE_DIFFUSION is true.
Kd_turb => NULL(), &!< The turbulent diapycnal diffusivity at the interfaces
!! between each layer, in m2 s-1.
Kv_turb => NULL(), &!< The turbulent vertical viscosity at the interfaces
!! between each layer, in m2 s-1.
Kd_shear => NULL(), &!< The shear-driven turbulent diapycnal diffusivity
!! at the interfaces between each layer, in m2 s-1.
Kv_shear => NULL(), &!< The shear-driven turbulent vertical viscosity
!! at the interfaces between each layer, in m2 s-1.
Kv_slow => NULL(), &!< The turbulent vertical viscosity component due to
!! "slow" processes (e.g., tidal, background,
!! convection etc).
TKE_turb => NULL() !< The turbulent kinetic energy per unit mass defined
!! at the interfaces between each layer, in m2 s-2.
end type vertvisc_type
Expand Down
6 changes: 3 additions & 3 deletions src/initialization/MOM_state_initialization.F90
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,7 @@ module MOM_state_initialization
use Rossby_front_2d_initialization, only : Rossby_front_initialize_temperature_salinity
use Rossby_front_2d_initialization, only : Rossby_front_initialize_velocity
use SCM_idealized_hurricane, only : SCM_idealized_hurricane_TS_init
use SCM_CVmix_tests, only: SCM_CVmix_tests_TS_init
use SCM_CVMix_tests, only: SCM_CVMix_tests_TS_init
use dyed_channel_initialization, only : dyed_channel_set_OBC_tracer_data
use dyed_obcs_initialization, only : dyed_obcs_set_OBC_data
use supercritical_initialization, only : supercritical_set_OBC_data
Expand Down Expand Up @@ -338,7 +338,7 @@ subroutine MOM_initialize_state(u, v, h, tv, Time, G, GV, PF, dirs, &
" \t dumbbell - sloshing channel ICs. \n"//&
" \t rossby_front - a mixed layer front in thermal wind balance.\n"//&
" \t SCM_ideal_hurr - used in the SCM idealized hurricane test.\n"//&
" \t SCM_CVmix_tests - used in the SCM CVmix tests.\n"//&
" \t SCM_CVMix_tests - used in the SCM CVMix tests.\n"//&
" \t USER - call a user modified routine.", &
fail_if_missing=new_sim, do_not_log=just_read)
! " \t baroclinic_zone - an analytic baroclinic zone. \n"//&
Expand Down Expand Up @@ -371,7 +371,7 @@ subroutine MOM_initialize_state(u, v, h, tv, Time, G, GV, PF, dirs, &
tv%S, h, G, GV, PF, eos, just_read_params=just_read)
case ("SCM_ideal_hurr"); call SCM_idealized_hurricane_TS_init ( tv%T, &
tv%S, h, G, GV, PF, just_read_params=just_read)
case ("SCM_CVmix_tests"); call SCM_CVmix_tests_TS_init (tv%T, &
case ("SCM_CVMix_tests"); call SCM_CVMix_tests_TS_init (tv%T, &
tv%S, h, G, GV, PF, just_read_params=just_read)
case ("dense"); call dense_water_initialize_TS(G, GV, PF, eos, tv%T, tv%S, &
h, just_read_params=just_read)
Expand Down
270 changes: 270 additions & 0 deletions src/parameterizations/vertical/MOM_CVMix_conv.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,270 @@
!> Interface to CVMix convection scheme.
module MOM_CVMix_conv

! This file is part of MOM6. See LICENSE.md for the license.

use MOM_diag_mediator, only : diag_ctrl, time_type, register_diag_field
use MOM_diag_mediator, only : post_data
use MOM_EOS, only : calculate_density
use MOM_variables, only : thermo_var_ptrs
use MOM_error_handler, only : MOM_error, is_root_pe, FATAL, WARNING, NOTE
use MOM_file_parser, only : openParameterBlock, closeParameterBlock
use MOM_debugging, only : hchksum
use MOM_grid, only : ocean_grid_type
use MOM_verticalGrid, only : verticalGrid_type
use MOM_file_parser, only : get_param, log_version, param_file_type
use CVMix_convection, only : CVMix_init_conv, CVMix_coeffs_conv
use CVMix_kpp, only : CVMix_kpp_compute_kOBL_depth

implicit none ; private

#include <MOM_memory.h>

public CVMix_conv_init, calculate_CVMix_conv, CVMix_conv_end, CVMix_conv_is_used

!> Control structure including parameters for CVMix convection.
type, public :: CVMix_conv_cs

! Parameters
real :: kd_conv_const !< diffusivity constant used in convective regime (m2/s)
real :: kv_conv_const !< viscosity constant used in convective regime (m2/s)
real :: bv_sqr_conv !< Threshold for squared buoyancy frequency
!! needed to trigger Brunt-Vaisala parameterization (1/s^2)
real :: min_thickness !< Minimum thickness allowed (m)
logical :: debug !< If true, turn on debugging

! Daignostic handles and pointers
type(diag_ctrl), pointer :: diag => NULL()
integer :: id_N2 = -1, id_kd_conv = -1, id_kv_conv = -1

! Diagnostics arrays
real, allocatable, dimension(:,:,:) :: N2 !< Squared Brunt-Vaisala frequency (1/s2)
real, allocatable, dimension(:,:,:) :: kd_conv !< Diffusivity added by convection (m2/s)
real, allocatable, dimension(:,:,:) :: kv_conv !< Viscosity added by convection (m2/s)

end type CVMix_conv_cs

character(len=40) :: mdl = "MOM_CVMix_conv" !< This module's name.

contains

!> Initialized the CVMix convection mixing routine.
logical function CVMix_conv_init(Time, G, GV, param_file, diag, CS)

type(time_type), intent(in) :: Time !< The current time.
type(ocean_grid_type), intent(in) :: G !< Grid structure.
type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure.
type(param_file_type), intent(in) :: param_file !< Run-time parameter file handle
type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure.
type(CVMix_conv_cs), pointer :: CS !< This module's control structure.

! Local variables
real :: prandtl_conv !< Turbulent Prandtl number used in convective instabilities.
logical :: useEPBL !< If True, use the ePBL boundary layer scheme.

! This include declares and sets the variable "version".
#include "version_variable.h"

if (associated(CS)) then
call MOM_error(WARNING, "CVMix_conv_init called with an associated "// &
"control structure.")
return
endif
allocate(CS)

! Read parameters
call log_version(param_file, mdl, version, &
"Parameterization of enhanced mixing due to convection via CVMix")
call get_param(param_file, mdl, "USE_CVMix_CONVECTION", CVMix_conv_init, &
"If true, turns on the enhanced mixing due to convection \n"// &
"via CVMix. This scheme increases diapycnal diffs./viscs. \n"// &
" at statically unstable interfaces. Relevant parameters are \n"// &
"contained in the CVMix_CONVECTION% parameter block.", &
default=.false.)

if (.not. CVMix_conv_init) return

call get_param(param_file, mdl, "ENERGETICS_SFC_PBL", useEPBL, default=.false., &
do_not_log=.true.)

! Warn user if EPBL is being used, since in this case mixing due to convection will
! be aplied in the boundary layer
if (useEPBL) then
call MOM_error(WARNING, 'MOM_CVMix_conv_init: '// &
'CVMix convection may not be properly applied when ENERGETICS_SFC_PBL = True'//&
'as convective mixing might occur in the boundary layer.')
endif

call get_param(param_file, mdl, 'DEBUG', CS%debug, default=.False., do_not_log=.True.)

call get_param(param_file, mdl, 'MIN_THICKNESS', CS%min_thickness, default=0.001, do_not_log=.True.)

call openParameterBlock(param_file,'CVMix_CONVECTION')

call get_param(param_file, mdl, "PRANDTL_CONV", prandtl_conv, &
"The turbulent Prandtl number applied to convective \n"//&
"instabilities (i.e., used to convert KD_CONV into KV_CONV)", &
units="nondim", default=1.0)

call get_param(param_file, mdl, 'KD_CONV', CS%kd_conv_const, &
"Diffusivity used in convective regime. Corresponding viscosity \n" // &
"(KV_CONV) will be set to KD_CONV * PRANDTL_TURB.", &
units='m2/s', default=1.00)

call get_param(param_file, mdl, 'BV_SQR_CONV', CS%bv_sqr_conv, &
"Threshold for squared buoyancy frequency needed to trigger \n" // &
"Brunt-Vaisala parameterization.", &
units='1/s^2', default=0.0)

call closeParameterBlock(param_file)

! set kv_conv_const based on kd_conv_const and prandtl_conv
CS%kv_conv_const = CS%kd_conv_const * prandtl_conv

! allocate arrays and set them to zero
allocate(CS%N2(SZI_(G), SZJ_(G), SZK_(G)+1)); CS%N2(:,:,:) = 0.
allocate(CS%kd_conv(SZI_(G), SZJ_(G), SZK_(G)+1)); CS%kd_conv(:,:,:) = 0.
allocate(CS%kv_conv(SZI_(G), SZJ_(G), SZK_(G)+1)); CS%kv_conv(:,:,:) = 0.

! Register diagnostics
CS%diag => diag
CS%id_N2 = register_diag_field('ocean_model', 'N2_conv', diag%axesTi, Time, &
'Square of Brunt-Vaisala frequency used by MOM_CVMix_conv module', '1/s2')
CS%id_kd_conv = register_diag_field('ocean_model', 'kd_conv', diag%axesTi, Time, &
'Additional diffusivity added by MOM_CVMix_conv module', 'm2/s')
CS%id_kv_conv = register_diag_field('ocean_model', 'kv_conv', diag%axesTi, Time, &
'Additional viscosity added by MOM_CVMix_conv module', 'm2/s')

call CVMix_init_conv(convect_diff=CS%kd_conv_const, &
convect_visc=CS%kv_conv_const, &
lBruntVaisala=.true., &
BVsqr_convect=CS%bv_sqr_conv)

end function CVMix_conv_init

!> Subroutine for calculating enhanced diffusivity/viscosity
!! due to convection via CVMix
subroutine calculate_CVMix_conv(h, tv, G, GV, CS, hbl)

type(ocean_grid_type), intent(in) :: G !< Grid structure.
type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure.
real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thickness, in m or kg m-2.
type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamics structure.
type(CVMix_conv_cs), pointer :: CS !< The control structure returned
!! by a previous call to CVMix_conv_init.
real, dimension(:,:), optional, pointer :: hbl!< Depth of ocean boundary layer (m)

! local variables
real, dimension(SZK_(G)) :: rho_lwr !< Adiabatic Water Density, this is a dummy
!! variable since here convection is always
!! computed based on Brunt Vaisala.
real, dimension(SZK_(G)) :: rho_1d !< water density in a column, this is also
!! a dummy variable, same reason as above.
real, dimension(SZK_(G)+1) :: iFaceHeight !< Height of interfaces (m)
real, dimension(SZK_(G)) :: cellHeight !< Height of cell centers (m)
integer :: kOBL !< level of OBL extent
real :: pref, g_o_rho0, rhok, rhokm1, dz, dh, hcorr
integer :: i, j, k

g_o_rho0 = GV%g_Earth / GV%Rho0

! initialize dummy variables
rho_lwr(:) = 0.0; rho_1d(:) = 0.0

if (.not. associated(hbl)) then
allocate(hbl(SZI_(G), SZJ_(G)));
hbl(:,:) = 0.0
endif

do j = G%jsc, G%jec
do i = G%isc, G%iec

! set N2 to zero at the top- and bottom-most interfaces
CS%N2(i,j,1) = 0.
CS%N2(i,j,G%ke+1) = 0.

! skip calling at land points
!if (G%mask2dT(i,j) == 0.) cycle

pRef = 0.
! Compute Brunt-Vaisala frequency (static stability) on interfaces
do k=2,G%ke

! pRef is pressure at interface between k and km1.
pRef = pRef + GV%H_to_Pa * h(i,j,k)
call calculate_density(tv%t(i,j,k), tv%s(i,j,k), pref, rhok, tv%eqn_of_state)
call calculate_density(tv%t(i,j,k-1), tv%s(i,j,k-1), pref, rhokm1, tv%eqn_of_state)

dz = ((0.5*(h(i,j,k-1) + h(i,j,k))+GV%H_subroundoff)*GV%H_to_m)
CS%N2(i,j,k) = g_o_rho0 * (rhok - rhokm1) / dz ! Can be negative

enddo

iFaceHeight(1) = 0.0 ! BBL is all relative to the surface
hcorr = 0.0
! compute heights at cell center and interfaces
do k=1,G%ke
dh = h(i,j,k) * GV%H_to_m ! Nominal thickness to use for increment
dh = dh + hcorr ! Take away the accumulated error (could temporarily make dh<0)
hcorr = min( dh - CS%min_thickness, 0. ) ! If inflating then hcorr<0
dh = max( dh, CS%min_thickness ) ! Limit increment dh>=min_thickness
cellHeight(k) = iFaceHeight(k) - 0.5 * dh
iFaceHeight(k+1) = iFaceHeight(k) - dh
enddo

kOBL = CVMix_kpp_compute_kOBL_depth(iFaceHeight, cellHeight,hbl(i,j))

call CVMix_coeffs_conv(Mdiff_out=CS%kv_conv(i,j,:), &
Tdiff_out=CS%kd_conv(i,j,:), &
Nsqr=CS%N2(i,j,:), &
dens=rho_1d(:), &
dens_lwr=rho_lwr(:), &
nlev=G%ke, &
max_nlev=G%ke, &
OBL_ind=kOBL)

! Do not apply mixing due to convection within the boundary layer
do k=1,kOBL
CS%kv_conv(i,j,k) = 0.0
CS%kd_conv(i,j,k) = 0.0
enddo

enddo
enddo

if (CS%debug) then
call hchksum(CS%N2, "MOM_CVMix_conv: N2",G%HI,haloshift=0)
call hchksum(CS%kd_conv, "MOM_CVMix_conv: kd_conv",G%HI,haloshift=0)
call hchksum(CS%kv_conv, "MOM_CVMix_conv: kv_conv",G%HI,haloshift=0)
endif

! send diagnostics to post_data
if (CS%id_N2 > 0) call post_data(CS%id_N2, CS%N2, CS%diag)
if (CS%id_kd_conv > 0) call post_data(CS%id_kd_conv, CS%kd_conv, CS%diag)
if (CS%id_kv_conv > 0) call post_data(CS%id_kv_conv, CS%kv_conv, CS%diag)

end subroutine calculate_CVMix_conv

!> Reads the parameter "USE_CVMix_CONVECTION" and returns state.
!! This function allows other modules to know whether this parameterization will
!! be used without needing to duplicate the log entry.
logical function CVMix_conv_is_used(param_file)
type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
call get_param(param_file, mdl, "USE_CVMix_CONVECTION", CVMix_conv_is_used, &
default=.false., do_not_log = .true.)

end function CVMix_conv_is_used

!> Clear pointers and dealocate memory
subroutine CVMix_conv_end(CS)
type(CVMix_conv_cs), pointer :: CS ! Control structure

deallocate(CS%N2)
deallocate(CS%kd_conv)
deallocate(CS%kv_conv)
deallocate(CS)

end subroutine CVMix_conv_end


end module MOM_CVMix_conv
Loading

0 comments on commit 3c70582

Please sign in to comment.