ROMS
Loading...
Searching...
No Matches
background_std.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#ifdef STD_MODEL
5!
6!!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! Standard deviation formulation used to model the 4D-Var Background !
14! Error Covariance matrix, B. It follows the method Mogensen et al. !
15! (2012) that assumes that the background errors are proportional to !
16! vertical derivatives of the state. That is, the difference between !
17! the background state value S_b and true value S_t is due to a !
18! vertical displacement of the profile: !
19! !
20! S_t(z) ~ S_b(z+deltaz) + [d(S_b)/d(z)]*deltaz !
21! !
22! where deltaz is the displacement. !
23! !
24! I assumes that the background and true profiles have similar shape !
25! and the true value of S lies somewhere in the background water !
26! column. The error is then [d(S_b)/d(z)]*deltaz. !
27! !
28! References: !
29! !
30! Kara A., P. Rochford, and E. Hulburt, 2000: An optimal definition !
31! for ocean mixed layer depth, J. Geophys. Res. vol 105, NoC7, !
32! pp 16, 803-16, 821. !
33! !
34! Mogensen, K., M.A. Balmaseda, and A.T. Weaver, 2012: The NEMOVAR !
35! ocean data assimilation system as implemented in the ECMWF ocean !
36! analysis for system 4. ECMWF Tech. Memorandum 668, 59. !
37! !
38! Moore, A., J. Zavala-Garay, H.G. Arango, C.A. Edwards, J. Anderson, !
39! and T. Hoar, 2020: Regional and basin scale applications of !
40! ensemble adjustement Kalman filter and 4D-Var ocean data !
41! assimilation systems, Progress in Oceanography, 189, 102450, !
42! https://doi.org/10.1016/j.pocean.2020.102450. !
43! !
44!=======================================================================
45!
46 USE mod_kinds
47 USE mod_param
48 USE mod_parallel
49# ifdef SOLVE3D
50 USE mod_coupling
51# endif
52 USE mod_grid
53 USE mod_iounits
54 USE mod_ncparam
55 USE mod_ocean
56 USE mod_scalars
57!
60# ifdef SOLVE3D
63# endif
64# ifdef DISTRIBUTE
66# endif
67# ifdef SOLVE3D
68 USE set_depth_mod, ONLY : set_depth
69# endif
70!
71 implicit none
72!
73 PRIVATE
74 PUBLIC :: background_std
75!
76 CONTAINS
77!
78!***********************************************************************
79 SUBROUTINE background_std (ng, tile, Lbck, Lstd)
80!***********************************************************************
81!
82! Imported variable declarations.
83!
84 integer, intent(in) :: ng, tile, Lbck, Lstd
85!
86! Local variable declarations.
87!
88# include "tile.h"
89!
90# ifdef SOLVE3D
91!
92! Compute background state thickness, depth arrays, thermal expansion,
93! and saline contraction coefficients.
94!
95 coupling(ng) % Zt_avg1 = 0.0_r8
96
97 CALL set_depth (ng, tile, inlm)
98!
99# endif
100 CALL background_std_tile (ng, tile, &
101 & lbi, ubi, lbj, ubj, &
102 & imins, imaxs, jmins, jmaxs, &
103 & lbck, lstd, &
104# ifdef SOLVE3D
105 & grid(ng) % Hz, &
106 & grid(ng) % z_r, &
107# endif
108# ifdef MASKING
109 & grid(ng) % rmask, &
110 & grid(ng) % umask, &
111 & grid(ng) % vmask, &
112# endif
113# ifdef SOLVE3D
114 & ocean(ng) % t, &
115 & ocean(ng) % u, &
116 & ocean(ng) % v, &
117 & ocean(ng) % e_t, &
118 & ocean(ng) % e_u, &
119 & ocean(ng) % e_v, &
120# endif
121 & ocean(ng) % e_ubar, &
122 & ocean(ng) % e_vbar, &
123 & ocean(ng) % e_zeta)
124!
125 IF (master) WRITE (stdout,10)
126 10 FORMAT (/,2x,'BACKGROUND_STD - computing standard deviation', &
127 & ' from prior.',/)
128!
129 RETURN
130 END SUBROUTINE background_std
131!
132!***********************************************************************
133 SUBROUTINE background_std_tile (ng, tile, &
134 & LBi, UBi, LBj, UBj, &
135 & IminS, ImaxS, JminS, JmaxS, &
136 & Lbck, Lstd, &
137# ifdef SOLVE3D
138 & Hz, z_r, &
139# endif
140# ifdef MASKING
141 & rmask, umask, vmask, &
142# endif
143# ifdef SOLVE3D
144 & t, u, v, &
145 & e_t, e_u, e_v, &
146# endif
147 & e_ubar, e_vbar, e_zeta)
148!***********************************************************************
149!
150! Imported variable declarations.
151!
152 integer, intent(in) :: ng, tile
153 integer, intent(in) :: LBi, UBi, LBj, UBj
154 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
155 integer, intent(in) :: Lbck, Lstd
156!
157# ifdef ASSUMED_SHAPE
158# ifdef SOLVE3D
159 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
160 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
161# endif
162# ifdef MASKING
163 real(r8), intent(in) :: rmask(LBi:,LBj:)
164 real(r8), intent(in) :: umask(LBi:,LBj:)
165 real(r8), intent(in) :: vmask(LBi:,LBj:)
166# endif
167# ifdef SOLVE3D
168 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
169 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
170 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
171 real(r8), intent(inout) :: e_t(LBi:,LBj:,:,:,:)
172 real(r8), intent(inout) :: e_u(LBi:,LBj:,:,:)
173 real(r8), intent(inout) :: e_v(LBi:,LBj:,:,:)
174# endif
175 real(r8), intent(inout) :: e_ubar(LBi:,LBj:,:)
176 real(r8), intent(inout) :: e_vbar(LBi:,LBj:,:)
177 real(r8), intent(inout) :: e_zeta(LBi:,LBj:,:)
178
179# else
180
181# ifdef SOLVE3D
182 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
183 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
184# endif
185# ifdef MASKING
186 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
189# endif
190# ifdef SOLVE3D
191 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
192 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,2,N(ng))
193 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,2,N(ng))
194 real(r8), intent(inout) :: e_t(LBi:UBi,LBj:UBj,N(ng),NSA,NT(ng))
195 real(r8), intent(inout) :: e_u(LBi:UBi,LBj:UBj,N(ng),NSA)
196 real(r8), intent(inout) :: e_v(LBi:UBi,LBj:UBj,N(ng),NSA)
197# endif
198 real(r8), intent(inout) :: e_ubar(LBi:UBi,LBj:UBj,NSA)
199 real(r8), intent(inout) :: e_vbar(LBi:UBi,LBj:UBj,NSA)
200 real(r8), intent(inout) :: e_zeta(LBi:UBi,LBj:UBj,NSA)
201
202# endif
203!
204! Local variable declarations.
205!
206 logical :: base_reached = .false.
207 logical :: ml_reached = .false.
208!
209 integer, parameter :: Norder = 2 ! Shapiro filter order
210 integer :: i, j, k, kref, khref, order
211!
212 real(r8) :: Temp_ref, T_dep, T_high, T_low, T_thvalue
213 real(r8) :: sigmabS, sigmabT, sigmabU, SigmabV
214 real(r8) :: cff, cff1, cff2, fac
215 real(r8) :: href
216!
217! Shapiro filter coefficients.
218!
219 real(r8), dimension(20) :: filter_coef = &
220 & (/ 2.500000E-1_r8, 6.250000E-2_r8, 1.562500E-2_r8, &
221 & 3.906250E-3_r8, 9.765625E-4_r8, 2.44140625E-4_r8, &
222 & 6.103515625E-5_r8, 1.5258789063E-5_r8, 3.814697E-6_r8, &
223 & 9.536743E-7_r8, 2.384186E-7_r8, 5.960464E-8_r8, &
224 & 1.490116E-8_r8, 3.725290E-9_r8, 9.313226E-10_r8, &
225 & 2.328306E-10_r8, 5.820766E-11_r8, 1.455192E-11_r8, &
226 & 3.637979E-12_r8, 9.094947E-13_r8 /)
227
228 real(r8), dimension(N(ng)) :: dSdT, dSdT_filter
229
230 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
231
232 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: dTdz
233 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: dUdz
234 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: dVdz
235 real(r8), dimension(N(ng)) :: dTdz_filter
236 real(r8), dimension(N(ng)) :: dUdz_filter, dVdz_filter
237# ifdef SALINITY
238 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: dSdz
239 real(r8), dimension(N(ng)) :: dSdz_filter
240# endif
241 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: mld
242
243# include "set_bounds.h"
244!
245!-----------------------------------------------------------------------
246! Compute temperature (dTdz) and salinity (dSdz) shears.
247!-----------------------------------------------------------------------
248!
249 DO j=jstrr,jendr
250 DO i=istrr,iendr
251 fc(i,0)=0.0_r8
252 dtdz(i,j,0)=0.0_r8
253# ifdef SALINITY
254 dsdz(i,j,0)=0.0_r8
255# endif
256 END DO
257 DO k=1,n(ng)-1
258 DO i=istrr,iendr
259 cff=1.0_r8/(2.0_r8*hz(i,j,k+1)+ &
260 & hz(i,j,k)*(2.0_r8-fc(i,k-1)))
261 fc(i,k)=cff*hz(i,j,k+1)
262 dtdz(i,j,k)=cff*(6.0_r8*(t(i,j,k+1,lbck,itemp)- &
263 & t(i,j,k ,lbck,itemp))- &
264 & hz(i,j,k)*dtdz(i,j,k-1))
265# ifdef SALINITY
266 dsdz(i,j,k)=cff*(6.0_r8*(t(i,j,k+1,lbck,isalt)- &
267 & t(i,j,k ,lbck,isalt))- &
268 & hz(i,j,k)*dsdz(i,j,k-1))
269# endif
270 END DO
271 END DO
272 DO i=istrr,iendr
273 dtdz(i,j,n(ng))=0.0_r8
274# ifdef SALINITY
275 dsdz(i,j,n(ng))=0.0_r8
276# endif
277 END DO
278 DO k=n(ng)-1,1,-1
279 DO i=istrr,iendr
280 dtdz(i,j,k)=dtdz(i,j,k)-fc(i,k)*dtdz(i,j,k+1)
281# ifdef SALINITY
282 dsdz(i,j,k)=dsdz(i,j,k)-fc(i,k)*dsdz(i,j,k+1)
283# endif
284 END DO
285 END DO
286!
287! Shapiro filter of the vertical derivatives.
288!
289 DO i=istrr,iendr
290 DO order=1,norder/2
291 IF (order.ne.norder/2) THEN
292# ifdef SALINITY
293 dsdz_filter(1)=2.0_r8*(dsdz(i,j,1)-dsdz(i,j,2))
294 dsdz_filter(n(ng))=2.0_r8*(dsdz(i,j,n(ng) )- &
295 & dsdz(i,j,n(ng)-1))
296# endif
297 dtdz_filter(1)=2.0_r8*(dtdz(i,j,1)-dtdz(i,j,2))
298 dtdz_filter(n(ng))=2.0_r8*(dtdz(i,j,n(ng) )- &
299 & dtdz(i,j,n(ng)-1))
300 ELSE
301# ifdef SALINITY
302 dsdz_filter(1)=0.0_r8
303 dsdz_filter(n(ng))=0.0_r8
304# endif
305 dtdz_filter(1)=0.0_r8
306 dtdz_filter(n(ng))=0.0_r8
307 END IF
308 DO k=2,n(ng)-1
309# ifdef SALINITY
310 dsdz_filter(k)=2.0_r8*dsdz(i,j,k)- &
311 & dsdz(i,j,k-1)-dsdz(i,j,k+1)
312# endif
313 dtdz_filter(k)=2.0_r8*dtdz(i,j,k)- &
314 & dtdz(i,j,k-1)-dtdz(i,j,k+1)
315 END DO
316 DO k=1,n(ng)
317# ifdef SALINITY
318 dsdz(i,j,k)=dsdz(i,j,k)- &
319 & filter_coef(norder/2)*dsdz_filter(k)
320# endif
321 dtdz(i,j,k)=dtdz(i,j,k)- &
322 & filter_coef(norder/2)*dtdz_filter(k)
323 END DO
324 END DO
325 END DO
326 END DO
327!
328! Compute velocity shears dUdz and dVdz.
329!
330 DO k=0,n(ng)
331 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
332 DO j=jstr-1,jend+1
333 DO i=istru-1,iend+1
334 dudz(i,j,k)=0.0_r8
335 END DO
336 END DO
337 DO j=jstrv-1,jend+1
338 DO i=istr-1,iend+1
339 dvdz(i,j,k)=0.0_r8
340 END DO
341 END DO
342 ELSE
343 DO j=jstr-1,jend+1
344 DO i=istru-1,iend+1
345 cff=1.0_r8/(0.5_r8*(z_r(i-1,j,k+1)-z_r(i-1,j,k)+ &
346 & z_r(i ,j,k+1)-z_r(i ,j,k)))
347 dudz(i,j,k)=cff*(u(i,j,k+1,lbck)- &
348 & u(i,j,k ,lbck))
349 END DO
350 END DO
351 DO j=jstrv-1,jend+1
352 DO i=istr-1,iend+1
353 cff=1.0_r8/(0.5_r8*(z_r(i,j-1,k+1)-z_r(i,j-1,k)+ &
354 & z_r(i,j ,k+1)-z_r(i,j ,k)))
355 dvdz(i,j,k)=cff*(v(i,j,k+1,lbck)- &
356 & v(i,j,k ,lbck))
357 END DO
358 END DO
359 END IF
360 END DO
361!
362! Shapiro filter of dUdz.
363!
364 DO j=jstr-1,jend+1
365 DO i=istru-1,iend+1
366 DO order=1,norder/2
367 IF (order.ne.norder/2) THEN
368 dudz_filter(1)=2.0_r8*(dudz(i,j,1)-dudz(i,j,2))
369 dudz_filter(n(ng))=2.0_r8*(dudz(i,j,n(ng) )- &
370 & dudz(i,j,n(ng)-1))
371 ELSE
372 dudz_filter(1)=0.0_r8
373 dudz_filter(n(ng))=0.0_r8
374 END IF
375 DO k=2,n(ng)-1
376 dudz_filter(k)=2.0_r8*dudz(i,j,k)- &
377 & dudz(i,j,k-1)-dudz(i,j,k+1)
378 END DO
379 DO k=1,n(ng)
380 dudz(i,j,k)=dudz(i,j,k)- &
381 & filter_coef(norder/2)*dudz_filter(k)
382 END DO
383 END DO
384 END DO
385 END DO
386!
387! Shapiro filter of dVdz.
388!
389 DO j=jstrv-1,jend+1
390 DO i=istr-1,iend+1
391 DO order=1,norder/2
392 IF (order.ne.norder/2) THEN
393 dvdz_filter(1)=2.0_r8*(dvdz(i,j,1)-dvdz(i,j,2))
394 dvdz_filter(n(ng))=2.0_r8*(dvdz(i,j,n(ng) )- &
395 & dvdz(i,j,n(ng)-1))
396 ELSE
397 dvdz_filter(1)=0.0_r8
398 dvdz_filter(n(ng))=0.0_r8
399 END IF
400 DO k=2,n(ng)-1
401 dvdz_filter(k)=2.0_r8*dvdz(i,j,k)- &
402 & dvdz(i,j,k-1)-dvdz(i,j,k+1)
403 END DO
404 DO k=1,n(ng)
405 dvdz(i,j,k)=dvdz(i,j,k)- &
406 & filter_coef(norder/2)*dvdz_filter(k)
407 END DO
408 END DO
409 END DO
410 END DO
411
412# ifdef COMPUTE_MLD
413!
414!-----------------------------------------------------------------------
415! Compute mixed layer depth (mld) according to criterion
416! on temperature from the work of Kara et al. (2000).
417!-----------------------------------------------------------------------
418!
419! Criteria on temperature : 0.8
420!
421! Start searching for the mixed layer depth at href.
422! CAUTION: href MUST be greater than the minimum depth in the model
423! domain (usually over land points).
424!
425 t_thvalue=0.8_r8
426! href = -10.0_r8
427 href = -9.0_r8
428!
429 DO j=jstrr,jendr
430 DO i=istrr,iendr
431!
432! Compute temperature at reference depth, href.
433!
434 DO k=n(ng)-1,1,-1
435 IF ((z_r(i,j,k+1).gt.href).and.(z_r(i,j,k).le.href)) THEN
436 khref = k
437 END IF
438 END DO
439 base_reached =.false.
440 kref= khref
441! Search for an uniform temperature region
442 DO k=khref-1,1,-1
443 t_dep = abs(t(i,j,k+1,lbck,itemp)-t(i,j,k,lbck,itemp))/ &
444 & abs(z_r(i,j,k+1)-z_r(i,j,k))
445 IF (.not.base_reached.and.(t_dep.lt.0.01_r8*t_thvalue)) THEN
446 kref = k
447 ELSE
448 base_reached = .true.
449 END IF
450 END DO
451!
452 IF (.not.base_reached) THEN
453 kref = khref
454 END IF
455!
456! Linearly interpolate temperature to href.
457!
458 IF (kref.eq.khref) THEN
459 DO k=1,n(ng)-1
460 IF ((href.ge.z_r(i,j,k)).and.(href.lt.z_r(i,j,k+1))) THEN
461 fac=(href-z_r(i,j,k))/(z_r(i,j,k+1)-z_r(i,j,k))
462 temp_ref=fac*t(i,j,k+1,lbck,itemp)+ &
463 (1.0_r8-fac)*t(i,j,k,lbck,itemp)
464 mld(i,j)=href
465 END IF
466 END DO
467 ELSE
468 temp_ref=t(i,j,kref,lbck,itemp)
469 mld(i,j)=z_r(i,j,kref)
470 END IF
471!
472 ml_reached=.false.
473!
474 DO k=kref-1,1,-1
475 IF (.not.ml_reached.and. &
476 & (abs(t(i,j,k,lbck,itemp)-temp_ref).gt.t_thvalue)) THEN
477 t_high=abs(t(i,j,k+1,lbck,itemp)-temp_ref)
478 t_low =abs(t(i,j,k ,lbck,itemp)-temp_ref)
479 mld(i,j)=((t_thvalue-t_high)*z_r(i,j,k)+ &
480 & (t_low-t_thvalue)*z_r(i,j,k+1))/(t_low-t_high)
481 ml_reached =.true.
482 END IF
483 END DO
484!
485 IF (.not.ml_reached) THEN
486 mld(i,j)=z_r(i,j,1)
487 END IF
488 END DO
489 END DO
490# else
491 DO j=jstrr,jendr
492 DO i=istrr,iendr
493 mld(i,j)=mld_uniform(ng)
494 END DO
495 END DO
496# endif
497!
498!-----------------------------------------------------------------------
499! Compute the background (prior) state vector standard deviations.
500!-----------------------------------------------------------------------
501!
502! Temperature and salinity.
503!
504 DO j=jstrr,jendr
505 DO i=istrr,iendr
506 DO k=1,n(ng)
507 cff1=0.5_r8*(dtdz(i,j,k-1)+dtdz(i,j,k))
508 sigmabt=min(abs(sigma_dz(istvar(itemp),ng)*cff1), &
509 & sigma_max(istvar(itemp),ng))
510 IF (z_r(i,j,k).ge.mld(i,j)) THEN
511 sigmabt=max(sigmabt, sigma_ml(istvar(itemp),ng))
512 ELSE
513 sigmabt=max(sigmabt, sigma_do(istvar(itemp),ng))
514 END IF
515 e_t(i,j,k,lstd,itemp)=sigmabt
516# ifdef MASKING
517 e_t(i,j,k,lstd,itemp)=e_t(i,j,k,lstd,itemp)*rmask(i,j)
518# endif
519# ifdef SALINITY
520 cff2=0.5_r8*(dsdz(i,j,k-1)+dsdz(i,j,k))
521 sigmabs=min(abs(sigma_dz(istvar(isalt),ng)*cff2), &
522 sigma_max(istvar(isalt),ng))
523 IF (z_r(i,j,k).ge.mld(i,j)) THEN
524 sigmabs=max(sigmabs, sigma_ml(istvar(isalt),ng))
525 ELSE
526 sigmabs=max(sigmabs, sigma_do(istvar(isalt),ng))
527 END IF
528 e_t(i,j,k,lstd,isalt)=sigmabs
529# ifdef MASKING
530 e_t(i,j,k,lstd,isalt)=e_t(i,j,k,lstd,isalt)*rmask(i,j)
531# endif
532# endif
533 END DO
534 END DO
535 END DO
536!
537! U-velocity component.
538!
539 DO j=jstr-1,jend+1
540 DO i=istru-1,iend+1
541 DO k=1,n(ng)
542 cff=0.5_r8*(dudz(i,j,k-1)+dudz(i,j,k))
543 sigmabu=min(abs(sigma_dz(isuvel,ng)*cff), &
544 & sigma_max(isuvel,ng))
545 IF (z_r(i,j,k).ge.mld(i,j)) THEN
546 sigmabu=max(sigmabu, sigma_ml(isuvel,ng))
547 ELSE
548 sigmabu=max(sigmabu, sigma_do(isuvel,ng))
549 END IF
550 e_u(i,j,k,lstd)=sigmabu
551# ifdef MASKING
552 e_u(i,j,k,lstd)=e_u(i,j,k,lstd)*umask(i,j)
553# endif
554 END DO
555 END DO
556 END DO
557!
558! V-velocity component.
559!
560 DO j=jstrv-1,jend+1
561 DO i=istr-1,iend+1
562 DO k=1,n(ng)
563 cff=0.5_r8*(dvdz(i,j,k-1)+dvdz(i,j,k))
564 sigmabv=min(abs(sigma_dz(isvvel,ng)*cff), &
565 & sigma_max(isvvel,ng))
566 IF (z_r(i,j,k).ge.mld(i,j)) THEN
567 sigmabv=max(sigmabv, sigma_ml(isvvel,ng))
568 ELSE
569 sigmabv=max(sigmabv, sigma_do(isvvel,ng))
570 END IF
571 e_v(i,j,k,lstd)=sigmabv
572# ifdef MASKING
573 e_v(i,j,k,lstd)=e_v(i,j,k,lstd)*vmask(i,j)
574# endif
575 END DO
576 END DO
577 END DO
578!
579! Vertically integrated velocity components. Not used, but needed
580! for I/O manipulations in the split schemes.
581!
582 DO j=jstrr,jendr
583 DO i=istr,iendr
584 e_ubar(i,j,lstd)=sigma_max(isubar,ng)
585# ifdef MASKING
586 e_ubar(i,j,lstd)=e_ubar(i,j,lstd)*umask(i,j)
587# endif
588 END DO
589 END DO
590
591 DO j=jstr,jendr
592 DO i=istrr,iendr
593 e_vbar(i,j,lstd)=sigma_max(isvbar,ng)
594# ifdef MASKING
595 e_vbar(i,j,lstd)=e_vbar(i,j,lstd)*vmask(i,j)
596# endif
597 END DO
598 END DO
599!
600! Free surface.
601!
602 DO j=jstrr,jendr
603 DO i=istrr,iendr
604 e_zeta(i,j,lstd)=sigma_max(isfsur,ng)
605# ifdef MASKING
606 e_zeta(i,j,lstd)=e_zeta(i,j,lstd)*rmask(i,j)
607# endif
608 END DO
609 END DO
610!
611! Exchange boundary information.
612
613 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
614 CALL exchange_r3d_tile (ng, tile, &
615 & lbi, ubi, lbj, ubj, 1, n(ng), &
616 & e_t(:,:,:,lstd,itemp))
617# ifdef SALINITY
618 CALL exchange_r3d_tile (ng, tile, &
619 & lbi, ubi, lbj, ubj, 1, n(ng), &
620 & e_t(:,:,:,lstd,isalt))
621# endif
622
623 CALL exchange_u3d_tile (ng, tile, &
624 & lbi, ubi, lbj, ubj, 1, n(ng), &
625 & e_u(:,:,:,lstd))
626 CALL exchange_v3d_tile (ng, tile, &
627 & lbi, ubi, lbj, ubj, 1, n(ng), &
628 & e_v(:,:,:,lstd))
629 CALL exchange_u2d_tile (ng, tile, &
630 & lbi, ubi, lbj, ubj, &
631 & e_ubar(:,:,lstd))
632 CALL exchange_v2d_tile (ng, tile, &
633 & lbi, ubi, lbj, ubj, &
634 & e_vbar(:,:,lstd))
635 CALL exchange_r2d_tile (ng, tile, &
636 & lbi, ubi, lbj, ubj, &
637 & e_zeta(:,:,lstd))
638 END IF
639
640# ifdef DISTRIBUTE
641!
642 CALL mp_exchange3d (ng, tile, inlm, 1, &
643 & lbi, ubi, lbj, ubj, 1, n(ng), &
644 & nghostpoints, &
645 & ewperiodic(ng), nsperiodic(ng), &
646 & e_t(:,:,:,lstd,itemp))
647# ifdef SALINITY
648 CALL mp_exchange3d (ng, tile, inlm, 1, &
649 & lbi, ubi, lbj, ubj, 1, n(ng), &
650 & nghostpoints, &
651 & ewperiodic(ng), nsperiodic(ng), &
652 & e_t(:,:,:,lstd,isalt))
653# endif
654 CALL mp_exchange3d (ng, tile, inlm, 2, &
655 & lbi, ubi, lbj, ubj, 1, n(ng), &
656 & nghostpoints, &
657 & ewperiodic(ng), nsperiodic(ng), &
658 & e_u(:,:,:,lstd), &
659 & e_v(:,:,:,lstd))
660 CALL mp_exchange2d (ng, tile, inlm, 3, &
661 & lbi, ubi, lbj, ubj, &
662 & nghostpoints, &
663 & ewperiodic(ng), nsperiodic(ng), &
664 & e_ubar(:,:,lstd), &
665 & e_vbar(:,:,lstd), &
666 & e_zeta(:,:,lstd))
667# endif
668!
669 RETURN
670 END SUBROUTINE background_std_tile
671#endif
672 END MODULE background_std_mod
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_v3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_coupling), dimension(:), allocatable coupling
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer stdout
integer isvvel
integer isvbar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer isubar
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer, parameter inlm
Definition mod_param.F:662
integer nghostpoints
Definition mod_param.F:710
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
integer isalt
integer itemp
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public set_depth(ng, tile, model)
Definition set_depth.F:34