mom_wave_interface module reference
Interface for surface waves.
Data Types
Container for all surface wave related parameters. |
Functions/Subroutines
Initializes parameters related to MOM_wave_interface. |
|
Set the parameters that are used to determine the averaged Stokes drift and Langmuir numbers. |
|
This interface provides the caller with information from the waves control structure. |
|
Subroutine that handles updating of surface wave/Stokes drift related properties. |
|
Constructs the Stokes Drift profile on the model grid based on desired coupling options. |
|
Return the value of (1 - exp(-x))/x [nondim], using an accurate expression for small values of x. |
|
Return the value of (1 - exp(-x)) [nondim], using an accurate expression for small values of x. |
|
A subroutine to fill the Stokes drift from a NetCDF file using the data_override procedures. |
|
Interface to get Langmuir number based on options stored in wave structure. |
|
function to return the wave method string set in the param file |
|
Get SL averaged Stokes drift from Li/FK 17 method. |
|
Get SL Averaged Stokes drift from a Stokes drift Profile. |
|
Get SL averaged Stokes drift from the banded Spectrum method. |
|
Compute the Stokes drift at a given depth. |
|
Explicit solver for Stokes mixing. |
|
Solver to add Coriolis-Stokes to model Still in development and not meant for general use. |
|
Computes tendency due to Stokes pressure gradient force anomaly including analytical integration of Stokes shear using multiple-exponential decay Stokes drift profile and vertical integration of the resulting pressure anomaly to the total pressure gradient force. |
|
Computes wind speed from ustar_air based on COARE 3.5 Cd relationship Probably doesn't belong in this module, but it is used here to estimate wind speed for wind-wave relationships. |
|
Clear pointers, deallocate memory. |
|
Register wave restart fields. |
Detailed Description
This module should be moved as wave coupling progresses and likely will should mirror the iceberg or sea-ice model set-up.
This module is meant to contain the routines to read in and interpret surface wave data for MOM6. In its original form, the capabilities include setting the Stokes drift in the model (from a variety of sources including prescribed, empirical, and input files). In short order, the plan is to also amend the subroutine to accept Stokes drift information from an external coupler. Eventually, it will be necessary to break this file apart so that general wave information may be stored in the control structure and the Stokes drift effect can be isolated from processes such as sea-state dependent momentum fluxes, gas fluxes, and other wave related air-sea interaction and boundary layer phenomenon.
The Stokes drift are stored on the C-grid with the conventional protocol to interpolate to the h-grid to compute Langmuir number, the primary quantity needed for Langmuir turbulence parameterizations in both the ePBL and KPP approach. This module also computes full 3d Stokes drift profiles, which will be useful if second-order type boundary layer parameterizations are implemented (perhaps via GOTM, work in progress).
Type Documentation
- type mom_wave_interface/wave_parameters_cs
Container for all surface wave related parameters.
- Type fields:
% id_pfu_stokes ::
integer, publicDiagnostic handles.% id_pfv_stokes ::
integer, publicDiagnostic handles.% id_3dstokes_x_from_ddt ::
integer, publicDiagnostic handles.% id_3dstokes_y_from_ddt ::
integer, publicDiagnostic handles.% id_p_deltastokes_l ::
integerDiagnostic handles.% id_p_deltastokes_i ::
integerDiagnostic handles.% id_surfacestokes_x ::
integerDiagnostic handles.% id_surfacestokes_y ::
integerDiagnostic handles.% id_3dstokes_x ::
integerDiagnostic handles.% id_3dstokes_y ::
integerDiagnostic handles.% id_ddt_3dstokes_x ::
integerDiagnostic handles.% id_ddt_3dstokes_y ::
integerDiagnostic handles.% id_la_turb ::
integerDiagnostic handles.% usewaves ::
logical, publicFlag to enable surface gravity wave feature.% stokes_vf ::
logical, publicTrue if Stokes vortex force is used.% passive_stokes_vf ::
logical, publicComputes Stokes VF, but doesn’t affect dynamics.% stokes_pgf ::
logical, publicTrue if Stokes shear pressure Gradient force is used.% robust_stokes_pgf ::
logical, publicIf true, use expressions to calculate the Stokes-induced pressure gradient anomalies that are more accurate in the limit of thin layers.% passive_stokes_pgf ::
logical, publicKeeps Stokes_PGF on, but doesn’t affect dynamics.% stokes_ddt ::
logical, publicDevelopmental: True if Stokes d/dt is used.% passive_stokes_ddt ::
logical, publicKeeps Stokes_DDT on, but doesn’t affect dynamics.% homogenize_surfbands ::
logicalTrue to homogenize surface band Stokes drift in the horizontal.% us_x ::
real, dimension(:,:,:), allocatable, public3d zonal Stokes drift profile [L T-1 ~> m s-1]% us_y ::
real, dimension(:,:,:), allocatable, public3d meridional Stokes drift profile [L T-1 ~> m s-1]% ddt_us_x ::
real, dimension(:,:,:), allocatable, public3d time tendency of zonal Stokes drift profile [L T-2 ~> m s-2]% ddt_us_y ::
real, dimension(:,:,:), allocatable, public3d time tendency of meridional Stokes drift profile [L T-2 ~> m s-2]% us_x_from_ddt ::
real, dimension(:,:,:), allocatable, publicCheck of 3d zonal Stokes drift profile [L T-1 ~> m s-1].% us_y_from_ddt ::
real, dimension(:,:,:), allocatable, publicCheck of 3d meridional Stokes drift profile [L T-1 ~> m s-1].% us_x_prev ::
real, dimension(:,:,:), allocatable, public3d zonal Stokes drift profile, previous dynamics call [L T-1 ~> m s-1]% us_y_prev ::
real, dimension(:,:,:), allocatable, public3d meridional Stokes drift profile, previous dynamics call [L T-1 ~> m s-1]% kvs ::
real, dimension(:,:,:), allocatable, publicViscosity for Stokes Drift shear [H Z T-1 ~> m2 s-1 or Pa s].% wavenum_cen ::
real, dimension(:), allocatable, publicWavenumber bands for read/coupled [Z-1 ~> m-1].% ustk_hb ::
real, dimension(:,:,:), allocatable, publicSurface Stokes Drift spectrum (zonal) [L T-1 ~> m s-1].% vstk_hb ::
real, dimension(:,:,:), allocatable, publicSurface Stokes Drift spectrum (meridional) [L T-1 ~> m s-1].% omega_w2x ::
real, dimension(:,:), allocatable, publicwind direction ccw from model x- axis [nondim radians]% numbands ::
integer, publicNumber of wavenumber/frequency partitions Must match the number of bands provided via either coupling or file.% wavemethod ::
integerOptions for including wave information Valid (tested) choices are: 0 - Test Profile 1 - Surface Stokes Drift Bands 2 - DHH85 3 - LF17 -99 - No waves computed, but empirical Langmuir number used.% lagrangianmixing ::
logicalThis feature is in development and not ready True if Stokes drift is present and mixing should be applied to Lagrangian current (mean current + Stokes drift). See Reichl et al., 2016 KPP-LT approach.% stokesmixing ::
logicalThis feature is in development and not ready. True if vertical mixing of momentum should be applied directly to Stokes current (with separate mixing parameter for Eulerian mixing contribution). See Harcourt 2013, 2015 Second-Moment approach.% coriolisstokes ::
logicalThis feature is in development and not ready.% stokes_min_thick_avg ::
realA layer thickness below which the cell-center Stokes drift is used instead of the cell average [Z ~> m]. This is only used if WAVE_INTERFACE_ANSWER_DATE < 20230101.% answer_date ::
integerThe vintage of the order of arithmetic and expressions in the surface wave calculations. Values below 20230101 recover the answers from the end of 2022, while higher values use updated and more robust forms of the same expressions.% partitionmode ::
integerMethod for partition mode (meant to check input) 0 - wavenumbers 1 - frequencies.% datasource ::
integerInteger that specifies where the model Looks for data Valid choices are: 1 - FMS DataOverride Routine 2 - Reserved For Coupler 3 - User input (fixed values, useful for 1d testing)% surfbandfilename ::
character(len=40)Filename if using DataOverride.% land_speed ::
realA large Stokes velocity that can be used to indicate land values in a data override file [L T-1 ~> m s-1]. Stokes drift components larger than this are set to zero in data override calls for the Stokes drift.% dataover_initialized ::
logicalFlag for DataOverride Initialization.% la_frachbl ::
realFraction of OSBL for averaging Langmuir number [nondim].% la_hbl_min ::
realMinimum boundary layer depth for averaging Langmuir number [Z ~> m].% la_misalignment ::
logicalFlag to use misalignment in Langmuir number.% la_misalign_bug ::
logicalFlag to use code with a sign error when calculating the misalignment between the shear and waves in the Langmuir number calculation.% g_earth ::
realThe gravitational acceleration, equivalent to GVg_Earth but with different dimensional rescaling appropriate for deep-water gravity waves [Z T-2 ~> m s-2].% i_g_earth ::
realThe inverse of the gravitational acceleration, with dimensional rescaling appropriate for deep-water gravity waves [T2 Z-1 ~> s2 m-1].% freq_cen ::
real, dimension(:), allocatableCentral frequency for wave bands, including a factor of 2*pi [T-1 ~> s-1].% prescribedsurfstkx ::
real, dimension(:), allocatableSurface Stokes drift if prescribed [L T-1 ~> m s-1].% prescribedsurfstky ::
real, dimension(:), allocatableSurface Stokes drift if prescribed [L T-1 ~> m s-1].% la_turb ::
real, dimension(:,:), allocatableAligned Turbulent Langmuir number [nondim].% us0_x ::
real, dimension(:,:), allocatableSurface Stokes Drift (zonal) [L T-1 ~> m s-1].% us0_y ::
real, dimension(:,:), allocatableSurface Stokes Drift (meridional) [L T-1 ~> m s-1].% stkx0 ::
real, dimension(:,:,:), allocatableStokes Drift spectrum (zonal) [L T-1 ~> m s-1].% stky0 ::
real, dimension(:,:,:), allocatableStokes Drift spectrum (meridional) [L T-1 ~> m s-1].% la_min ::
realAn arbitrary lower-bound on the Langmuir number [nondim]. Langmuir number is sqrt(u_star/u_stokes). When both are small but u_star is orders of magnitude smaller, the Langmuir number could have unintended consequences. Since both are small it can be safely capped to avoid such consequences.% la_stk_backgnd ::
realA small background Stokes velocity used in the denominator of some expressions for the Langmuir number [L T-1 ~> m s-1].% vonkar ::
realThe von Karman coefficient as used in the MOM_wave_interface module [nondim].% rho_air ::
realA typical density of air at sea level, as used in wave calculations [R ~> kg m-3].% nu_air ::
realThe viscosity of air, as used in wave calculations [Z2 T-1 ~> m2 s-1].% rho_ocn ::
realA typical surface density of seawater, as used in wave calculations in comparison with the density of air [R ~> kg m-3]. The default is RHO_0.% swh_from_u10sq ::
realA factor for converting the square of the 10 m wind speed to the significant wave height [Z T2 L-2 ~> s2 m-1].% charnock_min ::
realThe minimum value of the Charnock coefficient, which relates the square of the air friction velocity divided by the gravitational acceleration to the wave roughness length [nondim].% charnock_slope_u10 ::
realThe partial derivative of the Charnock coefficient with the 10 m wind speed [T L-1 ~> s m-1]. Note that in eq. 13 of the Edson et al. 2013 describing the COARE 3.5 bulk flux algorithm, this slope is given as 0.017. However, 0.0017 reproduces the curve in their figure 6, so that is the default value used in MOM6.% charnock_intercept ::
realThe intercept of the fit for the Charnock coefficient in the limit of no wind [nondim]. Note that this can be negative because CHARNOCK_MIN will keep the final value for the Charnock coefficient from being from being negative.% tp_stkx0 ::
realTest profile x-stokes drift amplitude [L T-1 ~> m s-1].% tp_stky0 ::
realTest profile y-stokes drift amplitude [L T-1 ~> m s-1].% tp_wvl ::
realTest profile wavelength [Z ~> m].% waveagepeakfreq ::
logicalFlag to use wave age to determine the peak frequency with DHH85.% staticwaves ::
logicalFlag to disable updating DHH85 Stokes drift.% dhh85_is_set ::
logicalThe if the wave properties have been set when WaveMethod = DHH85.% waveage ::
realThe fixed wave age used with the DHH85 spectrum [nondim].% wavewind ::
realWind speed for the DHH85 spectrum [L T-1 ~> m s-1].% omega_min ::
realMinimum wave frequency with the DHH85 spectrum [T-1 ~> s-1].% omega_max ::
realMaximum wave frequency with the DHH85 spectrum [T-1 ~> s-1].% time ::
type(time_type), pointerA pointer to the ocean model’s clock.% diag ::
type(diag_ctrl), pointerA structure that is used to regulate the timing of diagnostic output.
Function/Subroutine Documentation
- subroutine mom_wave_interface/mom_wave_interface_init(time, G, GV, US, param_file, CS, diag)
Initializes parameters related to MOM_wave_interface.
- Parameters:
time ::
time[in] Model timeg :: [inout] Grid structure
gv :: [in] Vertical grid structure
us :: [in] A dimensional unit scaling type
param_file ::
param_file[in] Input parameter structurecs :: Wave parameter control structure
diag ::
diag[inout] Diagnostic Pointer
- Call to:
couplerdataovrdhh85dhh85_stringefactorefactor_stringinputlf17lf17_stringmom_error_handler::mom_errornull_stringnull_wavemethodset_lf17_wave_paramssurfbandssurfbands_stringtestproftestprof_string- Called from:
mom6
- subroutine mom_wave_interface/set_lf17_wave_params(param_file, mdl, GV, US, CS)
Set the parameters that are used to determine the averaged Stokes drift and Langmuir numbers.
- Parameters:
param_file ::
param_file[in] Input parameter structuremdl ::
mdl[in] A module name to use in the get_param callsgv :: [in] Vertical grid structure
us :: [in] A dimensional unit scaling type
cs :: Wave parameter control structure
- Called from:
- subroutine mom_wave_interface/query_wave_properties(CS, NumBands, WaveNumbers, US)
This interface provides the caller with information from the waves control structure.
- Parameters:
cs :: Wave parameter Control structure
numbands :: [out] If present, this returns the number of < wavenumber partitions in the wave discretization
wavenumbers :: [out] If present this returns the characteristic wavenumbers of the wave discretization [m-1] or [Z-1 ~> m-1]
us :: [in] A dimensional unit scaling type that is used to undo the dimensional scaling of the output variables, if present
- Call to:
- subroutine mom_wave_interface/update_surface_waves(G, GV, US, Time_present, dt, CS, forces)
Subroutine that handles updating of surface wave/Stokes drift related properties.
- Parameters:
cs :: Wave parameter Control structure
g :: [inout] Grid structure
gv :: [in] Vertical grid structure
us :: [in] A dimensional unit scaling type
time_present :: [in] Model Time
dt ::
dt[in] Time increment as a time-typeforces ::
forces[in] MOM_forcing_type
- Call to:
couplerdataovrinputmom_error_handler::mom_errorsurface_bands_by_data_overridesurfbandstestprof- Called from:
- subroutine mom_wave_interface/update_stokes_drift(G, GV, US, CS, dz, ustar, dt, dynamics_step)
Constructs the Stokes Drift profile on the model grid based on desired coupling options.
- Parameters:
cs :: Wave parameter Control structure
g :: [inout] Grid structure
gv :: [in] Vertical grid structure
us :: [in] A dimensional unit scaling type
dz ::
dz[in] Thickness in height units [Z ~> m]ustar ::
ustar[in] Wind friction velocity [Z T-1 ~> m s-1].dt ::
dt[in] Time-step for computing Stokes-tendency [T ~> s]dynamics_step ::
dynamics_step[in] True if this call is on a dynamics step
- Call to:
dhh85dhh85_midefactorget_langmuir_numberone_minus_exp_xsurfbandstestprof- Called from:
- function mom_wave_interface/one_minus_exp_x(x)
Return the value of (1 - exp(-x))/x [nondim], using an accurate expression for small values of x.
- Parameters:
x ::
x[in] The argument of the function ((1 - exp(-x))/x) [nondim]- Called from:
- function mom_wave_interface/one_minus_exp(x)
Return the value of (1 - exp(-x)) [nondim], using an accurate expression for small values of x.
- Parameters:
x ::
x[in] The argument of the function ((1 - exp(-x))/x) [nondim]- Called from:
- subroutine mom_wave_interface/surface_bands_by_data_override(Time, G, GV, US, CS)
A subroutine to fill the Stokes drift from a NetCDF file using the data_override procedures.
- Parameters:
time :: [in] Time to get Stokes drift bands
cs :: Wave structure
g :: [inout] Grid structure
gv :: [in] Vertical grid structure
us :: [in] A dimensional unit scaling type
- Call to:
mom_spatial_means::global_area_meanmom_error_handler::mom_error- Called from:
- subroutine mom_wave_interface/get_langmuir_number(LA, G, GV, US, HBL, ustar, i, j, dz, Waves, U_H, V_H, Override_MA)
Interface to get Langmuir number based on options stored in wave structure.
Note this can be called with an unallocated Waves pointer, which is okay if we want the wind-speed only dependent Langmuir number. Therefore, we need to be careful about what we try to access here.
- Parameters:
g :: [in] Ocean grid structure
gv :: [in] Ocean vertical grid structure
la :: [out] Langmuir number [nondim]
us :: [in] A dimensional unit scaling type
hbl :: [in] (Positive) thickness of boundary layer [Z ~> m]
ustar ::
ustar[in] Friction velocity [Z T-1 ~> m s-1]i ::
i[in] Meridional index of h-pointj ::
j[in] Zonal index of h-pointdz ::
dz[in] Grid layer thickness [Z ~> m]waves :: Surface wave control structure.
u_h :: [in] Zonal velocity at H point [L T-1 ~> m s-1] or [m s-1]
v_h :: [in] Meridional velocity at H point [L T-1 ~> m s-1] or [m s-1]
override_ma :: [in] Override to use misalignment in LA calculation. This can be used if diagnostic LA outputs are desired that are different than those used by the dynamical model.
- Call to:
dhh85get_sl_average_bandget_sl_average_profget_stokessl_lifoxkemperlf17mom_error_handler::mom_errornull_wavemethodsurfbandstestprof- Called from:
mom_energetic_pbl::epbl_columnmom_cvmix_kpp::kpp_compute_bldupdate_stokes_drift
- function mom_wave_interface/get_wave_method(CS)
function to return the wave method string set in the param file
- Parameters:
cs :: Control structure
- Call to:
dhh85dhh85_stringefactorefactor_stringlf17lf17_stringnull_stringnull_wavemethodsurfbandssurfbands_stringtestproftestprof_string
- subroutine mom_wave_interface/get_stokessl_lifoxkemper(ustar, hbl, GV, US, CS, UStokes_SL, LA)
Get SL averaged Stokes drift from Li/FK 17 method.
Original description: * This function returns the enhancement factor, given the 10-meter wind [m s-1], friction velocity [m s-1] and the boundary layer depth [m].
Update (Jan/25): * Converted from function to subroutine, now returns Langmuir number.
Compute 10m wind internally, so only ustar and hbl need passed to subroutine.
Qing Li, 160606 * BGR port from CVMix to MOM6 Jan/25/2017
BGR change output to LA from Efactor
BGR remove u10 input
BGR note: fixed parameter values should be changed to “get_params”
- Parameters:
ustar ::
ustar[in] water-side surface friction velocity [Z T-1 ~> m s-1].hbl ::
hbl[in] boundary layer depth [Z ~> m].gv :: [in] Ocean vertical grid structure
us :: [in] A dimensional unit scaling type
cs :: Wave parameter Control structure
ustokes_sl :: [out] Surface layer averaged Stokes drift [L T-1 ~> m s-1]
la :: [out] Langmuir number [nondim]
- Call to:
- Called from:
- subroutine mom_wave_interface/get_sl_average_prof(GV, AvgDepth, dz, Profile, Average)
Get SL Averaged Stokes drift from a Stokes drift Profile.
- Parameters:
gv :: [in] Ocean vertical grid structure
avgdepth :: [in] Depth to average over (negative) [Z ~> m]
dz ::
dz[in] Grid thickness [Z ~> m]profile :: [in] Profile of quantity to be averaged in arbitrary units [A]
average :: [out] Output quantity averaged over depth AvgDepth [A] (used here for Stokes drift)
- Called from:
- subroutine mom_wave_interface/get_sl_average_band(GV, AvgDepth, NB, WaveNumbers, SurfStokes, Average)
Get SL averaged Stokes drift from the banded Spectrum method.
- Parameters:
gv :: [in] Ocean vertical grid
avgdepth :: [in] Depth to average over [Z ~> m].
nb :: [in] Number of bands used
wavenumbers :: [in] Wavenumber corresponding to each band [Z-1 ~> m-1]
surfstokes :: [in] Surface Stokes drift for each band [L T-1 ~> m s-1]
average :: [out] Output average Stokes drift over depth AvgDepth [L T-1 ~> m s-1]
- Called from:
- subroutine mom_wave_interface/dhh85_mid(GV, US, CS, zpt, UStokes)
Compute the Stokes drift at a given depth.
Taken from Qing Li (Brown) use for comparing MOM6 simulation to his LES computed at z mid point (I think) and not depth averaged. Should be fine to integrate in frequency from 0.1 to sqrt(-0.2*grav*2pi/dz
- Parameters:
gv :: [in] Ocean vertical grid
us :: [in] A dimensional unit scaling type
cs :: Wave parameter Control structure
zpt ::
zpt[in] Depth to get Stokes drift [Z ~> m].ustokes :: [out] Stokes drift [L T-1 ~> m s-1]
- Called from:
- subroutine mom_wave_interface/stokesmixing(G, GV, dt, h, dz, u, v, Waves)
Explicit solver for Stokes mixing. Still in development do not use.
- Parameters:
g :: [in] Ocean grid
gv :: [in] Ocean vertical grid
dt ::
dt[in] Time step of MOM6 [T ~> s] for explicit solverh ::
h[in] Layer thicknesses [H ~> m or kg m-2]dz ::
dz[in] Vertical distance between interfaces around a layer [Z ~> m]u ::
u[inout] Velocity i-component [L T-1 ~> m s-1]v ::
v[inout] Velocity j-component [L T-1 ~> m s-1]waves :: Surface wave related control structure.
- subroutine mom_wave_interface/coriolisstokes(G, GV, dt, h, u, v, Waves)
Solver to add Coriolis-Stokes to model Still in development and not meant for general use. Can be activated (with code intervention) for LES comparison CHECK THAT RIGHT TIMESTEP IS PASSED IF YOU USE THIS**.
Not accessed in the standard code.
- Parameters:
g :: [in] Ocean grid
gv :: [in] Ocean vertical grid
dt ::
dt[in] Time step of MOM6 [T ~> s]h ::
h[in] Layer thicknesses [H ~> m or kg m-2]u ::
u[inout] Velocity i-component [L T-1 ~> m s-1]v ::
v[inout] Velocity j-component [L T-1 ~> m s-1]waves :: Surface wave related control structure.
- subroutine mom_wave_interface/stokes_pgf(G, GV, US, dz, u, v, PFu_Stokes, PFv_Stokes, CS)
Computes tendency due to Stokes pressure gradient force anomaly including analytical integration of Stokes shear using multiple-exponential decay Stokes drift profile and vertical integration of the resulting pressure anomaly to the total pressure gradient force.
- Parameters:
g :: [in] Ocean grid
gv :: [in] Ocean vertical grid
us :: [in] A dimensional unit scaling type
dz ::
dz[in] Layer thicknesses in height units [Z ~> m]u ::
u[in] Lagrangian Velocity i-component [L T-1 ~> m s-1]v ::
v[in] Lagrangian Velocity j-component [L T-1 ~> m s-1]pfu_stokes :: [out] PGF Stokes-shear i-component [L T-2 ~> m s-2]
pfv_stokes :: [out] PGF Stokes-shear j-component [L T-2 ~> m s-2]
cs :: Surface wave related control structure.
- Call to:
- Called from:
mom_dynamics_split_rk2::step_mom_dyn_split_rk2mom_dynamics_split_rk2b::step_mom_dyn_split_rk2b
- subroutine mom_wave_interface/ust_2_u10_coare3p5(USTair, U10, GV, US, CS)
Computes wind speed from ustar_air based on COARE 3.5 Cd relationship Probably doesn’t belong in this module, but it is used here to estimate wind speed for wind-wave relationships. Should be a fine way to estimate the neutral wind-speed as written here.
- Parameters:
ustair :: [in] Wind friction velocity [Z T-1 ~> m s-1]
u10 :: [out] 10-m neutral wind speed [L T-1 ~> m s-1]
gv :: [in] vertical grid type
us :: [in] A dimensional unit scaling type
cs :: Wave parameter Control structure
- Call to:
- Called from:
- subroutine mom_wave_interface/waves_end(CS)
Clear pointers, deallocate memory.
- Parameters:
cs :: Control structure
- subroutine mom_wave_interface/waves_register_restarts(CS, HI, GV, US, param_file, restart_CSp)
Register wave restart fields. To be called before MOM_wave_interface_init.
- Parameters:
cs :: Wave parameter Control structure
hi :: [inout] Grid structure
gv :: [in] Vertical grid structure
us :: [in] A dimensional unit scaling type
param_file ::
param_file[in] Input parameter structurerestart_csp :: Restart structure, data intent(inout)
- Call to: