ROMS
Loading...
Searching...
No Matches
tl_omega.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined TANGENT && defined SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine computes S-coordinate vertical velocity (m^3/s), !
13! !
14! W=[Hz/(m*n)]*omega, !
15! !
16! diagnostically at horizontal RHO-points and vertical W-points. !
17! !
18! BASIC STATE variables needed: W, z_w. !
19! !
20! NOTE: We need to recompute basic state W in this routine since !
21! ---- intermediate values of W are needed by the tangent linear !
22! and adjoint routines. !
23! !
24!=======================================================================
25!
26 implicit none
27!
28 PRIVATE
29 PUBLIC :: tl_omega
30!
31 CONTAINS
32!
33!***********************************************************************
34 SUBROUTINE tl_omega (ng, tile, model)
35!***********************************************************************
36!
37 USE mod_param
38 USE mod_grid
39 USE mod_ocean
40# if defined SEDIMENT && defined SED_MORPH
41 USE mod_sedbed
42 USE mod_stepping
43# endif
44!
45! Imported variable declarations.
46!
47 integer, intent(in) :: ng, tile, model
48!
49! Local variable declarations.
50!
51 character (len=*), parameter :: myfile = &
52 & __FILE__
53!
54# include "tile.h"
55!
56# ifdef PROFILE
57 CALL wclock_on (ng, model, 13, __line__, myfile)
58# endif
59 CALL tl_omega_tile (ng, tile, model, &
60 & lbi, ubi, lbj, ubj, &
61 & imins, imaxs, jmins, jmaxs, &
62# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
63 & nstp(ng), nnew(ng), &
64 & grid(ng) % omn, &
65 & sedbed(ng) % bed_thick, &
66 & sedbed(ng) % tl_bed_thick, &
67# endif
68 & grid(ng) % Huon, &
69 & grid(ng) % Hvom, &
70 & grid(ng) % z_w, &
71 & grid(ng) % tl_Huon, &
72 & grid(ng) % tl_Hvom, &
73 & grid(ng) % tl_z_w, &
74 & ocean(ng) % W, &
75 & ocean(ng) % tl_W)
76# ifdef PROFILE
77 CALL wclock_off (ng, model, 13, __line__, myfile)
78# endif
79!
80 RETURN
81 END SUBROUTINE tl_omega
82!
83!***********************************************************************
84 SUBROUTINE tl_omega_tile (ng, tile, model, &
85 & LBi, UBi, LBj, UBj, &
86 & IminS, ImaxS, JminS, JmaxS, &
87# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
88 & nstp, nnew, &
89 & omn,
90 & bed_thick, tl_bed_thick, &
91# endif
92 & Huon, Hvom, z_w, &
93 & tl_Huon, tl_Hvom, tl_z_w, &
94 & W, tl_W)
95!***********************************************************************
96!
97 USE mod_param
98 USE mod_scalars
99 USE mod_sources
100!
101 USE bc_3d_mod, ONLY : bc_w3d_tile
102# ifdef DISTRIBUTE
104# endif
105!
106! Imported variable declarations.
107!
108 integer, intent(in) :: ng, tile, model
109 integer, intent(in) :: LBi, UBi, LBj, UBj
110 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
111# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
112 integer, intent(in) :: nstp, nnew
113# endif
114!
115# ifdef ASSUMED_SHAPE
116# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
117 real(r8), intent(in) :: omn(LBi:,LBj:)
118 real(r8), intent(in):: bed_thick(LBi:,LBj:,:)
119 real(r8), intent(in):: tl_bed_thick(LBi:,LBj:,:)
120# endif
121 real(r8), intent(in) :: Huon(LBi:,LBj:,:)
122 real(r8), intent(in) :: Hvom(LBi:,LBj:,:)
123 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
124 real(r8), intent(in) :: tl_Huon(LBi:,LBj:,:)
125 real(r8), intent(in) :: tl_Hvom(LBi:,LBj:,:)
126 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
127
128 real(r8), intent(out) :: W(LBi:,LBj:,0:)
129 real(r8), intent(out) :: tl_W(LBi:,LBj:,0:)
130
131# else
132
133# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
134 real(r8), intent(in) :: omn(LBi:UBi,LBj:UBj)
135 real(r8), intent(in):: bed_thick(LBi:UBi,LBj:UBj,2)
136 real(r8), intent(in):: tl_bed_thick(LBi:UBi,LBj:UBj,2)
137# endif
138 real(r8), intent(in) :: Huon(LBi:UBi,LBj:UBj,N(ng))
139 real(r8), intent(in) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
140 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
141 real(r8), intent(in) :: tl_Huon(LBi:UBi,LBj:UBj,N(ng))
142 real(r8), intent(in) :: tl_Hvom(LBi:UBi,LBj:UBj,N(ng))
143 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
144
145 real(r8), intent(out) :: W(LBi:UBi,LBj:UBj,0:N(ng))
146 real(r8), intent(out) :: tl_W(LBi:UBi,LBj:UBj,0:N(ng))
147# endif
148!
149! Local variable declarations.
150!
151 integer :: i, ii, is, j, jj, k
152 real(r8) :: cff, tl_cff
153# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
154 real(r8) :: cff1
155# endif
156 real(r8), dimension(IminS:ImaxS) :: wrk
157 real(r8), dimension(IminS:ImaxS) :: tl_wrk
158
159# include "set_bounds.h"
160!
161!------------------------------------------------------------------------
162! Vertically integrate horizontal mass flux divergence.
163!------------------------------------------------------------------------
164!
165! Starting with zero vertical velocity at the bottom, integrate
166! from the bottom (k=0) to the free-surface (k=N). The w(:,:,N(ng))
167! contains the vertical velocity at the free-surface, d(zeta)/d(t).
168! Notice that barotropic mass flux divergence is not used directly.
169!
170# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
171 cff1=1.0_r8/dt(ng)
172# endif
173 DO j=jstr,jend
174 DO i=istr,iend
175# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
176 w(i,j,0)=-cff1*(bed_thick(i,j,nstp)- &
177 & bed_thick(i,j,nnew))*omn(i,j)
178 tl_w(i,j,0)=-cff1*(tl_bed_thick(i,j,nstp)- &
179 & tl_bed_thick(i,j,nnew))*omn(i,j)
180# else
181 w(i,j,0)=0.0_r8
182 tl_w(i,j,0)=0.0_r8
183# endif
184 END DO
185!
186! Code added to clear tl_W to be consistent with adjoint.
187!
188 DO k=1,n(ng)
189 DO i=istr,iend
190 tl_w(i,j,k)=0.0_r8
191 END DO
192 END DO
193 DO k=1,n(ng)
194 DO i=istr,iend
195 w(i,j,k)=w(i,j,k-1)- &
196 & (huon(i+1,j,k)-huon(i,j,k)+ &
197 & hvom(i,j+1,k)-hvom(i,j,k))
198 tl_w(i,j,k)=tl_w(i,j,k-1)- &
199 & (tl_huon(i+1,j,k)-tl_huon(i,j,k)+ &
200 & tl_hvom(i,j+1,k)-tl_hvom(i,j,k))
201 END DO
202 END DO
203!
204! Apply mass point sources (volume vertical influx), if any.
205!
206! Overwrite W(Isrc,Jsrc,k) with the same divergence of Huon,Hvom as
207! above but add in point source Qsrc(k) and reaccumulate the vertical
208! sum to obtain the correct net Qbar given in user input - J. Levin
209! (Jupiter Intelligence Inc.) and J. Wilkin
210!
211! Dsrc(is) = 2, flow across grid cell w-face (positive or negative)
212!
213 IF (lwsrc(ng)) THEN
214 DO is=1,nsrc(ng)
215 IF (int(sources(ng)%Dsrc(is)).eq.2) THEN
216 ii=sources(ng)%Isrc(is)
217 jj=sources(ng)%Jsrc(is)
218 IF (((istrr.le.ii).and.(ii.le.iendr)).and. &
219 & ((jstrr.le.jj).and.(jj.le.jendr)).and. &
220 & (j.eq.jj)) THEN
221 DO k=1,n(ng)
222 w(ii,jj,k)=w(ii,jj,k-1)- &
223 & (huon(ii+1,jj,k)-huon(ii,jj,k)+ &
224 & hvom(ii,jj+1,k)-hvom(ii,jj,k))+ &
225 & sources(ng)%Qsrc(is,k)
226 tl_w(ii,jj,k)=tl_w(ii,jj,k-1)- &
227 & (tl_huon(ii+1,jj,k)-tl_huon(ii,jj,k)+ &
228 & tl_hvom(ii,jj+1,k)-tl_hvom(ii,jj,k))+ &
229 & sources(ng)%tl_Qsrc(is,k)
230 END DO
231 END IF
232 END IF
233 END DO
234 END IF
235!
236 DO i=istr,iend
237 cff=1.0_r8/(z_w(i,j,n(ng))-z_w(i,j,0))
238 tl_cff=-cff*cff*(tl_z_w(i,j,n(ng))-tl_z_w(i,j,0))
239 wrk(i)=cff*w(i,j,n(ng))
240 tl_wrk(i)=tl_cff*w(i,j,n(ng))+cff*tl_w(i,j,n(ng))
241 END DO
242!
243! In order to insure zero vertical velocity at the free-surface,
244! subtract the vertical velocities of the moving S-coordinates
245! isosurfaces. These isosurfaces are proportional to d(zeta)/d(t).
246! The proportionally coefficients are a linear function of the
247! S-coordinate with zero value at the bottom (k=0) and unity at
248! the free-surface (k=N).
249!
250 DO k=n(ng)-1,1,-1
251 DO i=istr,iend
252 w(i,j,k)=w(i,j,k)- &
253 & wrk(i)*(z_w(i,j,k)-z_w(i,j,0))
254 tl_w(i,j,k)=tl_w(i,j,k)- &
255 & tl_wrk(i)*(z_w(i,j,k)-z_w(i,j,0))- &
256 & wrk(i)*(tl_z_w(i,j,k)-tl_z_w(i,j,0))
257 END DO
258 END DO
259 DO i=istr,iend
260 w(i,j,n(ng))=0.0_r8
261 tl_w(i,j,n(ng))=0.0_r8
262 END DO
263 END DO
264!
265! Set lateral boundary conditions.
266!
267 CALL bc_w3d_tile (ng, tile, &
268 & lbi, ubi, lbj, ubj, 0, n(ng), &
269 & w)
270 CALL bc_w3d_tile (ng, tile, &
271 & lbi, ubi, lbj, ubj, 0, n(ng), &
272 & tl_w)
273
274# ifdef DISTRIBUTE
275 CALL mp_exchange3d (ng, tile, model, 1, &
276 & lbi, ubi, lbj, ubj, 0, n(ng), &
277 & nghostpoints, &
278 & ewperiodic(ng), nsperiodic(ng), &
279 & w)
280 CALL mp_exchange3d (ng, tile, model, 1, &
281 & lbi, ubi, lbj, ubj, 0, n(ng), &
282 & nghostpoints, &
283 & ewperiodic(ng), nsperiodic(ng), &
284 & tl_w)
285# endif
286!
287 RETURN
288 END SUBROUTINE tl_omega_tile
289#endif
290 END MODULE tl_omega_mod
subroutine bc_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
Definition bc_3d.F:591
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer nghostpoints
Definition mod_param.F:710
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable lwsrc
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
type(t_sources), dimension(:), allocatable sources
Definition mod_sources.F:90
integer, dimension(:), allocatable nsrc
Definition mod_sources.F:97
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine tl_omega_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, omn,
Definition tl_omega.F:90
subroutine, public tl_omega(ng, tile, model)
Definition tl_omega.F:35
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3