MOM6
MOM_wave_interface.F90
Go to the documentation of this file.
1 !> Interface for surface waves
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_diag_mediator, only : post_data, register_diag_field, safe_alloc_alloc
8 use mom_domains, only : pass_var, pass_vector, agrid
9 use mom_domains, only : to_south, to_west, to_all
10 use mom_error_handler, only : mom_error, fatal, warning
12 use mom_grid, only : ocean_grid_type
14 use mom_time_manager, only : time_type, operator(+), operator(/)
18 use data_override_mod, only : data_override_init, data_override
19 
20 implicit none ; private
21 
22 #include <MOM_memory.h>
23 
24 public mom_wave_interface_init ! Public interface to fully initialize the wave routines.
25 public mom_wave_interface_init_lite ! Public interface to quick initialize this module.
26 public update_surface_waves ! Public interface to update wave information at the
27  ! coupler/driver level.
28 public update_stokes_drift ! Public interface to update the Stokes drift profiles
29  ! called in step_mom.
30 public get_langmuir_number ! Public interface to compute Langmuir number called from
31  ! ePBL or KPP routines.
32 public stokesmixing ! NOT READY - Public interface to add down-Stokes gradient
33  ! momentum mixing (e.g. the approach of Harcourt 2013/2015)
34 public coriolisstokes ! NOT READY - Public interface to add Coriolis-Stokes acceleration
35  ! of the mean currents, needed for comparison with LES. It is
36  ! presently advised against implementing in non-1d settings without
37  ! serious consideration of the full 3d wave-averaged Navier-Stokes
38  ! CL2 effects.
39 public waves_end ! public interface to deallocate and free wave related memory.
40 
41 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
42 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
43 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
44 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
45 
46 !> Container for all surface wave related parameters
47 type, public :: wave_parameters_cs ; private
48 
49  !Main surface wave options
50  logical, public :: usewaves !< Flag to enable surface gravity wave feature
51  logical, public :: lagrangianmixing !< This feature is in development and not ready
52  !! True if Stokes drift is present and mixing
53  !! should be applied to Lagrangian current
54  !! (mean current + Stokes drift).
55  !! See Reichl et al., 2016 KPP-LT approach
56  logical, public :: stokesmixing !< This feature is in development and not ready.
57  !! True if vertical mixing of momentum
58  !! should be applied directly to Stokes current
59  !! (with separate mixing parameter for Eulerian
60  !! mixing contribution).
61  !! See Harcourt 2013, 2015 Second-Moment approach
62  logical, public :: coriolisstokes !< This feature is in development and not ready.
63  ! True if Coriolis-Stokes acceleration should be applied.
64  integer, public :: stklevelmode=1 !< Sets if Stokes drift is defined at mid-points
65  !! or layer averaged. Set to 0 if mid-point and set to
66  !! 1 if average value of Stokes drift over level.
67  !! If advecting with Stokes transport, 1 is the correct
68  !! approach.
69 
70  ! Surface Wave Dependent 1d/2d/3d vars
71  real, allocatable, dimension(:), public :: &
72  wavenum_cen !< Wavenumber bands for read/coupled [m-1]
73  real, allocatable, dimension(:), public :: &
74  freq_cen !< Frequency bands for read/coupled [s-1]
75  real, allocatable, dimension(:), public :: &
76  prescribedsurfstkx !< Surface Stokes drift if prescribed [m s-1]
77  real, allocatable, dimension(:), public :: &
78  prescribedsurfstky !< Surface Stokes drift if prescribed [m s-1]
79  real, allocatable, dimension(:,:,:), public :: &
80  us_x !< 3d zonal Stokes drift profile [m s-1]
81  !! Horizontal -> U points
82  !! Vertical -> Mid-points
83  real, allocatable, dimension(:,:,:), public :: &
84  us_y !< 3d meridional Stokes drift profile [m s-1]
85  !! Horizontal -> V points
86  !! Vertical -> Mid-points
87  real, allocatable, dimension(:,:), public :: &
88  la_sl,& !< SL Langmuir number (directionality factored later)
89  !! Horizontal -> H points
90  la_turb !< Aligned Turbulent Langmuir number
91  !! Horizontal -> H points
92  real, allocatable, dimension(:,:), public :: &
93  us0_x !< Surface Stokes Drift (zonal, m/s)
94  !! Horizontal -> U points
95  real, allocatable, dimension(:,:), public :: &
96  us0_y !< Surface Stokes Drift (meridional, m/s)
97  !! Horizontal -> V points
98  real, allocatable, dimension(:,:,:), public :: &
99  stkx0 !< Stokes Drift spectrum (zonal, m/s)
100  !! Horizontal -> U points
101  !! 3rd dimension -> Freq/Wavenumber
102  real, allocatable, dimension(:,:,:), public :: &
103  stky0 !< Stokes Drift spectrum (meridional, m/s)
104  !! Horizontal -> V points
105  !! 3rd dimension -> Freq/Wavenumber
106  real, allocatable, dimension(:,:,:), public :: &
107  kvs !< Viscosity for Stokes Drift shear [Z2 T-1 ~> m2 s-1]
108 
109  ! Pointers to auxiliary fields
110  type(time_type), pointer, public :: time !< A pointer to the ocean model's clock.
111  type(diag_ctrl), pointer, public :: diag !< A structure that is used to regulate the
112  !! timing of diagnostic output.
113 
114  !> An arbitrary lower-bound on the Langmuir number. Run-time parameter.
115  !! Langmuir number is sqrt(u_star/u_stokes). When both are small
116  !! but u_star is orders of magnitude smaller the Langmuir number could
117  !! have unintended consequences. Since both are small it can be safely capped
118  !! to avoid such consequences.
119  real :: la_min = 0.05
120 
121  !>@{ Diagnostic handles
122  integer, public :: id_surfacestokes_x = -1 , id_surfacestokes_y = -1
123  integer, public :: id_3dstokes_x = -1 , id_3dstokes_y = -1
124  integer, public :: id_la_turb = -1
125  !!@}
126 
127 end type wave_parameters_cs
128 
129 ! Options not needed outside of this module
130 
131 integer :: wavemethod=-99 !< Options for including wave information
132  !! Valid (tested) choices are:
133  !! 0 - Test Profile
134  !! 1 - Surface Stokes Drift Bands
135  !! 2 - DHH85
136  !! 3 - LF17
137  !! -99 - No waves computed, but empirical Langmuir number used.
138  !! \todo Module variable! Move into a control structure.
139 
140 ! Options if WaveMethod is Surface Stokes Drift Bands (1)
141 integer, public :: numbands =0 !< Number of wavenumber/frequency partitions to receive
142  !! This needs to match the number of bands provided
143  !! via either coupling or file.
144  !! \todo Module variable! Move into a control structure.
145 integer, public :: partitionmode !< Method for partition mode (meant to check input)
146  !! 0 - wavenumbers
147  !! 1 - frequencies
148  !! \todo Module variable! Move into a control structure.
149 integer :: datasource !< Integer that specifies where the Model Looks for Data
150  !! Valid choices are:
151  !! 1 - FMS DataOverride Routine
152  !! 2 - Reserved For Coupler
153  !! 3 - User input (fixed values, useful for 1d testing)
154  !! \todo Module variable! Move into a control structure.
155 
156 ! Options if using FMS DataOverride Routine
157 character(len=40) :: surfbandfilename !< Filename if using DataOverride
158  !! \todo Module variable! Move into a control structure.
159 logical :: dataoverrideisinitialized !< Flag for DataOverride Initialization
160  !! \todo Module variable! Move into a control structure.
161 
162 ! Options for computing Langmuir number
163 real :: la_frachbl !< Fraction of OSBL for averaging Langmuir number
164  !! \todo Module variable! Move into a control structure.
165 logical :: la_misalignment = .false. !< Flag to use misalignment in Langmuir number
166  !! \todo Module variable! Move into a control structure.
167 
168 ! This include declares and sets the variable "version".
169 #include "version_variable.h"
170 
171 character(len=40) :: mdl = "MOM_wave_interface" !< This module's name.
172 
173 !>@{ Undocumented parameters.
174 !! \todo These module variables need to be documented as static/private variables or moved
175 !! into a control structure.
176 ! Switches needed in import_stokes_drift
177 integer, parameter :: testprof = 0, surfbands = 1, &
178  dhh85 = 2, lf17 = 3, null_wavemethod=-99, &
179  dataovr = 1, coupler = 2, input = 3
180 
181 ! Options For Test Prof
183 logical :: waveagepeakfreq ! Flag to use W
186 real :: pi
187 !!@}
188 
189 contains
190 
191 !> Initializes parameters related to MOM_wave_interface
192 subroutine mom_wave_interface_init(time, G, GV, US, param_file, CS, diag )
193  type(time_type), target, intent(in) :: time !< Model time
194  type(ocean_grid_type), intent(inout) :: g !< Grid structure
195  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
196  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
197  type(param_file_type), intent(in) :: param_file !< Input parameter structure
198  type(wave_parameters_cs), pointer :: cs !< Wave parameter control structure
199  type(diag_ctrl), target, intent(inout) :: diag !< Diagnostic Pointer
200  ! Local variables
201  ! I/O
202  character*(13) :: tmpstring1,tmpstring2
203  character*(5), parameter :: null_string = "EMPTY"
204  character*(12), parameter :: testprof_string = "TEST_PROFILE"
205  character*(13), parameter :: surfbands_string = "SURFACE_BANDS"
206  character*(5), parameter :: dhh85_string = "DHH85"
207  character*(4), parameter :: lf17_string = "LF17"
208  character*(12), parameter :: dataovr_string = "DATAOVERRIDE"
209  character*(7), parameter :: coupler_string = "COUPLER"
210  character*(5), parameter :: input_string = "INPUT"
211 
212  ! Dummy Check
213  if (associated(cs)) then
214  call mom_error(fatal, "wave_interface_init called with an associated"//&
215  "control structure.")
216  return
217  endif
218 
219  pi=4.0*atan(1.0)
220 
221  ! Allocate CS and set pointers
222  allocate(cs)
223 
224  cs%diag => diag
225  cs%Time => time
226 
227  ! Add any initializations needed here
228  dataoverrideisinitialized = .false.
229 
230  ! The only way to get here is with UseWaves enabled.
231  cs%UseWaves=.true.
232 
233  call log_version(param_file, mdl, version)
234 
235  ! Wave modified physics
236  ! Presently these are all in research mode
237  call get_param(param_file, mdl, "LAGRANGIAN_MIXING", cs%LagrangianMixing, &
238  "Flag to use Lagrangian Mixing of momentum", units="", &
239  default=.false.)
240  if (cs%LagrangianMixing) then
241  ! Force Code Intervention
242  call mom_error(fatal,"Should you be enabling Lagrangian Mixing? Code not ready.")
243  endif
244  call get_param(param_file, mdl, "STOKES_MIXING", cs%StokesMixing, &
245  "Flag to use Stokes Mixing of momentum", units="", &
246  default=.false.)
247  if (cs%StokesMixing) then
248  ! Force Code Intervention
249  call mom_error(fatal,"Should you be enabling Stokes Mixing? Code not ready.")
250  endif
251  call get_param(param_file, mdl, "CORIOLIS_STOKES", cs%CoriolisStokes, &
252  "Flag to use Coriolis Stokes acceleration", units="", &
253  default=.false.)
254  if (cs%CoriolisStokes) then
255  ! Force Code Intervention
256  call mom_error(fatal,"Should you be enabling Coriolis-Stokes? Code not ready.")
257  endif
258 
259  ! Get Wave Method and write to integer WaveMethod
260  call get_param(param_file,mdl,"WAVE_METHOD",tmpstring1, &
261  "Choice of wave method, valid options include: \n"// &
262  " TEST_PROFILE - Prescribed from surface Stokes drift \n"// &
263  " and a decay wavelength.\n"// &
264  " SURFACE_BANDS - Computed from multiple surface values \n"// &
265  " and decay wavelengths.\n"// &
266  " DHH85 - Uses Donelan et al. 1985 empirical \n"// &
267  " wave spectrum with prescribed values. \n"// &
268  " LF17 - Infers Stokes drift profile from wind \n"// &
269  " speed following Li and Fox-Kemper 2017.\n", &
270  units='', default=null_string)
271  select case (trim(tmpstring1))
272  case (null_string)! No Waves
273  call mom_error(fatal, "wave_interface_init called with no specified "//&
274  "WAVE_METHOD.")
275  case (testprof_string)! Test Profile
277  call get_param(param_file,mdl,"TP_STKX_SURF",tp_stkx0,&
278  'Surface Stokes (x) for test profile',&
279  units='m/s',default=0.1)
280  call get_param(param_file,mdl,"TP_STKY_SURF",tp_stky0,&
281  'Surface Stokes (y) for test profile',&
282  units='m/s',default=0.0)
283  call get_param(param_file,mdl,"TP_WVL",tp_wvl,&
284  units='m', default=50.0, scale=us%m_to_Z)
285  case (surfbands_string)! Surface Stokes Drift Bands
287  call get_param(param_file, mdl, "SURFBAND_SOURCE",tmpstring2, &
288  "Choice of SURFACE_BANDS data mode, valid options include: \n"// &
289  " DATAOVERRIDE - Read from NetCDF using FMS DataOverride. \n"// &
290  " COUPLER - Look for variables from coupler pass \n"// &
291  " INPUT - Testing with fixed values.", &
292  units='', default=null_string)
293  select case (trim(tmpstring2))
294  case (null_string)! Default
295  call mom_error(fatal, "wave_interface_init called with SURFACE_BANDS"//&
296  " but no SURFBAND_SOURCE.")
297  case (dataovr_string)! Using Data Override
299  call get_param(param_file, mdl, "SURFBAND_FILENAME", surfbandfilename, &
300  "Filename of surface Stokes drift input band data.", default="StkSpec.nc")
301  case (coupler_string)! Reserved for coupling
303  case (input_string)! A method to input the Stokes band (globally uniform)
304  datasource = input
305  call get_param(param_file,mdl,"SURFBAND_NB",numbands, &
306  "Prescribe number of wavenumber bands for Stokes drift. "// &
307  "Make sure this is consistnet w/ WAVENUMBERS, STOKES_X, and "// &
308  "STOKES_Y, there are no safety checks in the code.", &
309  units='', default=1)
310  allocate( cs%WaveNum_Cen(1:numbands) )
311  cs%WaveNum_Cen(:) = 0.0
312  allocate( cs%PrescribedSurfStkX(1:numbands))
313  cs%PrescribedSurfStkX(:) = 0.0
314  allocate( cs%PrescribedSurfStkY(1:numbands))
315  cs%PrescribedSurfStkY(:) = 0.0
316  allocate( cs%STKx0(g%isdB:g%iedB,g%jsd:g%jed,1:numbands))
317  cs%STKx0(:,:,:) = 0.0
318  allocate( cs%STKy0(g%isd:g%ied,g%jsdB:g%jedB,1:numbands))
319  cs%STKy0(:,:,:) = 0.0
320  partitionmode=0
321  call get_param(param_file,mdl,"SURFBAND_WAVENUMBERS",cs%WaveNum_Cen, &
322  "Central wavenumbers for surface Stokes drift bands.",units='rad/m', &
323  default=0.12566)
324  call get_param(param_file,mdl,"SURFBAND_STOKES_X",cs%PrescribedSurfStkX, &
325  "X-direction surface Stokes drift for bands.",units='m/s', &
326  default=0.15)
327  call get_param(param_file,mdl,"SURFBAND_STOKES_Y",cs%PrescribedSurfStkY, &
328  "Y-direction surface Stokes drift for bands.",units='m/s', &
329  default=0.0)
330  case default! No method provided
331  call mom_error(fatal,'Check WAVE_METHOD.')
332  end select
333 
334  case (dhh85_string)!Donelan et al., 1985 spectrum
335  wavemethod = dhh85
336  call mom_error(warning,"DHH85 only ever set-up for uniform cases w/"//&
337  " Stokes drift in x-direction.")
338  call get_param(param_file,mdl,"DHH85_AGE_FP",waveagepeakfreq, &
339  "Choose true to use waveage in peak frequency.", &
340  units='', default=.false.)
341  call get_param(param_file,mdl,"DHH85_AGE",waveage, &
342  "Wave Age for DHH85 spectrum.", &
343  units='', default=1.2)
344  call get_param(param_file,mdl,"DHH85_WIND",wavewind, &
345  "Wind speed for DHH85 spectrum.", &
346  units='', default=10.0)
347  call get_param(param_file,mdl,"STATIC_DHH85",staticwaves, &
348  "Flag to disable updating DHH85 Stokes drift.", &
349  default=.false.)
350  case (lf17_string)!Li and Fox-Kemper 17 wind-sea Langmuir number
351  wavemethod = lf17
352  case default
353  call mom_error(fatal,'Check WAVE_METHOD.')
354  end select
355 
356  ! Langmuir number Options
357  call get_param(param_file, mdl, "LA_DEPTH_RATIO", la_frachbl, &
358  "The depth (normalized by BLD) to average Stokes drift over in "//&
359  "Langmuir number calculation, where La = sqrt(ust/Stokes).", &
360  units="nondim",default=0.04)
361  call get_param(param_file, mdl, "LA_MISALIGNMENT", la_misalignment, &
362  "Flag (logical) if using misalignment bt shear and waves in LA",&
363  default=.false.)
364  call get_param(param_file, mdl, "MIN_LANGMUIR", cs%La_min, &
365  "A minimum value for all Langmuir numbers that is not physical, "//&
366  "but is likely only encountered when the wind is very small and "//&
367  "therefore its effects should be mostly benign.",units="nondim",&
368  default=0.05)
369 
370  ! Allocate and initialize
371  ! a. Stokes driftProfiles
372  allocate(cs%Us_x(g%isdB:g%IedB,g%jsd:g%jed,g%ke))
373  cs%Us_x(:,:,:) = 0.0
374  allocate(cs%Us_y(g%isd:g%Ied,g%jsdB:g%jedB,g%ke))
375  cs%Us_y(:,:,:) = 0.0
376  ! b. Surface Values
377  allocate(cs%US0_x(g%isdB:g%iedB,g%jsd:g%jed))
378  cs%US0_x(:,:) = 0.0
379  allocate(cs%US0_y(g%isd:g%ied,g%jsdB:g%jedB))
380  cs%US0_y(:,:) = 0.0
381  ! c. Langmuir number
382  allocate(cs%La_SL(g%isc:g%iec,g%jsc:g%jec))
383  allocate(cs%La_turb(g%isc:g%iec,g%jsc:g%jec))
384  cs%La_SL(:,:) = 0.0
385  cs%La_turb (:,:) = 0.0
386  ! d. Viscosity for Stokes drift
387  if (cs%StokesMixing) then
388  allocate(cs%KvS(g%isd:g%Ied,g%jsd:g%jed,g%ke))
389  cs%KvS(:,:,:) = 0.0
390  endif
391 
392  ! Initialize Wave related outputs
393  cs%id_surfacestokes_y = register_diag_field('ocean_model','surface_stokes_y', &
394  cs%diag%axesCu1,time,'Surface Stokes drift (y)','m s-1')
395  cs%id_surfacestokes_x = register_diag_field('ocean_model','surface_stokes_x', &
396  cs%diag%axesCv1,time,'Surface Stokes drift (x)','m s-1')
397  cs%id_3dstokes_y = register_diag_field('ocean_model','3d_stokes_y', &
398  cs%diag%axesCvL,time,'3d Stokes drift (y)','m s-1')
399  cs%id_3dstokes_x = register_diag_field('ocean_model','3d_stokes_x', &
400  cs%diag%axesCuL,time,'3d Stokes drift (y)','m s-1')
401  cs%id_La_turb = register_diag_field('ocean_model','La_turbulent',&
402  cs%diag%axesT1,time,'Surface (turbulent) Langmuir number','nondim')
403 
404  return
405 end subroutine mom_wave_interface_init
406 
407 !> A 'lite' init subroutine to initialize a few inputs needed if using wave information
408 !! with the wind-speed dependent Stokes drift formulation of LF17
409 subroutine mom_wave_interface_init_lite(param_file)
410  type(param_file_type), intent(in) :: param_file !< Input parameter structure
411  character*(5), parameter :: null_string = "EMPTY"
412  character*(4), parameter :: lf17_string = "LF17"
413  character*(13) :: tmpstring1
414  logical :: statisticalwaves
415 
416  ! Langmuir number Options
417  call get_param(param_file, mdl, "LA_DEPTH_RATIO", la_frachbl, &
418  "The depth (normalized by BLD) to average Stokes drift over in "//&
419  "Langmuir number calculation, where La = sqrt(ust/Stokes).", &
420  units="nondim",default=0.04)
421 
422  ! Check if using LA_LI2016
423  call get_param(param_file,mdl,"USE_LA_LI2016",statisticalwaves, &
424  do_not_log=.true.,default=.false.)
425  if (statisticalwaves) then
426  wavemethod = lf17
427  pi=4.0*atan(1.0)
428  else
430  end if
431 
432  return
433 end subroutine mom_wave_interface_init_lite
434 
435 !> Subroutine that handles updating of surface wave/Stokes drift related properties
436 subroutine update_surface_waves(G, GV, US, Day, dt, CS)
437  type(wave_parameters_cs), pointer :: cs !< Wave parameter Control structure
438  type(ocean_grid_type), intent(inout) :: g !< Grid structure
439  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
440  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
441  type(time_type), intent(in) :: day !< Current model time
442  type(time_type), intent(in) :: dt !< Timestep as a time-type
443  ! Local variables
444  integer :: ii, jj, kk, b
445  type(time_type) :: day_center
446 
447  ! Computing central time of time step
448  day_center = day + dt/2
449 
450  if (wavemethod == testprof) then
451  ! Do nothing
452  elseif (wavemethod==surfbands) then
453  if (datasource==dataovr) then
454  call surface_bands_by_data_override(day_center, g, gv, us, cs)
455  elseif (datasource==coupler) then
456  ! Reserve for coupler hooks
457  elseif (datasource==input) then
458  do b=1,numbands
459  do ii=g%isdB,g%iedB
460  do jj=g%jsd,g%jed
461  cs%STKx0(ii,jj,b) = cs%PrescribedSurfStkX(b)
462  enddo
463  enddo
464  do ii=g%isd,g%ied
465  do jj=g%jsdB, g%jedB
466  cs%STKY0(ii,jj,b) = cs%PrescribedSurfStkY(b)
467  enddo
468  enddo
469  enddo
470  endif
471  endif
472 
473  return
474 end subroutine update_surface_waves
475 
476 !> Constructs the Stokes Drift profile on the model grid based on
477 !! desired coupling options
478 subroutine update_stokes_drift(G, GV, US, CS, h, ustar)
479  type(wave_parameters_cs), pointer :: cs !< Wave parameter Control structure
480  type(ocean_grid_type), intent(inout) :: g !< Grid structure
481  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
482  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
483  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
484  intent(in) :: h !< Thickness [H ~> m or kg m-2]
485  real, dimension(SZI_(G),SZJ_(G)), &
486  intent(in) :: ustar !< Wind friction velocity [Z T-1 ~> m s-1].
487  ! Local Variables
488  real :: top, midpoint, bottom, one_cm
489  real :: decayscale
490  real :: cmn_fac, wn, ustokes
491  real :: la
492  integer :: ii, jj, kk, b, iim1, jjm1
493 
494  one_cm = 0.01*us%m_to_Z
495 
496  ! 1. If Test Profile Option is chosen
497  ! Computing mid-point value from surface value and decay wavelength
498  if (wavemethod==testprof) then
499  decayscale = 4.*pi / tp_wvl !4pi
500  do ii = g%isdB,g%iedB
501  do jj = g%jsd,g%jed
502  iim1 = max(1,ii-1)
503  bottom = 0.0
504  midpoint = 0.0
505  do kk = 1,g%ke
506  top = bottom
507  midpoint = bottom - gv%H_to_Z*0.25*(h(ii,jj,kk)+h(iim1,jj,kk))
508  bottom = bottom - gv%H_to_Z*0.5*(h(ii,jj,kk)+h(iim1,jj,kk))
509  cs%Us_x(ii,jj,kk) = tp_stkx0*exp(midpoint*decayscale)
510  enddo
511  enddo
512  enddo
513  do ii = g%isd,g%ied
514  do jj = g%jsdB,g%jedB
515  jjm1 = max(1,jj-1)
516  bottom = 0.0
517  midpoint = 0.0
518  do kk = 1,g%ke
519  top = bottom
520  midpoint = bottom - gv%H_to_Z*0.25*(h(ii,jj,kk)+h(ii,jjm1,kk))
521  bottom = bottom - gv%H_to_Z*0.5*(h(ii,jj,kk)+h(ii,jjm1,kk))
522  cs%Us_y(ii,jj,kk) = tp_stky0*exp(midpoint*decayscale)
523  enddo
524  enddo
525  enddo
526  ! 2. If Surface Bands is chosen
527  ! In wavenumber mode compute integral for layer averaged Stokes drift.
528  ! In frequency mode compuate value at midpoint.
529  elseif (wavemethod==surfbands) then
530  cs%Us_x(:,:,:) = 0.0
531  cs%Us_y(:,:,:) = 0.0
532  cs%Us0_x(:,:) = 0.0
533  cs%Us0_y(:,:) = 0.0
534  ! Computing X direction Stokes drift
535  do ii = g%isdB,g%iedB
536  do jj = g%jsd,g%jed
537  ! 1. First compute the surface Stokes drift
538  ! by integrating over the partitionas.
539  do b = 1,numbands
540  if (partitionmode==0) then
541  ! In wavenumber we are averaging over (small) level
542  cmn_fac = (1.0-exp(-one_cm*2.*cs%WaveNum_Cen(b))) / &
543  (one_cm*2.*cs%WaveNum_Cen(b))
544  elseif (partitionmode==1) then
545  ! In frequency we are not averaging over level and taking top
546  cmn_fac = 1.0
547  endif
548  cs%US0_x(ii,jj) = cs%US0_x(ii,jj) + cs%STKx0(ii,jj,b)*cmn_fac
549  enddo
550  ! 2. Second compute the level averaged Stokes drift
551  bottom = 0.0
552  do kk = 1,g%ke
553  top = bottom
554  iim1 = max(ii-1,1)
555  midpoint = bottom - gv%H_to_Z*0.25*(h(ii,jj,kk)+h(iim1,jj,kk))
556  bottom = bottom - gv%H_to_Z*0.5*(h(ii,jj,kk)+h(iim1,jj,kk))
557  do b = 1,numbands
558  if (partitionmode==0) then
559  ! In wavenumber we are averaging over level
560  cmn_fac = (exp(top*2.*cs%WaveNum_Cen(b))-exp(bottom*2.*cs%WaveNum_Cen(b)))&
561  / ((top-bottom)*(2.*cs%WaveNum_Cen(b)))
562  elseif (partitionmode==1) then
563  if (cs%StkLevelMode==0) then
564  ! Take the value at the midpoint
565  cmn_fac = exp(midpoint*2.*(2.*pi*cs%Freq_Cen(b)*us%T_to_s)**2/(us%L_to_Z**2*gv%g_Earth))
566  elseif (cs%StkLevelMode==1) then
567  ! Use a numerical integration and then
568  ! divide by layer thickness
569  wn = (2.*pi*cs%Freq_Cen(b)*us%T_to_s)**2 / (us%L_to_Z**2*gv%g_Earth) !bgr bug-fix missing g
570  cmn_fac = (exp(2.*wn*top)-exp(2.*wn*bottom)) / (2.*wn*(top-bottom))
571  endif
572  endif
573  cs%US_x(ii,jj,kk) = cs%US_x(ii,jj,kk) + cs%STKx0(ii,jj,b)*cmn_fac
574  enddo
575  enddo
576  enddo
577  enddo
578  ! Computing Y direction Stokes drift
579  do ii = g%isd,g%ied
580  do jj = g%jsdB,g%jedB
581  ! Compute the surface values.
582  do b = 1,numbands
583  if (partitionmode==0) then
584  ! In wavenumber we are averaging over (small) level
585  cmn_fac = (1.0-exp(-one_cm*2.*cs%WaveNum_Cen(b))) / &
586  (one_cm*2.*cs%WaveNum_Cen(b))
587  elseif (partitionmode==1) then
588  ! In frequency we are not averaging over level and taking top
589  cmn_fac = 1.0
590  endif
591  cs%US0_y(ii,jj) = cs%US0_y(ii,jj) + cs%STKy0(ii,jj,b)*cmn_fac
592  enddo
593  ! Compute the level averages.
594  bottom = 0.0
595  do kk = 1,g%ke
596  top = bottom
597  jjm1 = max(jj-1,1)
598  midpoint = bottom - gv%H_to_Z*0.25*(h(ii,jj,kk)+h(ii,jjm1,kk))
599  bottom = bottom - gv%H_to_Z*0.5*(h(ii,jj,kk)+h(ii,jjm1,kk))
600  do b = 1,numbands
601  if (partitionmode==0) then
602  ! In wavenumber we are averaging over level
603  cmn_fac = (exp(top*2.*cs%WaveNum_Cen(b)) - &
604  exp(bottom*2.*cs%WaveNum_Cen(b))) / &
605  ((top-bottom)*(2.*cs%WaveNum_Cen(b)))
606  elseif (partitionmode==1) then
607  if (cs%StkLevelMode==0) then
608  ! Take the value at the midpoint
609  cmn_fac = exp(midpoint*2.*(2.*pi*cs%Freq_Cen(b)*us%T_to_s)**2/(us%L_to_Z**2*gv%g_Earth))
610  elseif (cs%StkLevelMode==1) then
611  ! Use a numerical integration and then
612  ! divide by layer thickness
613  wn = (2.*pi*cs%Freq_Cen(b)*us%T_to_s)**2 / (us%L_to_Z**2*gv%g_Earth)
614  cmn_fac = (exp(2.*wn*top)-exp(2.*wn*bottom)) / (2.*wn*(top-bottom))
615  endif
616  endif
617  cs%US_y(ii,jj,kk) = cs%US_y(ii,jj,kk) + cs%STKy0(ii,jj,b)*cmn_fac
618  enddo
619  enddo
620  enddo
621  enddo
622  elseif (wavemethod==dhh85) then
623  if (.not.(staticwaves .and. dhh85_is_set)) then
624  do ii = g%isdB,g%iedB
625  do jj = g%jsd,g%jed
626  bottom = 0.0
627  do kk = 1,g%ke
628  top = bottom
629  iim1 = max(ii-1,1)
630  midpoint = bottom - gv%H_to_Z*0.25*(h(ii,jj,kk)+h(iim1,jj,kk))
631  bottom = bottom - gv%H_to_Z*0.5*(h(ii,jj,kk)+h(iim1,jj,kk))
632  !bgr note that this is using a u-point ii on h-point ustar
633  ! this code has only been previous used for uniform
634  ! grid cases. This needs fixed if DHH85 is used for non
635  ! uniform cases.
636  call dhh85_mid(gv, us, midpoint, ustokes)
637  ! Putting into x-direction (no option for direction
638  cs%US_x(ii,jj,kk) = ustokes
639  enddo
640  enddo
641  enddo
642  do ii = g%isd,g%ied
643  do jj = g%jsdB,g%jedB
644  bottom = 0.0
645  do kk=1, g%ke
646  top = bottom
647  jjm1 = max(jj-1,1)
648  midpoint = bottom - gv%H_to_Z*0.25*(h(ii,jj,kk)+h(ii,jjm1,kk))
649  bottom = bottom - gv%H_to_Z*0.5*(h(ii,jj,kk)+h(ii,jjm1,kk))
650  !bgr note that this is using a v-point jj on h-point ustar
651  ! this code has only been previous used for uniform
652  ! grid cases. This needs fixed if DHH85 is used for non
653  ! uniform cases.
654  ! call DHH85_mid(GV, US, Midpoint, UStokes)
655  ! Putting into x-direction, so setting y direction to 0
656  cs%US_y(ii,jj,kk) = 0.0 !### Note that =0 should be =US - RWH
657  ! bgr - see note above, but this is true
658  ! if this is used for anything
659  ! other than simple LES comparison
660  enddo
661  enddo
662  enddo
663  dhh85_is_set = .true.
664  endif
665  else! Keep this else, fallback to 0 Stokes drift
666  do kk= 1,g%ke
667  do ii = g%isdB,g%iedB
668  do jj = g%jsd,g%jed
669  cs%Us_x(ii,jj,kk) = 0.
670  enddo
671  enddo
672  do ii = g%isd,g%ied
673  do jj = g%jsdB,g%jedB
674  cs%Us_y(ii,jj,kk) = 0.
675  enddo
676  enddo
677  enddo
678  endif
679 
680  ! Turbulent Langmuir number is computed here and available to use anywhere.
681  ! SL Langmuir number requires mixing layer depth, and therefore is computed
682  ! in the routine it is needed by (e.g. KPP or ePBL).
683  do ii = g%isc,g%iec
684  do jj = g%jsc, g%jec
685  top = h(ii,jj,1)*gv%H_to_Z
686  call get_langmuir_number( la, g, gv, us, top, ustar(ii,jj), ii, jj, &
687  h(ii,jj,:),override_ma=.false.,waves=cs)
688  cs%La_turb(ii,jj) = la
689  enddo
690  enddo
691 
692  ! Output any desired quantities
693  if (cs%id_surfacestokes_y>0) &
694  call post_data(cs%id_surfacestokes_y, cs%us0_y, cs%diag)
695  if (cs%id_surfacestokes_x>0) &
696  call post_data(cs%id_surfacestokes_x, cs%us0_x, cs%diag)
697  if (cs%id_3dstokes_y>0) &
698  call post_data(cs%id_3dstokes_y, cs%us_y, cs%diag)
699  if (cs%id_3dstokes_x>0) &
700  call post_data(cs%id_3dstokes_x, cs%us_x, cs%diag)
701  if (cs%id_La_turb>0) &
702  call post_data(cs%id_La_turb, cs%La_turb, cs%diag)
703 
704 end subroutine update_stokes_drift
705 
706 !> A subroutine to fill the Stokes drift from a NetCDF file
707 !! using the data_override procedures.
708 subroutine surface_bands_by_data_override(day_center, G, GV, US, CS)
709  use netcdf
710  type(time_type), intent(in) :: day_center !< Center of timestep
711  type(wave_parameters_cs), pointer :: CS !< Wave structure
712  type(ocean_grid_type), intent(inout) :: G !< Grid structure
713  type(verticalgrid_type), intent(in) :: GV !< Vertical grid structure
714  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
715  ! Local variables
716  real :: temp_x(SZI_(G),SZJ_(G)) ! Pseudo-zonal Stokes drift of band at h-points [m s-1]
717  real :: temp_y(SZI_(G),SZJ_(G)) ! Psuedo-meridional Stokes drift of band at h-points [m s-1]
718  real :: Top, MidPoint
719  integer :: b
720  integer :: i, j
721  integer, dimension(4) :: start, counter, dims, dim_id
722  character(len=12) :: dim_name(4)
723  character(20) :: varname, varread1, varread2
724  integer :: rcode_fr, rcode_wn, ncid, varid_fr, varid_wn, id, ndims
725 
726  if (.not.dataoverrideisinitialized) then
727  call data_override_init(ocean_domain_in=g%Domain%mpp_domain)
729 
730  ! Read in number of wavenumber bands in file to set number to be read in
731  ! Hardcoded filename/variables
732  varread1 = 'wavenumber' !Old method gives wavenumber
733  varread2 = 'frequency' !New method gives frequency
734  rcode_wn = nf90_open(trim(surfbandfilename), nf90_nowrite, ncid)
735  if (rcode_wn /= 0) then
736  call mom_error(fatal,"error opening file "//trim(surfbandfilename)//&
737  " in MOM_wave_interface.")
738  endif
739 
740  ! Check if rcode_wn or rcode_fr is 0 (checks if input has wavenumber or frequency)
741  rcode_wn = nf90_inq_varid(ncid, varread1, varid_wn)
742  rcode_fr = nf90_inq_varid(ncid, varread2, varid_fr)
743 
744  if (rcode_wn /= 0 .and. rcode_fr /= 0) then
745  call mom_error(fatal,"error finding variable "//trim(varread1)//&
746  " or "//trim(varread2)//" in file "//trim(surfbandfilename)//" in MOM_wave_interface.")
747 
748  elseif (rcode_wn == 0) then
749  ! wavenumbers found:
750  partitionmode = 0
751  rcode_wn = nf90_inquire_variable(ncid, varid_wn, ndims=ndims, &
752  dimids=dims)
753  if (rcode_wn /= 0) then
754  call mom_error(fatal, &
755  'error inquiring dimensions MOM_wave_interface.')
756  endif
757  rcode_wn = nf90_inquire_dimension(ncid, dims(1), dim_name(1), len=id)
758  if (rcode_wn /= 0) then
759  call mom_error(fatal,"error reading dimension 1 data for "// &
760  trim(varread1)//" in file "// trim(surfbandfilename)// &
761  " in MOM_wave_interface.")
762  endif
763  rcode_wn = nf90_inq_varid(ncid, dim_name(1), dim_id(1))
764  if (rcode_wn /= 0) then
765  call mom_error(fatal,"error finding variable "//trim(dim_name(1))//&
766  " in file "//trim(surfbandfilename)//" in MOM_wave_interace.")
767  endif
768  ! Allocating size of wavenumber bins
769  allocate( cs%WaveNum_Cen(1:id) )
770  cs%WaveNum_Cen(:) = 0.0
771  elseif (rcode_fr == 0) then
772  ! frequencies found:
773  partitionmode = 1
774  rcode_fr = nf90_inquire_variable(ncid, varid_fr, ndims=ndims, &
775  dimids=dims)
776  if (rcode_fr /= 0) then
777  call mom_error(fatal,&
778  'error inquiring dimensions MOM_wave_interface.')
779  endif
780  rcode_fr = nf90_inquire_dimension(ncid, dims(1), dim_name(1), len=id)
781  if (rcode_fr /= 0) then
782  call mom_error(fatal,"error reading dimension 1 data for "// &
783  trim(varread2)//" in file "// trim(surfbandfilename)// &
784  " in MOM_wave_interface.")
785  endif
786  rcode_fr = nf90_inq_varid(ncid, dim_name(1), dim_id(1))
787  if (rcode_fr /= 0) then
788  call mom_error(fatal,"error finding variable "//trim(dim_name(1))//&
789  " in file "//trim(surfbandfilename)//" in MOM_wave_interace.")
790  endif
791  ! Allocating size of frequency bins
792  allocate( cs%Freq_Cen(1:id) )
793  cs%Freq_Cen(:) = 0.0
794  ! Allocating size of wavenumber bins
795  allocate( cs%WaveNum_Cen(1:id) )
796  cs%WaveNum_Cen(:) = 0.0
797  allocate( cs%STKx0(g%isdB:g%iedB,g%jsd:g%jed,1:id))
798  cs%STKx0(:,:,:) = 0.0
799  allocate( cs%STKy0(g%isd:g%ied,g%jsdB:g%jedB,1:id))
800  cs%STKy0(:,:,:) = 0.0
801  endif
802 
803  ! Reading wavenumber bins/Frequencies
804  start(:) = 1 ! Set all start to 1
805  counter(:) = 1 ! Set all counter to 1
806  counter(1) = id ! Set counter(1) to id (number of frequency bins)
807  if (partitionmode==0) then
808  rcode_wn = nf90_get_var(ncid, dim_id(1), cs%WaveNum_Cen, start, counter)
809  if (rcode_wn /= 0) then
810  call mom_error(fatal,&
811  "error reading dimension 1 values for var_name "// &
812  trim(varread1)//",dim_name "//trim(dim_name(1))// &
813  " in file "// trim(surfbandfilename)//" in MOM_wave_interface")
814  endif
815  numbands = id
816  do b = 1,numbands ; cs%WaveNum_Cen(b) = us%Z_to_m*cs%WaveNum_Cen(b) ; enddo
817  elseif (partitionmode==1) then
818  rcode_fr = nf90_get_var(ncid, dim_id(1), cs%Freq_Cen, start, counter)
819  if (rcode_fr /= 0) then
820  call mom_error(fatal,&
821  "error reading dimension 1 values for var_name "// &
822  trim(varread2)//",dim_name "//trim(dim_name(1))// &
823  " in file "// trim(surfbandfilename)//" in MOM_wave_interface")
824  endif
825  numbands = id
826  do b = 1,numbands
827  cs%WaveNum_Cen(b) = (2.*pi*cs%Freq_Cen(b)*us%T_to_s)**2 / (us%L_to_Z**2*gv%g_Earth)
828  enddo
829  endif
830 
831  endif
832 
833  do b = 1,numbands
834  temp_x(:,:) = 0.0
835  temp_y(:,:) = 0.0
836  varname = ' '
837  write(varname,"(A3,I0)")'Usx',b
838  call data_override('OCN',trim(varname), temp_x, day_center)
839  varname = ' '
840  write(varname,'(A3,I0)')'Usy',b
841  call data_override('OCN',trim(varname), temp_y, day_center)
842  ! Disperse into halo on h-grid
843  call pass_vector(temp_x, temp_y, g%Domain, to_all, agrid)
844  !Filter land values
845  do j = g%jsd,g%jed
846  do i = g%Isd,g%Ied
847  if (abs(temp_x(i,j)) > 10. .or. abs(temp_y(i,j)) > 10.) then
848  ! Assume land-mask and zero out
849  temp_x(i,j) = 0.0
850  temp_y(i,j) = 0.0
851  endif
852  enddo
853  enddo
854 
855  ! Interpolate to u/v grids
856  do j = g%jsc,g%jec
857  do i = g%IscB,g%IecB
858  cs%STKx0(i,j,b) = 0.5 * (temp_x(i,j) + temp_x(i+1,j))
859  enddo
860  enddo
861  do j = g%JscB,g%JecB
862  do i = g%isc,g%iec
863  cs%STKy0(i,j,b) = 0.5 * (temp_y(i,j) + temp_y(i,j+1))
864  enddo
865  enddo
866  ! Disperse into halo on u/v grids
867  call pass_vector(cs%STKx0(:,:,b),cs%STKy0(:,:,b), g%Domain, to_all)
868  enddo !Closes b-loop
869 
870 end subroutine surface_bands_by_data_override
871 
872 !> Interface to get Langmuir number based on options stored in wave structure
873 !!
874 !! Note this can be called with an unallocated Waves pointer, which is okay if we
875 !! want the wind-speed only dependent Langmuir number. Therefore, we need to be
876 !! careful about what we try to access here.
877 subroutine get_langmuir_number( LA, G, GV, US, HBL, ustar, i, j, &
878  H, U_H, V_H, Override_MA, Waves )
879  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
880  type(verticalgrid_type), intent(in) :: gv !< Ocean vertical grid structure
881  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
882  integer, intent(in) :: i !< Meridional index of h-point
883  integer, intent(in) :: j !< Zonal index of h-point
884  real, intent(in) :: ustar !< Friction velocity [Z T-1 ~> m s-1].
885  real, intent(in) :: hbl !< (Positive) thickness of boundary layer [Z ~> m].
886  logical, optional, intent(in) :: override_ma !< Override to use misalignment in LA
887  !! calculation. This can be used if diagnostic
888  !! LA outputs are desired that are different than
889  !! those used by the dynamical model.
890  real, dimension(SZK_(GV)), optional, &
891  intent(in) :: h !< Grid layer thickness [H ~> m or kg m-2]
892  real, dimension(SZK_(GV)), optional, &
893  intent(in) :: u_h !< Zonal velocity at H point [m s-1]
894  real, dimension(SZK_(GV)), optional, &
895  intent(in) :: v_h !< Meridional velocity at H point [m s-1]
896  type(wave_parameters_cs), &
897  pointer :: waves !< Surface wave control structure.
898 
899  real, intent(out) :: la !< Langmuir number
900 
901 !Local Variables
902  real :: top, bottom, midpoint
903  real :: dpt_lasl, sheardirection, wavedirection
904  real :: la_stkx, la_stky, la_stk ! Stokes velocities in [m s-1]
905  logical :: continueloop, use_ma
906  real, dimension(SZK_(G)) :: us_h, vs_h
907  real, dimension(NumBands) :: stkband_x, stkband_y
908  integer :: kk, bb
909 
910  ! Compute averaging depth for Stokes drift (negative)
911  dpt_lasl = min(-0.1*us%m_to_Z, -la_frachbl*hbl)
912 
913  use_ma = la_misalignment
914  if (present(override_ma)) use_ma = override_ma
915 
916  ! If requesting to use misalignment in the Langmuir number compute the Shear Direction
917  if (use_ma) then
918  if (.not.(present(h).and.present(u_h).and.present(v_h))) then
919  call mom_error(fatal,'Get_LA_waves requested to consider misalignment.')
920  endif
921  continueloop = .true.
922  bottom = 0.0
923  do kk = 1,g%ke
924  top = bottom
925  midpoint = bottom + gv%H_to_Z*0.5*h(kk)
926  bottom = bottom + gv%H_to_Z*h(kk)
927  if (midpoint > dpt_lasl .and. kk > 1 .and. continueloop) then
928  sheardirection = atan2(v_h(1)-v_h(kk),u_h(1)-u_h(kk))
929  continueloop = .false.
930  endif
931  enddo
932  endif
933 
934  if (wavemethod==testprof) then
935  do kk = 1,g%ke
936  us_h(kk) = 0.5*(waves%US_X(i,j,kk)+waves%US_X(i-1,j,kk))
937  vs_h(kk) = 0.5*(waves%US_Y(i,j,kk)+waves%US_Y(i,j-1,kk))
938  enddo
939  call get_sl_average_prof( gv, dpt_lasl, h, us_h, la_stkx)
940  call get_sl_average_prof( gv, dpt_lasl, h, vs_h, la_stky)
941  la_stk = sqrt(la_stkx*la_stkx+la_stky*la_stky)
942  elseif (wavemethod==surfbands) then
943  do bb = 1,numbands
944  stkband_x(bb) = 0.5*(waves%STKx0(i,j,bb)+waves%STKx0(i-1,j,bb))
945  stkband_y(bb) = 0.5*(waves%STKy0(i,j,bb)+waves%STKy0(i,j-1,bb))
946  enddo
947  call get_sl_average_band(gv, dpt_lasl, numbands, waves%WaveNum_Cen, stkband_x, la_stkx )
948  call get_sl_average_band(gv, dpt_lasl, numbands, waves%WaveNum_Cen, stkband_y, la_stky )
949  la_stk = sqrt(la_stkx**2 + la_stky**2)
950  elseif (wavemethod==dhh85) then
951  ! Temporarily integrating profile rather than spectrum for simplicity
952  do kk = 1,gv%ke
953  us_h(kk) = 0.5*(waves%US_X(i,j,kk)+waves%US_X(i-1,j,kk))
954  vs_h(kk) = 0.5*(waves%US_Y(i,j,kk)+waves%US_Y(i,j-1,kk))
955  enddo
956  call get_sl_average_prof( gv, dpt_lasl, h, us_h, la_stkx)
957  call get_sl_average_prof( gv, dpt_lasl, h, vs_h, la_stky)
958  la_stk = sqrt(la_stkx**2 + la_stky**2)
959  elseif (wavemethod==lf17) then
960  call get_stokessl_lifoxkemper(ustar, hbl*la_frachbl, gv, us, la_stk, la)
961  elseif (wavemethod==null_wavemethod) then
962  call mom_error(fatal, "Get_Langmuir_number called without defining a WaveMethod. "//&
963  "Suggest to make sure USE_LT is set/overridden to False or "//&
964  "choose a wave method (or set USE_LA_LI2016 to use statistical "//&
965  "waves.")
966  endif
967 
968  if (.not.(wavemethod==lf17)) then
969  ! This is an arbitrary lower bound on Langmuir number.
970  ! We shouldn't expect values lower than this, but
971  ! there is also no good reason to cap it here other then
972  ! to prevent large enhancements in unconstrained parts of
973  ! the curve fit parameterizations.
974  ! Note the dimensional constant background Stokes velocity of 10^-10 m s-1.
975  la = max(waves%La_min, sqrt(us%Z_to_m*us%s_to_T*ustar / (la_stk+1.e-10)))
976  endif
977 
978  if (use_ma) then
979  wavedirection = atan2(la_stky, la_stkx)
980  la = la / sqrt(max(1.e-8, cos( wavedirection - sheardirection)))
981  endif
982 
983  return
984 end subroutine get_langmuir_number
985 
986 !> Get SL averaged Stokes drift from Li/FK 17 method
987 !!
988 !! Original description:
989 !! - This function returns the enhancement factor, given the 10-meter
990 !! wind [m s-1], friction velocity [m s-1] and the boundary layer depth [m].
991 !!
992 !! Update (Jan/25):
993 !! - Converted from function to subroutine, now returns Langmuir number.
994 !! - Computs 10m wind internally, so only ustar and hbl need passed to
995 !! subroutine.
996 !!
997 !! Qing Li, 160606
998 !! - BGR port from CVMix to MOM6 Jan/25/2017
999 !! - BGR change output to LA from Efactor
1000 !! - BGR remove u10 input
1001 !! - BGR note: fixed parameter values should be changed to "get_params"
1002 subroutine get_stokessl_lifoxkemper(ustar, hbl, GV, US, UStokes_SL, LA)
1003  real, intent(in) :: ustar !< water-side surface friction velocity [Z T-1 ~> m s-1].
1004  real, intent(in) :: hbl !< boundary layer depth [Z ~> m].
1005  type(verticalgrid_type), intent(in) :: GV !< Ocean vertical grid structure
1006  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1007  real, intent(out) :: UStokes_SL !< Surface layer averaged Stokes drift [m s-1]
1008  real, intent(out) :: LA !< Langmuir number
1009  ! Local variables
1010  ! parameters
1011  real, parameter :: &
1012  ! ratio of U19.5 to U10 (Holthuijsen, 2007)
1013  u19p5_to_u10 = 1.075, &
1014  ! ratio of mean frequency to peak frequency for
1015  ! Pierson-Moskowitz spectrum (Webb, 2011)
1016  fm_into_fp = 1.296, &
1017  ! ratio of surface Stokes drift to U10
1018  us_to_u10 = 0.0162, &
1019  ! loss ratio of Stokes transport
1020  r_loss = 0.667
1021  real :: UStokes, hm0, fm, fp, vstokes, kphil, kstar
1022  real :: z0, z0i, r1, r2, r3, r4, tmp, lasl_sqr_i
1023  real :: u10
1024 
1025  if (ustar > 0.0) then
1026  ! Computing u10 based on u_star and COARE 3.5 relationships
1027  call ust_2_u10_coare3p5(us%Z_to_m*us%s_to_T*ustar*sqrt(us%R_to_kg_m3*gv%Rho0/1.225), u10, gv, us)
1028  ! surface Stokes drift
1029  ustokes = us_to_u10*u10
1030  !
1031  ! significant wave height from Pierson-Moskowitz
1032  ! spectrum (Bouws, 1998)
1033  hm0 = 0.0246 *u10**2
1034  !
1035  ! peak frequency (PM, Bouws, 1998)
1036  tmp = 2.0 * pi * u19p5_to_u10 * u10
1037  fp = 0.877 * gv%mks_g_Earth / tmp
1038  !
1039  ! mean frequency
1040  fm = fm_into_fp * fp
1041  !
1042  ! total Stokes transport (a factor r_loss is applied to account
1043  ! for the effect of directional spreading, multidirectional waves
1044  ! and the use of PM peak frequency and PM significant wave height
1045  ! on estimating the Stokes transport)
1046  vstokes = 0.125 * pi * r_loss * fm * hm0**2
1047  !
1048  ! the general peak wavenumber for Phillips' spectrum
1049  ! (Breivik et al., 2016) with correction of directional spreading
1050  kphil = 0.176 * ustokes / vstokes
1051  !
1052  ! surface layer averaged Stokes dirft with Stokes drift profile
1053  ! estimated from Phillips' spectrum (Breivik et al., 2016)
1054  ! the directional spreading effect from Webb and Fox-Kemper, 2015
1055  ! is also included
1056  kstar = kphil * 2.56
1057  ! surface layer
1058  z0 = abs(us%Z_to_m*hbl)
1059  z0i = 1.0 / z0
1060  ! term 1 to 4
1061  r1 = ( 0.151 / kphil * z0i -0.84 ) * &
1062  ( 1.0 - exp(-2.0 * kphil * z0) )
1063  r2 = -( 0.84 + 0.0591 / kphil * z0i ) * &
1064  sqrt( 2.0 * pi * kphil * z0 ) * &
1065  erfc( sqrt( 2.0 * kphil * z0 ) )
1066  r3 = ( 0.0632 / kstar * z0i + 0.125 ) * &
1067  (1.0 - exp(-2.0 * kstar * z0) )
1068  r4 = ( 0.125 + 0.0946 / kstar * z0i ) * &
1069  sqrt( 2.0 * pi *kstar * z0) * &
1070  erfc( sqrt( 2.0 * kstar * z0 ) )
1071  ustokes_sl = ustokes * (0.715 + r1 + r2 + r3 + r4)
1072  la = sqrt(us%Z_to_m*us%s_to_T*ustar / ustokes_sl)
1073  else
1074  ustokes_sl = 0.0
1075  la=1.e8
1076  endif
1077 
1078 end subroutine get_stokessl_lifoxkemper
1079 
1080 !> Get SL Averaged Stokes drift from a Stokes drift Profile
1081 subroutine get_sl_average_prof( GV, AvgDepth, H, Profile, Average )
1083  intent(in) :: GV !< Ocean vertical grid structure
1084  real, intent(in) :: AvgDepth !< Depth to average over [Z ~> m].
1085  real, dimension(SZK_(GV)), &
1086  intent(in) :: H !< Grid thickness [H ~> m or kg m-2]
1087  real, dimension(SZK_(GV)), &
1088  intent(in) :: Profile !< Profile of quantity to be averaged [arbitrary]
1089  !! (used here for Stokes drift)
1090  real, intent(out) :: Average !< Output quantity averaged over depth AvgDepth [arbitrary]
1091  !! (used here for Stokes drift)
1092  !Local variables
1093  real :: top, midpoint, bottom ! Depths [Z ~> m].
1094  real :: Sum
1095  integer :: kk
1096 
1097  ! Initializing sum
1098  sum = 0.0
1099 
1100  ! Integrate
1101  bottom = 0.0
1102  do kk = 1, gv%ke
1103  top = bottom
1104  midpoint = bottom - gv%H_to_Z * 0.5*h(kk)
1105  bottom = bottom - gv%H_to_Z * h(kk)
1106  if (avgdepth < bottom) then !Whole cell within H_LA
1107  sum = sum + profile(kk) * (gv%H_to_Z * h(kk))
1108  elseif (avgdepth < top) then !partial cell within H_LA
1109  sum = sum + profile(kk) * (top-avgdepth)
1110  endif
1111  enddo
1112 
1113  ! Divide by AvgDepth !### Consider dividing by the depth in the column if that is smaller. -RWH
1114  average = sum / abs(avgdepth)
1115 
1116  return
1117 end subroutine get_sl_average_prof
1118 
1119 !> Get SL averaged Stokes drift from the banded Spectrum method
1120 subroutine get_sl_average_band( GV, AvgDepth, NB, WaveNumbers, SurfStokes, Average )
1122  intent(in) :: GV !< Ocean vertical grid
1123  real, intent(in) :: AvgDepth !< Depth to average over [Z ~> m].
1124  integer, intent(in) :: NB !< Number of bands used
1125  real, dimension(NB), &
1126  intent(in) :: WaveNumbers !< Wavenumber corresponding to each band [Z-1 ~> m-1]
1127  real, dimension(NB), &
1128  intent(in) :: SurfStokes !< Surface Stokes drift for each band [m s-1]
1129  real, intent(out) :: Average !< Output average Stokes drift over depth AvgDepth [m s-1]
1130 
1131  ! Local variables
1132  integer :: bb
1133 
1134  ! Loop over bands
1135  average = 0.0
1136  do bb = 1, nb
1137  ! Factor includes analytical integration of e(2kz)
1138  ! - divided by (-H_LA) to get average from integral.
1139  average = average + surfstokes(bb) * &
1140  (1.-exp(-abs(avgdepth * 2.0 * wavenumbers(bb)))) / &
1141  abs(avgdepth * 2.0 * wavenumbers(bb))
1142  enddo
1143 
1144  return
1145 end subroutine get_sl_average_band
1146 
1147 !> Compute the Stokes drift at a given depth
1148 !!
1149 !! Taken from Qing Li (Brown)
1150 !! use for comparing MOM6 simulation to his LES
1151 !! computed at z mid point (I think) and not depth averaged.
1152 !! Should be fine to integrate in frequency from 0.1 to sqrt(-0.2*grav*2pi/dz
1153 subroutine dhh85_mid(GV, US, zpt, UStokes)
1154  type(verticalgrid_type), intent(in) :: GV !< Ocean vertical grid
1155  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1156  real, intent(in) :: ZPT !< Depth to get Stokes drift [Z ~> m]. !### THIS IS NOT USED YET.
1157  real, intent(out) :: UStokes !< Stokes drift [m s-1]
1158  !
1159  real :: ann, Bnn, Snn, Cnn, Dnn
1160  real :: omega_peak, omega, u10, WA, domega
1161  real :: omega_min, omega_max, wavespec, Stokes
1162  integer :: Nomega, OI
1163 
1164  wa = waveage
1165  u10 = wavewind
1166 
1167  !/
1168  omega_min = 0.1 ! Hz
1169  ! Cut off at 30cm for now...
1170  omega_max = 10. ! ~sqrt(0.2*GV%mks_g_Earth*2*pi/0.3)
1171  nomega = 1000
1172  domega = (omega_max-omega_min)/real(nomega)
1173 
1174  !
1175  if (waveagepeakfreq) then
1176  omega_peak = gv%mks_g_Earth / (wa * u10)
1177  else
1178  omega_peak = 2. * pi * 0.13 * gv%mks_g_Earth / u10
1179  endif
1180  !/
1181  ann = 0.006 * waveage**(-0.55)
1182  bnn = 1.0
1183  snn = 0.08 * (1.0 + 4.0 * waveage**3)
1184  cnn = 1.7
1185  if (wa < 1.) then
1186  cnn = cnn - 6.0*log10(wa)
1187  endif
1188  !/
1189  ustokes = 0.0
1190  omega = omega_min + 0.5*domega
1191  do oi = 1,nomega-1
1192  dnn = exp( -0.5 * (omega-omega_peak)**2 / (snn**2 * omega_peak**2) )
1193  ! wavespec units = m2s
1194  wavespec = (ann * gv%mks_g_Earth**2 / (omega_peak*omega**4 ) ) * &
1195  exp(-bnn*(omega_peak/omega)**4)*cnn**dnn
1196  ! Stokes units m (multiply by frequency range for units of m/s)
1197  stokes = 2.0 * wavespec * omega**3 * &
1198  exp( 2.0 * omega**2 * zpt / gv%mks_g_Earth) / gv%mks_g_Earth
1199  ustokes = ustokes + stokes*domega
1200  omega = omega + domega
1201  enddo
1202 
1203  return
1204 end subroutine dhh85_mid
1205 
1206 !> Explicit solver for Stokes mixing.
1207 !! Still in development do not use.
1208 subroutine stokesmixing(G, GV, dt, h, u, v, Waves )
1210  intent(in) :: g !< Ocean grid
1211  type(verticalgrid_type), &
1212  intent(in) :: gv !< Ocean vertical grid
1213  real, intent(in) :: dt !< Time step of MOM6 [T ~> s] for explicit solver
1214  real, dimension(SZI_(G),SZJ_(G),SZK_(G)),&
1215  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
1216  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
1217  intent(inout) :: u !< Velocity i-component [m s-1]
1218  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
1219  intent(inout) :: v !< Velocity j-component [m s-1]
1220  type(wave_parameters_cs), &
1221  pointer :: waves !< Surface wave related control structure.
1222  ! Local variables
1223  real :: dtauup, dtaudn ! Vertical momentum fluxes [Z T-1 m s-1]
1224  real :: h_lay ! The layer thickness at a velocity point [Z ~> m].
1225  integer :: i,j,k
1226 
1227 ! This is a template to think about down-Stokes mixing.
1228 ! This is not ready for use...
1229 
1230  do k = 1, g%ke
1231  do j = g%jsc, g%jec
1232  do i = g%iscB, g%iecB
1233  h_lay = gv%H_to_Z*0.5*(h(i,j,k)+h(i+1,j,k))
1234  dtauup = 0.0
1235  if (k > 1) &
1236  dtauup = 0.5*(waves%Kvs(i,j,k)+waves%Kvs(i+1,j,k)) * &
1237  (waves%us_x(i,j,k-1)-waves%us_x(i,j,k)) / &
1238  (0.5*(h_lay + gv%H_to_Z*0.5*(h(i,j,k-1)+h(i+1,j,k-1)) ))
1239  dtaudn = 0.0
1240  if (k < g%ke-1) &
1241  dtaudn = 0.5*(waves%Kvs(i,j,k+1)+waves%Kvs(i+1,j,k+1)) * &
1242  (waves%us_x(i,j,k)-waves%us_x(i,j,k+1)) / &
1243  (0.5*(h_lay + gv%H_to_Z*0.5*(h(i,j,k+1)+h(i+1,j,k+1)) ))
1244  u(i,j,k) = u(i,j,k) + dt * (dtauup-dtaudn) / h_lay
1245  enddo
1246  enddo
1247  enddo
1248 
1249  do k = 1, g%ke
1250  do j = g%jscB, g%jecB
1251  do i = g%isc, g%iec
1252  h_lay = gv%H_to_Z*0.5*(h(i,j,k)+h(i,j+1,k))
1253  dtauup = 0.
1254  if (k > 1) &
1255  dtauup = 0.5*(waves%Kvs(i,j,k)+waves%Kvs(i,j+1,k)) * &
1256  (waves%us_y(i,j,k-1)-waves%us_y(i,j,k)) / &
1257  (0.5*(h_lay + gv%H_to_Z*0.5*(h(i,j,k-1)+h(i,j+1,k-1)) ))
1258  dtaudn = 0.0
1259  if (k < g%ke-1) &
1260  dtaudn =0.5*(waves%Kvs(i,j,k+1)+waves%Kvs(i,j+1,k+1)) * &
1261  (waves%us_y(i,j,k)-waves%us_y(i,j,k+1)) / &
1262  (0.5*(h_lay + gv%H_to_Z*0.5*(h(i,j,k+1)+h(i,j+1,k+1)) ))
1263  v(i,j,k) = v(i,j,k) + dt * (dtauup-dtaudn) / h_lay
1264  enddo
1265  enddo
1266  enddo
1267 
1268 end subroutine stokesmixing
1269 
1270 !> Solver to add Coriolis-Stokes to model
1271 !! Still in development and not meant for general use.
1272 !! Can be activated (with code intervention) for LES comparison
1273 !! CHECK THAT RIGHT TIMESTEP IS PASSED IF YOU USE THIS**
1274 !!
1275 !! Not accessed in the standard code.
1276 subroutine coriolisstokes(G, GV, DT, h, u, v, WAVES, US)
1278  intent(in) :: g !< Ocean grid
1279  type(verticalgrid_type), &
1280  intent(in) :: gv !< Ocean vertical grid
1281  real, intent(in) :: dt !< Time step of MOM6 [s] CHECK IF PASSING RIGHT TIMESTEP
1282  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
1283  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
1284  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
1285  intent(inout) :: u !< Velocity i-component [m s-1]
1286  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
1287  intent(inout) :: v !< Velocity j-component [m s-1]
1288  type(wave_parameters_cs), &
1289  pointer :: waves !< Surface wave related control structure.
1290  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
1291  ! Local variables
1292  real :: dvel ! A rescaled velocity change [m s-1 T-1 ~> m s-2]
1293  integer :: i,j,k
1294 
1295  do k = 1, g%ke
1296  do j = g%jsc, g%jec
1297  do i = g%iscB, g%iecB
1298  dvel = 0.25*(waves%us_y(i,j+1,k)+waves%us_y(i-1,j+1,k))*g%CoriolisBu(i,j+1) + &
1299  0.25*(waves%us_y(i,j,k)+waves%us_y(i-1,j,k))*g%CoriolisBu(i,j)
1300  u(i,j,k) = u(i,j,k) + dvel*us%s_to_T*dt
1301  enddo
1302  enddo
1303  enddo
1304 
1305  do k = 1, g%ke
1306  do j = g%jscB, g%jecB
1307  do i = g%isc, g%iec
1308  dvel = 0.25*(waves%us_x(i+1,j,k)+waves%us_x(i+1,j-1,k))*g%CoriolisBu(i+1,j) + &
1309  0.25*(waves%us_x(i,j,k)+waves%us_x(i,j-1,k))*g%CoriolisBu(i,j)
1310  v(i,j,k) = v(i,j,k) - dvel*us%s_to_T*dt
1311  enddo
1312  enddo
1313  enddo
1314 end subroutine coriolisstokes
1315 
1316 !> Computes wind speed from ustar_air based on COARE 3.5 Cd relationship
1317 !! Probably doesn't belong in this module, but it is used here to estimate
1318 !! wind speed for wind-wave relationships. Should be a fine way to estimate
1319 !! the neutral wind-speed as written here.
1320 subroutine ust_2_u10_coare3p5(USTair, U10, GV, US)
1321  real, intent(in) :: USTair !< Wind friction velocity [m s-1]
1322  real, intent(out) :: U10 !< 10-m neutral wind speed [m s-1]
1323  type(verticalgrid_type), intent(in) :: GV !< vertical grid type
1324  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1325 
1326  ! Local variables
1327  real, parameter :: vonkar = 0.4 ! Should access a get_param von karman
1328  real, parameter :: nu=1e-6 ! Should access a get_param air-viscosity
1329  real :: z0sm, z0, z0rough, u10a, alpha, CD
1330  integer :: CT
1331 
1332  ! Uses empirical formula for z0 to convert ustar_air to u10 based on the
1333  ! COARE 3.5 paper (Edson et al., 2013)
1334  ! alpha=m*U10+b
1335  ! Note in Edson et al. 2013, eq. 13 m is given as 0.017. However,
1336  ! m=0.0017 reproduces the curve in their figure 6.
1337 
1338  z0sm = 0.11 * nu * us%m_to_Z / ustair !Compute z0smooth from ustar guess
1339  u10 = ustair/sqrt(0.001) !Guess for u10
1340  u10a = 1000
1341 
1342  ct=0
1343  do while (abs(u10a/u10-1.) > 0.001)
1344  ct=ct+1
1345  u10a = u10
1346  alpha = min(0.028, 0.0017 * u10 - 0.005)
1347  z0rough = alpha * (us%m_s_to_L_T*ustair)**2 / gv%g_Earth ! Compute z0rough from ustar guess
1348  z0 = z0sm + z0rough
1349  cd = ( vonkar / log(10.*us%m_to_Z / z0) )**2 ! Compute CD from derived roughness
1350  u10 = ustair/sqrt(cd) ! Compute new u10 from derived CD, while loop
1351  ! ends and checks for convergence...CT counter
1352  ! makes sure loop doesn't run away if function
1353  ! doesn't converge. This code was produced offline
1354  ! and converged rapidly (e.g. 2 cycles)
1355  ! for ustar=0.0001:0.0001:10.
1356  if (ct>20) then
1357  u10 = ustair/sqrt(0.0015) ! I don't expect to get here, but just
1358  ! in case it will output a reasonable value.
1359  exit
1360  endif
1361  enddo
1362  return
1363 end subroutine ust_2_u10_coare3p5
1364 
1365 !> Clear pointers, deallocate memory
1366 subroutine waves_end(CS)
1367  type(wave_parameters_cs), pointer :: cs !< Control structure
1368 
1369  if (allocated(cs%WaveNum_Cen)) deallocate( cs%WaveNum_Cen )
1370  if (allocated(cs%Freq_Cen)) deallocate( cs%Freq_Cen )
1371  if (allocated(cs%Us_x)) deallocate( cs%Us_x )
1372  if (allocated(cs%Us_y)) deallocate( cs%Us_y )
1373  if (allocated(cs%La_SL)) deallocate( cs%La_SL )
1374  if (allocated(cs%La_turb)) deallocate( cs%La_turb )
1375  if (allocated(cs%STKx0)) deallocate( cs%STKx0 )
1376  if (allocated(cs%STKy0)) deallocate( cs%STKy0 )
1377  if (allocated(cs%KvS)) deallocate( cs%KvS )
1378  if (allocated(cs%Us0_y)) deallocate( cs%Us0_y )
1379  if (allocated(cs%Us0_x)) deallocate( cs%Us0_x )
1380 
1381  deallocate( cs )
1382 
1383  return
1384 end subroutine waves_end
1385 
1386 !> \namespace mom_wave_interface
1387 !!
1388 !! \author Brandon Reichl, 2018.
1389 !!
1390 !! This module should be moved as wave coupling progresses and
1391 !! likely will should mirror the iceberg or sea-ice model set-up.
1392 !!
1393 !! This module is meant to contain the routines to read in and
1394 !! interpret surface wave data for MOM6. In its original form, the
1395 !! capabilities include setting the Stokes drift in the model (from a
1396 !! variety of sources including prescribed, empirical, and input
1397 !! files). In short order, the plan is to also ammend the subroutine
1398 !! to accept Stokes drift information from an external coupler.
1399 !! Eventually, it will be necessary to break this file apart so that
1400 !! general wave information may be stored in the control structure
1401 !! and the Stokes drift effect can be isolated from processes such as
1402 !! sea-state dependent momentum fluxes, gas fluxes, and other wave
1403 !! related air-sea interaction and boundary layer phenomenon.
1404 !!
1405 !! The Stokes drift are stored on the C-grid with the conventional
1406 !! protocol to interpolate to the h-grid to compute Langmuir number,
1407 !! the primary quantity needed for Langmuir turbulence
1408 !! parameterizations in both the ePBL and KPP approach. This module
1409 !! also computes full 3d Stokes drift profiles, which will be useful
1410 !! if second-order type boundary layer parameterizations are
1411 !! implemented (perhaps via GOTM, work in progress).
1412 
1413 end module mom_wave_interface
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_wave_interface::waves_end
subroutine, public waves_end(CS)
Clear pointers, deallocate memory.
Definition: MOM_wave_interface.F90:1367
mom_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
mom_wave_interface::testprof
integer, parameter testprof
Undocumented parameters.
Definition: MOM_wave_interface.F90:177
mom_wave_interface::stokesmixing
subroutine, public stokesmixing(G, GV, dt, h, u, v, Waves)
Explicit solver for Stokes mixing. Still in development do not use.
Definition: MOM_wave_interface.F90:1209
mom_wave_interface::waveagepeakfreq
logical waveagepeakfreq
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:183
mom_safe_alloc
Convenience functions for safely allocating memory without accidentally reallocating pointer and caus...
Definition: MOM_safe_alloc.F90:3
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_wave_interface::dhh85_mid
subroutine dhh85_mid(GV, US, zpt, UStokes)
Compute the Stokes drift at a given depth.
Definition: MOM_wave_interface.F90:1154
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_wave_interface::wavewind
real wavewind
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:185
mom_wave_interface::coriolisstokes
subroutine, public coriolisstokes(G, GV, DT, h, u, v, WAVES, US)
Solver to add Coriolis-Stokes to model Still in development and not meant for general use....
Definition: MOM_wave_interface.F90:1277
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_diag_mediator::register_diag_field
integer function, public register_diag_field(module_name, field_name, axes_in, init_time, long_name, units, missing_value, range, mask_variant, standard_name, verbose, do_not_log, err_msg, interp_method, tile_count, cmor_field_name, cmor_long_name, cmor_units, cmor_standard_name, cell_methods, x_cell_method, y_cell_method, v_cell_method, conversion, v_extensive)
Returns the "diag_mediator" handle for a group (native, CMOR, z-coord, ...) of diagnostics derived fr...
Definition: MOM_diag_mediator.F90:1878
mom_wave_interface::get_stokessl_lifoxkemper
subroutine get_stokessl_lifoxkemper(ustar, hbl, GV, US, UStokes_SL, LA)
Get SL averaged Stokes drift from Li/FK 17 method.
Definition: MOM_wave_interface.F90:1003
mom_wave_interface::update_surface_waves
subroutine, public update_surface_waves(G, GV, US, Day, dt, CS)
Subroutine that handles updating of surface wave/Stokes drift related properties.
Definition: MOM_wave_interface.F90:437
mom_wave_interface::tp_stky0
real tp_stky0
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:182
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_wave_interface::wavemethod
integer wavemethod
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:131
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_wave_interface::numbands
integer, public numbands
Number of wavenumber/frequency partitions to receive This needs to match the number of bands provided...
Definition: MOM_wave_interface.F90:141
mom_wave_interface::pi
real pi
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:186
mom_wave_interface::dataovr
integer, parameter dataovr
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_wave_interface::partitionmode
integer, public partitionmode
Method for partition mode (meant to check input) 0 - wavenumbers 1 - frequencies.
Definition: MOM_wave_interface.F90:145
mom_wave_interface::mdl
character(len=40) mdl
This module's name.
Definition: MOM_wave_interface.F90:171
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_wave_interface::coupler
integer, parameter coupler
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_wave_interface::la_frachbl
real la_frachbl
Fraction of OSBL for averaging Langmuir number.
Definition: MOM_wave_interface.F90:163
mom_wave_interface::mom_wave_interface_init
subroutine, public mom_wave_interface_init(time, G, GV, US, param_file, CS, diag)
Initializes parameters related to MOM_wave_interface.
Definition: MOM_wave_interface.F90:193
mom_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_wave_interface
Interface for surface waves.
Definition: MOM_wave_interface.F90:2
mom_wave_interface::wave_parameters_cs
Container for all surface wave related parameters.
Definition: MOM_wave_interface.F90:47
mom_wave_interface::get_sl_average_band
subroutine get_sl_average_band(GV, AvgDepth, NB, WaveNumbers, SurfStokes, Average)
Get SL averaged Stokes drift from the banded Spectrum method.
Definition: MOM_wave_interface.F90:1121
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_wave_interface::dhh85
integer, parameter dhh85
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_wave_interface::staticwaves
logical staticwaves
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:184
mom_wave_interface::lf17
integer, parameter lf17
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_wave_interface::surface_bands_by_data_override
subroutine surface_bands_by_data_override(day_center, G, GV, US, CS)
A subroutine to fill the Stokes drift from a NetCDF file using the data_override procedures.
Definition: MOM_wave_interface.F90:709
mom_wave_interface::null_wavemethod
integer, parameter null_wavemethod
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_wave_interface::mom_wave_interface_init_lite
subroutine, public mom_wave_interface_init_lite(param_file)
A 'lite' init subroutine to initialize a few inputs needed if using wave information with the wind-sp...
Definition: MOM_wave_interface.F90:410
mom_wave_interface::tp_wvl
real tp_wvl
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:182
mom_wave_interface::surfbandfilename
character(len=40) surfbandfilename
Filename if using DataOverride.
Definition: MOM_wave_interface.F90:157
mom_wave_interface::la_misalignment
logical la_misalignment
Flag to use misalignment in Langmuir number.
Definition: MOM_wave_interface.F90:165
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_wave_interface::input
integer, parameter input
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_wave_interface::waveage
real waveage
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:185
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_wave_interface::update_stokes_drift
subroutine, public update_stokes_drift(G, GV, US, CS, h, ustar)
Constructs the Stokes Drift profile on the model grid based on desired coupling options.
Definition: MOM_wave_interface.F90:479
mom_wave_interface::tp_stkx0
real tp_stkx0
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:182
mom_wave_interface::get_sl_average_prof
subroutine get_sl_average_prof(GV, AvgDepth, H, Profile, Average)
Get SL Averaged Stokes drift from a Stokes drift Profile.
Definition: MOM_wave_interface.F90:1082
mom_wave_interface::dhh85_is_set
logical dhh85_is_set
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:184
mom_safe_alloc::safe_alloc_ptr
Allocate a pointer to a 1-d, 2-d or 3-d array.
Definition: MOM_safe_alloc.F90:12
mom_wave_interface::datasource
integer datasource
Integer that specifies where the Model Looks for Data Valid choices are: 1 - FMS DataOverride Routine...
Definition: MOM_wave_interface.F90:149
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239
mom_wave_interface::get_langmuir_number
subroutine, public get_langmuir_number(LA, G, GV, US, HBL, ustar, i, j, H, U_H, V_H, Override_MA, Waves)
Interface to get Langmuir number based on options stored in wave structure.
Definition: MOM_wave_interface.F90:879
mom_wave_interface::ust_2_u10_coare3p5
subroutine ust_2_u10_coare3p5(USTair, U10, GV, US)
Computes wind speed from ustar_air based on COARE 3.5 Cd relationship Probably doesn't belong in this...
Definition: MOM_wave_interface.F90:1321
mom_wave_interface::surfbands
integer, parameter surfbands
Options for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stoke...
Definition: MOM_wave_interface.F90:177
mom_wave_interface::dataoverrideisinitialized
logical dataoverrideisinitialized
Flag for DataOverride Initialization.
Definition: MOM_wave_interface.F90:159