From 1705c024e02e0a1928de5c8901c9c9bbf19c8d46 Mon Sep 17 00:00:00 2001 From: Luke Van Roekel Date: Wed, 10 Jul 2024 15:32:48 -0500 Subject: [PATCH 1/5] Fixes to KPP interface Two portions of the KPP interface had bugs 1. the computation of surface buoyancy forcing and friction velocity were incorrect. The surfaceBuoyancyForcing subtracted the wrong values from the heat fluxes from thickness fluxes. The surfaceFrictionVelocity had incorrect interpolation to centers 2. the KPP interface for matchBoth is fixed and the indexing is fixed for matching interior to boundary layer diffusivities. --- .../src/shared/mpas_ocn_diagnostics.F | 26 +++---------------- .../src/shared/mpas_ocn_vmix_cvmix.F | 13 +++++----- 2 files changed, 11 insertions(+), 28 deletions(-) diff --git a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F index 52ef08267fc1..1b5df485a20d 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F @@ -3303,8 +3303,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & evapTemperatureFlux, & icebergTemperatureFlux, & seaIceTemperatureFlux, & - surfaceStress, & - surfaceStressMagnitude + sfcStressMag, & real (kind=RKIND), dimension(:,:), pointer :: & layerThickness, &! layer thickness @@ -3381,10 +3380,8 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & surfaceThicknessFluxRunoff) call mpas_pool_get_array(forcingPool, 'penetrativeTemperatureFlux', & penetrativeTemperatureFlux) - call mpas_pool_get_array(forcingPool, 'surfaceStress', & - surfaceStress) call mpas_pool_get_array(forcingPool, 'surfaceStressMagnitude', & - surfaceStressMagnitude) + sfcStressMag) call mpas_pool_get_array(forcingPool, 'rainTemperatureFlux', & rainTemperatureFlux) call mpas_pool_get_array(forcingPool, 'evapTemperatureFlux', & @@ -3477,9 +3474,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & - fracAbsorbed* (rainTemperatureFlux(iCell) + & evapTemperatureFlux(iCell) + & seaIceTemperatureFlux(iCell) + & - icebergTemperatureFlux(iCell)) & - - fracAbsorbedRunoff* & - activeTracersSurfaceFluxRunoff(indexTempFlux,iCell) + icebergTemperatureFlux(iCell)) nonLocalSurfaceTracerFlux(indexSaltFlux,iCell) = & activeTracersSurfaceFlux(indexSaltFlux,iCell) & @@ -3500,21 +3495,8 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & surfaceBuoyancyForcing(iCell) = surfaceBuoyancyForcing(iCell)* & gravity - ! compute magnitude of surface stress - sumSurfaceStressSquared = 0.0_RKIND - do i = 1, nEdgesOnCell(iCell) - iEdge = edgesOnCell(i, iCell) - sumSurfaceStressSquared = sumSurfaceStressSquared & - + edgeAreaFractionOfCell(i,iCell)* & - surfaceStress(iEdge)**2 - enddo - - ! NOTE that the factor of 2 is from averaging dot products - ! to cell centers on a C-grid - surfaceStressMagnitude(iCell) = & - sqrt(2.0_RKIND * sumSurfaceStressSquared) surfaceFrictionVelocity(iCell) = & - sqrt(surfaceStressMagnitude(iCell) / rho_sw) + sqrt(sfcStressMag(iCell) / rho_sw) enddo !$omp end do diff --git a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F index fcede65e420c..c9a9f52c6f0a 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F @@ -65,7 +65,7 @@ module ocn_vmix_cvmix type(cvmix_tidal_params_type) :: cvmix_tidal_params logical :: cvmixOn, cvmixConvectionOn, cvmixKPPOn - real (kind=RKIND) :: backgroundVisc, backgroundDiff + real (kind=RKIND) :: multVal, backgroundVisc, backgroundDiff integer :: cvmixBackgroundChoice ! user choice of cvmix background scheme @@ -725,9 +725,9 @@ subroutine ocn_vmix_coefs_cvmix_build(meshPool, statePool, forcingPool, err, tim ! intent out of BoundaryLayerDepth is boundary layer depth measured in meters and vertical index indexBoundaryLayerDepth(iCell) = cvmix_variables % kOBL_depth - do k = minLevelCell(iCell), maxLevelCell(iCell) + 1 - vertViscTopOfCell(k, iCell) = cvmix_variables % Mdiff_iface(k) - vertDiffTopOfCell(k, iCell) = cvmix_variables % Tdiff_iface(k) + do k = minLevelCell(iCell), floor(indexBoundaryLayerDepth(iCell)) + 1 + vertViscTopOfCell(k, iCell) = multVal*vertViscTopOfCell(k,iCell) + cvmix_variables % Mdiff_iface(k) + vertDiffTopOfCell(k, iCell) = multVal*vertDiffTopOfCell(k,iCell) + cvmix_variables % Tdiff_iface(k) end do ! store non-local flux terms @@ -751,7 +751,7 @@ subroutine ocn_vmix_coefs_cvmix_build(meshPool, statePool, forcingPool, err, tim ! add convective mixing to vertical viscosity/diffusivity ! if using KPP, then do not apply convective mixing within the ocean boundary layer if(config_use_cvmix_kpp) then - do k = ceiling(indexBoundaryLayerDepth(iCell)) + 1, maxLevelCell(iCell) + do k = floor(indexBoundaryLayerDepth(iCell)) + 1, maxLevelCell(iCell) vertViscTopOfCell(k,iCell) = vertViscTopOfCell(k,iCell) + cvmix_variables % Mdiff_iface(k) vertDiffTopOfCell(k,iCell) = vertDiffTopOfCell(k,iCell) + cvmix_variables % Tdiff_iface(k) enddo @@ -1081,13 +1081,14 @@ subroutine ocn_vmix_cvmix_init(domain,err)!{{{ call cvmix_init_kpp ( & ri_crit = config_cvmix_kpp_criticalBulkRichardsonNumber, & interp_type = config_cvmix_kpp_interpolationOMLType, & - interp_type2 = config_cvmix_kpp_interpolationOMLType, & + interp_type2 = 'LMD94', & lEkman = config_cvmix_kpp_EkmanOBL, & lMonOb = config_cvmix_kpp_MonObOBL, & MatchTechnique = config_cvmix_kpp_matching, & surf_layer_ext = config_cvmix_kpp_surface_layer_extent, & langmuir_mixing_str = config_cvmix_kpp_langmuir_mixing_opt, & langmuir_entrainment_str = config_cvmix_kpp_langmuir_entrainment_opt, & + lnoDGat1 = .true., & lenhanced_diff = config_cvmix_kpp_use_enhanced_diff) endif From 8ad3fa20730ec1b15329ad89a1cd14dfc3773242 Mon Sep 17 00:00:00 2001 From: Luke Van Roekel Date: Thu, 11 Jul 2024 15:47:40 -0500 Subject: [PATCH 2/5] updates to the KPP buoyancy flux calculation --- .../mpas-ocean/src/shared/mpas_ocn_diagnostics.F | 13 +------------ 1 file changed, 1 insertion(+), 12 deletions(-) diff --git a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F index 1b5df485a20d..5489f3ad0422 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F @@ -3382,14 +3382,6 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & penetrativeTemperatureFlux) call mpas_pool_get_array(forcingPool, 'surfaceStressMagnitude', & sfcStressMag) - call mpas_pool_get_array(forcingPool, 'rainTemperatureFlux', & - rainTemperatureFlux) - call mpas_pool_get_array(forcingPool, 'evapTemperatureFlux', & - evapTemperatureFlux) - call mpas_pool_get_array(forcingPool, 'seaIceTemperatureFlux', & - seaIceTemperatureFlux) - call mpas_pool_get_array(forcingPool, 'icebergTemperatureFlux', & - icebergTemperatureFlux) ! allocate scratch space displaced density computation ncells = nCellsAll @@ -3471,10 +3463,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & activeTracersSurfaceFlux(indexTempFlux,iCell) & + penetrativeTemperatureFlux(iCell) & - penetrativeTemperatureFluxOBL(iCell) & - - fracAbsorbed* (rainTemperatureFlux(iCell) + & - evapTemperatureFlux(iCell) + & - seaIceTemperatureFlux(iCell) + & - icebergTemperatureFlux(iCell)) + - fracAbsorbed*activeTracers(indexTmpFlux,iCell) nonLocalSurfaceTracerFlux(indexSaltFlux,iCell) = & activeTracersSurfaceFlux(indexSaltFlux,iCell) & From 295480d6b4a52685911d3d14526a13dda17e081e Mon Sep 17 00:00:00 2001 From: vanroekel Date: Thu, 11 Jul 2024 20:21:49 -0700 Subject: [PATCH 3/5] Fixup to mpas_vmix_cvmix --- components/mpas-ocean/src/shared/Makefile | 2 +- .../mpas-ocean/src/shared/mpas_ocn_tendency.F | 2 -- components/mpas-ocean/src/shared/mpas_ocn_vmix.F | 6 ++++++ .../mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F | 15 ++++++++++----- 4 files changed, 17 insertions(+), 8 deletions(-) diff --git a/components/mpas-ocean/src/shared/Makefile b/components/mpas-ocean/src/shared/Makefile index 35017551c754..9199d2469d4a 100644 --- a/components/mpas-ocean/src/shared/Makefile +++ b/components/mpas-ocean/src/shared/Makefile @@ -157,7 +157,7 @@ mpas_ocn_tracer_short_wave_absorption_variable.o: mpas_ocn_constants.o mpas_ocn_ mpas_ocn_tracer_short_wave_absorption_jerlov.o: mpas_ocn_constants.o mpas_ocn_config.o -mpas_ocn_vmix.o: mpas_ocn_vmix_cvmix.o mpas_ocn_vmix_coefs_redi.o mpas_ocn_constants.o mpas_ocn_config.o mpas_ocn_diagnostics_variables.o mpas_ocn_vmix_gotm.o +mpas_ocn_vmix.o: mpas_ocn_diagnostics.o mpas_ocn_vmix_cvmix.o mpas_ocn_vmix_coefs_redi.o mpas_ocn_constants.o mpas_ocn_config.o mpas_ocn_diagnostics_variables.o mpas_ocn_vmix_gotm.o mpas_ocn_vmix_cvmix.o: mpas_ocn_constants.o mpas_ocn_config.o mpas_ocn_diagnostics_variables.o mpas_ocn_mesh.o mpas_ocn_stokes_drift.o diff --git a/components/mpas-ocean/src/shared/mpas_ocn_tendency.F b/components/mpas-ocean/src/shared/mpas_ocn_tendency.F index a16e52f14499..b9138a1ce88b 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_tendency.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_tendency.F @@ -1246,8 +1246,6 @@ subroutine ocn_tend_tracer(tendPool, statePool, forcingPool, & if (config_use_cvmix_kpp) then call mpas_timer_start("non-local flux from KPP") if (.not. config_cvmix_kpp_nonlocal_with_implicit_mix) then - call ocn_compute_KPP_input_fields(statePool, forcingPool,& - meshPool, timeLevel) if (computeBudgets) then !$omp parallel diff --git a/components/mpas-ocean/src/shared/mpas_ocn_vmix.F b/components/mpas-ocean/src/shared/mpas_ocn_vmix.F index 68ef9d82a6bf..2d27fbd2f439 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_vmix.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_vmix.F @@ -35,6 +35,7 @@ module ocn_vmix use ocn_vmix_coefs_redi use ocn_diagnostics_variables use ocn_surface_land_ice_fluxes + use ocn_diagnostics implicit none private @@ -206,6 +207,11 @@ subroutine ocn_vmix_coefs(meshPool, statePool, forcingPool, scratchPool, err, ti !$acc update host(vertViscTopOfEdge, vertDiffTopOfCell) #endif + if(config_use_cvmix_kpp .or. config_use_gotm) then + call ocn_compute_mixing_input_fields(statePool, forcingPool,& + meshPool, timeLevel) + endif + call ocn_vmix_coefs_cvmix_build(meshPool, statePool, forcingPool, err1, timeLevel) call ocn_vmix_coefs_redi_build(meshPool, statePool, err2, timeLevel) call ocn_vmix_coefs_gotm_build(statePool, forcingPool, err3, timeLevel) diff --git a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F index c9a9f52c6f0a..8e7694849ad9 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F @@ -1042,11 +1042,7 @@ subroutine ocn_vmix_cvmix_init(domain,err)!{{{ ! initialize KPP boundary layer scheme ! if (config_use_cvmix_kpp) then - if(config_cvmix_kpp_matching.eq."MatchBoth") then - call mpas_log_write( & - "Use of option MatchBoth is discouraged, use SimpleShapes instead", & - MPAS_LOG_WARN) - elseif(.not. config_cvmix_kpp_matching.eq."SimpleShapes") then + if(.not. config_cvmix_kpp_matching.eq."SimpleShapes" .or. .not. config_cvmix_kpp_matching.eq."MatchBoth") then call mpas_log_write( & "Unknown value for config_cvmix_kpp_matching., supported values are:" // & " SimpleShapes or MatchBoth", & @@ -1055,6 +1051,15 @@ subroutine ocn_vmix_cvmix_init(domain,err)!{{{ return endif + !multVal sets the parameter for adding diffusivity/viscosity in the + !boundary layer + + if (config_cvmix_kpp_matching.eq."SimpleShapes") then + multVal = 1.0_RKIND + else + multVal = 0.0_RKIND + endif + if (trim(config_cvmix_kpp_langmuir_mixing_opt) .ne. "NONE" .and. & trim(config_cvmix_kpp_langmuir_mixing_opt) .ne. "LWF16" .and. & trim(config_cvmix_kpp_langmuir_mixing_opt) .ne. "RWHGK16") then From 0e861cfea5ee66162369998bb1d76aa916850a5f Mon Sep 17 00:00:00 2001 From: vanroekel Date: Fri, 12 Jul 2024 07:05:55 -0700 Subject: [PATCH 4/5] compile time fixes --- .../src/shared/mpas_ocn_diagnostics.F | 22 +++++++++---------- .../src/shared/mpas_ocn_vmix_cvmix.F | 2 +- 2 files changed, 12 insertions(+), 12 deletions(-) diff --git a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F index 5489f3ad0422..74922640630c 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F @@ -64,7 +64,7 @@ module ocn_diagnostics ocn_fuperp, & ocn_filter_btr_mode_tend_vel, & ocn_reconstruct_eddy_vectors, & - ocn_compute_kpp_input_fields, & + ocn_compute_mixing_input_fields, & ocn_validate_state, & ocn_build_log_filename, & ocn_diagnostics_init @@ -3241,12 +3241,12 @@ end subroutine ocn_filter_btr_mode_tend_vel!}}} !*********************************************************************** ! -! routine ocn_compute_KPP_input_fields +! routine ocn_compute_mixing_input_fields ! !> \brief -!> Compute fields necessary to drive the CVMix KPP module -!> \author Todd Ringler -!> \date 20 August 2013 +!> Compute fields necessary to drive the CVMix KPP and gotm modules +!> \author Todd Ringler, Luke Van Roekel +!> \date 11 July 2024 !> \details !> CVMix/KPP requires the following fields as intent(in): !> surfaceBuoyancyForcing @@ -3255,7 +3255,7 @@ end subroutine ocn_filter_btr_mode_tend_vel!}}} ! !----------------------------------------------------------------------- - subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & + subroutine ocn_compute_mixing_input_fields(statePool, forcingPool, & meshPool, timeLevelIn)!{{{ !----------------------------------------------------------------- @@ -3303,7 +3303,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & evapTemperatureFlux, & icebergTemperatureFlux, & seaIceTemperatureFlux, & - sfcStressMag, & + sfcStressMag real (kind=RKIND), dimension(:,:), pointer :: & layerThickness, &! layer thickness @@ -3335,7 +3335,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & !----------------------------------------------------------------- ! Begin code - call mpas_timer_start('KPP input fields') + call mpas_timer_start('Mixing input fields') if (present(timeLevelIn)) then timeLevel = timeLevelIn @@ -3463,7 +3463,7 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & activeTracersSurfaceFlux(indexTempFlux,iCell) & + penetrativeTemperatureFlux(iCell) & - penetrativeTemperatureFluxOBL(iCell) & - - fracAbsorbed*activeTracers(indexTmpFlux,iCell) + - fracAbsorbed*activeTracers(indexTempFlux,kmin,iCell) nonLocalSurfaceTracerFlux(indexSaltFlux,iCell) = & activeTracersSurfaceFlux(indexSaltFlux,iCell) & @@ -3496,11 +3496,11 @@ subroutine ocn_compute_KPP_input_fields(statePool, forcingPool, & salineContractionCoeff, & densitySurfaceDisplaced) - call mpas_timer_stop('KPP input fields') + call mpas_timer_stop('Mixing input fields') !------------------------------------------------------------------- - end subroutine ocn_compute_KPP_input_fields!}}} + end subroutine ocn_compute_mixing_input_fields!}}} !*********************************************************************** ! diff --git a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F index 8e7694849ad9..127de4ee87a5 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F @@ -1042,7 +1042,7 @@ subroutine ocn_vmix_cvmix_init(domain,err)!{{{ ! initialize KPP boundary layer scheme ! if (config_use_cvmix_kpp) then - if(.not. config_cvmix_kpp_matching.eq."SimpleShapes" .or. .not. config_cvmix_kpp_matching.eq."MatchBoth") then + if(.not. config_cvmix_kpp_matching.eq."SimpleShapes" .and. .not. config_cvmix_kpp_matching.eq."MatchBoth") then call mpas_log_write( & "Unknown value for config_cvmix_kpp_matching., supported values are:" // & " SimpleShapes or MatchBoth", & From 0e9cd30e17dc88d81945cc78891d3a875fe1a52f Mon Sep 17 00:00:00 2001 From: vanroekel Date: Fri, 12 Jul 2024 14:15:38 -0700 Subject: [PATCH 5/5] runtime fixes --- components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F | 3 ++- components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F | 6 +++--- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F index 74922640630c..26cdab567510 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_diagnostics.F @@ -3463,7 +3463,8 @@ subroutine ocn_compute_mixing_input_fields(statePool, forcingPool, & activeTracersSurfaceFlux(indexTempFlux,iCell) & + penetrativeTemperatureFlux(iCell) & - penetrativeTemperatureFluxOBL(iCell) & - - fracAbsorbed*activeTracers(indexTempFlux,kmin,iCell) + - fracAbsorbed*surfaceThicknessFlux(iCell) * & + activeTracers(indexTempFlux,kmin,iCell) nonLocalSurfaceTracerFlux(indexSaltFlux,iCell) = & activeTracersSurfaceFlux(indexSaltFlux,iCell) & diff --git a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F index 127de4ee87a5..7d1a477b7b9d 100644 --- a/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F +++ b/components/mpas-ocean/src/shared/mpas_ocn_vmix_cvmix.F @@ -725,9 +725,9 @@ subroutine ocn_vmix_coefs_cvmix_build(meshPool, statePool, forcingPool, err, tim ! intent out of BoundaryLayerDepth is boundary layer depth measured in meters and vertical index indexBoundaryLayerDepth(iCell) = cvmix_variables % kOBL_depth - do k = minLevelCell(iCell), floor(indexBoundaryLayerDepth(iCell)) + 1 - vertViscTopOfCell(k, iCell) = multVal*vertViscTopOfCell(k,iCell) + cvmix_variables % Mdiff_iface(k) - vertDiffTopOfCell(k, iCell) = multVal*vertDiffTopOfCell(k,iCell) + cvmix_variables % Tdiff_iface(k) + do k = minLevelCell(iCell), maxLevelCell(iCell) !floor(indexBoundaryLayerDepth(iCell)) + 1 + vertViscTopOfCell(k, iCell) = cvmix_variables % Mdiff_iface(k) + vertDiffTopOfCell(k, iCell) = cvmix_variables % Tdiff_iface(k) end do ! store non-local flux terms