Line data Source code
1 : module gw_front
2 :
3 : !
4 : ! This module handles gravity waves from frontal sources, and was extracted
5 : ! from gw_drag in May 2013.
6 : !
7 :
8 : use gw_utils, only: r8
9 : use gw_common, only: GWBand, pi, pver
10 :
11 : implicit none
12 : private
13 : save
14 :
15 : public :: CMSourceDesc
16 : public :: flat_cm_desc
17 : public :: gaussian_cm_desc
18 : public :: gw_cm_src
19 :
20 : ! Tuneable settings.
21 : type CMSourceDesc
22 : ! Source level.
23 : integer :: ksrc
24 : ! Level at which to check whether the frontogenesis function is above
25 : ! the critical threshold.
26 : integer :: kfront
27 : ! Frontogenesis function critical threshold.
28 : real(r8) :: frontgfc = huge(1._r8)
29 : ! The stress spectrum to produce at the source level.
30 : real(r8), allocatable :: src_tau(:)
31 : end type CMSourceDesc
32 :
33 : contains
34 :
35 : ! Create a flat profile to be launched (all wavenumbers have the same
36 : ! source strength, except that l=0 is excluded).
37 0 : function flat_cm_desc(band, ksrc, kfront, frontgfc, taubgnd) result(desc)
38 : ! Wavelengths triggered by frontogenesis.
39 : type(GWBand), intent(in) :: band
40 : ! The following are used to set the corresponding object components.
41 : integer, intent(in) :: ksrc
42 : integer, intent(in) :: kfront
43 : real(r8), intent(in) :: frontgfc
44 : ! Amount of stress to launch from each wavelength.
45 : real(r8), intent(in) :: taubgnd
46 :
47 : type(CMSourceDesc) :: desc
48 :
49 0 : desc%ksrc = ksrc
50 0 : desc%kfront = kfront
51 0 : desc%frontgfc = frontgfc
52 :
53 0 : allocate(desc%src_tau(-band%ngwv:band%ngwv))
54 0 : desc%src_tau = taubgnd
55 :
56 : ! Prohibit wavenumber 0.
57 0 : desc%src_tau(0) = 0._r8
58 :
59 0 : end function flat_cm_desc
60 :
61 : ! Create a source tau profile that is a gaussian over wavenumbers (l=0 is
62 : ! excluded).
63 1536 : function gaussian_cm_desc(band, ksrc, kfront, frontgfc, height, width) &
64 : result(desc)
65 :
66 : use shr_spfn_mod, only: erfc => shr_spfn_erfc
67 :
68 : ! Wavelengths triggered by frontogenesis.
69 : type(GWBand), intent(in) :: band
70 : ! The following are used to set the corresponding object components.
71 : integer, intent(in) :: ksrc
72 : integer, intent(in) :: kfront
73 : real(r8), intent(in) :: frontgfc
74 : ! Parameters of gaussian.
75 : real(r8), intent(in) :: height
76 : real(r8), intent(in) :: width
77 :
78 : type(CMSourceDesc) :: desc
79 :
80 : ! Bounds used to average bins of the gaussian.
81 3072 : real(r8) :: gaussian_bounds(2*band%ngwv+2)
82 :
83 : ! Wavenumber index.
84 : integer :: l
85 :
86 1536 : desc%ksrc = ksrc
87 1536 : desc%kfront = kfront
88 1536 : desc%frontgfc = frontgfc
89 :
90 4608 : allocate(desc%src_tau(-band%ngwv:band%ngwv))
91 :
92 : ! Find the boundaries of each bin.
93 101376 : gaussian_bounds(:2*band%ngwv+1) = band%cref - 0.5_r8*band%dc
94 1536 : gaussian_bounds(2*band%ngwv+2) = band%cref(band%ngwv) + 0.5_r8*band%dc
95 :
96 : ! Integral of the gaussian at bin interfaces (from the point to
97 : ! positive infinity).
98 : gaussian_bounds = &
99 101376 : [( erfc(gaussian_bounds(l)/width)*height*width*sqrt(pi)/2._r8, &
100 204288 : l = 1, 2*band%ngwv+2 )]
101 :
102 : ! Get average in each bin using integral from right to left side.
103 : desc%src_tau = &
104 102912 : (gaussian_bounds(:2*band%ngwv+1) - gaussian_bounds(2:)) / band%dc
105 :
106 : ! Prohibit wavenumber 0.
107 1536 : desc%src_tau(0) = 0._r8
108 :
109 1536 : end function gaussian_cm_desc
110 :
111 : !==========================================================================
112 2978352 : subroutine gw_cm_src(ncol, band, desc, u, v, frontgf, &
113 1489176 : src_level, tend_level, tau, ubm, ubi, xv, yv, c)
114 : use gw_utils, only: get_unit_vector, dot_2d, midpoint_interp
115 : !-----------------------------------------------------------------------
116 : ! Driver for multiple gravity wave drag parameterization.
117 : !
118 : ! The parameterization is assumed to operate only where water vapor
119 : ! concentrations are negligible in determining the density.
120 : !-----------------------------------------------------------------------
121 :
122 : !------------------------------Arguments--------------------------------
123 : ! Column dimension.
124 : integer, intent(in) :: ncol
125 :
126 : ! Wavelengths triggered by frontogenesis.
127 : type(GWBand), intent(in) :: band
128 :
129 : ! Bandification of how to produce the gravity wave bandtrum.
130 : type(CMSourceDesc), intent(in) :: desc
131 :
132 : ! Midpoint zonal/meridional winds.
133 : real(r8), intent(in) :: u(ncol,pver), v(ncol,pver)
134 : ! Frontogenesis function.
135 : real(r8), intent(in) :: frontgf(:,:)
136 :
137 : ! Indices of top gravity wave source level and lowest level where wind
138 : ! tendencies are allowed.
139 : integer, intent(out) :: src_level(ncol)
140 : integer, intent(out) :: tend_level(ncol)
141 :
142 : ! Wave Reynolds stress.
143 : real(r8), intent(out) :: tau(ncol,-band%ngwv:band%ngwv,pver+1)
144 : ! Projection of wind at midpoints and interfaces.
145 : real(r8), intent(out) :: ubm(ncol,pver), ubi(ncol,pver+1)
146 : ! Unit vectors of source wind (zonal and meridional components).
147 : real(r8), intent(out) :: xv(ncol), yv(ncol)
148 : ! Phase speeds.
149 : real(r8), intent(out) :: c(ncol,-band%ngwv:band%ngwv)
150 :
151 : !---------------------------Local Storage-------------------------------
152 : ! Column and wavenumber indices.
153 : integer :: k, l
154 :
155 : ! Whether or not to launch waves in this column.
156 2978352 : logical :: launch_wave(ncol)
157 :
158 : ! Zonal/meridional wind averaged over source region.
159 2978352 : real(r8) :: usrc(ncol), vsrc(ncol)
160 :
161 : !------------------------------------------------------------------------
162 : ! Determine the source layer wind and unit vectors, then project winds.
163 : !------------------------------------------------------------------------
164 :
165 : ! Just use the source level interface values for the source wind speed
166 : ! and direction (unit vector).
167 24865776 : src_level = desc%ksrc
168 24865776 : tend_level = desc%ksrc
169 24865776 : usrc = 0.5_r8*(u(:,desc%ksrc+1)+u(:,desc%ksrc))
170 24865776 : vsrc = 0.5_r8*(v(:,desc%ksrc+1)+v(:,desc%ksrc))
171 :
172 : ! Get the unit vector components and magnitude at the surface.
173 1489176 : call get_unit_vector(usrc, vsrc, xv, yv, ubi(:,desc%ksrc+1))
174 :
175 : ! Project the local wind at midpoints onto the source wind.
176 102753144 : do k = 1, desc%ksrc
177 102753144 : ubm(:,k) = dot_2d(u(:,k), v(:,k), xv, yv)
178 : end do
179 :
180 : ! Compute the interface wind projection by averaging the midpoint winds.
181 : ! Use the top level wind at the top interface.
182 26354952 : ubi(:,1) = ubm(:,1)
183 :
184 1489176 : ubi(:,2:desc%ksrc) = midpoint_interp(ubm(:,1:desc%ksrc))
185 :
186 : !-----------------------------------------------------------------------
187 : ! Gravity wave sources
188 : !-----------------------------------------------------------------------
189 :
190 >15207*10^7 : tau = 0._r8
191 :
192 : ! GW generation depends on frontogenesis at specified level (may be below
193 : ! actual source level).
194 24865776 : launch_wave = (frontgf(:,desc%kfront) > desc%frontgfc)
195 :
196 98285616 : do l = -band%ngwv, band%ngwv
197 1617764616 : where (launch_wave)
198 96796440 : tau(:,l,desc%ksrc+1) = desc%src_tau(l)
199 : end where
200 : end do
201 :
202 : ! Set phase speeds as reference speeds plus the wind speed at the source
203 : ! level.
204 : c = spread(band%cref, 1, ncol) + &
205 1617764616 : spread(ubi(:,desc%ksrc+1), 2, 2*band%ngwv+1)
206 :
207 1489176 : end subroutine gw_cm_src
208 :
209 0 : end module gw_front
|
