ROMS
Loading...
Searching...
No Matches
rp_uv3drelax_mod Module Reference

Functions/Subroutines

subroutine, public rp_uv3drelax (ng, tile)
 
subroutine rp_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, u, v, tl_u, tl_v)
 

Function/Subroutine Documentation

◆ rp_uv3drelax()

subroutine, public rp_uv3drelax_mod::rp_uv3drelax ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 28 of file rp_uv3drelax.F.

29!***********************************************************************
30!
31 USE mod_param
32 USE mod_coupling
33 USE mod_grid
34 USE mod_ocean
35 USE mod_stepping
36!
37! Imported variable declarations.
38!
39 integer, intent(in) :: ng, tile
40!
41! Local variable declarations.
42!
43 character (len=*), parameter :: MyFile = &
44 & __FILE__
45!
46# include "tile.h"
47!
48# ifdef PROFILE
49 CALL wclock_on (ng, irpm, 30, __line__, myfile)
50# endif
51 CALL rp_uv3drelax_tile (ng, tile, &
52 & lbi, ubi, lbj, ubj, &
53 & imins, imaxs, jmins, jmaxs, &
54 & nrhs(ng), nnew(ng), &
55# ifdef MASKING
56 & grid(ng) % pmask, &
57# endif
58 & grid(ng) % Hz, &
59 & grid(ng) % pm, &
60 & grid(ng) % pmon_p, &
61 & grid(ng) % pmon_r, &
62 & grid(ng) % pn, &
63 & grid(ng) % pnom_p, &
64 & grid(ng) % pnom_r, &
65 & ocean(ng) % u, &
66 & ocean(ng) % v, &
67 & ocean(ng) % tl_u, &
68 & ocean(ng) % tl_v)
69# ifdef PROFILE
70 CALL wclock_off (ng, irpm, 30, __line__, myfile)
71# endif
72!
73 RETURN
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter irpm
Definition mod_param.F:664
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

References mod_grid::grid, mod_param::irpm, mod_stepping::nnew, mod_stepping::nrhs, mod_ocean::ocean, rp_uv3drelax_tile(), wclock_off(), and wclock_on().

Referenced by rp_rhs3d_mod::rp_rhs3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ rp_uv3drelax_tile()

subroutine rp_uv3drelax_mod::rp_uv3drelax_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
integer, intent(in) nnew,
real(r8), dimension(lbi:,lbj:), intent(in) pmask,
real(r8), dimension(lbi:,lbj:,:), intent(in) hz,
real(r8), dimension(lbi:,lbj:), intent(in) pm,
real(r8), dimension(lbi:,lbj:), intent(in) pmon_p,
real(r8), dimension(lbi:,lbj:), intent(in) pmon_r,
real(r8), dimension(lbi:,lbj:), intent(in) pn,
real(r8), dimension(lbi:,lbj:), intent(in) pnom_p,
real(r8), dimension(lbi:,lbj:), intent(in) pnom_r,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) v,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_v )
private

Definition at line 78 of file rp_uv3drelax.F.

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(in) :: u(LBi:,LBj:,:,:)
115 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
116
117 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
118 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
119# else
120# ifdef MASKING
121 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
122# endif
123 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
124 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
125 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
126 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
127 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
128 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
129 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
130
131 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
132 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
133
134 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
135 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
136# endif
137!
138! Local variable declarations.
139!
140 integer :: i, j, k
141
142 real(r8) :: cff
143
144 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
145 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
146 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
147 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
148
149# include "set_bounds.h"
150!
151!-----------------------------------------------------------------------
152! Compute horizontal diffusion relaxation of 3D momentum between
153! current and previous representer tangent linear Picard iteration
154! trajectory (basic state).
155!-----------------------------------------------------------------------
156!
157 IF (tl_m3diff(ng).gt.0.0_r8) THEN
158!
159 k_loop : DO k=1,n(ng)
160!
161! Compute flux-components of the diffusive relaxation (m4/s2) in XI-
162! and ETA-directions.
163!
164 DO j=jstr,jend
165 DO i=istru-1,iend
166 ufx(i,j)=tl_m3diff(ng)*pmon_r(i,j)* &
167 & hz(i,j,k)* &
168 & (tl_u(i+1,j,k,nrhs)-u(i+1,j,k,nrhs)- &
169 & tl_u(i ,j,k,nrhs)+u(i ,j,k,nrhs))
170 END DO
171 END DO
172 DO j=jstr,jend+1
173 DO i=istru,iend
174 ufe(i,j)=tl_m3diff(ng)*pnom_p(i,j)* &
175 & 0.25_r8*(hz(i,j ,k)+hz(i-1,j ,k)+ &
176 & hz(i,j-1,k)+hz(i-1,j-1,k))* &
177 & (tl_u(i,j ,k,nrhs)-u(i,j ,k,nrhs)- &
178 & tl_u(i,j-1,k,nrhs)+u(i,j-1,k,nrhs))
179# ifdef MASKING
180 ufe(i,j)=ufe(i,j)*pmask(i,j)
181# endif
182 END DO
183 END DO
184 DO j=jstrv,jend
185 DO i=istr,iend+1
186 vfx(i,j)=tl_m3diff(ng)*pmon_p(i,j)* &
187 & 0.25_r8*(hz(i,j ,k)+hz(i-1,j ,k)+ &
188 & hz(i,j-1,k)+hz(i-1,j-1,k))* &
189 & (tl_v(i ,j,k,nrhs)-v(i ,j,k,nrhs)- &
190 & tl_v(i-1,j,k,nrhs)+v(i-1,j,k,nrhs))
191# ifdef MASKING
192 vfx(i,j)=vfx(i,j)*pmask(i,j)
193# endif
194 END DO
195 END DO
196 DO j=jstrv-1,jend
197 DO i=istr,iend
198 vfe(i,j)=tl_m3diff(ng)*pnom_r(i,j)* &
199 & hz(i,j,k)* &
200 & (tl_v(i,j+1,k,nrhs)-v(i,j+1,k,nrhs)- &
201 & tl_v(i,j ,k,nrhs)+v(i,j ,k,nrhs))
202 END DO
203 END DO
204!
205! Time-step diffusive relaxation term. Notice that momentum at this
206! stage is HzU and HzV and has m2/s units.
207!
208 cff=dt(ng)*0.25_r8
209 DO j=jstr,jend
210 DO i=istru,iend
211 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)+ &
212 & cff*(pm(i-1,j)+pm(i,j))* &
213 & (pn(i-1,j)+pn(i,j))* &
214 & (ufx(i,j)-ufx(i-1,j)+ &
215 & ufe(i,j+1)-ufe(i,j))
216
217 END DO
218 END DO
219 DO j=jstrv,jend
220 DO i=istr,iend
221 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)+ &
222 & cff*(pm(i,j)+pm(i,j-1))* &
223 & (pn(i,j)+pn(i,j-1))* &
224 & (vfx(i+1,j)-vfx(i,j)+ &
225 & vfe(i,j)-vfe(i,j-1))
226 END DO
227 END DO
228 END DO k_loop
229 END IF
230!
231 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
real(dp), dimension(:), allocatable dt
real(r8), dimension(:), allocatable tl_m3diff

References mod_scalars::dt, and mod_scalars::tl_m3diff.

Referenced by rp_uv3drelax().

Here is the caller graph for this function: