From 28268d3e632ff1a406c997c9da582a88bc21ea39 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Thu, 23 Dec 2021 09:02:02 -0700 Subject: [PATCH 1/2] fixes in MOM_wave_interface array indices --- src/user/MOM_wave_interface.F90 | 14 +++++--------- 1 file changed, 5 insertions(+), 9 deletions(-) diff --git a/src/user/MOM_wave_interface.F90 b/src/user/MOM_wave_interface.F90 index f89d26451d..340322379c 100644 --- a/src/user/MOM_wave_interface.F90 +++ b/src/user/MOM_wave_interface.F90 @@ -373,7 +373,7 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag, restar "or the model will fail.",units='', default=1) allocate( CS%WaveNum_Cen(CS%NumBands), source=0.0 ) allocate( CS%STKx0(G%isdB:G%iedB,G%jsd:G%jed,CS%NumBands), source=0.0 ) - allocate( CS%STKy0(G%isdB:G%iedB,G%jsd:G%jed,CS%NumBands), source=0.0 ) + allocate( CS%STKy0(G%isd:G%ied,G%jsdB:G%jedB,CS%NumBands), source=0.0 ) CS%PartitionMode = 0 call get_param(param_file, mdl, "SURFBAND_WAVENUMBERS", CS%WaveNum_Cen, & "Central wavenumbers for surface Stokes drift bands.", & @@ -442,9 +442,7 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag, restar ! a. Stokes driftProfiles if (CS%Stokes_DDT) then allocate(CS%ddt_Us_x(G%isdB:G%IedB,G%jsd:G%jed,G%ke), source=0.0) - CS%ddt_Us_x(:,:,:) = 0.0 allocate(CS%ddt_Us_y(G%isd:G%Ied,G%jsdB:G%jedB,G%ke), source=0.0) - CS%ddt_Us_y(:,:,:) = 0.0 endif ! b. Surface Values allocate(CS%US0_x(G%isdB:G%iedB,G%jsd:G%jed), source=0.0) @@ -825,7 +823,7 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar, dt) enddo CS%DHH85_is_set = .true. endif - call pass_vector(CS%Us_x(:,:),CS%Us_y(:,:), G%Domain) + call pass_vector(CS%Us_x(:,:,:),CS%Us_y(:,:,:), G%Domain) else! Keep this else, fallback to 0 Stokes drift do kk= 1,GV%ke do jj = G%jsd,G%jed @@ -1520,7 +1518,7 @@ subroutine CoriolisStokes(G, GV, dt, h, u, v, WAVES, US) integer :: i,j,k do k = 1, GV%ke -do j = G%jsc, G%jec + do j = G%jsc, G%jec do I = G%iscB, G%iecB DVel = 0.25*(WAVES%us_y(i,j+1,k)+WAVES%us_y(i-1,j+1,k))*G%CoriolisBu(i,j+1) + & 0.25*(WAVES%us_y(i,j,k)+WAVES%us_y(i-1,j,k))*G%CoriolisBu(i,j) @@ -1788,10 +1786,8 @@ subroutine waves_register_restarts(CS, HI, GV, param_file, restart_CSp) if (.not.(use_waves .or. StatisticalWaves)) return ! Allocate wave fields needed for restart file - allocate(CS%Us_x(HI%isdB:HI%IedB,HI%jsd:HI%jed,GV%ke)) - CS%Us_x(:,:,:) = 0.0 - allocate(CS%Us_y(HI%isd:HI%Ied,HI%jsdB:HI%jedB,GV%ke)) - CS%Us_y(:,:,:) = 0.0 + allocate(CS%Us_x(HI%isdB:HI%IedB,HI%jsd:HI%jed,GV%ke), source=0.0) + allocate(CS%Us_y(HI%isd:HI%Ied,HI%jsdB:HI%jedB,GV%ke), source=0.0) call get_param(param_file,mdl,"WAVE_METHOD",wave_method_str, do_not_log=.true., default=NULL_STRING) From 4f592f2988a168152fc8372bda09ff28250b776e Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Thu, 23 Dec 2021 13:24:01 -0700 Subject: [PATCH 2/2] More changes for surfbands wave coupling: - for ustkb and vstkb, call pass_var, instead of pass_vector because they are on h grd. -remove unused ustk0 and vstk0 arrays. - correct the order of wave-related get param & allocate_forcing calls in nuopc cap. - add lamult flag to allocate_forcing_by_group because lamult should only be allocated if wave_method=="EFACTOR". --- config_src/drivers/nuopc_cap/mom_cap.F90 | 12 +++------- .../nuopc_cap/mom_ocean_model_nuopc.F90 | 7 +++--- .../nuopc_cap/mom_surface_forcing_nuopc.F90 | 10 ++------ src/core/MOM_forcing_type.F90 | 23 +++++-------------- 4 files changed, 15 insertions(+), 37 deletions(-) diff --git a/config_src/drivers/nuopc_cap/mom_cap.F90 b/config_src/drivers/nuopc_cap/mom_cap.F90 index 409f3bd166..8d75c814f7 100644 --- a/config_src/drivers/nuopc_cap/mom_cap.F90 +++ b/config_src/drivers/nuopc_cap/mom_cap.F90 @@ -714,16 +714,10 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) Ice_ocean_boundary%lamult = 0.0 else if (wave_method == "SURFACE_BANDS") then call query_ocean_state(ocean_state, NumWaveBands=Ice_ocean_boundary%num_stk_bands) - 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 + allocate(Ice_ocean_boundary%ustkb(isc:iec,jsc:jec,Ice_ocean_boundary%num_stk_bands), source=0.0) + allocate(Ice_ocean_boundary%vstkb(isc:iec,jsc:jec,Ice_ocean_boundary%num_stk_bands), source=0.0) + allocate(Ice_ocean_boundary%stk_wavenumbers(Ice_ocean_boundary%num_stk_bands), source=0.0) call query_ocean_state(ocean_state, WaveNumbers=Ice_ocean_boundary%stk_wavenumbers, unscale=.true.) - Ice_ocean_boundary%ustkb = 0.0 - Ice_ocean_boundary%vstkb = 0.0 else call MOM_error(FATAL, "Unsupported WAVE_METHOD encountered in NUOPC cap.") endif 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 5e0242be13..c5ac1d44b4 100644 --- a/config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90 +++ b/config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90 @@ -373,6 +373,8 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i call get_param(param_file, mdl, "USE_CFC_CAP", use_CFC, & default=.false., do_not_log=.true.) + call get_param(param_file, mdl, "USE_WAVES", OS%Use_Waves, & + "If true, enables surface wave modules.", default=.false.) ! Consider using a run-time flag to determine whether to do the diagnostic ! vertical integrals, since the related 3-d sums are not negligible in cost. @@ -393,12 +395,11 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i call allocate_forcing_type(OS%grid, OS%fluxes, shelf=.true.) endif - call allocate_forcing_type(OS%grid, OS%fluxes, waves=.true.) - 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.) endif + call allocate_forcing_type(OS%grid, OS%fluxes, waves=.true., lamult=(trim(OS%wave_method)=="EFACTOR")) + ! MOM_wave_interface_init is called regardless of the value of USE_WAVES because ! it also initializes statistical waves. call MOM_wave_interface_init(OS%Time, OS%grid, OS%GV, OS%US, param_file, OS%Waves, OS%diag, OS%restart_CSp) diff --git a/config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90 b/config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90 index c704214930..421ada487f 100644 --- a/config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90 +++ b/config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90 @@ -193,8 +193,6 @@ module MOM_surface_forcing_nuopc !! for divergence damping, as determined !! outside of the ocean model in [m3/s] real, pointer, dimension(:,:) :: lamult => NULL() !< Langmuir enhancement factor [nondim] - real, pointer, dimension(:,:) :: ustk0 => NULL() !< Surface Stokes drift, zonal [m/s] - real, pointer, dimension(:,:) :: vstk0 => NULL() !< Surface Stokes drift, meridional [m/s] real, pointer, dimension(:) :: stk_wavenumbers => NULL() !< The central wave number of Stokes bands [rad/m] real, pointer, dimension(:,:,:) :: ustkb => NULL() !< Stokes Drift spectrum, zonal [m/s] !! Horizontal - u points @@ -893,17 +891,13 @@ subroutine convert_IOB_to_forces(IOB, forces, index_bounds, Time, G, US, CS) if ( associated(IOB%ustkb) ) then forces%stk_wavenumbers(:) = IOB%stk_wavenumbers - do j=js,je; do i=is,ie - forces%ustk0(i,j) = IOB%ustk0(i-I0,j-J0) ! How to be careful here that the domains are right? - forces%vstk0(i,j) = IOB%vstk0(i-I0,j-J0) - enddo ; enddo - call pass_vector(forces%ustk0,forces%vstk0, G%domain ) do istk = 1,IOB%num_stk_bands do j=js,je; do i=is,ie forces%ustkb(i,j,istk) = IOB%ustkb(i-I0,j-J0,istk) forces%vstkb(i,j,istk) = IOB%vstkb(i-I0,j-J0,istk) enddo; enddo - call pass_vector(forces%ustkb(:,:,istk),forces%vstkb(:,:,istk), G%domain ) + call pass_var(forces%ustkb(:,:,istk), G%domain ) + call pass_var(forces%vstkb(:,:,istk), G%domain ) enddo endif diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90 index 2429ce9d2d..58041dbdac 100644 --- a/src/core/MOM_forcing_type.F90 +++ b/src/core/MOM_forcing_type.F90 @@ -256,9 +256,6 @@ module MOM_forcing_type logical :: accumulate_rigidity = .false. !< If true, the rigidity due to various types of !! ice needs to be accumulated, and the rigidity explicitly !! reset to zero at the driver level when appropriate. - real, pointer, dimension(:,:) :: & - ustk0 => NULL(), & !< Surface Stokes drift, zonal [m/s] - vstk0 => NULL() !< Surface Stokes drift, meridional [m/s] real, pointer, dimension(:) :: & stk_wavenumbers => NULL() !< The central wave number of Stokes bands [rad/m] real, pointer, dimension(:,:,:) :: & @@ -2953,7 +2950,7 @@ end subroutine forcing_diagnostics !> Conditionally allocate fields within the forcing type subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, & - shelf, iceberg, salt, fix_accum_bug, cfc, waves) + shelf, iceberg, salt, fix_accum_bug, cfc, waves, lamult) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields logical, optional, intent(in) :: water !< If present and true, allocate water fluxes @@ -2967,6 +2964,7 @@ subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, & !! accumulation of ustar_gustless logical, optional, intent(in) :: cfc !< If present and true, allocate cfc fluxes logical, optional, intent(in) :: waves !< If present and true, allocate wave fields + logical, optional, intent(in) :: lamult !< If present and true, allocate langmuir enhancement factor ! Local variables integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB @@ -3030,7 +3028,7 @@ subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, & !These fields should only on allocated when wave coupling is activated. call myAlloc(fluxes%ice_fraction,isd,ied,jsd,jed, waves) - call myAlloc(fluxes%lamult,isd,ied,jsd,jed, waves) + call myAlloc(fluxes%lamult,isd,ied,jsd,jed, lamult) if (present(fix_accum_bug)) fluxes%gustless_accum_bug = .not.fix_accum_bug end subroutine allocate_forcing_by_group @@ -3125,8 +3123,6 @@ subroutine allocate_mech_forcing_by_group(G, forces, stress, ustar, shelf, & call myAlloc(forces%mass_berg,isd,ied,jsd,jed, iceberg) !These fields should only be allocated when waves - call myAlloc(forces%ustk0,isd,ied,jsd,jed, waves) - call myAlloc(forces%vstk0,isd,ied,jsd,jed, waves) if (present(waves)) then; if (waves) then; if (.not. present(num_stk_bands)) then call MOM_error(FATAL,"Requested to & @@ -3134,20 +3130,13 @@ subroutine allocate_mech_forcing_by_group(G, forces, stress, ustar, shelf, & endif if (num_stk_bands > 0) then if (.not.associated(forces%ustkb)) then - allocate(forces%stk_wavenumbers(num_stk_bands)) - forces%stk_wavenumbers(:) = 0.0 - allocate(forces%ustkb(isd:ied,jsd:jed,num_stk_bands)) - forces%ustkb(isd:ied,jsd:jed,:) = 0.0 + allocate(forces%stk_wavenumbers(num_stk_bands), source=0.0) + allocate(forces%ustkb(isd:ied,jsd:jed,num_stk_bands), source=0.0) + allocate(forces%vstkb(isd:ied,jsd:jed,num_stk_bands), source=0.0) endif endif endif ; endif - - if (present(waves)) then; if (waves) then; if (.not.associated(forces%vstkb)) then - allocate(forces%vstkb(isd:ied,jsd:jed,num_stk_bands)) - forces%vstkb(isd:ied,jsd:jed,:) = 0.0 - endif ; endif ; endif - end subroutine allocate_mech_forcing_by_group