From 0d591560846c96bb3b667e2f59bc7324bfe106ed Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Sun, 13 Jun 2021 22:56:47 -0600 Subject: [PATCH] introduce a new wave coupling method: VR12-MA (not fully implemented yet) --- config_src/drivers/nuopc_cap/mom_cap.F90 | 49 ++++++++++--------- .../nuopc_cap/mom_ocean_model_nuopc.F90 | 16 ++++-- src/user/MOM_wave_interface.F90 | 14 +++++- 3 files changed, 51 insertions(+), 28 deletions(-) diff --git a/config_src/drivers/nuopc_cap/mom_cap.F90 b/config_src/drivers/nuopc_cap/mom_cap.F90 index 2d79674606..8d383c3d3f 100644 --- a/config_src/drivers/nuopc_cap/mom_cap.F90 +++ b/config_src/drivers/nuopc_cap/mom_cap.F90 @@ -696,17 +696,21 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) Ice_ocean_boundary%frunoff = 0.0 if (ocean_state%use_waves) then - Ice_ocean_boundary%num_stk_bands=ocean_state%Waves%NumBands - allocate ( Ice_ocean_boundary% ustk0 (isc:iec,jsc:jec), & - Ice_ocean_boundary% vstk0 (isc:iec,jsc:jec), & - Ice_ocean_boundary% ustkb (isc:iec,jsc:jec,Ice_ocean_boundary%num_stk_bands), & - Ice_ocean_boundary% vstkb (isc:iec,jsc:jec,Ice_ocean_boundary%num_stk_bands), & - Ice_ocean_boundary%stk_wavenumbers (Ice_ocean_boundary%num_stk_bands)) - Ice_ocean_boundary%ustk0 = 0.0 - Ice_ocean_boundary%vstk0 = 0.0 - Ice_ocean_boundary%stk_wavenumbers = ocean_state%Waves%WaveNum_Cen - Ice_ocean_boundary%ustkb = 0.0 - Ice_ocean_boundary%vstkb = 0.0 + if (ocean_state%wave_method == "VR12-MA") then + continue + else + Ice_ocean_boundary%num_stk_bands=ocean_state%Waves%NumBands + allocate ( Ice_ocean_boundary% ustk0 (isc:iec,jsc:jec), & + Ice_ocean_boundary% vstk0 (isc:iec,jsc:jec), & + Ice_ocean_boundary% ustkb (isc:iec,jsc:jec,Ice_ocean_boundary%num_stk_bands), & + Ice_ocean_boundary% vstkb (isc:iec,jsc:jec,Ice_ocean_boundary%num_stk_bands), & + Ice_ocean_boundary%stk_wavenumbers (Ice_ocean_boundary%num_stk_bands)) + Ice_ocean_boundary%ustk0 = 0.0 + Ice_ocean_boundary%vstk0 = 0.0 + Ice_ocean_boundary%stk_wavenumbers = ocean_state%Waves%WaveNum_Cen + Ice_ocean_boundary%ustkb = 0.0 + Ice_ocean_boundary%vstkb = 0.0 + endif endif ocean_internalstate%ptr%ocean_state_type_ptr => ocean_state @@ -723,9 +727,6 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) !call fld_list_add(fldsToOcn_num, fldsToOcn, "Sw_ustokes" , "will provide") !call fld_list_add(fldsToOcn_num, fldsToOcn, "Sw_vstokes" , "will provide") !call fld_list_add(fldsToOcn_num, fldsToOcn, "Sw_hstokes" , "will provide") - !call fld_list_add(fldsToOcn_num, fldsToOcn, "Fioi_melth" , "will provide") - !call fld_list_add(fldsToOcn_num, fldsToOcn, "Fioi_meltw" , "will provide") - !call fld_list_add(fldsFrOcn_num, fldsFrOcn, "So_fswpen" , "will provide") else !call fld_list_add(fldsToOcn_num, fldsToOcn, "mass_of_overlying_sea_ice" , "will provide") !call fld_list_add(fldsFrOcn_num, fldsFrOcn, "sea_lev" , "will provide") @@ -753,15 +754,19 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) !call fld_list_add(fldsToOcn_num, fldsToOcn, "mean_runoff_heat_flx" , "will provide") !call fld_list_add(fldsToOcn_num, fldsToOcn, "mean_calving_heat_flx" , "will provide") if (ocean_state%use_waves) then - if (Ice_ocean_boundary%num_stk_bands > 3) then - call MOM_error(FATAL, "Number of Stokes Bands > 3, NUOPC cap not set up for this") + if (ocean_state%wave_method == "VR12-MA") then + continue + else + if (Ice_ocean_boundary%num_stk_bands > 3) then + call MOM_error(FATAL, "Number of Stokes Bands > 3, NUOPC cap not set up for this") + endif + call fld_list_add(fldsToOcn_num, fldsToOcn, "eastward_partitioned_stokes_drift_1" , "will provide") + call fld_list_add(fldsToOcn_num, fldsToOcn, "northward_partitioned_stokes_drift_1", "will provide") + call fld_list_add(fldsToOcn_num, fldsToOcn, "eastward_partitioned_stokes_drift_2" , "will provide") + call fld_list_add(fldsToOcn_num, fldsToOcn, "northward_partitioned_stokes_drift_2", "will provide") + call fld_list_add(fldsToOcn_num, fldsToOcn, "eastward_partitioned_stokes_drift_3" , "will provide") + call fld_list_add(fldsToOcn_num, fldsToOcn, "northward_partitioned_stokes_drift_3", "will provide") endif - call fld_list_add(fldsToOcn_num, fldsToOcn, "eastward_partitioned_stokes_drift_1" , "will provide") - call fld_list_add(fldsToOcn_num, fldsToOcn, "northward_partitioned_stokes_drift_1", "will provide") - call fld_list_add(fldsToOcn_num, fldsToOcn, "eastward_partitioned_stokes_drift_2" , "will provide") - call fld_list_add(fldsToOcn_num, fldsToOcn, "northward_partitioned_stokes_drift_2", "will provide") - call fld_list_add(fldsToOcn_num, fldsToOcn, "eastward_partitioned_stokes_drift_3" , "will provide") - call fld_list_add(fldsToOcn_num, fldsToOcn, "northward_partitioned_stokes_drift_3", "will provide") endif !--------- export fields ------------- diff --git a/config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90 b/config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90 index 493762f4bc..7e7716d386 100644 --- a/config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90 +++ b/config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90 @@ -145,6 +145,7 @@ module MOM_ocean_model_nuopc integer :: nstep = 0 !< The number of calls to update_ocean. logical :: use_ice_shelf !< If true, the ice shelf model is enabled. logical,public :: use_waves !< If true use wave coupling. + character(len=40),public :: wave_method !< Wave coupling method. logical :: icebergs_alter_ocean !< If true, the icebergs can change ocean the !! ocean dynamics and forcing fluxes. @@ -387,10 +388,15 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i call get_param(param_file, mdl, "USE_WAVES", OS%Use_Waves, & "If true, enables surface wave modules.", default=.false.) if (OS%use_waves) then + call get_param(param_file, mdl, "WAVE_METHOD", OS%wave_method, default="EMPTY", do_not_log=.true.) call MOM_wave_interface_init(OS%Time, OS%grid, OS%GV, OS%US, param_file, OS%Waves, OS%diag) - call get_param(param_file,mdl,"SURFBAND_WAVENUMBERS",OS%Waves%WaveNum_Cen, & - "Central wavenumbers for surface Stokes drift bands.",units='rad/m', & - default=0.12566) + if (OS%wave_method == "VR12-MA") then + continue + else + call get_param(param_file,mdl,"SURFBAND_WAVENUMBERS",OS%Waves%WaveNum_Cen, & + "Central wavenumbers for surface Stokes drift bands.",units='rad/m', & + default=0.12566) + endif else call MOM_wave_interface_init_lite(param_file) endif @@ -575,7 +581,9 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & call set_net_mass_forcing(OS%fluxes, OS%forces, OS%grid, OS%US) if (OS%use_waves) then - call Update_Surface_Waves(OS%grid, OS%GV, OS%US, OS%time, ocean_coupling_time_step, OS%waves, OS%forces) + if (OS%wave_method /= "VR12-MA") then + call Update_Surface_Waves(OS%grid, OS%GV, OS%US, OS%time, ocean_coupling_time_step, OS%waves, OS%forces) + endif endif if (OS%nstep==0) then diff --git a/src/user/MOM_wave_interface.F90 b/src/user/MOM_wave_interface.F90 index a8e3e207a6..1ec1e9314c 100644 --- a/src/user/MOM_wave_interface.F90 +++ b/src/user/MOM_wave_interface.F90 @@ -141,6 +141,7 @@ module MOM_wave_interface !! 1 - Surface Stokes Drift Bands !! 2 - DHH85 !! 3 - LF17 + !! 4 - VR12MA !! -99 - No waves computed, but empirical Langmuir number used. !! \todo Module variable! Move into a control structure. @@ -178,7 +179,7 @@ module MOM_wave_interface !! into a control structure. ! Switches needed in import_stokes_drift integer, parameter :: TESTPROF = 0, SURFBANDS = 1, & - DHH85 = 2, LF17 = 3, NULL_WaveMethod=-99, & + DHH85 = 2, LF17 = 3, VR12MA = 4, NULL_WaveMethod=-99, & DATAOVR = 1, COUPLER = 2, INPUT = 3 ! Options For Test Prof @@ -208,6 +209,7 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag ) character*(13), parameter :: SURFBANDS_STRING = "SURFACE_BANDS" character*(5), parameter :: DHH85_STRING = "DHH85" character*(4), parameter :: LF17_STRING = "LF17" + character*(7), parameter :: VR12MA_STRING = "VR12-MA" character*(12), parameter :: DATAOVR_STRING = "DATAOVERRIDE" character*(7), parameter :: COUPLER_STRING = "COUPLER" character*(5), parameter :: INPUT_STRING = "INPUT" @@ -268,7 +270,11 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag ) " DHH85 - Uses Donelan et al. 1985 empirical \n"// & " wave spectrum with prescribed values. \n"// & " LF17 - Infers Stokes drift profile from wind \n"// & - " speed following Li and Fox-Kemper 2017.\n", & + " speed following Li and Fox-Kemper 2017.\n"// & + " VR12-MA - Applies an enhancement factor to the KPP\n"//& + " turbulent velocity scale received from \n"// & + " WW3 and is based on the surface layer \n"// & + " and projected Langmuir number. (Li 2016)\n", & units='', default=NULL_STRING) select case (TRIM(TMPSTRING1)) case (NULL_STRING)! No Waves @@ -363,6 +369,8 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag ) default=.false.) case (LF17_STRING)!Li and Fox-Kemper 17 wind-sea Langmuir number WaveMethod = LF17 + case (VR12MA_STRING)!Li and Fox-Kemper 16 + WaveMethod = VR12MA case default call MOM_error(FATAL,'Check WAVE_METHOD.') end select @@ -732,6 +740,8 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar) enddo DHH85_is_set = .true. endif + elseif (WaveMethod==VR12MA) then + return ! todo else! Keep this else, fallback to 0 Stokes drift do kk= 1,GV%ke do jj = G%jsd,G%jed