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

Functions/Subroutines

subroutine, public ini_fields (ng, tile, model)
 
subroutine ini_fields_tile (ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kstp, krhs, knew, nstp, nnew, rmask, umask, vmask, rmask_wet, umask_wet, vmask_wet, ru, rv, rubar, rvbar, rzeta, bottom, bed, bed_frac, bed_mass, hz, t, u, v, ubar, vbar, zeta)
 
subroutine, public ini_zeta (ng, tile, model)
 
subroutine ini_zeta_tile (ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kstp, krhs, knew, rmask, rmask_wet, h, bed, bed_thick0, bed_thick, zeta)
 
subroutine, public set_zeta_timeavg (ng, tile, model)
 
subroutine set_zeta_timeavg_tile (ng, tile, model, lbi, ubi, lbj, ubj, kstp, zeta, zt_avg1)
 

Function/Subroutine Documentation

◆ ini_fields()

subroutine, public ini_fields_mod::ini_fields ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )

Definition at line 66 of file ini_fields.F.

67!***********************************************************************
68!
69 USE mod_stepping
70!
71! Imported variable declarations.
72!
73 integer, intent(in) :: ng, tile, model
74!
75! Local variable declarations.
76!
77 character (len=*), parameter :: MyFile = &
78 & __FILE__
79!
80# include "tile.h"
81!
82# ifdef PROFILE
83 CALL wclock_on (ng, model, 2, __line__, myfile)
84# endif
85 CALL ini_fields_tile (ng, tile, model, &
86 & lbi, ubi, lbj, ubj, &
87 & imins, imaxs, jmins, jmaxs, &
88 & kstp(ng), krhs(ng), knew(ng), &
89# ifdef SOLVE3D
90 & nstp(ng), nnew(ng), &
91# endif
92# ifdef MASKING
93 & grid(ng) % rmask, &
94 & grid(ng) % umask, &
95 & grid(ng) % vmask, &
96# endif
97# ifdef WET_DRY
98 & grid(ng) % rmask_wet, &
99 & grid(ng) % umask_wet, &
100 & grid(ng) % vmask_wet, &
101# endif
102# ifdef PERFECT_RESTART
103# ifdef SOLVE3D
104 & ocean(ng) % ru, &
105 & ocean(ng) % rv, &
106# endif
107 & ocean(ng) % rubar, &
108 & ocean(ng) % rvbar, &
109 & ocean(ng) % rzeta, &
110# endif
111# ifdef SOLVE3D
112# if defined SEDIMENT || defined BBL_MODEL
113 & sedbed(ng) % bottom, &
114# endif
115# if defined SEDIMENT
116 & sedbed(ng) % bed, &
117 & sedbed(ng) % bed_frac, &
118 & sedbed(ng) % bed_mass, &
119# endif
120 & grid(ng) % Hz, &
121 & ocean(ng) % t, &
122 & ocean(ng) % u, &
123 & ocean(ng) % v, &
124# endif
125 & ocean(ng) % ubar, &
126 & ocean(ng) % vbar, &
127 & ocean(ng) % zeta)
128# ifdef PROFILE
129 CALL wclock_off (ng, model, 2, __line__, myfile)
130# endif
131!
132 RETURN
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable krhs
integer, dimension(:), allocatable nstp
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, ini_fields_tile(), mod_stepping::knew, mod_stepping::krhs, mod_stepping::kstp, mod_stepping::nnew, mod_stepping::nstp, mod_ocean::ocean, mod_sedbed::sedbed, wclock_off(), and wclock_on().

Referenced by post_initial_mod::post_initial().

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

◆ ini_fields_tile()

subroutine ini_fields_mod::ini_fields_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
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) kstp,
integer, intent(in) krhs,
integer, intent(in) knew,
integer, intent(in) nstp,
integer, intent(in) nnew,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:), intent(in) rmask_wet,
real(r8), dimension(lbi:,lbj:), intent(inout) umask_wet,
real(r8), dimension(lbi:,lbj:), intent(inout) vmask_wet,
real(r8), dimension(lbi:,lbj:,0:,:), intent(inout) ru,
real(r8), dimension(lbi:,lbj:,0:,:), intent(inout) rv,
real(r8), dimension(lbi:,lbj:,:), intent(inout) rubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) rvbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) rzeta,
real(r8), dimension(lbi:,lbj:,:), intent(inout) bottom,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) bed,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) bed_frac,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) bed_mass,
real(r8), dimension(lbi:,lbj:,:), intent(in) hz,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) vbar,
real(r8), dimension(lbi:,lbj:,:), intent(in) zeta )
private

Definition at line 136 of file ini_fields.F.

166!***********************************************************************
167!
168! Imported variable declarations.
169!
170 integer, intent(in) :: ng, tile, model
171 integer, intent(in) :: LBi, UBi, LBj, UBj
172 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
173 integer, intent(in) :: kstp, krhs, knew
174# ifdef SOLVE3D
175 integer, intent(in) :: nstp, nnew
176# endif
177!
178# ifdef ASSUMED_SHAPE
179# ifdef MASKING
180 real(r8), intent(in) :: rmask(LBi:,LBj:)
181 real(r8), intent(in) :: umask(LBi:,LBj:)
182 real(r8), intent(in) :: vmask(LBi:,LBj:)
183# endif
184# ifdef WET_DRY
185 real(r8), intent(in) :: rmask_wet(LBi:,LBj:)
186 real(r8), intent(inout) :: umask_wet(LBi:,LBj:)
187 real(r8), intent(inout) :: vmask_wet(LBi:,LBj:)
188# endif
189# ifdef SOLVE3D
190 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
191# endif
192 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
193# ifdef SOLVE3D
194# if defined SEDIMENT || defined BBL_MODEL
195 real(r8), intent(inout) :: bottom(LBi:,LBj:,:)
196# endif
197# if defined SEDIMENT
198 real(r8), intent(inout) :: bed(LBi:,LBj:,:,:)
199 real(r8), intent(inout) :: bed_frac(LBi:,LBj:,:,:)
200 real(r8), intent(inout) :: bed_mass(LBi:,LBj:,:,:,:)
201# endif
202# endif
203# ifdef PERFECT_RESTART
204# ifdef SOLVE3D
205 real(r8), intent(inout) :: ru(LBi:,LBj:,0:,:)
206 real(r8), intent(inout) :: rv(LBi:,LBj:,0:,:)
207# endif
208 real(r8), intent(inout) :: rubar(LBi:,LBj:,:)
209 real(r8), intent(inout) :: rvbar(LBi:,LBj:,:)
210 real(r8), intent(inout) :: rzeta(LBi:,LBj:,:)
211# endif
212# ifdef SOLVE3D
213 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
214 real(r8), intent(inout) :: u(LBi:,LBj:,:,:)
215 real(r8), intent(inout) :: v(LBi:,LBj:,:,:)
216# endif
217 real(r8), intent(inout) :: ubar(LBi:,LBj:,:)
218 real(r8), intent(inout) :: vbar(LBi:,LBj:,:)
219
220# else
221
222# ifdef MASKING
223 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
224 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
225 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
226# endif
227# ifdef WET_DRY
228 real(r8), intent(in) :: rmask_wet(LBi:UBi,LBj:UBj)
229 real(r8), intent(inout) :: umask_wet(LBi:UBi,LBj:UBj)
230 real(r8), intent(inout) :: vmask_wet(LBi:UBi,LBj:UBj)
231# endif
232# ifdef SOLVE3D
233 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
234# endif
235 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
236# ifdef SOLVE3D
237# if defined SEDIMENT || defined BBL_MODEL
238 real(r8), intent(inout) :: bottom(LBi:UBi,LBj:UBj,MBOTP)
239# endif
240# if defined SEDIMENT
241 real(r8), intent(inout) :: bed(LBi:UBi,LBj:UBj,Nbed,MBEDP)
242 real(r8), intent(inout) :: bed_frac(LBi:UBi,LBj:UBj,Nbed,NST)
243 real(r8), intent(inout) :: bed_mass(LBi:UBi,LBj:UBj,Nbed,2,NST)
244# endif
245# endif
246# ifdef PERFECT_RESTART
247# ifdef SOLVE3D
248 real(r8), intent(inout) :: ru(LBi:UBi,LBj:UBj,0:N(ng),2)
249 real(r8), intent(inout) :: rv(LBi:UBi,LBj:UBj,0:N(ng),2)
250# endif
251 real(r8), intent(inout) :: rubar(LBi:UBi,LBj:UBj,2)
252 real(r8), intent(inout) :: rvbar(LBi:UBi,LBj:UBj,2)
253 real(r8), intent(inout) :: rzeta(LBi:UBi,LBj:UBj,2)
254# endif
255# ifdef SOLVE3D
256 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
257 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
258 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
259# endif
260 real(r8), intent(inout) :: ubar(LBi:UBi,LBj:UBj,:)
261 real(r8), intent(inout) :: vbar(LBi:UBi,LBj:UBj,:)
262# endif
263!
264! Local variable declarations.
265!
266 integer :: i, ic, itrc, j, k
267# if defined SEDIMENT || defined BBL_MODEL
268 integer :: ised
269# endif
270
271 real(r8) :: cff1, cff2
272# ifdef SOLVE3D
273 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
274 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
275# endif
276
277# include "set_bounds.h"
278
279# ifdef SOLVE3D
280!
281!-----------------------------------------------------------------------
282! If not perfect restart, initialize other time levels for 3D momentum.
283!-----------------------------------------------------------------------
284!
285 IF (.not.perfectrst(ng)) THEN
286 DO j=jstrb,jendb
287 DO k=1,n(ng)
288 DO i=istrm,iendb
289 cff1=u(i,j,k,nstp)
290# ifdef MASKING
291 cff1=cff1*umask(i,j)
292# endif
293# ifdef WET_DRY
294 cff1=cff1*umask_wet(i,j)
295# endif
296 u(i,j,k,nstp)=cff1
297 END DO
298 END DO
299!
300 IF (j.ge.jstrm) THEN
301 DO k=1,n(ng)
302 DO i=istrb,iendb
303 cff2=v(i,j,k,nstp)
304# ifdef MASKING
305 cff2=cff2*vmask(i,j)
306# endif
307# ifdef WET_DRY
308 cff2=cff2*vmask_wet(i,j)
309# endif
310 v(i,j,k,nstp)=cff2
311 END DO
312 END DO
313 END IF
314 END DO
315!
316! Apply boundary conditions.
317!
318 CALL u3dbc_tile (ng, tile, &
319 & lbi, ubi, lbj, ubj, n(ng), &
320 & imins, imaxs, jmins, jmaxs, &
321 & nstp, nstp, &
322 & u)
323 CALL v3dbc_tile (ng, tile, &
324 & lbi, ubi, lbj, ubj, n(ng), &
325 & imins, imaxs, jmins, jmaxs, &
326 & nstp, nstp, &
327 & v)
328 END IF
329!
330 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
331 CALL exchange_u3d_tile (ng, tile, &
332 & lbi, ubi, lbj, ubj, 1, n(ng), &
333 & u(:,:,:,nstp))
334 CALL exchange_v3d_tile (ng, tile, &
335 & lbi, ubi, lbj, ubj, 1, n(ng), &
336 & v(:,:,:,nstp))
337 END IF
338
339# ifdef DISTRIBUTE
340!
341 CALL mp_exchange3d (ng, tile, model, 2, &
342 & lbi, ubi, lbj, ubj, 1, n(ng), &
343 & nghostpoints, &
344 & ewperiodic(ng), nsperiodic(ng), &
345 & u(:,:,:,nstp), v(:,:,:,nstp))
346# endif
347# endif
348
349# ifdef SOLVE3D
350!
351!-----------------------------------------------------------------------
352! If not perfect restart, compute vertically-integrated momentum
353! (ubar, vbar) from initial 3D momentum (u, v).
354!-----------------------------------------------------------------------
355!
356! Here DC(i,1:N) are the grid cell thicknesses, DC(i,0) is the total
357! depth of the water column, and CF(i,0) is the vertical integral.
358!
359 IF (.not.perfectrst(ng)) THEN
360 DO j=jstrb,jendb
361 DO i=istrm,iendb
362 dc(i,0)=0.0_r8
363 cf(i,0)=0.0_r8
364 END DO
365 DO k=1,n(ng)
366 DO i=istrm,iendb
367 dc(i,k)=0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
368 dc(i,0)=dc(i,0)+dc(i,k)
369 cf(i,0)=cf(i,0)+dc(i,k)*u(i,j,k,nstp)
370 END DO
371 END DO
372 DO i=istrm,iendb
373 cff1=1.0_r8/dc(i,0)
374 cff2=cf(i,0)*cff1
375# ifdef MASKING
376 cff2=cff2*umask(i,j)
377# endif
378# ifdef WET_DRY
379 cff2=cff2*umask_wet(i,j)
380# endif
381 ubar(i,j,kstp)=cff2
382 END DO
383!
384 IF (j.ge.jstrm) THEN
385 DO i=istrb,iendb
386 dc(i,0)=0.0_r8
387 cf(i,0)=0.0_r8
388 END DO
389 DO k=1,n(ng)
390 DO i=istrb,iendb
391 dc(i,k)=0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
392 dc(i,0)=dc(i,0)+dc(i,k)
393 cf(i,0)=cf(i,0)+dc(i,k)*v(i,j,k,nstp)
394 END DO
395 END DO
396 DO i=istrb,iendb
397 cff1=1.0_r8/dc(i,0)
398 cff2=cf(i,0)*cff1
399# ifdef MASKING
400 cff2=cff2*vmask(i,j)
401# endif
402# ifdef WET_DRY
403 cff2=cff2*vmask_wet(i,j)
404# endif
405 vbar(i,j,kstp)=cff2
406 END DO
407 END IF
408 END DO
409!
410! Apply boundary conditions.
411!
412 IF (.not.(any(lbc(:,isubar,ng)%radiation).or. &
413 & any(lbc(:,isvbar,ng)%radiation).or. &
414 & any(lbc(:,isubar,ng)%Flather).or. &
415 & any(lbc(:,isvbar,ng)%Flather))) THEN
416 CALL u2dbc_tile (ng, tile, &
417 & lbi, ubi, lbj, ubj, &
418 & imins, imaxs, jmins, jmaxs, &
419 & krhs, kstp, kstp, &
420 & ubar, vbar, zeta)
421 CALL v2dbc_tile (ng, tile, &
422 & lbi, ubi, lbj, ubj, &
423 & imins, imaxs, jmins, jmaxs, &
424 & krhs, kstp, kstp, &
425 & ubar, vbar, zeta)
426 END IF
427 END IF
428!
429 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
430 CALL exchange_u2d_tile (ng, tile, &
431 & lbi, ubi, lbj, ubj, &
432 & ubar(:,:,kstp))
433 CALL exchange_v2d_tile (ng, tile, &
434 & lbi, ubi, lbj, ubj, &
435 & vbar(:,:,kstp))
436 END IF
437
438# ifdef DISTRIBUTE
439!
440 CALL mp_exchange2d (ng, tile, model, 2, &
441 & lbi, ubi, lbj, ubj, &
442 & nghostpoints, &
443 & ewperiodic(ng), nsperiodic(ng), &
444 & ubar(:,:,kstp), vbar(:,:,kstp))
445# endif
446
447# else
448!
449!-----------------------------------------------------------------------
450! If not perfect restart, initialize other time levels for 2D momentum
451! (shallow-water model).
452!-----------------------------------------------------------------------
453!
454 IF (.not.perfectrst(ng)) THEN
455 DO j=jstrb,jendb
456 DO i=istrm,iendb
457 cff1=ubar(i,j,kstp)
458# ifdef MASKING
459 cff1=cff1*umask(i,j)
460# endif
461# ifdef WET_DRY
462 cff1=cff1*umask_wet(i,j)
463# endif
464 ubar(i,j,kstp)=cff1
465 END DO
466!
467 IF (j.ge.jstrm) THEN
468 DO i=istrb,iendb
469 cff2=vbar(i,j,kstp)
470# ifdef MASKING
471 cff2=cff2*vmask(i,j)
472# endif
473# ifdef WET_DRY
474 cff2=cff2*vmask_wet(i,j)
475# endif
476 vbar(i,j,kstp)=cff2
477 END DO
478 END IF
479 END DO
480!
481! Apply boundary conditions.
482!
483 IF (.not.(any(lbc(:,isubar,ng)%radiation).or. &
484 & any(lbc(:,isvbar,ng)%radiation).or. &
485 & any(lbc(:,isubar,ng)%Flather).or. &
486 & any(lbc(:,isvbar,ng)%Flather))) THEN
487 CALL u2dbc_tile (ng, tile, &
488 & lbi, ubi, lbj, ubj, &
489 & imins, imaxs, jmins, jmaxs, &
490 & krhs, kstp, kstp, &
491 & ubar, vbar, zeta)
492 CALL v2dbc_tile (ng, tile, &
493 & lbi, ubi, lbj, ubj, &
494 & imins, imaxs, jmins, jmaxs, &
495 & krhs, kstp, kstp, &
496 & ubar, vbar, zeta)
497 END IF
498 END IF
499!
500 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
501 CALL exchange_u2d_tile (ng, tile, &
502 & lbi, ubi, lbj, ubj, &
503 & ubar(:,:,kstp))
504 CALL exchange_v2d_tile (ng, tile, &
505 & lbi, ubi, lbj, ubj, &
506 & vbar(:,:,kstp))
507 IF (perfectrst(ng)) THEN
508 CALL exchange_u2d_tile (ng, tile, &
509 & lbi, ubi, lbj, ubj, &
510 & ubar(:,:,knew))
511 CALL exchange_v2d_tile (ng, tile, &
512 & lbi, ubi, lbj, ubj, &
513 & vbar(:,:,knew))
514 END IF
515 END IF
516
517# ifdef DISTRIBUTE
518!
519 CALL mp_exchange2d (ng, tile, model, 2, &
520 & lbi, ubi, lbj, ubj, &
521 & nghostpoints, &
522 & ewperiodic(ng), nsperiodic(ng), &
523 & ubar(:,:,kstp), vbar(:,:,kstp))
524 IF (perfectrst(ng)) THEN
525 CALL mp_exchange2d (ng, tile, model, 2, &
526 & lbi, ubi, lbj, ubj, &
527 & nghostpoints, &
528 & ewperiodic(ng), nsperiodic(ng), &
529 & ubar(:,:,knew), vbar(:,:,knew))
530 END IF
531# endif
532# endif
533
534# ifdef SOLVE3D
535!
536!-----------------------------------------------------------------------
537! If not perfect restart, initialize other time levels for tracers.
538!-----------------------------------------------------------------------
539!
540 ic=0
541 IF (.not.perfectrst(ng)) THEN
542 DO itrc=1,nt(ng)
543 IF (ltracerclm(itrc,ng).and.lnudgetclm(itrc,ng)) THEN
544 ic=ic+1 ! OBC nudging coefficient index
545 END IF
546 DO k=1,n(ng)
547 DO j=jstrb,jendb
548 DO i=istrb,iendb
549 cff1=t(i,j,k,nstp,itrc)
550# ifdef MASKING
551 cff1=cff1*rmask(i,j)
552# endif
553 t(i,j,k,nstp,itrc)=cff1
554 END DO
555 END DO
556 END DO
557!
558! Apply boundary conditions.
559!
560 CALL t3dbc_tile (ng, tile, itrc, ic, &
561 & lbi, ubi, lbj, ubj, n(ng), nt(ng), &
562 & imins, imaxs, jmins, jmaxs, &
563 & nstp, nstp, &
564 & t)
565 END DO
566 END IF
567!
568 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
569 DO itrc=1,nt(ng)
570 CALL exchange_r3d_tile (ng, tile, &
571 & lbi, ubi, lbj, ubj, 1, n(ng), &
572 & t(:,:,:,nstp,itrc))
573 END DO
574 END IF
575
576# ifdef DISTRIBUTE
577!
578 CALL mp_exchange4d (ng, tile, model, 1, &
579 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
580 & nghostpoints, &
581 & ewperiodic(ng), nsperiodic(ng), &
582 & t(:,:,:,nstp,:))
583# endif
584
585# if defined BBL_MODEL || defined SEDIMENT
586!
587!-----------------------------------------------------------------------
588! Initialize sediment bed properties.
589!-----------------------------------------------------------------------
590!
591# ifdef SEDIMENT
592 IF (.not.perfectrst(ng)) THEN
593 DO k=1,nbed
594 DO j=jstrt,jendt
595 DO i=istrt,iendt
596 DO ised=1,nst
597 cff1=bed_mass(i,j,k,1,ised)
598 bed_mass(i,j,k,2,ised)=cff1
599 END DO
600 END DO
601 END DO
602 END DO
603 END IF
604# endif
605!
606 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
607# ifdef SEDIMENT
608 DO ised=1,nst
609 CALL exchange_r3d_tile (ng, tile, &
610 & lbi, ubi, lbj, ubj, 1, nbed, &
611 & bed_frac(:,:,:,ised))
612 CALL exchange_r3d_tile (ng, tile, &
613 & lbi, ubi, lbj, ubj, 1, nbed, &
614 & bed_mass(:,:,:,1,ised))
615 CALL exchange_r3d_tile (ng, tile, &
616 & lbi, ubi, lbj, ubj, 1, nbed, &
617 & bed_mass(:,:,:,2,ised))
618 END DO
619 DO ised=1,mbedp
620 CALL exchange_r3d_tile (ng, tile, &
621 & lbi, ubi, lbj, ubj, 1, nbed, &
622 & bed(:,:,:,ised))
623 END DO
624# endif
625 CALL exchange_r3d_tile (ng, tile, &
626 & lbi, ubi, lbj, ubj, 1, mbotp, &
627 & bottom)
628 END IF
629
630# ifdef DISTRIBUTE
631# ifdef SEDIMENT
632 CALL mp_exchange4d (ng, tile, model, 1, &
633 & lbi, ubi, lbj, ubj, 1, nbed, 1, nst, &
634 & nghostpoints, &
635 & ewperiodic(ng), nsperiodic(ng), &
636 & bed_frac)
637 CALL mp_exchange4d (ng, tile, model, 2, &
638 & lbi, ubi, lbj, ubj, 1, nbed, 1, nst, &
639 & nghostpoints, &
640 & ewperiodic(ng), nsperiodic(ng), &
641 & bed_mass(:,:,:,1,:),bed_mass(:,:,:,2,:))
642 CALL mp_exchange4d (ng, tile, model, 1, &
643 & lbi, ubi, lbj, ubj, 1, nbed, 1, mbedp, &
644 & nghostpoints, &
645 & ewperiodic(ng), nsperiodic(ng), &
646 & bed)
647# endif
648 CALL mp_exchange3d (ng, tile, model, 1, &
649 & lbi, ubi, lbj, ubj, 1, mbotp, &
650 & nghostpoints, &
651 & ewperiodic(ng), nsperiodic(ng), &
652 & bottom)
653# endif
654# endif
655# endif
656
657# ifdef PERFECT_RESTART
658!
659!-----------------------------------------------------------------------
660! If perfect restart, apply periodic boundary condition to
661! right-hand-side terms.
662!-----------------------------------------------------------------------
663!
664 IF (perfectrst(ng)) THEN
665 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
666 DO i=1,2
667 CALL exchange_u2d_tile (ng, tile, &
668 & lbi, ubi, lbj, ubj, &
669 & rubar(:,:,i))
670 CALL exchange_v2d_tile (ng, tile, &
671 & lbi, ubi, lbj, ubj, &
672 & rvbar(:,:,i))
673 CALL exchange_r2d_tile (ng, tile, &
674 & lbi, ubi, lbj, ubj, &
675 & rzeta(:,:,i))
676# ifdef SOLVE3D
677 CALL exchange_u3d_tile (ng, tile, &
678 & lbi, ubi, lbj, ubj, 0, n(ng), &
679 & ru(:,:,:,i))
680 CALL exchange_v3d_tile (ng, tile, &
681 & lbi, ubi, lbj, ubj, 0, n(ng), &
682 & rv(:,:,:,i))
683# endif
684 END DO
685 END IF
686
687# ifdef DISTRIBUTE
688 CALL mp_exchange3d (ng, tile, model, 3, &
689 & lbi, ubi, lbj, ubj, 1, 2, &
690 & nghostpoints, &
691 & ewperiodic(ng), nsperiodic(ng), &
692 & rubar, rvbar, rzeta)
693# ifdef SOLVE3D
694 CALL mp_exchange4d (ng, tile, model, 2, &
695 & lbi, ubi, lbj, ubj, 0, n(ng), 1, 2, &
696 & nghostpoints, &
697 & ewperiodic(ng), nsperiodic(ng), &
698 & ru, rv)
699# endif
700# endif
701 END IF
702# endif
703!
704 RETURN

References mod_scalars::ewperiodic, exchange_2d_mod::exchange_r2d_tile(), exchange_3d_mod::exchange_r3d_tile(), exchange_2d_mod::exchange_u2d_tile(), exchange_3d_mod::exchange_u3d_tile(), exchange_2d_mod::exchange_v2d_tile(), exchange_3d_mod::exchange_v3d_tile(), mod_ncparam::isubar, mod_ncparam::isvbar, mod_param::lbc, mod_scalars::lnudgetclm, mod_scalars::ltracerclm, mp_exchange_mod::mp_exchange2d(), mp_exchange_mod::mp_exchange3d(), mp_exchange_mod::mp_exchange4d(), mod_param::nghostpoints, mod_scalars::nsperiodic, mod_scalars::perfectrst, t3dbc_mod::t3dbc_tile(), u2dbc_mod::u2dbc_tile(), u3dbc_mod::u3dbc_tile(), v2dbc_mod::v2dbc_tile(), and v3dbc_mod::v3dbc_tile().

Referenced by ini_fields().

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

◆ ini_zeta()

subroutine, public ini_fields_mod::ini_zeta ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )

Definition at line 708 of file ini_fields.F.

709!***********************************************************************
710!
711 USE mod_stepping
712!
713! Imported variable declarations.
714!
715 integer, intent(in) :: ng, tile, model
716!
717! Local variable declarations.
718!
719 character (len=*), parameter :: MyFile = &
720 & __FILE__//", ini_zeta"
721!
722# include "tile.h"
723!
724# ifdef PROFILE
725 CALL wclock_on (ng, model, 2, __line__, myfile)
726# endif
727 CALL ini_zeta_tile (ng, tile, model, &
728 & lbi, ubi, lbj, ubj, &
729 & imins, imaxs, jmins, jmaxs, &
730 & kstp(ng), krhs(ng), knew(ng), &
731# ifdef MASKING
732 & grid(ng) % rmask, &
733# endif
734# ifdef WET_DRY
735 & grid(ng) % rmask_wet, &
736 & grid(ng) % h, &
737# endif
738# ifdef SOLVE3D
739# if defined SEDIMENT && defined SED_MORPH
740 & sedbed(ng) % bed, &
741 & sedbed(ng) % bed_thick0, &
742 & sedbed(ng) % bed_thick, &
743# endif
744# endif
745 & ocean(ng) % zeta)
746# ifdef PROFILE
747 CALL wclock_off (ng, model, 2, __line__, myfile)
748# endif
749!
750 RETURN

References mod_grid::grid, ini_zeta_tile(), mod_stepping::knew, mod_stepping::krhs, mod_stepping::kstp, mod_ocean::ocean, mod_sedbed::sedbed, wclock_off(), and wclock_on().

Referenced by post_initial_mod::post_initial().

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

◆ ini_zeta_tile()

subroutine ini_fields_mod::ini_zeta_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
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) kstp,
integer, intent(in) krhs,
integer, intent(in) knew,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(inout) rmask_wet,
real(r8), dimension(lbi:,lbj:), intent(in) h,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) bed,
real(r8), dimension(lbi:,lbj:), intent(inout) bed_thick0,
real(r8), dimension(lbi:,lbj:,:), intent(inout) bed_thick,
real(r8), dimension(lbi:,lbj:,:), intent(inout) zeta )
private

Definition at line 754 of file ini_fields.F.

771!***********************************************************************
772!
773! Imported variable declarations.
774!
775 integer, intent(in) :: ng, tile, model
776 integer, intent(in) :: LBi, UBi, LBj, UBj
777 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
778 integer, intent(in) :: kstp, krhs, knew
779!
780# ifdef ASSUMED_SHAPE
781# ifdef MASKING
782 real(r8), intent(in) :: rmask(LBi:,LBj:)
783# endif
784# ifdef WET_DRY
785 real(r8), intent(in) :: h(LBi:,LBj:)
786# endif
787# if defined SOLVE3D && defined SEDIMENT && defined SED_MORPH
788 real(r8), intent(in) :: bed(LBi:,LBj:,:,:)
789 real(r8), intent(inout) :: bed_thick0(LBi:,LBj:)
790 real(r8), intent(inout) :: bed_thick(LBi:,LBj:,:)
791# endif
792# ifdef WET_DRY
793 real(r8), intent(inout) :: rmask_wet(LBi:,LBj:)
794# endif
795 real(r8), intent(inout) :: zeta(LBi:,LBj:,:)
796
797# else
798
799# ifdef MASKING
800 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
801# endif
802# ifdef WET_DRY
803 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
804# endif
805# if defined SOLVE3D && defined SEDIMENT && defined SED_MORPH
806 real(r8), intent(in) :: bed(LBi:UBi,LBj:UBj,Nbed,MBEDP)
807 real(r8), intent(inout) :: bed_thick0(LBi:UBi,LBj:UBj)
808 real(r8), intent(inout) :: bed_thick(LBi:UBi,LBj:UBj,3)
809# endif
810# ifdef WET_DRY
811 real(r8), intent(inout) :: rmask_wet(LBi:UBi,LBj:UBj)
812# endif
813 real(r8), intent(inout) :: zeta(LBi:UBi,LBj:UBj,:)
814# endif
815!
816! Local variable declarations.
817!
818 integer :: Imin, Imax, Jmin, Jmax
819 integer :: i, j, kbed
820
821 real(r8) :: cff1
822
823# include "set_bounds.h"
824!
825!-----------------------------------------------------------------------
826! Initialize other time levels for free-surface.
827!-----------------------------------------------------------------------
828!
829 IF (.not.perfectrst(ng)) THEN
830 IF (.not.(any(lbc(:,isfsur,ng)%radiation).or. &
831 & any(lbc(:,isfsur,ng)%Chapman_explicit).or. &
832 & any(lbc(:,isfsur,ng)%Chapman_implicit))) THEN
833 imin=istrb
834 imax=iendb
835 jmin=jstrb
836 jmax=jendb
837 ELSE
838 imin=istrt
839 imax=iendt
840 jmin=jstrt
841 jmax=jendt
842 END IF
843 DO j=jmin,jmax
844 DO i=imin,imax
845 cff1=zeta(i,j,kstp)
846# ifdef MASKING
847 cff1=cff1*rmask(i,j)
848# endif
849# ifdef WET_DRY
850 IF (cff1.le.(dcrit(ng)-h(i,j))) THEN
851 cff1=dcrit(ng)-h(i,j)
852!! rmask_wet(i,j)=0.0_r8
853!! ELSE
854!! rmask_wet(i,j)=1.0_r8*rmask(i,j)
855 END IF
856# endif
857 zeta(i,j,kstp)=cff1
858 END DO
859 END DO
860!
861! Apply boundary conditions.
862!
863 IF (.not.(any(lbc(:,isfsur,ng)%radiation).or. &
864 & any(lbc(:,isfsur,ng)%Chapman_explicit).or. &
865 & any(lbc(:,isfsur,ng)%Chapman_implicit))) THEN
866 CALL zetabc_tile (ng, tile, &
867 & lbi, ubi, lbj, ubj, &
868 & imins, imaxs, jmins, jmaxs, &
869 & krhs, kstp, kstp, &
870 & zeta)
871 END IF
872 END IF
873!
874 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
875 CALL exchange_r2d_tile (ng, tile, &
876 & lbi, ubi, lbj, ubj, &
877 & zeta(:,:,kstp))
878 IF (perfectrst(ng)) THEN
879# ifdef SOLVE3D
880 CALL exchange_r2d_tile (ng, tile, &
881 & lbi, ubi, lbj, ubj, &
882 & zeta(:,:,krhs))
883# else
884 CALL exchange_r2d_tile (ng, tile, &
885 & lbi, ubi, lbj, ubj, &
886 & zeta(:,:,knew))
887# endif
888 END IF
889
890# ifdef WET_DRY
891 CALL exchange_r2d_tile (ng, tile, &
892 & lbi, ubi, lbj, ubj, &
893 & rmask_wet)
894# endif
895 END IF
896
897# ifdef DISTRIBUTE
898 CALL mp_exchange2d (ng, tile, model, 1, &
899 & lbi, ubi, lbj, ubj, &
900 & nghostpoints, &
901 & ewperiodic(ng), nsperiodic(ng), &
902 & zeta(:,:,kstp))
903 IF (perfectrst(ng)) THEN
904# ifdef SOLVE3D
905 CALL mp_exchange2d (ng, tile, model, 1, &
906 & lbi, ubi, lbj, ubj, &
907 & nghostpoints, &
908 & ewperiodic(ng), nsperiodic(ng), &
909 & zeta(:,:,krhs))
910# else
911 CALL mp_exchange2d (ng, tile, model, 1, &
912 & lbi, ubi, lbj, ubj, &
913 & nghostpoints, &
914 & ewperiodic(ng), nsperiodic(ng), &
915 & zeta(:,:,knew))
916# endif
917 END IF
918
919# ifdef WET_DRY
920 CALL mp_exchange2d (ng, tile, model, 1, &
921 & lbi, ubi, lbj, ubj, &
922 & nghostpoints, &
923 & ewperiodic(ng), nsperiodic(ng), &
924 & rmask_wet)
925# endif
926# endif
927
928# ifdef SOLVE3D
929!
930!-----------------------------------------------------------------------
931! Initialize fast-time averaged free-surface (Zt_avg1) with the inital
932! free-surface
933!-----------------------------------------------------------------------
934!
935 CALL set_zeta_timeavg (ng, tile, model)
936
937# if defined SEDIMENT && defined SED_MORPH
938!
939!-----------------------------------------------------------------------
940! Compute initial total thickness for all sediment bed layers.
941!-----------------------------------------------------------------------
942!
943 DO j=jstrt,jendt
944 DO i=istrt,iendt
945 bed_thick0(i,j)=0.0_r8
946 DO kbed=1,nbed
947 bed_thick0(i,j)=bed_thick0(i,j)+bed(i,j,kbed,ithck)
948 END DO
949 bed_thick(i,j,1)=bed_thick0(i,j)
950 bed_thick(i,j,2)=bed_thick0(i,j)
951 END DO
952 END DO
953!
954 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
955 CALL exchange_r2d_tile (ng, tile, &
956 & lbi, ubi, lbj, ubj, &
957 & bed_thick0)
958 CALL exchange_r2d_tile (ng, tile, &
959 & lbi, ubi, lbj, ubj, &
960 & bed_thick(:,:,1))
961 CALL exchange_r2d_tile (ng, tile, &
962 & lbi, ubi, lbj, ubj, &
963 & bed_thick(:,:,2))
964 END IF
965
966# ifdef DISTRIBUTE
967 CALL mp_exchange2d (ng, tile, model, 3, &
968 & lbi, ubi, lbj, ubj, &
969 & nghostpoints, &
970 & ewperiodic(ng), nsperiodic(ng), &
971 & bed_thick0, &
972 & bed_thick(:,:,1), bed_thick(:,:,2))
973# endif
974# endif
975# endif
976!
977 RETURN
978!

References mod_scalars::dcrit, mod_scalars::ewperiodic, exchange_2d_mod::exchange_r2d_tile(), mod_ncparam::isfsur, mod_sediment::ithck, mod_param::lbc, mp_exchange_mod::mp_exchange2d(), mod_param::nghostpoints, mod_scalars::nsperiodic, mod_scalars::perfectrst, set_zeta_timeavg(), and zetabc_mod::zetabc_tile().

Referenced by ini_zeta().

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

◆ set_zeta_timeavg()

subroutine, public ini_fields_mod::set_zeta_timeavg ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )

Definition at line 985 of file ini_fields.F.

986!***********************************************************************
987!
988 USE mod_stepping
989!
990! Imported variable declarations.
991!
992 integer, intent(in) :: ng, tile, model
993!
994! Local variable declarations.
995!
996 character (len=*), parameter :: MyFile = &
997 & __FILE__//", set_zeta_timeavg"
998!
999# include "tile.h"
1000!
1001# ifdef PROFILE
1002 CALL wclock_on (ng, model, 2, __line__, myfile)
1003# endif
1004 CALL set_zeta_timeavg_tile (ng, tile, model, &
1005 & lbi, ubi, lbj, ubj, &
1006 & kstp(ng), &
1007 & ocean(ng) % zeta, &
1008 & coupling(ng) % Zt_avg1)
1009# ifdef PROFILE
1010 CALL wclock_off (ng, model, 2, __line__, myfile)
1011# endif
1012!
1013 RETURN

References mod_coupling::coupling, mod_stepping::kstp, mod_ocean::ocean, set_zeta_timeavg_tile(), wclock_off(), and wclock_on().

Referenced by ad_initial(), ad_post_initial_mod::ad_post_initial(), ini_zeta_tile(), initial(), tl_initial(), and tl_post_initial_mod::tl_post_initial().

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

◆ set_zeta_timeavg_tile()

subroutine ini_fields_mod::set_zeta_timeavg_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) kstp,
real(r8), dimension(lbi:,lbj:,:), intent(in) zeta,
real(r8), dimension(lbi:,lbj:), intent(inout) zt_avg1 )
private

Definition at line 1017 of file ini_fields.F.

1022!***********************************************************************
1023!
1024! Imported variable declarations.
1025!
1026 integer, intent(in) :: ng, tile, model
1027 integer, intent(in) :: LBi, UBi, LBj, UBj
1028 integer, intent(in) :: kstp
1029!
1030# ifdef ASSUMED_SHAPE
1031 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
1032 real(r8), intent(inout) :: Zt_avg1(LBi:,LBj:)
1033# else
1034 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
1035 real(r8), intent(inout) :: Zt_avg1(LBi:UBi,LBj:UBj)
1036# endif
1037!
1038! Local variable declarations.
1039!
1040 integer :: i, j
1041!
1042# include "set_bounds.h"
1043!
1044!-----------------------------------------------------------------------
1045! Initialize fast-time averaged free-surface (Zt_avg1) with the inital
1046! free-surface.
1047!-----------------------------------------------------------------------
1048!
1049 DO j=jstrt,jendt
1050 DO i=istrt,iendt
1051 zt_avg1(i,j)=zeta(i,j,kstp)
1052 END DO
1053 END DO
1054!
1055 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1056 CALL exchange_r2d_tile (ng, tile, &
1057 & lbi, ubi, lbj, ubj, &
1058 & zt_avg1)
1059 END IF
1060
1061# ifdef DISTRIBUTE
1062!
1063 CALL mp_exchange2d (ng, tile, model, 1, &
1064 & lbi, ubi, lbj, ubj, &
1065 & nghostpoints, &
1066 & ewperiodic(ng), nsperiodic(ng), &
1067 & zt_avg1)
1068# endif
1069!
1070 RETURN

References mod_scalars::ewperiodic, exchange_2d_mod::exchange_r2d_tile(), mp_exchange_mod::mp_exchange2d(), mod_param::nghostpoints, and mod_scalars::nsperiodic.

Referenced by set_zeta_timeavg().

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