ROMS
Loading...
Searching...
No Matches
frc_iau.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#if defined RBL4DVAR && defined RPCG
5!
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2019 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.txt !
10!=======================================================================
11! !
12! This routine is used to compute the IAU forcing for the NLM. !
13! !
14! The method is described in: !
15! Bloom et al., 1996: Data Assimilation Using Incremental Analysis !
16! Update. Mon. Wea. Rev., 124, 1256-1271. !
17! !
18!=======================================================================
19!
20 USE mod_param
21 USE mod_forces
22 USE mod_fourdvar
23 USE mod_ocean
24 USE mod_scalars
25 USE mod_stepping
26
27
28!
29 implicit none
30!
31 PRIVATE
32 PUBLIC :: frc_iau
33 PUBLIC :: frc_iau_ini
34!
35 CONTAINS
36!
37 SUBROUTINE frc_iau_ini (ng, tile, irec)
38!
39!=======================================================================
40! !
41! This subroutine computes the incremental analysis update (IAU) !
42! that is applied as a forcing term in the non-linear model during !
43! the first timeIAU seconds of the assimilation window. !
44! !
45! The adjoint arrays ad_var(irec) are used as temporary storage for !
46! the IAU. !
47! !
48! On Input: !
49! !
50! ng Nested grid number. !
51! tile Domain partition. !
52! irec ad_var record containing the fields on input. !
53! !
54!=======================================================================
55!
56! Imported variable declarations.
57!
58 integer, intent(in) :: ng, tile, irec
59!
60! Local variable declarations.
61!
62# include "tile.h"
63!
64# ifdef PROFILE
65 CALL wclock_on (ng, iadm, 7, __line__, __file__)
66# endif
67 CALL frc_iau_ini_tile (ng, tile, &
68 & lbi, ubi, lbj, ubj, &
69 & imins, imaxs, jmins, jmaxs, &
70 & irec, &
71 & ocean(ng) % ad_zeta, &
72# ifdef SOLVE3D
73 & ocean(ng) % ad_u, &
74 & ocean(ng) % ad_v, &
75 & ocean(ng) % ad_t, &
76# else
77 & ocean(ng) % ad_ubar, &
78 & ocean(ng) % ad_vbar, &
79# endif
80 & ocean(ng) % f_zeta, &
81# ifdef SOLVE3D
82 & ocean(ng) % f_u, &
83 & ocean(ng) % f_v, &
84 & ocean(ng) % f_t)
85# else
86 & ocean(ng) % f_ubar, &
87 & ocean(ng) % f_vbar)
88# endif
89# ifdef PROFILE
90 CALL wclock_off (ng, iadm, 7, __line__, __file__)
91# endif
92!
93 RETURN
94 END SUBROUTINE frc_iau_ini
95!
96!***********************************************************************
97 SUBROUTINE frc_iau_ini_tile (ng, tile, &
98 & LBi, UBi, LBj, UBj, &
99 & IminS, ImaxS, JminS, JmaxS, &
100 & irec, &
101 & ad_zeta, &
102# ifdef SOLVE3D
103 & ad_u, ad_v, ad_t, &
104# else
105 & ad_ubar, ad_vbar, &
106# endif
107 & f_zeta, &
108# ifdef SOLVE3D
109 & f_u, f_v, f_t)
110# else
111 & f_ubar, f_vbar)
112# endif
113!***********************************************************************
114!
115! Imported variable declarations.
116!
117 integer, intent(in) :: ng, tile, irec
118 integer, intent(in) :: LBi, UBi, LBj, UBj
119 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
120!
121# ifdef ASSUMED_SHAPE
122 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
123 real(r8), intent(inout) :: f_zeta(LBi:,LBj:)
124# ifdef SOLVE3D
125 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
126 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
127 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
128 real(r8), intent(inout) :: f_u(LBi:,LBj:,:)
129 real(r8), intent(inout) :: f_v(LBi:,LBj:,:)
130 real(r8), intent(inout) :: f_t(LBi:,LBj:,:,:)
131# else
132 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
133 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
134 real(r8), intent(inout) :: f_ubar(LBi:,LBj:)
135 real(r8), intent(inout) :: f_vbar(LBi:,LBj:)
136# endif
137# else
138 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,3)
139 real(r8), intent(inout) :: f_zeta(LBi:UBi,LBj:UBj)
140# ifdef SOLVE3D
141 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
142 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
143 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),2,NT(ng))
144 real(r8), intent(inout) :: f_u(LBi:UBi,LBj:UBj,N(ng))
145 real(r8), intent(inout) :: f_v(LBi:UBi,LBj:UBj,N(ng))
146 real(r8), intent(inout) :: f_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
147# else
148 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,3)
149 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,3)
150 real(r8), intent(inout) :: f_ubar(LBi:UBi,LBj:UBj)
151 real(r8), intent(inout) :: f_vbar(LBi:UBi,LBj:UBj)
152# endif
153# endif
154!
155! Local variable declarations.
156!
157 integer :: i, it1, it2, j, k
158# ifdef SOLVE3D
159 integer :: itrc
160# endif
161 real(r8) :: fac1
162
163# include "set_bounds.h"
164!
165!-----------------------------------------------------------------------
166! Compute IAU forcing terms used in the nonlinear model.
167! Use ad_var arrays as temporary storage for incremental analysis
168! update forcing - needed when weak constraing is also activated
169! and used in in the routine frc_iau.
170!-----------------------------------------------------------------------
171!
172! Set uniform weights in time.
173!
174 fac1=real(timeiau(ng)/dt(ng),r8)
175 fac1=1.0_r8/fac1
176!
177 DO j=jstrr,jendr
178 DO i=istrr,iendr
179 f_zeta(i,j)=fac1*ad_zeta(i,j,irec)
180 ad_zeta(i,j,irec)=f_zeta(i,j)
181 END DO
182 END DO
183
184# ifndef SOLVE3D
185!
186! Compute 2D-momentum iau forcing terms.
187!
188 DO j=jstrr,jendr
189 DO i=istr,iendr
190 f_ubar(i,j)=fac1*ad_ubar(i,j,irec)
191 ad_ubar(i,j,irec)=f_ubar(i,j)
192 END DO
193 END DO
194 DO j=jstr,jendr
195 DO i=istrr,iendr
196 f_vbar(i,j)=fac1*ad_vbar(i,j,irec)
197 ad_vbar(i,j,irec)=f_vbar(i,j)
198 END DO
199 END DO
200# endif
201# ifdef SOLVE3D
202!
203! Compute 3D-momentum iau forcing terms.
204!
205 DO k=1,n(ng)
206 DO j=jstrr,jendr
207 DO i=istr,iendr
208 f_u(i,j,k)=fac1*ad_u(i,j,k,irec)
209 ad_u(i,j,k,irec)=f_u(i,j,k)
210 END DO
211 END DO
212 END DO
213 DO k=1,n(ng)
214 DO j=jstr,jendr
215 DO i=istrr,iendr
216 f_v(i,j,k)=fac1*ad_v(i,j,k,irec)
217 ad_v(i,j,k,irec)=f_v(i,j,k)
218 END DO
219 END DO
220 END DO
221!
222! Compute tracer iau forcing terms.
223!
224 DO itrc=1,nt(ng)
225 DO k=1,n(ng)
226 DO j=jstrr,jendr
227 DO i=istrr,iendr
228 f_t(i,j,k,itrc)=fac1*ad_t(i,j,k,irec,itrc)
229 ad_t(i,j,k,irec,itrc)=f_t(i,j,k,itrc)
230 END DO
231 END DO
232 END DO
233 END DO
234# endif
235!
236 RETURN
237 END SUBROUTINE frc_iau_ini_tile
238!
239 SUBROUTINE frc_iau (ng, tile, irec)
240!
241!=======================================================================
242! !
243! This subroutine computes the combined weak constraint forcing !
244! and the incremental analysis update (IAU) when both weak constraint !
245! data assimilation and the IAU are activated. !
246! The input adjoint arrays ad_var(irec) hold the IAU forcing. !
247! !
248! On Input: !
249! !
250! ng Nested grid number. !
251! tile Domain partition. !
252! irec ad_var record containing the fields on input. !
253! !
254!=======================================================================
255!
256! Imported variable declarations.
257!
258 integer, intent(in) :: ng, tile, irec
259!
260! Local variable declarations.
261!
262# include "tile.h"
263!
264# ifdef PROFILE
265 CALL wclock_on (ng, iadm, 7, __line__, __file__)
266# endif
267 CALL frc_iau_tile (ng, tile, &
268 & lbi, ubi, lbj, ubj, &
269 & imins, imaxs, jmins, jmaxs, &
270 & irec, &
271 & ocean(ng) % ad_zeta, &
272# ifdef SOLVE3D
273 & ocean(ng) % ad_u, &
274 & ocean(ng) % ad_v, &
275 & ocean(ng) % ad_t, &
276# else
277 & ocean(ng) % ad_ubar, &
278 & ocean(ng) % ad_vbar, &
279# endif
280 & ocean(ng) % f_zeta, &
281# ifdef SOLVE3D
282 & ocean(ng) % f_u, &
283 & ocean(ng) % f_v, &
284 & ocean(ng) % f_t)
285# else
286 & ocean(ng) % f_ubar, &
287 & ocean(ng) % f_vbar)
288# endif
289# ifdef PROFILE
290 CALL wclock_off (ng, iadm, 7, __line__, __file__)
291# endif
292!
293 RETURN
294 END SUBROUTINE frc_iau
295!
296!***********************************************************************
297 SUBROUTINE frc_iau_tile (ng, tile, &
298 & LBi, UBi, LBj, UBj, &
299 & IminS, ImaxS, JminS, JmaxS, &
300 & irec, &
301 & ad_zeta, &
302# ifdef SOLVE3D
303 & ad_u, ad_v, ad_t, &
304# else
305 & ad_ubar, ad_vbar, &
306# endif
307 & f_zeta, &
308# ifdef SOLVE3D
309 & f_u, f_v, f_t)
310# else
311 & f_ubar, f_vbar)
312# endif
313!***********************************************************************
314!
315! Imported variable declarations.
316!
317 integer, intent(in) :: ng, tile, irec
318 integer, intent(in) :: LBi, UBi, LBj, UBj
319 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
320!
321# ifdef ASSUMED_SHAPE
322 real(r8), intent(in) :: ad_zeta(LBi:,LBj:,:)
323 real(r8), intent(inout) :: f_zeta(LBi:,LBj:)
324# ifdef SOLVE3D
325 real(r8), intent(in) :: ad_u(LBi:,LBj:,:,:)
326 real(r8), intent(in) :: ad_v(LBi:,LBj:,:,:)
327 real(r8), intent(in) :: ad_t(LBi:,LBj:,:,:,:)
328 real(r8), intent(inout) :: f_u(LBi:,LBj:,:)
329 real(r8), intent(inout) :: f_v(LBi:,LBj:,:)
330 real(r8), intent(inout) :: f_t(LBi:,LBj:,:,:)
331# else
332 real(r8), intent(in) :: ad_ubar(LBi:,LBj:,:)
333 real(r8), intent(in) :: ad_vbar(LBi:,LBj:,:)
334 real(r8), intent(inout) :: f_ubar(LBi:,LBj:)
335 real(r8), intent(inout) :: f_vbar(LBi:,LBj:)
336# endif
337# else
338 real(r8), intent(in) :: ad_zeta(LBi:UBi,LBj:UBj,3)
339 real(r8), intent(inout) :: f_zeta(LBi:UBi,LBj:UBj)
340# ifdef SOLVE3D
341 real(r8), intent(in) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
342 real(r8), intent(in) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
343 real(r8), intent(in) :: ad_t(LBi:UBi,LBj:UBj,N(ng),2,NT(ng))
344 real(r8), intent(inout) :: f_u(LBi:UBi,LBj:UBj,N(ng))
345 real(r8), intent(inout) :: f_v(LBi:UBi,LBj:UBj,N(ng))
346 real(r8), intent(inout) :: f_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
347# else
348 real(r8), intent(in) :: ad_ubar(LBi:UBi,LBj:UBj,3)
349 real(r8), intent(in) :: ad_vbar(LBi:UBi,LBj:UBj,3)
350 real(r8), intent(inout) :: f_ubar(LBi:UBi,LBj:UBj)
351 real(r8), intent(inout) :: f_vbar(LBi:UBi,LBj:UBj)
352# endif
353# endif
354!
355! Local variable declarations.
356!
357 integer :: i, it1, it2, j, k
358# ifdef SOLVE3D
359 integer :: itrc
360# endif
361
362# include "set_bounds.h"
363!
364!-----------------------------------------------------------------------
365! Update the weak constraint forcing arrays with the IAU forcing
366! which is stored in the ad_var(irec) arrays.
367!-----------------------------------------------------------------------
368!
369 DO j=jstrr,jendr
370 DO i=istrr,iendr
371 f_zeta(i,j)=f_zeta(i,j)+ad_zeta(i,j,irec)
372 END DO
373 END DO
374
375# ifndef SOLVE3D
376!
377! Add 2D-momentum iau forcing terms.
378!
379 DO j=jstrr,jendr
380 DO i=istr,iendr
381 f_ubar(i,j)=f_ubar(i,j)+ad_ubar(i,j,irec)
382 END DO
383 END DO
384 DO j=jstr,jendr
385 DO i=istrr,iendr
386 f_vbar(i,j)=f_vbar(i,j)+ad_vbar(i,j,irec)
387 END DO
388 END DO
389# endif
390# ifdef SOLVE3D
391!
392! Add 3D-momentum iau forcing terms.
393!
394 DO k=1,n(ng)
395 DO j=jstrr,jendr
396 DO i=istr,iendr
397 f_u(i,j,k)=f_u(i,j,k)+ad_u(i,j,k,irec)
398 END DO
399 END DO
400 END DO
401 DO k=1,n(ng)
402 DO j=jstr,jendr
403 DO i=istrr,iendr
404 f_v(i,j,k)=f_v(i,j,k)+ad_v(i,j,k,irec)
405 END DO
406 END DO
407 END DO
408!
409! Add tracer iau forcing terms.
410!
411 DO itrc=1,nt(ng)
412 DO k=1,n(ng)
413 DO j=jstrr,jendr
414 DO i=istrr,iendr
415 f_t(i,j,k,itrc)=f_t(i,j,k,itrc)+ad_t(i,j,k,irec,itrc)
416 END DO
417 END DO
418 END DO
419 END DO
420# endif
421!
422 RETURN
423 END SUBROUTINE frc_iau_tile
424#endif
425 END MODULE frc_iau_mod
subroutine, public frc_iau_ini(ng, tile, irec)
Definition frc_iau.F:38
subroutine frc_iau_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, irec, ad_zeta, ad_u, ad_v, ad_t, f_zeta, f_u, f_v, f_t)
Definition frc_iau.F:310
subroutine, public frc_iau(ng, tile, irec)
Definition frc_iau.F:240
subroutine frc_iau_ini_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, irec, ad_zeta, ad_u, ad_v, ad_t, f_zeta, f_u, f_v, f_t)
Definition frc_iau.F:110
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter iadm
Definition mod_param.F:665
real(dp), dimension(:), allocatable dt
real(dp), dimension(:), allocatable timeiau
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