ROMS
Loading...
Searching...
No Matches
tl_uv3drelax.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#if defined TL_IOMS && \
5 defined rpm_relaxation && \
6 defined r4dvar_ana_sensitivity && \
7 defined solve3d
8!
9!git $Id$
10!=================================================== Andrew M. Moore ===
11! Copyright (c) 2002-2025 The ROMS Group Hernan G. Arango !
12! Licensed under a MIT/X style license !
13! See License_ROMS.md !
14!=======================================================================
15! !
16! This routine relaxes current representer tangent linear 3D !
17! momentum to previous Picard iteration solution (basic state) !
18! to improve stability and convergence. !
19! !
20!=======================================================================
21!
22 implicit none
23
24 PRIVATE
25
26 PUBLIC tl_uv3drelax
27
28 CONTAINS
29!
30!***********************************************************************
31 SUBROUTINE tl_uv3drelax (ng, tile)
32!***********************************************************************
33!
34 USE mod_param
35 USE mod_coupling
36 USE mod_grid
37 USE mod_ocean
38 USE mod_stepping
39!
40! Imported variable declarations.
41!
42 integer, intent(in) :: ng, tile
43!
44! Local variable declarations.
45!
46 character (len=*), parameter :: myfile = &
47 & __FILE__
48!
49# include "tile.h"
50!
51# ifdef PROFILE
52 CALL wclock_on (ng, itlm, 30, __line__, myfile)
53# endif
54 CALL tl_uv3drelax_tile (ng, tile, &
55 & lbi, ubi, lbj, ubj, &
56 & imins, imaxs, jmins, jmaxs, &
57 & nrhs(ng), nnew(ng), &
58# ifdef MASKING
59 & grid(ng) % pmask, &
60# endif
61 & grid(ng) % Hz, &
62 & grid(ng) % pm, &
63 & grid(ng) % pmon_p, &
64 & grid(ng) % pmon_r, &
65 & grid(ng) % pn, &
66 & grid(ng) % pnom_p, &
67 & grid(ng) % pnom_r, &
68 & ocean(ng) % tl_u, &
69 & ocean(ng) % tl_v)
70# ifdef PROFILE
71 CALL wclock_off (ng, itlm, 30, __line__, myfile)
72# endif
73!
74 RETURN
75 END SUBROUTINE tl_uv3drelax
76
77!
78!***********************************************************************
79 SUBROUTINE tl_uv3drelax_tile (ng, tile, &
80 & LBi, UBi, LBj, UBj, &
81 & IminS, ImaxS, JminS, JmaxS, &
82 & nrhs, nnew, &
83# ifdef MASKING
84 & pmask, &
85# endif
86 & Hz, &
87 & pm, pmon_p, pmon_r, &
88 & pn, pnom_p, pnom_r, &
89 & tl_u, tl_v)
90!***********************************************************************
91!
92 USE mod_param
93 USE mod_scalars
94!
95! Imported variable declarations.
96!
97 integer, intent(in) :: ng, tile
98 integer, intent(in) :: LBi, UBi, LBj, UBj
99 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
100 integer, intent(in) :: nrhs, nnew
101
102# ifdef ASSUMED_SHAPE
103# ifdef MASKING
104 real(r8), intent(in) :: pmask(LBi:,LBj:)
105# endif
106 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
107 real(r8), intent(in) :: pm(LBi:,LBj:)
108 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
109 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
110 real(r8), intent(in) :: pn(LBi:,LBj:)
111 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
112 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
113
114 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
115 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
116# else
117# ifdef MASKING
118 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
119# endif
120 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
121 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
122 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
123 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
124 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
125 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
126 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
127
128 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
129 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
130# endif
131!
132! Local variable declarations.
133!
134 integer :: i, j, k
135
136 real(r8) :: cff
137
138 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
139 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
140 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
141 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
142
143# include "set_bounds.h"
144!
145!-----------------------------------------------------------------------
146! Compute horizontal diffusion relaxation of 3D momentum between
147! current and previous representer tangent linear Picard iteration
148! trajectory (basic state).
149!-----------------------------------------------------------------------
150!
151! This is the tangent linear of the relaxation terms that appear in
152! the representer model. Since the representer model is linear, we
153! DO NOT need to include the tangent linear of the Hz terms below
154! because Hz is computed from the background trajectory (AMM).
155!
156 IF (tl_m3diff(ng).gt.0.0_r8) THEN
157!
158 k_loop : DO k=1,n(ng)
159!
160! Compute flux-components of the diffusive relaxation (m4/s2) in XI-
161! and ETA-directions.
162!
163 DO j=jstr,jend
164 DO i=istru-1,iend
165 ufx(i,j)=tl_m3diff(ng)*pmon_r(i,j)* &
166 & hz(i,j,k)* &
167 & (tl_u(i+1,j,k,nrhs)- &
168 & tl_u(i ,j,k,nrhs))
169 END DO
170 END DO
171 DO j=jstr,jend+1
172 DO i=istru,iend
173 ufe(i,j)=tl_m3diff(ng)*pnom_p(i,j)* &
174 & 0.25_r8*(hz(i,j ,k)+hz(i-1,j ,k)+ &
175 & hz(i,j-1,k)+hz(i-1,j-1,k))* &
176 & (tl_u(i,j ,k,nrhs)- &
177 & tl_u(i,j-1,k,nrhs))
178# ifdef MASKING
179 ufe(i,j)=ufe(i,j)*pmask(i,j)
180# endif
181 END DO
182 END DO
183 DO j=jstrv,jend
184 DO i=istr,iend+1
185 vfx(i,j)=tl_m3diff(ng)*pmon_p(i,j)* &
186 & 0.25_r8*(hz(i,j ,k)+hz(i-1,j ,k)+ &
187 & hz(i,j-1,k)+hz(i-1,j-1,k))* &
188 & (tl_v(i ,j,k,nrhs)- &
189 & tl_v(i-1,j,k,nrhs))
190# ifdef MASKING
191 vfx(i,j)=vfx(i,j)*pmask(i,j)
192# endif
193 END DO
194 END DO
195 DO j=jstrv-1,jend
196 DO i=istr,iend
197 vfe(i,j)=tl_m3diff(ng)*pnom_r(i,j)* &
198 & hz(i,j,k)* &
199 & (tl_v(i,j+1,k,nrhs)- &
200 & tl_v(i,j ,k,nrhs))
201 END DO
202 END DO
203!
204! Time-step diffusive relaxation term. Notice that momentum at this
205! stage is HzU and HzV and has m2/s units.
206!
207 cff=dt(ng)*0.25_r8
208 DO j=jstr,jend
209 DO i=istru,iend
210 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)+ &
211 & cff*(pm(i-1,j)+pm(i,j))* &
212 & (pn(i-1,j)+pn(i,j))* &
213 & (ufx(i,j)-ufx(i-1,j)+ &
214 & ufe(i,j+1)-ufe(i,j))
215
216 END DO
217 END DO
218 DO j=jstrv,jend
219 DO i=istr,iend
220 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)+ &
221 & cff*(pm(i,j)+pm(i,j-1))* &
222 & (pn(i,j)+pn(i,j-1))* &
223 & (vfx(i+1,j)-vfx(i,j)+ &
224 & vfe(i,j)-vfe(i,j-1))
225 END DO
226 END DO
227 END DO k_loop
228 END IF
229!
230 RETURN
231 END SUBROUTINE tl_uv3drelax_tile
232#endif
233 END MODULE tl_uv3drelax_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter itlm
Definition mod_param.F:663
real(dp), dimension(:), allocatable dt
real(r8), dimension(:), allocatable tl_m3diff
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
subroutine, public tl_uv3drelax(ng, tile)
subroutine tl_uv3drelax_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, pmask, hz, pm, pmon_p, pmon_r, pn, pnom_p, pnom_r, tl_u, tl_v)
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