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

Functions/Subroutines

subroutine, public ad_forcing (ng, tile, kfrc, nfrc)
 
subroutine ad_forcing_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kfrc, nfrc, f_t, f_u, f_v, f_ubar, f_vbar, f_zeta, ad_t, ad_u, ad_v, ad_ubar, ad_vbar, ad_zt_avg1, ad_zeta, ad_t_sol, ad_u_sol, ad_v_sol, ad_zeta_sol)
 

Function/Subroutine Documentation

◆ ad_forcing()

subroutine, public ad_forcing_mod::ad_forcing ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) kfrc,
integer, intent(in) nfrc )

Definition at line 25 of file ad_forcing.F.

26!***********************************************************************
27!
28 USE mod_param
29 USE mod_ocean
30# ifdef SOLVE3D
31 USE mod_coupling
32# endif
33!
34! Imported variable declarations.
35!
36 integer, intent(in) :: ng, tile, Kfrc, Nfrc
37!
38! Local variable declarations.
39!
40# include "tile.h"
41!
42 CALL ad_forcing_tile (ng, tile, &
43 & lbi, ubi, lbj, ubj, &
44 & imins, imaxs, jmins, jmaxs, &
45 & kfrc, nfrc, &
46# ifdef SOLVE3D
47 & ocean(ng) % f_t, &
48 & ocean(ng) % f_u, &
49 & ocean(ng) % f_v, &
50# ifdef FORCING_SV
51 & ocean(ng) % f_ubar, &
52 & ocean(ng) % f_vbar, &
53# endif
54# else
55 & ocean(ng) % f_ubar, &
56 & ocean(ng) % f_vbar, &
57# endif
58 & ocean(ng) % f_zeta, &
59# ifdef SOLVE3D
60 & ocean(ng) % ad_t, &
61 & ocean(ng) % ad_u, &
62 & ocean(ng) % ad_v, &
63# ifdef FORCING_SV
64 & ocean(ng) % ad_ubar, &
65 & ocean(ng) % ad_vbar, &
66# endif
67# else
68 & ocean(ng) % ad_ubar, &
69 & ocean(ng) % ad_vbar, &
70# endif
71# ifdef SOLVE3D
72 & coupling(ng) % ad_Zt_avg1, &
73# endif
74 & ocean(ng) % ad_zeta, &
75# ifdef SOLVE3D
76 & ocean(ng) % ad_t_sol, &
77 & ocean(ng) % ad_u_sol, &
78 & ocean(ng) % ad_v_sol, &
79# else
80 & ocean(ng) % ad_ubar_sol, &
81 & ocean(ng) % ad_vbar_sol, &
82# endif
83 & ocean(ng) % ad_zeta_sol)
84
85 RETURN
type(t_coupling), dimension(:), allocatable coupling
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351

References ad_forcing_tile(), mod_coupling::coupling, and mod_ocean::ocean.

Referenced by ad_main3d().

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

◆ ad_forcing_tile()

subroutine ad_forcing_mod::ad_forcing_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) kfrc,
integer, intent(in) nfrc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) f_t,
real(r8), dimension(lbi:,lbj:,:), intent(inout) f_u,
real(r8), dimension(lbi:,lbj:,:), intent(inout) f_v,
real(r8), dimension(lbi:,lbj:), intent(inout) f_ubar,
real(r8), dimension(lbi:,lbj:), intent(inout) f_vbar,
real(r8), dimension(lbi:,lbj:), intent(inout) f_zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_vbar,
real(r8), dimension(lbi:,lbj:), intent(inout) ad_zt_avg1,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_zeta,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_t_sol,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_u_sol,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_v_sol,
real(r8), dimension(lbi:,lbj:), intent(inout) ad_zeta_sol )
private

Definition at line 89 of file ad_forcing.F.

120!***********************************************************************
121!
122 USE mod_param
123 USE mod_scalars
124!
125! Imported variable declarations.
126!
127 integer, intent(in) :: ng, tile
128 integer, intent(in) :: LBi, UBi, LBj, UBj
129 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
130 integer, intent(in) :: Kfrc
131 integer, intent(in) :: Nfrc
132!
133# ifdef ASSUMED_SHAPE
134# ifdef SOLVE3D
135 real(r8), intent(inout) :: f_t(LBi:,LBj:,:,:)
136 real(r8), intent(inout) :: f_u(LBi:,LBj:,:)
137 real(r8), intent(inout) :: f_v(LBi:,LBj:,:)
138# ifdef FORCING_SV
139 real(r8), intent(inout) :: f_ubar(LBi:,LBj:)
140 real(r8), intent(inout) :: f_vbar(LBi:,LBj:)
141# endif
142# else
143 real(r8), intent(inout) :: f_ubar(LBi:,LBj:)
144 real(r8), intent(inout) :: f_vbar(LBi:,LBj:)
145# endif
146 real(r8), intent(inout) :: f_zeta(LBi:,LBj:)
147# ifdef SOLVE3D
148 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
149 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
150 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
151# ifdef FORCING_SV
152 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
153 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
154# endif
155# else
156 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
157 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
158# endif
159# ifdef SOLVE3D
160 real(r8), intent(inout) :: ad_Zt_avg1(LBi:,LBj:)
161# endif
162 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
163# ifdef SOLVE3D
164 real(r8), intent(inout) :: ad_t_sol(LBi:,LBj:,:,:)
165 real(r8), intent(inout) :: ad_u_sol(LBi:,LBj:,:)
166 real(r8), intent(inout) :: ad_v_sol(LBi:,LBj:,:)
167# else
168 real(r8), intent(inout) :: ad_ubar_sol(LBi:,LBj:)
169 real(r8), intent(inout) :: ad_vbar_sol(LBi:,LBj:)
170# endif
171 real(r8), intent(inout) :: ad_zeta_sol(LBi:,LBj:)
172# else
173# ifdef SOLVE3D
174 real(r8), intent(inout) :: f_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
175 real(r8), intent(inout) :: f_u(LBi:UBi,LBj:UBj,N(ng))
176 real(r8), intent(inout) :: f_v(LBi:UBi,LBj:UBj,N(ng))
177# ifdef FORCING_SV
178 real(r8), intent(inout) :: f_ubar(LBi:UBi,LBj:UBj)
179 real(r8), intent(inout) :: f_vbar(LBi:UBi,LBj:UBj)
180# endif
181# else
182 real(r8), intent(inout) :: f_ubar(LBi:UBi,LBj:UBj)
183 real(r8), intent(inout) :: f_vbar(LBi:UBi,LBj:UBj)
184# endif
185 real(r8), intent(inout) :: f_zeta(LBi:UBi,LBj:UBj)
186# ifdef SOLVE3D
187 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
188 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
189 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
190# ifdef FORCING_SV
191 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
192 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
193# endif
194# else
195 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
196 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
197# endif
198# ifdef SOLVE3D
199 real(r8), intent(inout) :: ad_Zt_avg1(LBi:UBi,LBj:UBj)
200# endif
201 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,:)
202# ifdef SOLVE3D
203 real(r8), intent(inout) :: ad_t_sol(LBi:UBi,LBj:UBj,N(ng),NT(ng))
204 real(r8), intent(inout) :: ad_u_sol(LBi:UBi,LBj:UBj,N(ng))
205 real(r8), intent(inout) :: ad_v_sol(LBi:UBi,LBj:UBj,N(ng))
206# else
207 real(r8), intent(inout) :: ad_ubar_sol(LBi:UBi,LBj:UBj)
208 real(r8), intent(inout) :: ad_vbar_sol(LBi:UBi,LBj:UBj)
209# endif
210 real(r8), intent(inout) :: ad_zeta_sol(LBi:UBi,LBj:UBj)
211# endif
212!
213! Local variable declarations.
214!
215 integer :: i, j
216# ifdef SOLVE3D
217 integer :: itrc, k
218# endif
219
220# include "set_bounds.h"
221!
222!
223! Adjoint linear free-surface. The two different cases in the
224! case of SOLVE3D are due to the fact that tl_ini_fields is
225! also called on the first timestep. tl_forcing MUST be called
226! before tl_ini_fields.
227!
228# ifdef SOLVE3D
229 DO j=jstrr,jendr
230 DO i=istrr,iendr
231 f_zeta(i,j)=f_zeta(i,j)+ad_zt_avg1(i,j)
232 END DO
233 END DO
234# else
235 DO j=jstrr,jendr
236 DO i=istrr,iendr
237 f_zeta(i,j)=f_zeta(i,j)+ad_zeta(i,j,kfrc)
238 END DO
239 END DO
240# endif
241
242# ifndef SOLVE3D
243!
244! Adjoint linear 2D momentum.
245!
246 DO j=jstrr,jendr
247 DO i=istr,iendr
248 f_ubar(i,j)=f_ubar(i,j)+ad_ubar(i,j,kfrc)
249 END DO
250 END DO
251!
252 DO j=jstr,jendr
253 DO i=istrr,iendr
254 f_vbar(i,j)=f_vbar(i,j)+ad_vbar(i,j,kfrc)
255 END DO
256 END DO
257
258# else
259# ifdef FORCING_SV
260!
261! Adjoint linear 2D momentum.
262!
263 DO j=jstrr,jendr
264 DO i=istr,iendr
265 f_ubar(i,j)=f_ubar(i,j)+ad_ubar(i,j,1)+ad_ubar(i,j,2)
266 END DO
267 END DO
268!
269 DO j=jstr,jendr
270 DO i=istrr,iendr
271 f_vbar(i,j)=f_vbar(i,j)+ad_vbar(i,j,1)+ad_vbar(i,j,2)
272 END DO
273 END DO
274# endif
275!
276! Adjoint linear 3D momentum.
277!
278 DO k=1,n(ng)
279 DO j=jstrr,jendr
280 DO i=istr,iendr
281 f_u(i,j,k)=f_u(i,j,k)+ad_u(i,j,k,nfrc)
282 END DO
283 END DO
284 DO j=jstr,jendr
285 DO i=istrr,iendr
286 f_v(i,j,k)=f_v(i,j,k)+ad_v(i,j,k,nfrc)
287 END DO
288 END DO
289 END DO
290!
291! Adjoint linear tracers.
292!
293 DO itrc=1,nt(ng)
294 DO k=1,n(ng)
295 DO j=jstrr,jendr
296 DO i=istrr,iendr
297 f_t(i,j,k,itrc)=f_t(i,j,k,itrc)+ &
298 & ad_t(i,j,k,nfrc,itrc)
299 END DO
300 END DO
301 END DO
302 END DO
303# endif
304
305 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable nt
Definition mod_param.F:489

Referenced by ad_forcing().

Here is the caller graph for this function: