refactoring and changes to write stokes drift

profile to restart file when surfbands wave coupling mode is on.

config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90

@@ -279,7 +279,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
   call initialize_MOM(OS%Time, Time_init, param_file, OS%dirs, OS%MOM_CSp, &
                       OS%restart_CSp, Time_in, offline_tracer_mode=OS%offline_tracer_mode, &
                       input_restart_file=input_restart_file, &
-                      diag_ptr=OS%diag, count_calls=.true.)
+                      diag_ptr=OS%diag, count_calls=.true., waves_CSp=OS%Waves)
   call get_MOM_state_elements(OS%MOM_CSp, G=OS%grid, GV=OS%GV, US=OS%US, C_p=OS%C_p, &
                               C_p_scaled=OS%fluxes%C_p, use_temp=use_temperature)
 
@@ -401,7 +401,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
   endif
   ! 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)
+  call MOM_wave_interface_init(OS%Time, OS%grid, OS%GV, OS%US, param_file, OS%Waves, OS%diag, OS%restart_CSp)
 
   if (associated(OS%grid%Domain%maskmap)) then
     call initialize_ocean_public_type(OS%grid%Domain%mpp_domain, Ocean_sfc, &

src/core/MOM.F90

@@ -133,7 +133,7 @@ module MOM
 use MOM_verticalGrid,          only : verticalGrid_type, verticalGridInit, verticalGridEnd
 use MOM_verticalGrid,          only : fix_restart_scaling
 use MOM_verticalGrid,          only : get_thickness_units, get_flux_units, get_tr_flux_units
-use MOM_wave_interface,        only : wave_parameters_CS, waves_end
+use MOM_wave_interface,        only : wave_parameters_CS, waves_end, waves_register_restarts
 use MOM_wave_interface,        only : Update_Stokes_Drift
 
 ! ODA modules
@@ -1634,7 +1634,7 @@ end subroutine step_offline
 !! initializing the ocean state variables, and initializing subsidiary modules
 subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
                           Time_in, offline_tracer_mode, input_restart_file, diag_ptr, &
-                          count_calls, tracer_flow_CSp,  ice_shelf_CSp)
+                          count_calls, tracer_flow_CSp,  ice_shelf_CSp, waves_CSp)
   type(time_type), target,   intent(inout) :: Time        !< model time, set in this routine
   type(time_type),           intent(in)    :: Time_init   !< The start time for the coupled model's calendar
   type(param_file_type),     intent(out)   :: param_file  !< structure indicating parameter file to parse
@@ -1656,6 +1656,8 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
                                                           !! calls to step_MOM instead of the number of
                                                           !! dynamics timesteps.
   type(ice_shelf_CS), optional,     pointer :: ice_shelf_CSp !< A pointer to an ice shelf control structure
+  type(Wave_parameters_CS), &
+                   optional, pointer       :: Waves_CSp       !< An optional pointer to a wave property CS
   ! local variables
   type(ocean_grid_type),  pointer :: G => NULL()    ! A pointer to the metric grid use for the run
   type(ocean_grid_type),  pointer :: G_in => NULL() ! Pointer to the input grid
@@ -2348,6 +2350,9 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
   if (associated(CS%OBC)) &
     call open_boundary_register_restarts(dg%HI, GV, CS%OBC, CS%tracer_Reg, &
                           param_file, restart_CSp, use_temperature)
+  if (present(waves_CSp)) then
+    call waves_register_restarts(waves_CSp, dG%HI, GV, param_file, restart_CSp)
+  endif
 
   call callTree_waypoint("restart registration complete (initialize_MOM)")
 

src/user/MOM_wave_interface.F90

@@ -12,12 +12,15 @@ module MOM_wave_interface
 use MOM_file_parser,   only : get_param, log_version, param_file_type
 use MOM_forcing_type,  only : mech_forcing
 use MOM_grid,          only : ocean_grid_type
+use MOM_hor_index,     only : hor_index_type
 use MOM_io,            only : file_exists, get_var_sizes, read_variable
+use MOM_io,            only : vardesc, var_desc
 use MOM_safe_alloc,    only : safe_alloc_ptr
 use MOM_time_manager,  only : time_type, operator(+), operator(/)
 use MOM_unit_scaling,  only : unit_scale_type
 use MOM_variables,     only : thermo_var_ptrs, surface
 use MOM_verticalgrid,  only : verticalGrid_type
+use MOM_restart,       only : register_restart_field, MOM_restart_CS, query_initialized
 
 implicit none ; private
 
@@ -42,6 +45,7 @@ module MOM_wave_interface
                       ! CL2 effects.
 public Waves_end ! public interface to deallocate and free wave related memory.
 public get_wave_method ! public interface to obtain the wave method string
+public waves_register_restarts ! public interface to register wave restart fields
 
 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
@@ -210,14 +214,15 @@ module MOM_wave_interface
 contains
 
 !> Initializes parameters related to MOM_wave_interface
-subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag )
+subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag, restart_CSp)
   type(time_type), target, intent(in)    :: Time       !< Model time
   type(ocean_grid_type),   intent(inout) :: G          !< Grid structure
   type(verticalGrid_type), intent(in)    :: GV         !< Vertical grid structure
   type(unit_scale_type),   intent(in)    :: US         !< A dimensional unit scaling type
   type(param_file_type),   intent(in)    :: param_file !< Input parameter structure
   type(wave_parameters_CS), pointer      :: CS         !< Wave parameter control structure
   type(diag_ctrl), target, intent(inout) :: diag       !< Diagnostic Pointer
+  type(MOM_restart_CS), optional, pointer:: restart_CSp!< Restart control structure
 
   ! Local variables
   character(len=40)  :: mdl = "MOM_wave_interface" !< This module's name.
@@ -231,8 +236,8 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag )
   logical :: StatisticalWaves
 
   ! Dummy Check
-  if (associated(CS)) then
-    call MOM_error(FATAL, "wave_interface_init called with an associated control structure.")
+  if (.not. associated(CS)) then
+    call MOM_error(FATAL, "wave_interface_init called without an associated control structure.")
     return
   endif
 
@@ -245,9 +250,6 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag )
 
   if (.not.(use_waves .or. StatisticalWaves)) return
 
-  ! Allocate CS and set pointers
-  allocate(CS)
-
   CS%UseWaves = use_waves
   CS%diag => diag
   CS%Time => Time
@@ -441,10 +443,6 @@ subroutine MOM_wave_interface_init(time, G, GV, US, param_file, CS, diag )
 
   ! Allocate and initialize
   ! a. Stokes driftProfiles
-  allocate(CS%Us_x(G%isdB:G%IedB,G%jsd:G%jed,GV%ke))
-  CS%Us_x(:,:,:) = 0.0
-  allocate(CS%Us_y(G%isd:G%Ied,G%jsdB:G%jedB,GV%ke))
-  CS%Us_y(:,:,:) = 0.0
   if (CS%Stokes_DDT) then
     allocate(CS%ddt_Us_x(G%isdB:G%IedB,G%jsd:G%jed,G%ke))
     CS%ddt_Us_x(:,:,:) = 0.0
@@ -612,6 +610,10 @@ subroutine Update_Stokes_Drift(G, GV, US, CS, h, ustar, dt)
   min_level_thick_avg = 1.e-3*US%m_to_Z
   idt = 1.0/dt
 
+  if (allocated(CS%US_x) .and. allocated(CS%US_y)) then
+    call pass_vector(CS%US_x(:,:,:),CS%US_y(:,:,:), G%Domain)
+  endif
+
   ! Getting Stokes drift profile from previous step
   CS%ddt_us_x(:,:,:) = CS%US_x(:,:,:)
   CS%ddt_us_y(:,:,:) = CS%US_y(:,:,:)
@@ -1757,6 +1759,51 @@ subroutine Waves_end(CS)
 
 end subroutine Waves_end
 
+!> Register wave restart fields. To be called before MOM_wave_interface_init
+subroutine waves_register_restarts(CS, HI, GV, param_file, restart_CSp)
+  type(wave_parameters_CS), pointer       :: CS           !< Wave parameter Control structure
+  type(hor_index_type),     intent(inout) :: HI           !< Grid structure
+  type(verticalGrid_type),  intent(in)    :: GV           !< Vertical grid structure
+  type(param_file_type),    intent(in)    :: param_file   !< Input parameter structure
+  type(MOM_restart_CS),     pointer       :: restart_CSp  !< Restart structure, data intent(inout)
+  ! Local variables
+  type(vardesc) :: vd(2)
+  logical :: use_waves
+  logical :: StatisticalWaves
+  character*(13) :: wave_method_str
+  character(len=40)  :: mdl = "MOM_wave_interface" !< This module's name.
+
+  if (associated(CS)) then
+    call MOM_error(FATAL, "waves_register_restarts: Called with initialized waves control structure")
+  endif 
+  allocate(CS)
+
+  call get_param(param_file, mdl, "USE_WAVES", use_waves, &
+       "If true, enables surface wave modules.", do_not_log=.true., default=.false.)
+
+  ! Check if using LA_LI2016
+  call get_param(param_file,mdl,"USE_LA_LI2016",StatisticalWaves,     &
+                 do_not_log=.true.,default=.false.)
+
+  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
+
+  call get_param(param_file,mdl,"WAVE_METHOD",wave_method_str, do_not_log=.true., default=NULL_STRING)
+
+  if (trim(wave_method_str)== trim(SURFBANDS_STRING)) then
+    vd(1) = var_desc("US_x", "m s-1", "3d zonal Stokes drift profile")
+    vd(2) = var_desc("US_y", "m s-1", "3d meridional Stokes drift profile")
+    call register_restart_field(CS%US_x(:,:,:), vd(1), .true., restart_CSp)
+    call register_restart_field(CS%US_y(:,:,:), vd(2), .false., restart_CSp)
+  endif
+
+end subroutine waves_register_restarts
+
 !> \namespace  mom_wave_interface
 !!
 !! \author Brandon Reichl, 2018.