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

Functions/Subroutines

subroutine, public ad_variability (ng, tile, linp, lweak)
 
subroutine ad_variability_tile (ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, linp, lweak, t_obc_std, u_obc_std, v_obc_std, ubar_obc_std, vbar_obc_std, zeta_obc_std, sustr_std, svstr_std, stflx_std, t_std, u_std, v_std, ubar_std, vbar_std, zeta_std, ad_t_obc, ad_u_obc, ad_v_obc, ad_ubar_obc, ad_vbar_obc, ad_zeta_obc, ad_ustr, ad_vstr, ad_tflux, ad_t, ad_u, ad_v, ad_ubar, ad_vbar, ad_zeta)
 

Function/Subroutine Documentation

◆ ad_variability()

subroutine, public ad_variability_mod::ad_variability ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
logical, intent(in) lweak )

Definition at line 55 of file ad_variability.F.

56!***********************************************************************
57!
58 USE mod_param
59# ifdef ADJUST_BOUNDARY
60 USE mod_boundary
61# endif
62# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
63 USE mod_forces
64# endif
65 USE mod_ocean
66!
67! Imported variable declarations.
68!
69 logical, intent(in) :: Lweak
70
71 integer, intent(in) :: ng, tile, Linp
72!
73! Local variable declarations.
74!
75# include "tile.h"
76!
77 CALL ad_variability_tile (ng, tile, &
78 & lbi, ubi, lbj, ubj, lbij, ubij, &
79 & imins, imaxs, jmins, jmaxs, &
80 & linp, lweak, &
81# ifdef ADJUST_BOUNDARY
82# ifdef SOLVE3D
83 & boundary(ng) % e_t_obc, &
84 & boundary(ng) % e_u_obc, &
85 & boundary(ng) % e_v_obc, &
86# endif
87 & boundary(ng) % e_ubar_obc, &
88 & boundary(ng) % e_vbar_obc, &
89 & boundary(ng) % e_zeta_obc, &
90# endif
91# ifdef ADJUST_WSTRESS
92 & forces(ng) % e_sustr, &
93 & forces(ng) % e_svstr, &
94# endif
95# if defined ADJUST_STFLUX && defined SOLVE3D
96 & forces(ng) % e_stflx, &
97# endif
98# ifdef SOLVE3D
99 & ocean(ng) % e_t, &
100 & ocean(ng) % e_u, &
101 & ocean(ng) % e_v, &
102# if defined WEAK_CONSTRAINT && defined TIME_CONV
103 & ocean(ng) % e_ubar, &
104 & ocean(ng) % e_vbar, &
105# endif
106# else
107 & ocean(ng) % e_ubar, &
108 & ocean(ng) % e_vbar, &
109# endif
110 & ocean(ng) % e_zeta, &
111# ifdef ADJUST_BOUNDARY
112# ifdef SOLVE3D
113 & boundary(ng) % ad_t_obc, &
114 & boundary(ng) % ad_u_obc, &
115 & boundary(ng) % ad_v_obc, &
116# endif
117 & boundary(ng) % ad_ubar_obc, &
118 & boundary(ng) % ad_vbar_obc, &
119 & boundary(ng) % ad_zeta_obc, &
120# endif
121# ifdef ADJUST_WSTRESS
122 & forces(ng) % ad_ustr, &
123 & forces(ng) % ad_vstr, &
124# endif
125# if defined ADJUST_STFLUX && defined SOLVE3D
126 & forces(ng) % ad_tflux, &
127# endif
128# ifdef SOLVE3D
129 & ocean(ng) % ad_t, &
130 & ocean(ng) % ad_u, &
131 & ocean(ng) % ad_v, &
132# if defined WEAK_CONSTRAINT && defined TIME_CONV
133 & ocean(ng) % ad_ubar, &
134 & ocean(ng) % ad_vbar, &
135# endif
136# else
137 & ocean(ng) % ad_ubar, &
138 & ocean(ng) % ad_vbar, &
139# endif
140 & ocean(ng) % ad_zeta)
141
142 RETURN
type(t_boundary), dimension(:), allocatable boundary
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351

References ad_variability_tile(), mod_boundary::boundary, mod_forces::forces, and mod_ocean::ocean.

Referenced by convolve_mod::convolve(), convolve_mod::error_covariance(), i4dvar_mod::increment(), and convolve_mod::saddlec().

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

◆ ad_variability_tile()

subroutine ad_variability_mod::ad_variability_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) lbij,
integer, intent(in) ubij,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) linp,
logical, intent(in) lweak,
real(r8), dimension(lbij:,:,:,:), intent(in) t_obc_std,
real(r8), dimension(lbij:,:,:), intent(in) u_obc_std,
real(r8), dimension(lbij:,:,:), intent(in) v_obc_std,
real(r8), dimension(lbij:,:), intent(in) ubar_obc_std,
real(r8), dimension(lbij:,:), intent(in) vbar_obc_std,
real(r8), dimension(lbij:,:), intent(in) zeta_obc_std,
real(r8), dimension(lbi:,lbj:), intent(in) sustr_std,
real(r8), dimension(lbi:,lbj:), intent(in) svstr_std,
real(r8), dimension(lbi:,lbj:,:), intent(in) stflx_std,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) t_std,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) u_std,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) v_std,
real(r8), dimension(lbi:,lbj:,:), intent(in) ubar_std,
real(r8), dimension(lbi:,lbj:,:), intent(in) vbar_std,
real(r8), dimension(lbi:,lbj:,:), intent(in) zeta_std,
real(r8), dimension(lbij:,:,:,:,:,:), intent(inout) ad_t_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) ad_u_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) ad_v_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_ubar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_vbar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_zeta_obc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_vstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_tflux,
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_zeta )
private

Definition at line 146 of file ad_variability.F.

194!***********************************************************************
195!
196 USE mod_param
197 USE mod_ncparam
198 USE mod_scalars
199!
200# ifdef DISTRIBUTE
202# endif
203!
204! Imported variable declarations.
205!
206 logical, intent(in) :: Lweak
207
208 integer, intent(in) :: ng, tile
209 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
210 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
211 integer, intent(in) :: Linp
212!
213# ifdef ASSUMED_SHAPE
214# ifdef ADJUST_BOUNDARY
215# ifdef SOLVE3D
216 real(r8), intent(in) :: t_obc_std(LBij:,:,:,:)
217 real(r8), intent(in) :: u_obc_std(LBij:,:,:)
218 real(r8), intent(in) :: v_obc_std(LBij:,:,:)
219# endif
220 real(r8), intent(in) :: ubar_obc_std(LBij:,:)
221 real(r8), intent(in) :: vbar_obc_std(LBij:,:)
222 real(r8), intent(in) :: zeta_obc_std(LBij:,:)
223# endif
224# ifdef ADJUST_WSTRESS
225 real(r8), intent(in) :: sustr_std(LBi:,LBj:)
226 real(r8), intent(in) :: svstr_std(LBi:,LBj:)
227# endif
228# if defined ADJUST_STFLUX && defined SOLVE3D
229 real(r8), intent(in) :: stflx_std(LBi:,LBj:,:)
230# endif
231# ifdef SOLVE3D
232 real(r8), intent(in) :: t_std(LBi:,LBj:,:,:,:)
233 real(r8), intent(in) :: u_std(LBi:,LBj:,:,:)
234 real(r8), intent(in) :: v_std(LBi:,LBj:,:,:)
235# if defined WEAK_CONSTRAINT && defined TIME_CONV
236 real(r8), intent(in) :: ubar_std(LBi:,LBj:,:)
237 real(r8), intent(in) :: vbar_std(LBi:,LBj:,:)
238# endif
239# else
240 real(r8), intent(in) :: ubar_std(LBi:,LBj:,:)
241 real(r8), intent(in) :: vbar_std(LBi:,LBj:,:)
242# endif
243 real(r8), intent(in) :: zeta_std(LBi:,LBj:,:)
244# ifdef ADJUST_BOUNDARY
245# ifdef SOLVE3D
246 real(r8), intent(inout) :: ad_t_obc(LBij:,:,:,:,:,:)
247 real(r8), intent(inout) :: ad_u_obc(LBij:,:,:,:,:)
248 real(r8), intent(inout) :: ad_v_obc(LBij:,:,:,:,:)
249# endif
250 real(r8), intent(inout) :: ad_ubar_obc(LBij:,:,:,:)
251 real(r8), intent(inout) :: ad_vbar_obc(LBij:,:,:,:)
252 real(r8), intent(inout) :: ad_zeta_obc(LBij:,:,:,:)
253# endif
254# ifdef ADJUST_WSTRESS
255 real(r8), intent(inout) :: ad_ustr(LBi:,LBj:,:,:)
256 real(r8), intent(inout) :: ad_vstr(LBi:,LBj:,:,:)
257# endif
258# if defined ADJUST_STFLUX && defined SOLVE3D
259 real(r8), intent(inout) :: ad_tflux(LBi:,LBj:,:,:,:)
260# endif
261# ifdef SOLVE3D
262 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
263 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
264 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
265# if defined WEAK_CONSTRAINT && defined TIME_CONV
266 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
267 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
268# endif
269# else
270 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
271 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
272# endif
273 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
274
275# else
276
277# ifdef ADJUST_BOUNDARY
278# ifdef SOLVE3D
279 real(r8), intent(in) :: t_obc_std(LBij:UBij,N(ng),4,NT(ng))
280 real(r8), intent(in) :: u_obc_std(LBij:UBij,N(ng),4)
281 real(r8), intent(in) :: v_obc_std(LBij:UBij,N(ng),4)
282# endif
283 real(r8), intent(in) :: ubar_obc_std(LBij:UBij,4)
284 real(r8), intent(in) :: vbar_obc_std(LBij:UBij,4)
285 real(r8), intent(in) :: zeta_obc_std(LBij:UBij,4)
286# endif
287# ifdef ADJUST_WSTRESS
288 real(r8), intent(in) :: sustr_std(LBi:,LBj:)
289 real(r8), intent(in) :: svstr_std(LBi:,LBj:)
290# endif
291# if defined ADJUST_STFLUX && defined SOLVE3D
292 real(r8), intent(in) :: stflx_std(LBi:UBi,LBj:UBj,NT(ng))
293# endif
294# ifdef SOLVE3D
295 real(r8), intent(in) :: t_std(LBi:UBi,LBj:UBj,N(ng),NSA,NT(ng))
296 real(r8), intent(in) :: u_std(LBi:UBi,LBj:UBj,N(ng),NSA)
297 real(r8), intent(in) :: v_std(LBi:UBi,LBj:UBj,N(ng),NSA)
298# if defined WEAK_CONSTRAINT && defined TIME_CONV
299 real(r8), intent(in) :: ubar_std(LBi:UBi,LBj:UBj,NSA)
300 real(r8), intent(in) :: vbar_std(LBi:UBi,LBj:UBj,NSA)
301# endif
302# else
303 real(r8), intent(in) :: ubar_std(LBi:UBi,LBj:UBj,NSA)
304 real(r8), intent(in) :: vbar_std(LBi:UBi,LBj:UBj,NSA)
305# endif
306 real(r8), intent(in) :: zeta_std(LBi:UBi,LBj:UBj,NSA)
307# ifdef ADJUST_BOUNDARY
308# ifdef SOLVE3D
309 real(r8), intent(inout) :: ad_t_obc(LBij:UBij,N(ng),4, &
310 & Nbrec(ng),2,NT(ng))
311 real(r8), intent(inout) :: ad_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
312 real(r8), intent(inout) :: ad_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
313# endif
314 real(r8), intent(inout) :: ad_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
315 real(r8), intent(inout) :: ad_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
316 real(r8), intent(inout) :: ad_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
317# endif
318# ifdef ADJUST_WSTRESS
319 real(r8), intent(inout) :: ad_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
320 real(r8), intent(inout) :: ad_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
321# endif
322# if defined ADJUST_STFLUX && defined SOLVE3D
323 real(r8), intent(inout) :: ad_tflux(LBi:UBi,LBj:UBj, &
324 & Nfrec(ng),2,NT(ng))
325# endif
326# ifdef SOLVE3D
327 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
328 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
329 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
330# if defined WEAK_CONSTRAINT && defined TIME_CONV
331 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
332 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
333# endif
334# else
335 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
336 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
337# endif
338 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,:)
339# endif
340!
341! Local variable declarations.
342!
343 integer :: i, ib, ir, j, rec
344# ifdef SOLVE3D
345 integer :: it, k
346# endif
347
348# include "set_bounds.h"
349!
350!-----------------------------------------------------------------------
351! Determine error covariance standard deviation factors to use.
352!-----------------------------------------------------------------------
353!
354 IF (lweak) THEN
355 rec=2 ! weak constraint
356 ELSE
357 rec=1 ! strong constraint
358 END IF
359!
360!-----------------------------------------------------------------------
361! Initial conditions and model error covariance: Multiply adjoint state
362! by its corresponding standard deviation.
363!-----------------------------------------------------------------------
364!
365! Free-surface.
366!
367# ifdef DISTRIBUTE
368 CALL ad_mp_exchange2d (ng, tile, iadm, 1, &
369 & lbi, ubi, lbj, ubj, &
370 & nghostpoints, &
371 & ewperiodic(ng), nsperiodic(ng), &
372 & ad_zeta(:,:,linp))
373# endif
374 DO j=jstrt,jendt
375 DO i=istrt,iendt
376 ad_zeta(i,j,linp)=ad_zeta(i,j,linp)*zeta_std(i,j,rec)
377 END DO
378 END DO
379# if !defined SOLVE3D || (defined WEAK_CONSTRAINT && defined TIME_CONV)
380!
381! 2D U-momentum.
382!
383# ifdef DISTRIBUTE
384 CALL ad_mp_exchange2d (ng, tile, iadm, 2, &
385 & lbi, ubi, lbj, ubj, &
386 & nghostpoints, &
387 & ewperiodic(ng), nsperiodic(ng), &
388 & ad_ubar(:,:,linp), &
389 & ad_vbar(:,:,linp))
390# endif
391 DO j=jstrt,jendt
392 DO i=istrp,iendt
393 ad_ubar(i,j,linp)=ad_ubar(i,j,linp)*ubar_std(i,j,rec)
394 END DO
395 END DO
396 DO j=jstrp,jendt
397 DO i=istrt,iendt
398 ad_vbar(i,j,linp)=ad_vbar(i,j,linp)*vbar_std(i,j,rec)
399 END DO
400 END DO
401# endif
402# ifdef SOLVE3D
403!
404! 3D U-momentum.
405!
406# ifdef DISTRIBUTE
407 CALL ad_mp_exchange3d (ng, tile, iadm, 2, &
408 & lbi, ubi, lbj, ubj, 1, n(ng), &
409 & nghostpoints, &
410 & ewperiodic(ng), nsperiodic(ng), &
411 & ad_u(:,:,:,linp), &
412 & ad_v(:,:,:,linp))
413# endif
414 DO k=1,n(ng)
415 DO j=jstrt,jendt
416 DO i=istrp,iendt
417 ad_u(i,j,k,linp)=ad_u(i,j,k,linp)*u_std(i,j,k,rec)
418 END DO
419 END DO
420 DO j=jstrp,jendt
421 DO i=istrt,iendt
422 ad_v(i,j,k,linp)=ad_v(i,j,k,linp)*v_std(i,j,k,rec)
423 END DO
424 END DO
425 END DO
426!
427! Tracers.
428!
429# ifdef DISTRIBUTE
430 CALL ad_mp_exchange4d (ng, tile, iadm, 1, &
431 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
432 & nghostpoints, &
433 & ewperiodic(ng), nsperiodic(ng), &
434 & ad_t(:,:,:,linp,:))
435# endif
436 DO it=1,nt(ng)
437 DO k=1,n(ng)
438 DO j=jstrt,jendt
439 DO i=istrt,iendt
440 ad_t(i,j,k,linp,it)=ad_t(i,j,k,linp,it)* &
441 & t_std(i,j,k,rec,it)
442 END DO
443 END DO
444 END DO
445 END DO
446# endif
447
448# ifdef ADJUST_BOUNDARY
449!
450!-----------------------------------------------------------------------
451! Boundary conditions error covariance: Multiply adjoint state by its
452! corresponding standard deviation.
453!-----------------------------------------------------------------------
454!
455! Free-surface open boundaries.
456!
457 IF (.not.lweak.and.(any(lobc(:,isfsur,ng)))) THEN
458 DO ir=1,nbrec(ng)
459# ifdef DISTRIBUTE
460 DO ib=1,4
461 IF (lobc(ib,isfsur,ng)) THEN
462 CALL ad_mp_exchange2d_bry (ng, tile, iadm, 1, ib, &
463 & lbij, ubij, &
464 & nghostpoints, &
465 & ewperiodic(ng), nsperiodic(ng),&
466 & ad_zeta_obc(:,ib,ir,linp))
467 END IF
468 END DO
469# endif
470 IF ((lobc(iwest,isfsur,ng)).and. &
471 & domain(ng)%Western_Edge(tile)) THEN
472 ib=iwest
473 DO j=jstr,jend
474 ad_zeta_obc(j,ib,ir,linp)=ad_zeta_obc(j,ib,ir,linp)* &
475 & zeta_obc_std(j,ib)
476 END DO
477 END IF
478 IF ((lobc(ieast,isfsur,ng)).and. &
479 & domain(ng)%Eastern_Edge(tile)) THEN
480 ib=ieast
481 DO j=jstr,jend
482 ad_zeta_obc(j,ib,ir,linp)=ad_zeta_obc(j,ib,ir,linp)* &
483 & zeta_obc_std(j,ib)
484 END DO
485 END IF
486 IF ((lobc(isouth,isfsur,ng)).and. &
487 & domain(ng)%Southern_Edge(tile)) THEN
488 ib=isouth
489 DO i=istr,iend
490 ad_zeta_obc(i,ib,ir,linp)=ad_zeta_obc(i,ib,ir,linp)* &
491 & zeta_obc_std(i,ib)
492 END DO
493 END IF
494 IF ((lobc(inorth,isfsur,ng)).and. &
495 & domain(ng)%Northern_Edge(tile)) THEN
496 ib=inorth
497 DO i=istr,iend
498 ad_zeta_obc(i,ib,ir,linp)=ad_zeta_obc(i,ib,ir,linp)* &
499 & zeta_obc_std(i,ib)
500 END DO
501 END IF
502 END DO
503 END IF
504!
505! 2D U-momentum open boundaries.
506!
507 IF (.not.lweak.and.(any(lobc(:,isubar,ng)))) THEN
508 DO ir=1,nbrec(ng)
509# ifdef DISTRIBUTE
510 DO ib=1,4
511 IF (lobc(ib,isubar,ng)) THEN
512 CALL ad_mp_exchange2d_bry (ng, tile, iadm, 1, ib, &
513 & lbij, ubij, &
514 & nghostpoints, &
515 & ewperiodic(ng), nsperiodic(ng),&
516 & ad_ubar_obc(:,ib,ir,linp))
517 END IF
518 END DO
519# endif
520 IF ((lobc(iwest,isubar,ng)).and. &
521 & domain(ng)%Western_Edge(tile)) THEN
522 ib=iwest
523 DO j=jstr,jend
524 ad_ubar_obc(j,ib,ir,linp)=ad_ubar_obc(j,ib,ir,linp)* &
525 & ubar_obc_std(j,ib)
526 END DO
527 END IF
528 IF ((lobc(ieast,isubar,ng)).and. &
529 & domain(ng)%Eastern_Edge(tile)) THEN
530 ib=ieast
531 DO j=jstr,jend
532 ad_ubar_obc(j,ib,ir,linp)=ad_ubar_obc(j,ib,ir,linp)* &
533 & ubar_obc_std(j,ib)
534 END DO
535 END IF
536 IF ((lobc(isouth,isubar,ng)).and. &
537 & domain(ng)%Southern_Edge(tile)) THEN
538 ib=isouth
539 DO i=istru,iend
540 ad_ubar_obc(i,ib,ir,linp)=ad_ubar_obc(i,ib,ir,linp)* &
541 & ubar_obc_std(i,ib)
542 END DO
543 END IF
544 IF ((lobc(inorth,isubar,ng)).and. &
545 & domain(ng)%Northern_Edge(tile)) THEN
546 ib=inorth
547 DO i=istru,iend
548 ad_ubar_obc(i,ib,ir,linp)=ad_ubar_obc(i,ib,ir,linp)* &
549 & ubar_obc_std(i,ib)
550 END DO
551 END IF
552 END DO
553 END IF
554!
555! 2D V-momentum open boundaries.
556!
557 IF (.not.lweak.and.(any(lobc(:,isvbar,ng)))) THEN
558 DO ir=1,nbrec(ng)
559# ifdef DISTRIBUTE
560 DO ib=1,4
561 IF (lobc(ib,isvbar,ng)) THEN
562 CALL ad_mp_exchange2d_bry (ng, tile, iadm, 1, ib, &
563 & lbij, ubij, &
564 & nghostpoints, &
565 & ewperiodic(ng), nsperiodic(ng),&
566 & ad_vbar_obc(:,ib,ir,linp))
567 END IF
568 END DO
569# endif
570 IF ((lobc(iwest,isvbar,ng)).and. &
571 & domain(ng)%Western_Edge(tile)) THEN
572 ib=iwest
573 DO j=jstrv,jend
574 ad_vbar_obc(j,ib,ir,linp)=ad_vbar_obc(j,ib,ir,linp)* &
575 & vbar_obc_std(j,ib)
576 END DO
577 END IF
578 IF ((lobc(ieast,isvbar,ng)).and. &
579 & domain(ng)%Eastern_Edge(tile)) THEN
580 ib=ieast
581 DO j=jstrv,jend
582 ad_vbar_obc(j,ib,ir,linp)=ad_vbar_obc(j,ib,ir,linp)* &
583 & vbar_obc_std(j,ib)
584 END DO
585 END IF
586 IF ((lobc(isouth,isvbar,ng)).and. &
587 & domain(ng)%Southern_Edge(tile)) THEN
588 ib=isouth
589 DO i=istr,iend
590 ad_vbar_obc(i,ib,ir,linp)=ad_vbar_obc(i,ib,ir,linp)* &
591 & vbar_obc_std(i,ib)
592 END DO
593 END IF
594 IF ((lobc(inorth,isvbar,ng)).and. &
595 & domain(ng)%Northern_Edge(tile)) THEN
596 ib=inorth
597 DO i=istr,iend
598 ad_vbar_obc(i,ib,ir,linp)=ad_vbar_obc(i,ib,ir,linp)* &
599 & vbar_obc_std(i,ib)
600 END DO
601 END IF
602 END DO
603 END IF
604
605# ifdef SOLVE3D
606!
607! 3D U-momentum open boundaries.
608!
609 IF (.not.lweak.and.(any(lobc(:,isuvel,ng)))) THEN
610 DO ir=1,nbrec(ng)
611# ifdef DISTRIBUTE
612 DO ib=1,4
613 IF (lobc(ib,isuvel,ng)) THEN
614 CALL ad_mp_exchange3d_bry (ng, tile, iadm, 1, ib, &
615 & lbij, ubij, 1, n(ng), &
616 & nghostpoints, &
617 & ewperiodic(ng), nsperiodic(ng),&
618 & ad_u_obc(:,:,ib,ir,linp))
619 END IF
620 END DO
621# endif
622 IF ((lobc(iwest,isuvel,ng)).and. &
623 & domain(ng)%Western_Edge(tile)) THEN
624 ib=iwest
625 DO k=1,n(ng)
626 DO j=jstr,jend
627 ad_u_obc(j,k,ib,ir,linp)=ad_u_obc(j,k,ib,ir,linp)* &
628 & u_obc_std(j,k,ib)
629 END DO
630 END DO
631 END IF
632 IF ((lobc(ieast,isuvel,ng)).and. &
633 & domain(ng)%Eastern_Edge(tile)) THEN
634 ib=ieast
635 DO k=1,n(ng)
636 DO j=jstr,jend
637 ad_u_obc(j,k,ib,ir,linp)=ad_u_obc(j,k,ib,ir,linp)* &
638 & u_obc_std(j,k,ib)
639 END DO
640 END DO
641 END IF
642 IF ((lobc(isouth,isuvel,ng)).and. &
643 & domain(ng)%Southern_Edge(tile)) THEN
644 ib=isouth
645 DO k=1,n(ng)
646 DO i=istru,iend
647 ad_u_obc(i,k,ib,ir,linp)=ad_u_obc(i,k,ib,ir,linp)* &
648 & u_obc_std(i,k,ib)
649 END DO
650 END DO
651 END IF
652 IF ((lobc(inorth,isuvel,ng)).and. &
653 & domain(ng)%Northern_Edge(tile)) THEN
654 ib=inorth
655 DO k=1,n(ng)
656 DO i=istru,iend
657 ad_u_obc(i,k,ib,ir,linp)=ad_u_obc(i,k,ib,ir,linp)* &
658 & u_obc_std(i,k,ib)
659 END DO
660 END DO
661 END IF
662 END DO
663 END IF
664!
665! 3D V-momentum open boundaries.
666!
667 IF (.not.lweak.and.(any(lobc(:,isvvel,ng)))) THEN
668 DO ir=1,nbrec(ng)
669# ifdef DISTRIBUTE
670 DO ib=1,4
671 IF (lobc(ib,isvvel,ng)) THEN
672 CALL ad_mp_exchange3d_bry (ng, tile, iadm, 1, ib, &
673 & lbij, ubij, 1, n(ng), &
674 & nghostpoints, &
675 & ewperiodic(ng), nsperiodic(ng),&
676 & ad_v_obc(:,:,ib,ir,linp))
677 END IF
678 END DO
679# endif
680 IF ((lobc(iwest,isvvel,ng)).and. &
681 & domain(ng)%Western_Edge(tile)) THEN
682 ib=iwest
683 DO k=1,n(ng)
684 DO j=jstrv,jend
685 ad_v_obc(j,k,ib,ir,linp)=ad_v_obc(j,k,ib,ir,linp)* &
686 & v_obc_std(j,k,ib)
687 END DO
688 END DO
689 END IF
690 IF ((lobc(ieast,isvvel,ng)).and. &
691 & domain(ng)%Eastern_Edge(tile)) THEN
692 ib=ieast
693 DO k=1,n(ng)
694 DO j=jstrv,jend
695 ad_v_obc(j,k,ib,ir,linp)=ad_v_obc(j,k,ib,ir,linp)* &
696 & v_obc_std(j,k,ib)
697 END DO
698 END DO
699 END IF
700 IF ((lobc(isouth,isvvel,ng)).and. &
701 & domain(ng)%Southern_Edge(tile)) THEN
702 ib=isouth
703 DO k=1,n(ng)
704 DO i=istr,iend
705 ad_v_obc(i,k,ib,ir,linp)=ad_v_obc(i,k,ib,ir,linp)* &
706 & v_obc_std(i,k,ib)
707 END DO
708 END DO
709 END IF
710 IF ((lobc(inorth,isvvel,ng)).and. &
711 & domain(ng)%Northern_Edge(tile)) THEN
712 ib=inorth
713 DO k=1,n(ng)
714 DO i=istr,iend
715 ad_v_obc(i,k,ib,ir,linp)=ad_v_obc(i,k,ib,ir,linp)* &
716 & v_obc_std(i,k,ib)
717 END DO
718 END DO
719 END IF
720 END DO
721 END IF
722!
723! Tracers open boundaries.
724!
725 DO it=1,nt(ng)
726 IF (.not.lweak.and.(any(lobc(:,istvar(it),ng)))) THEN
727 DO ir=1,nbrec(ng)
728# ifdef DISTRIBUTE
729 DO ib=1,4
730 IF (lobc(ib,istvar(it),ng)) THEN
731 CALL ad_mp_exchange3d_bry (ng, tile, iadm, 1, ib, &
732 & lbij, ubij, 1, n(ng), &
733 & nghostpoints, &
734 & ewperiodic(ng), &
735 & nsperiodic(ng), &
736 & ad_t_obc(:,:,ib,ir,linp,it))
737 END IF
738 END DO
739# endif
740 IF ((lobc(iwest,istvar(it),ng)).and. &
741 & domain(ng)%Western_Edge(tile)) THEN
742 ib=iwest
743 DO k=1,n(ng)
744 DO j=jstr,jend
745 ad_t_obc(j,k,ib,ir,linp,it)= &
746 & ad_t_obc(j,k,ib,ir,linp,it)* &
747 & t_obc_std(j,k,ib,it)
748 END DO
749 END DO
750 END IF
751 IF ((lobc(ieast,istvar(it),ng)).and. &
752 & domain(ng)%Eastern_Edge(tile)) THEN
753 ib=ieast
754 DO k=1,n(ng)
755 DO j=jstr,jend
756 ad_t_obc(j,k,ib,ir,linp,it)= &
757 & ad_t_obc(j,k,ib,ir,linp,it)* &
758 & t_obc_std(j,k,ib,it)
759 END DO
760 END DO
761 END IF
762 IF ((lobc(isouth,istvar(it),ng)).and. &
763 & domain(ng)%Southern_Edge(tile)) THEN
764 ib=isouth
765 DO k=1,n(ng)
766 DO i=istr,iend
767 ad_t_obc(i,k,ib,ir,linp,it)= &
768 & ad_t_obc(i,k,ib,ir,linp,it)* &
769 & t_obc_std(i,k,ib,it)
770 END DO
771 END DO
772 END IF
773 IF ((lobc(inorth,istvar(it),ng)).and. &
774 & domain(ng)%Northern_Edge(tile)) THEN
775 ib=inorth
776 DO k=1,n(ng)
777 DO i=istr,iend
778 ad_t_obc(i,k,ib,ir,linp,it)= &
779 & ad_t_obc(i,k,ib,ir,linp,it)* &
780 & t_obc_std(i,k,ib,it)
781 END DO
782 END DO
783 END IF
784 END DO
785 END IF
786 END DO
787# endif
788# endif
789
790# if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
791!
792!-----------------------------------------------------------------------
793! Surface forcing error covariance: Multiply adjoint state by its
794! corresponding standard deviation.
795!-----------------------------------------------------------------------
796
797# ifdef ADJUST_WSTRESS
798!
799! Adjoint surface momentum stress.
800!
801 IF (.not.lweak) THEN
802# ifdef DISTRIBUTE
803 CALL ad_mp_exchange3d (ng, tile, iadm, 2, &
804 & lbi, ubi, lbj, ubj, 1, nfrec(ng), &
805 & nghostpoints, &
806 & ewperiodic(ng), nsperiodic(ng), &
807 & ad_ustr(:,:,:,linp), &
808 & ad_vstr(:,:,:,linp))
809# endif
810 DO ir=1,nfrec(ng)
811 DO j=jstrt,jendt
812 DO i=istrp,iendt
813 ad_ustr(i,j,ir,linp)=ad_ustr(i,j,ir,linp)*sustr_std(i,j)
814 END DO
815 END DO
816 DO j=jstrp,jendt
817 DO i=istrt,iendt
818 ad_vstr(i,j,ir,linp)=ad_vstr(i,j,ir,linp)*svstr_std(i,j)
819 END DO
820 END DO
821 END DO
822 END IF
823# endif
824# if defined ADJUST_STFLUX && defined SOLVE3D
825!
826! Surface tracer flux.
827!
828 IF (.not.lweak) THEN
829 DO it=1,nt(ng)
830 IF (lstflux(it,ng)) THEN
831# ifdef DISTRIBUTE
832 CALL ad_mp_exchange3d (ng, tile, iadm, 1, &
833 & lbi, ubi, lbj, ubj, 1,nfrec(ng), &
834 & nghostpoints, &
835 & ewperiodic(ng), nsperiodic(ng), &
836 & ad_tflux(:,:,:,linp,it))
837# endif
838 DO ir=1,nfrec(ng)
839 DO j=jstrt,jendt
840 DO i=istrt,iendt
841 ad_tflux(i,j,ir,linp,it)=ad_tflux(i,j,ir,linp,it)* &
842 & stflx_std(i,j,it)
843 END DO
844 END DO
845 END DO
846 END IF
847 END DO
848 END IF
849# endif
850# endif
851
852 RETURN
integer isvvel
integer isvbar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer isubar
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer nghostpoints
Definition mod_param.F:710
integer, parameter iadm
Definition mod_param.F:665
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable nt
Definition mod_param.F:489
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nbrec
subroutine ad_mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c, ad_d)
subroutine ad_mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c)
subroutine ad_mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c, ad_d)
subroutine ad_mp_exchange3d_bry(ng, tile, model, nvar, boundary, lbij, ubij, lbk, ubk, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c, ad_d)
subroutine ad_mp_exchange2d_bry(ng, tile, model, nvar, boundary, lbij, ubij, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c, ad_d)

References mp_exchange_mod::ad_mp_exchange2d(), mp_exchange_mod::ad_mp_exchange2d_bry(), mp_exchange_mod::ad_mp_exchange3d(), mp_exchange_mod::ad_mp_exchange3d_bry(), mp_exchange_mod::ad_mp_exchange4d(), mod_param::domain, mod_scalars::ewperiodic, mod_param::iadm, mod_scalars::ieast, mod_scalars::inorth, mod_ncparam::isfsur, mod_scalars::isouth, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::iwest, mod_scalars::lobc, mod_scalars::lstflux, mod_param::nghostpoints, and mod_scalars::nsperiodic.

Referenced by ad_variability().

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