MOM_wave_drag.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> Frequency-dependent linear wave drag
6
7module mom_wave_drag
8
9use mom_domains, only : pass_vector, to_all, scalar_pair
10use mom_error_handler, only : mom_error, note
11use mom_file_parser, only : get_param, log_param, param_file_type
12use mom_grid, only : ocean_grid_type
13use mom_io, only : mom_read_data, slasher, east_face, north_face
16
17implicit none ; private
18
20
21#include <MOM_memory.h>
22
23!> Control structure for the MOM_wave_drag module
24type, public :: wave_drag_cs ; private
25 integer :: nf !< Number of filters to be used in the simulation
26 real, allocatable, dimension(:,:,:) :: coef_u !< frequency-dependent drag coefficients [H T-1 ~> m s-1]
27 real, allocatable, dimension(:,:,:) :: coef_v !< frequency-dependent drag coefficients [H T-1 ~> m s-1]
28 real, allocatable, dimension(:,:,:) :: coef_uv !< frequency-dependent drag coefficients [H T-1 ~> m s-1]
29 real, allocatable, dimension(:,:,:) :: coef_vu !< frequency-dependent drag coefficients [H T-1 ~> m s-1]
30 logical :: tensor_drag !< If true, include the off-diagonal components of the
31 !! wave drag tensor for computing the wave drag
32end type wave_drag_cs
33
34contains
35
36!> This subroutine reads drag coefficients from file.
37subroutine wave_drag_init(param_file, wave_drag_file, G, GV, US, CS)
38 type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
39 character(len=*), intent(in) :: wave_drag_file !< The file from which to read drag coefficients
40 type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure
41 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
42 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
43 type(wave_drag_cs), intent(out) :: cs !< Control structure of MOM_wave_drag
44
45 ! Local variables
46 character(len=40) :: mdl = "MOM_wave_drag" !< This module's name
47 character(len=50) :: filter_name_str !< List of drag coefficients to be used
48 character(len=2), allocatable, dimension(:) :: filter_names !< Names of drag coefficients
49 character(len=80) :: var_names(4) !< Names of variables in wave_drag_file
50 character(len=200) :: mesg
51 real :: var_scale !< Scaling factors of drag coefficients [nondim]
52 integer :: c
53
54 ! The number and names of drag coefficients should match those of the streaming filters.
55 call get_param(param_file, mdl, "N_FILTERS", cs%nf, &
56 "Number of streaming band-pass filters to be used in the simulation.", &
57 default=0, do_not_log=.true.)
58 call get_param(param_file, mdl, "FILTER_NAMES", filter_name_str, &
59 "Names of streaming band-pass filters to be used in the simulation.", &
60 do_not_log=.true.)
61
62 allocate(cs%coef_u(g%IsdB:g%IedB,g%jsd:g%jed,cs%nf)) ; cs%coef_u(:,:,:) = 0.0
63 allocate(cs%coef_v(g%isd:g%ied,g%JsdB:g%JedB,cs%nf)) ; cs%coef_v(:,:,:) = 0.0
64 allocate(cs%coef_uv(g%IsdB:g%IedB,g%jsd:g%jed,cs%nf)) ; cs%coef_uv(:,:,:) = 0.0
65 allocate(cs%coef_vu(g%isd:g%ied,g%JsdB:g%JedB,cs%nf)) ; cs%coef_vu(:,:,:) = 0.0
66 allocate(filter_names(cs%nf)) ; read(filter_name_str, *) filter_names
67
68 cs%tensor_drag = .false.
69
70 if (len_trim(wave_drag_file) > 0) then
71 do c=1,cs%nf
72 call get_param(param_file, mdl, "BT_"//trim(filter_names(c))//"_DRAG_U", &
73 var_names(1), "The name of the variable in BT_WAVE_DRAG_FILE "//&
74 "for the drag coefficient of the "//trim(filter_names(c))//&
75 " frequency at u points.", default="")
76 call get_param(param_file, mdl, "BT_"//trim(filter_names(c))//"_DRAG_V", &
77 var_names(2), "The name of the variable in BT_WAVE_DRAG_FILE "//&
78 "for the drag coefficient of the "//trim(filter_names(c))//&
79 " frequency at v points.", default="")
80 call get_param(param_file, mdl, "BT_"//trim(filter_names(c))//"_DRAG_UV", &
81 var_names(3), "The name of the variable in BT_WAVE_DRAG_FILE "//&
82 "for the drag coefficient of the "//trim(filter_names(c))//&
83 " frequency at u points, corresponding to the off-diagonal "//&
84 "component of the wave drag tensor.", default="")
85 call get_param(param_file, mdl, "BT_"//trim(filter_names(c))//"_DRAG_VU", &
86 var_names(4), "The name of the variable in BT_WAVE_DRAG_FILE "//&
87 "for the drag coefficient of the "//trim(filter_names(c))//&
88 " frequency at v points, corresponding to the off-diagonal "//&
89 "component of the wave drag tensor.", default="")
90 call get_param(param_file, mdl, "BT_"//trim(filter_names(c))//"_DRAG_SCALE", &
91 var_scale, "A scaling factor for the drag coefficient of the "//&
92 trim(filter_names(c))//" frequency.", default=1.0, units="nondim")
93
94 if (len_trim(var_names(1))>0 .and. len_trim(var_names(2))>0 .and. var_scale>0.0) then
95 call mom_read_data(wave_drag_file, trim(var_names(1)), cs%coef_u(:,:,c), g%Domain, &
96 position=east_face, scale=var_scale*gv%m_to_H*us%T_to_s)
97 call mom_read_data(wave_drag_file, trim(var_names(2)), cs%coef_v(:,:,c), g%Domain, &
98 position=north_face, scale=var_scale*gv%m_to_H*us%T_to_s)
99 call pass_vector(cs%coef_u(:,:,c), cs%coef_v(:,:,c), g%domain, &
100 direction=to_all+scalar_pair)
101
102 if (len_trim(var_names(3))>0 .and. len_trim(var_names(4))>0) then
103 cs%tensor_drag = .true.
104
105 call mom_read_data(wave_drag_file, trim(var_names(3)), cs%coef_uv(:,:,c), g%Domain, &
106 position=east_face, scale=var_scale*gv%m_to_H*us%T_to_s)
107 call mom_read_data(wave_drag_file, trim(var_names(4)), cs%coef_vu(:,:,c), g%Domain, &
108 position=north_face, scale=var_scale*gv%m_to_H*us%T_to_s)
109 call pass_vector(cs%coef_uv(:,:,c), cs%coef_vu(:,:,c), g%domain, &
110 direction=to_all+scalar_pair)
111 endif
112
113 write(mesg, *) "MOM_wave_drag: ", trim(filter_names(c)), &
114 " coefficients read from file, scaling factor = ", var_scale
115 call mom_error(note, trim(mesg))
116 endif ! (len_trim(var_names(1))+len_trim(var_names(2))>0 .and. var_scale>0.0)
117 enddo ! k=1,CS%nf
118 endif ! (len_trim(wave_drag_file) > 0)
119
120end subroutine wave_drag_init
121
122!> This subroutine calculates the sum of the products of the tidal velocities and the scaled
123!! frequency-dependent drag for each tidal constituent specified in MOM_input.
124subroutine wave_drag_calc(u, v, drag_u, drag_v, G, CS)
125 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
126 type(wave_drag_cs), intent(in) :: cs !< Control structure of MOM_wave_drag
127 real, dimension(:,:,:), pointer, intent(in) :: u !< Zonal velocity from the output of
128 !! streaming band-pass filters [L T-1 ~> m s-1]
129 real, dimension(:,:,:), pointer, intent(in) :: v !< Meridional velocity from the output of
130 !! streaming band-pass filters [L T-1 ~> m s-1]
131 real, dimension(G%IsdB:G%IedB,G%jsd:G%jed), intent(out) :: drag_u !< Sum of products of filtered velocities
132 !! and scaled frequency-dependent drag [L2 T-2 ~> m2 s-2]
133 real, dimension(G%isd:G%ied,G%JsdB:G%JedB), intent(out) :: drag_v !< Sum of products of filtered velocities
134 !! and scaled frequency-dependent drag [L2 T-2 ~> m2 s-2]
135
136 ! Local variables
137 integer :: is, ie, js, je, i, j, c
138
139 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
140
141 drag_u(:,:) = 0.0 ; drag_v(:,:) = 0.0
142
143 if (cs%tensor_drag) then
144 call pass_vector(u(:,:,1:cs%nf), v(:,:,1:cs%nf), g%domain, direction=to_all+scalar_pair)
145 !$OMP do
146 do j=js,je ; do i=is-1,ie ; do c=1,cs%nf ; if (g%mask2dCu(i,j) * cs%coef_u(i,j,c) > 0.0) then
147 drag_u(i,j) = drag_u(i,j) + (u(i,j,c) * cs%coef_u(i,j,c) + &
148 0.25 * ((v(i+1,j,c) + v(i,j-1,c)) + (v(i,j,c) + v(i+1,j-1,c))) * cs%coef_uv(i,j,c))
149 endif ; enddo ; enddo ; enddo
150 !$OMP do
151 do j=js-1,je ; do i=is,ie ; do c=1,cs%nf ; if (g%mask2dCv(i,j) * cs%coef_v(i,j,c) > 0.0) then
152 drag_v(i,j) = drag_v(i,j) + (v(i,j,c) * cs%coef_v(i,j,c) + &
153 0.25 * ((u(i-1,j,c) + u(i,j+1,c)) + (u(i,j,c) + u(i-1,j+1,c))) * cs%coef_vu(i,j,c))
154 endif ; enddo ; enddo ; enddo
155 else ! (.not.CS%tensor_drag)
156 !$OMP do
157 do j=js,je ; do i=is-1,ie ; do c=1,cs%nf ; if (g%mask2dCu(i,j) * cs%coef_u(i,j,c) > 0.0) then
158 drag_u(i,j) = drag_u(i,j) + u(i,j,c) * cs%coef_u(i,j,c)
159 endif ; enddo ; enddo ; enddo
160 !$OMP do
161 do j=js-1,je ; do i=is,ie ; do c=1,cs%nf ; if (g%mask2dCv(i,j) * cs%coef_v(i,j,c) > 0.0) then
162 drag_v(i,j) = drag_v(i,j) + v(i,j,c) * cs%coef_v(i,j,c)
163 endif ; enddo ; enddo ; enddo
164 endif ! (CS%tensor_drag)
165
166end subroutine wave_drag_calc
167
168!> \namespace mom_wave_drag
169!!
170!! By Chengzhu Xu (chengzhu.xu@oregonstate.edu) and Edward D. Zaron
171!!
172!! This module calculates the net effects of the frequency-dependent internal wave drag applied to
173!! the tidal velocities, and returns the sum of products of frequency-dependent drag coefficients
174!! and tidal velocities for each constituent to the MOM_barotropic module for further calculations.
175!! It relies on the use of MOM_streaming_filter for determining the tidal velocities. Furthermore,
176!! the number of drag coefficients cannot exceed that of the streaming filters, and the names of
177!! drag coefficients should match those of the streaming filters. The frequency-dependent drag
178!! coefficients are read from the same file for the linear drag coefficients in MOM_barotropic.
179!!
180!! Reference: Xu, C., & Zaron, E. D. (2025). Parameterization of frequency-dependent internal wave drag.
181!! Journal of Advances in Modeling Earth Systems, 17, e2025MS005126. https://doi.org/10.1029/2025MS005126
182
183end module mom_wave_drag
184