ROMS
Loading...
Searching...
No Matches
ad_uv3dmix2_s.h
Go to the documentation of this file.
1 MODULE ad_uv3dmix2_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This routine computes adjoint harmonic mixing of momentum, along !
11! constant S-surfaces, from the horizontal divergence of the stress !
12! tensor. A transverse isotropy is assumed so the stress tensor is !
13! splitted into vertical and horizontal subtensors. !
14! !
15! Reference: !
16! !
17! Wajsowicz, R.C, 1993: A consistent formulation of the !
18! anisotropic stress tensor for use in models of the !
19! large-scale ocean circulation, JCP, 105, 333-338. !
20! !
21! Sadourny, R. and K. Maynard, 1997: Formulations of !
22! lateral diffusion in geophysical fluid dynamics !
23! models, In Numerical Methods of Atmospheric and !
24! Oceanic Modelling. Lin, Laprise, and Ritchie, !
25! Eds., NRC Research Press, 547-556. !
26! !
27! Griffies, S.M. and R.W. Hallberg, 2000: Biharmonic !
28! friction with a Smagorinsky-like viscosity for !
29! use in large-scale eddy-permitting ocean models, !
30! Monthly Weather Rev., 128, 8, 2935-2946. !
31! !
32! BASIC STATE variables needed: visc2, u, v, Hz !
33! !
34!=======================================================================
35!
36 implicit none
37!
38 PRIVATE
39 PUBLIC ad_uv3dmix2
40!
41 CONTAINS
42!
43!***********************************************************************
44 SUBROUTINE ad_uv3dmix2 (ng, tile)
45!***********************************************************************
46!
47 USE mod_param
48 USE mod_coupling
49!!#ifdef DIAGNOSTICS_UV
50!! USE mod_diags
51!!#endif
52 USE mod_grid
53 USE mod_mixing
54 USE mod_ocean
55 USE mod_stepping
56!
57! Imported variable declarations.
58!
59 integer, intent(in) :: ng, tile
60!
61! Local variable declarations.
62!
63 character (len=*), parameter :: MyFile = &
64 & __FILE__
65!
66#include "tile.h"
67!
68#ifdef PROFILE
69 CALL wclock_on (ng, iadm, 30, __line__, myfile)
70#endif
71 CALL ad_uv3dmix2_s_tile (ng, tile, &
72 & lbi, ubi, lbj, ubj, &
73 & imins, imaxs, jmins, jmaxs, &
74 & nrhs(ng), nnew(ng), &
75#ifdef MASKING
76 & grid(ng) % pmask, &
77#endif
78 & grid(ng) % Hz, &
79 & grid(ng) % ad_Hz, &
80 & grid(ng) % om_p, &
81 & grid(ng) % om_r, &
82 & grid(ng) % on_p, &
83 & grid(ng) % on_r, &
84 & grid(ng) % pm, &
85 & grid(ng) % pmon_p, &
86 & grid(ng) % pmon_r, &
87 & grid(ng) % pn, &
88 & grid(ng) % pnom_p, &
89 & grid(ng) % pnom_r, &
90 & mixing(ng) % visc2_p, &
91 & mixing(ng) % visc2_r, &
92!!#ifdef DIAGNOSTICS_UV
93!! & DIAGS(ng) % DiaRUfrc, &
94!! & DIAGS(ng) % DiaRVfrc, &
95!! & DIAGS(ng) % DiaU3wrk, &
96!! & DIAGS(ng) % DiaV3wrk, &
97!!#endif
98 & ocean(ng) % u, &
99 & ocean(ng) % v, &
100 & coupling(ng) % ad_rufrc, &
101 & coupling(ng) % ad_rvfrc, &
102 & ocean(ng) % ad_u, &
103 & ocean(ng) % ad_v)
104#ifdef PROFILE
105 CALL wclock_off (ng, iadm, 30, __line__, myfile)
106#endif
107!
108 RETURN
109 END SUBROUTINE ad_uv3dmix2
110
111!
112!***********************************************************************
113 SUBROUTINE ad_uv3dmix2_s_tile (ng, tile, &
114 & LBi, UBi, LBj, UBj, &
115 & IminS, ImaxS, JminS, JmaxS, &
116 & nrhs, nnew, &
117#ifdef MASKING
118 & pmask, &
119#endif
120 & Hz, ad_Hz, &
121 & om_p, om_r, on_p, on_r, &
122 & pm, pmon_p, pmon_r, &
123 & pn, pnom_p, pnom_r, &
124 & visc2_p, visc2_r, &
125!!#ifdef DIAGNOSTICS_UV
126!! & DiaRUfrc, DiaRVfrc, &
127!! & DiaU3wrk, DiaV3wrk, &
128!!#endif
129 & u, v, &
130 & ad_rufrc, ad_rvfrc, &
131 & ad_u, ad_v)
132!***********************************************************************
133!
134 USE mod_param
135 USE mod_scalars
136!
137! Imported variable declarations.
138!
139 integer, intent(in) :: ng, tile
140 integer, intent(in) :: LBi, UBi, LBj, UBj
141 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
142 integer, intent(in) :: nrhs, nnew
143
144#ifdef ASSUMED_SHAPE
145#ifdef MASKING
146 real(r8), intent(in) :: pmask(LBi:,LBj:)
147#endif
148 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
149 real(r8), intent(in) :: om_p(LBi:,LBj:)
150 real(r8), intent(in) :: om_r(LBi:,LBj:)
151 real(r8), intent(in) :: on_p(LBi:,LBj:)
152 real(r8), intent(in) :: on_r(LBi:,LBj:)
153 real(r8), intent(in) :: pm(LBi:,LBj:)
154 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
155 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
156 real(r8), intent(in) :: pn(LBi:,LBj:)
157 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
158 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
159 real(r8), intent(in) :: visc2_p(LBi:,LBj:)
160 real(r8), intent(in) :: visc2_r(LBi:,LBj:)
161
162 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
163 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
164
165 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
166 real(r8), intent(inout) :: ad_rufrc(LBi:,LBj:)
167 real(r8), intent(inout) :: ad_rvfrc(LBi:,LBj:)
168 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
169 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
170#else
171#ifdef MASKING
172 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
173#endif
174 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
175 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
178 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
179 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
180 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
181 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
182 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
183 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
184 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
185 real(r8), intent(in) :: visc2_p(LBi:UBi,LBj:UBj)
186 real(r8), intent(in) :: visc2_r(LBi:UBi,LBj:UBj)
187
188 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
189 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
190
191 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
192 real(r8), intent(inout) :: ad_rufrc(LBi:UBi,LBj:UBj)
193 real(r8), intent(inout) :: ad_rvfrc(LBi:UBi,LBj:UBj)
194 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
195 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
196#endif
197!
198! Local variable declarations.
199!
200 integer :: i, j, k
201
202 real(r8) :: cff, ad_cff, ad_cff1, ad_cff2
203 real(r8) :: adfac, adfac1, adfac2, adfac3, adfac4
204
205 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_UFe
206 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_VFe
207 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_UFx
208 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_VFx
209
210#include "set_bounds.h"
211!
212!-----------------------------------------------------------------------
213! Initialize adjoint private variables.
214!-----------------------------------------------------------------------
215!
216 ad_cff=0.0_r8
217 ad_cff1=0.0_r8
218 ad_cff2=0.0_r8
219
220 ad_ufe(imins:imaxs,jmins:jmaxs)=0.0_r8
221 ad_vfe(imins:imaxs,jmins:jmaxs)=0.0_r8
222 ad_ufx(imins:imaxs,jmins:jmaxs)=0.0_r8
223 ad_vfx(imins:imaxs,jmins:jmaxs)=0.0_r8
224!
225!-----------------------------------------------------------------------
226! Compute horizontal harmonic viscosity along constant S-surfaces.
227!-----------------------------------------------------------------------
228!
229 k_loop : DO k=1,n(ng)
230!
231! Time-step harmonic, S-surfaces viscosity term. Notice that momentum
232! at this stage is HzU and HzV and has m2/s units. Add contribution for
233! barotropic forcing terms.
234!
235 DO j=jstrv,jend
236 DO i=istr,iend
237 cff=0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
238!^ tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)+tl_cff2
239!^
240 ad_cff2=ad_cff2+ad_v(i,j,k,nnew)
241!^ tl_rvfrc(i,j)=tl_rvfrc(i,j)+tl_cff1
242!^
243 ad_cff1=ad_cff1+ad_rvfrc(i,j)
244!^ tl_cff2=dt(ng)*cff*tl_cff1
245!^
246 ad_cff1=ad_cff1+dt(ng)*cff*ad_cff2
247 ad_cff2=0.0_r8
248!^ tl_cff1=0.5_r8*((pn(i,j-1)+pn(i,j))* &
249!^ & (tl_VFx(i+1,j)-tl_VFx(i,j ))- &
250!^ & (pm(i,j-1)+pm(i,j))* &
251!^ & (tl_VFe(i ,j)-tl_VFe(i,j-1)))
252!^
253 adfac=0.5_r8*ad_cff1
254 adfac1=adfac*(pn(i,j-1)+pn(i,j))
255 adfac2=adfac*(pm(i,j-1)+pm(i,j))
256 ad_vfx(i ,j)=ad_vfx(i ,j)-adfac1
257 ad_vfx(i+1,j)=ad_vfx(i+1,j)+adfac1
258 ad_vfe(i,j-1)=ad_vfe(i,j-1)+adfac2
259 ad_vfe(i,j )=ad_vfe(i,j )-adfac2
260 ad_cff1=0.0_r8
261 END DO
262 END DO
263 DO j=jstr,jend
264 DO i=istru,iend
265 cff=0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
266!^ tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)+tl_cff2
267!^
268 ad_cff2=ad_cff2+ad_u(i,j,k,nnew)
269!^ tl_rufrc(i,j)=tl_rufrc(i,j)+tl_cff1
270!^
271 ad_cff1=ad_cff1+ad_rufrc(i,j)
272!^ tl_cff2=dt(ng)*cff*tl_cff1
273!^
274 ad_cff1=ad_cff1+dt(ng)*cff*ad_cff2
275 ad_cff2=0.0_r8
276!^ cff1=0.5_r8*((pn(i-1,j)+pn(i,j))* &
277!^ & (UFx(i,j )-UFx(i-1,j))+ &
278!^ & (pm(i-1,j)+pm(i,j))* &
279!^ & (UFe(i,j+1)-UFe(i ,j)))
280!^
281!^
282 adfac=0.5_r8*ad_cff1
283 adfac1=adfac*(pn(i-1,j)+pn(i,j))
284 adfac2=adfac*(pm(i-1,j)+pm(i,j))
285 ad_ufx(i-1,j)=ad_ufx(i-1,j)-adfac1
286 ad_ufx(i ,j)=ad_ufx(i ,j)+adfac1
287 ad_ufe(i,j )=ad_ufe(i,j )-adfac2
288 ad_ufe(i,j+1)=ad_ufe(i,j+1)+adfac2
289 ad_cff1=0.0_r8
290 END DO
291 END DO
292!
293! Compute flux-components of the horizontal divergence of the stress
294! tensor (m5/s2) in XI- and ETA-directions.
295!
296 DO j=jstr,jend+1
297 DO i=istr,iend+1
298!^ tl_VFx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
299!^ tl_UFe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
300!^
301 ad_cff=ad_cff+ &
302 & on_p(i,j)*on_p(i,j)*ad_vfx(i,j)+ &
303 & om_p(i,j)*om_p(i,j)*ad_ufe(i,j)
304 ad_vfx(i,j)=0.0_r8
305 ad_ufe(i,j)=0.0_r8
306#ifdef MASKING
307!^ tl_cff=tl_cff*pmask(i,j)
308!^
309 ad_cff=ad_cff*pmask(i,j)
310#endif
311!^ tl_cff=visc2_p(i,j)*0.125_r8* &
312!^ & ((tl_Hz(i-1,j ,k)+tl_Hz(i,j ,k)+ &
313!^ & tl_Hz(i-1,j-1,k)+tl_Hz(i,j-1,k))* &
314!^ & (pmon_p(i,j)* &
315!^ & ((pn(i ,j-1)+pn(i ,j))*v(i ,j,k,nrhs)- &
316!^ & (pn(i-1,j-1)+pn(i-1,j))*v(i-1,j,k,nrhs))+ &
317!^ & pnom_p(i,j)* &
318!^ & ((pm(i-1,j )+pm(i,j ))*u(i,j ,k,nrhs)- &
319!^ & (pm(i-1,j-1)+pm(i,j-1))*u(i,j-1,k,nrhs)))+ &
320!^ & (Hz(i-1,j ,k)+Hz(i,j ,k)+ &
321!^ & Hz(i-1,j-1,k)+Hz(i,j-1,k))* &
322!^ & (pmon_p(i,j)* &
323!^ & ((pn(i ,j-1)+pn(i ,j))*tl_v(i ,j,k,nrhs)- &
324!^ & (pn(i-1,j-1)+pn(i-1,j))*tl_v(i-1,j,k,nrhs))+ &
325!^ & pnom_p(i,j)* &
326!^ & ((pm(i-1,j )+pm(i,j ))*tl_u(i,j ,k,nrhs)- &
327!^ & (pm(i-1,j-1)+pm(i,j-1))*tl_u(i,j-1,k,nrhs))))
328!^
329
330 adfac=visc2_p(i,j)*0.125_r8*ad_cff
331 adfac1=adfac* &
332 & (pmon_p(i,j)* &
333 & ((pn(i ,j-1)+pn(i ,j))*v(i ,j,k,nrhs)- &
334 & (pn(i-1,j-1)+pn(i-1,j))*v(i-1,j,k,nrhs))+ &
335 & pnom_p(i,j)* &
336 & ((pm(i-1,j )+pm(i,j ))*u(i,j ,k,nrhs)- &
337 & (pm(i-1,j-1)+pm(i,j-1))*u(i,j-1,k,nrhs)))
338 adfac2=adfac*(hz(i-1,j ,k)+hz(i,j ,k)+ &
339 & hz(i-1,j-1,k)+hz(i,j-1,k))
340 adfac3=adfac2*pmon_p(i,j)
341 adfac4=adfac2*pnom_p(i,j)
342 ad_hz(i-1,j-1,k)=ad_hz(i-1,j-1,k)+adfac1
343 ad_hz(i ,j-1,k)=ad_hz(i ,j-1,k)+adfac1
344 ad_hz(i-1,j ,k)=ad_hz(i-1,j ,k)+adfac1
345 ad_hz(i ,j ,k)=ad_hz(i ,j ,k)+adfac1
346 ad_v(i-1,j,k,nrhs)=ad_v(i-1,j,k,nrhs)- &
347 & (pn(i-1,j-1)+pn(i-1,j))*adfac3
348 ad_v(i ,j,k,nrhs)=ad_v(i ,j,k,nrhs)+ &
349 & (pn(i ,j-1)+pn(i ,j))*adfac3
350 ad_u(i,j-1,k,nrhs)=ad_u(i,j-1,k,nrhs)- &
351 & (pm(i-1,j-1)+pm(i,j-1))*adfac4
352 ad_u(i,j ,k,nrhs)=ad_u(i,j ,k,nrhs)+ &
353 & (pm(i-1,j )+pm(i,j ))*adfac4
354 ad_cff=0.0_r8
355 END DO
356 END DO
357 DO j=jstrv-1,jend
358 DO i=istru-1,iend
359!^ tl_VFe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
360!^ tl_UFx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
361!^
362 ad_cff=ad_cff+ &
363 & om_r(i,j)*om_r(i,j)*ad_vfe(i,j)+ &
364 & on_r(i,j)*on_r(i,j)*ad_ufx(i,j)
365 ad_vfe(i,j)=0.0_r8
366 ad_ufx(i,j)=0.0_r8
367!^ tl_cff=visc2_r(i,j)*0.5_r8* &
368!^ & (tl_Hz(i,j,k)* &
369!^ & (pmon_r(i,j)* &
370!^ & ((pn(i ,j)+pn(i+1,j))*u(i+1,j,k,nrhs)- &
371!^ & (pn(i-1,j)+pn(i ,j))*u(i ,j,k,nrhs))- &
372!^ & pnom_r(i,j)* &
373!^ & ((pm(i,j )+pm(i,j+1))*v(i,j+1,k,nrhs)- &
374!^ & (pm(i,j-1)+pm(i,j ))*v(i,j ,k,nrhs)))+ &
375!^ & Hz(i,j,k)* &
376!^ & (pmon_r(i,j)* &
377!^ & ((pn(i ,j)+pn(i+1,j))*tl_u(i+1,j,k,nrhs)- &
378!^ & (pn(i-1,j)+pn(i ,j))*tl_u(i ,j,k,nrhs))- &
379!^ & pnom_r(i,j)* &
380!^ & ((pm(i,j )+pm(i,j+1))*tl_v(i,j+1,k,nrhs)- &
381!^ & (pm(i,j-1)+pm(i,j ))*tl_v(i,j ,k,nrhs))))
382!^
383 adfac=visc2_r(i,j)*0.5_r8*ad_cff
384 adfac1=adfac*hz(i,j,k)
385 adfac2=adfac1*pmon_r(i,j)
386 adfac3=adfac1*pnom_r(i,j)
387 ad_hz(i,j,k)=ad_hz(i,j,k)+ &
388 & (pmon_r(i,j)* &
389 & ((pn(i ,j)+pn(i+1,j))*u(i+1,j,k,nrhs)- &
390 & (pn(i-1,j)+pn(i ,j))*u(i ,j,k,nrhs))- &
391 & pnom_r(i,j)* &
392 & ((pm(i,j )+pm(i,j+1))*v(i,j+1,k,nrhs)- &
393 & (pm(i,j-1)+pm(i,j ))*v(i,j ,k,nrhs)))* &
394 & adfac
395 ad_u(i ,j,k,nrhs)=ad_u(i ,j,k,nrhs)- &
396 & (pn(i-1,j)+pn(i ,j))*adfac2
397 ad_u(i+1,j,k,nrhs)=ad_u(i+1,j,k,nrhs)+ &
398 & (pn(i ,j)+pn(i+1,j))*adfac2
399 ad_v(i,j ,k,nrhs)=ad_v(i,j ,k,nrhs)+ &
400 & (pm(i,j-1)+pm(i,j ))*adfac3
401 ad_v(i,j+1,k,nrhs)=ad_v(i,j+1,k,nrhs)- &
402 & (pm(i,j )+pm(i,j+1))*adfac3
403 ad_cff=0.0_r8
404 END DO
405 END DO
406 END DO k_loop
407!
408 RETURN
409 END SUBROUTINE ad_uv3dmix2_s_tile
410
411 END MODULE ad_uv3dmix2_mod
subroutine, public ad_uv3dmix2(ng, tile)
subroutine ad_uv3dmix2_s_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, pmask, hz, ad_hz, om_p, om_r, on_p, on_r, pm, pmon_p, pmon_r, pn, pnom_p, pnom_r, visc2_p, visc2_r, u, v, ad_rufrc, ad_rvfrc, ad_u, ad_v)
type(t_coupling), dimension(:), allocatable coupling
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter iadm
Definition mod_param.F:665
real(dp), dimension(:), allocatable dt
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
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