3***************************************************** john c. warner ***
4** copyright(c) 2002-2025 the roms group hernan g. arango **
5** licensed under a mit/x style license **
6** see license_roms.md **
7************************************************************************
9** these routines are
use couple roms to swan wave model using **
10** the model coupling toolkit(mct). **
12************************************************************************
15 SUBROUTINE initialize_ocn2wav_coupling (ng, tile)
34 integer,
intent(in) :: ng, tile
38 integer :: Istr, Iend, Jstr, Jend
39 integer :: IstrR, IendR, JstrR, JendR, IstrU, JstrV
40 integer :: Asize, Jsize, MyError
41 integer :: j, jc, nprocs
43 integer,
allocatable :: length(:)
44 integer,
allocatable :: start(:)
58 IF (
domain(ng)%Western_Edge(tile))
THEN
59 istrr=
bounds(ng)%Istr(tile)-1
61 istrr=
bounds(ng)%Istr(tile)
63 IF (
domain(ng)%Eastern_Edge(tile))
THEN
64 iendr=
bounds(ng)%Iend(tile)+1
66 iendr=
bounds(ng)%Iend(tile)
68 IF (
domain(ng)%Southern_Edge(tile))
THEN
69 jstrr=
bounds(ng)%Jstr(tile)-1
71 jstrr=
bounds(ng)%Jstr(tile)
73 IF (
domain(ng)%Northern_Edge(tile))
THEN
74 jendr=
bounds(ng)%Jend(tile)+1
76 jendr=
bounds(ng)%Jend(tile)
96 IF (.not.
allocated(start))
THEN
97 allocate ( start(jsize) )
99 IF (.not.
allocated(length))
THEN
100 allocate ( length(jsize) )
105 start(jc)=j*(
lm(ng)+2)+istrr+1
106 length(jc)=(iendr-istrr+1)
108 CALL globalsegmap_init (gsmaproms, start, length, 0, &
124 CALL attrvect_zero (ocn2wav_av)
132 IF (
allocated(start))
THEN
135 IF (
allocated(length))
THEN
140 END SUBROUTINE initialize_ocn2wav_coupling
142 SUBROUTINE ocn2wav_coupling (ng, tile)
179 integer,
intent(in) :: ng, tile
183 character (len=*),
parameter :: MyFile = &
191 CALL ocn2wav_coupling_tile (ng, tile, &
192 & lbi, ubi, lbj, ubj)
198 END SUBROUTINE ocn2wav_coupling
201 SUBROUTINE ocn2wav_coupling_tile (ng, tile, &
202 & LBi, UBi, LBj, UBj)
220 USE roms_import_mod,
ONLY : roms_import2d
221 USE roms_export_mod,
ONLY : roms_export2d
227 integer,
intent(in) :: ng, tile
228 integer,
intent(in) :: LBi, UBi, LBj, UBj
232 integer :: Istr, Iend, Jstr, Jend
233 integer :: IstrR, IendR, JstrR, JendR, IstrU, JstrV
234 integer :: Asize, Iimport, Iexport, MyError
235 integer :: gtype, i, id, ifield, ij, j, status
237 real(r8),
parameter :: Lwave_max=500.0_r8
239 real(r8) :: add_offset, scale
240 real(r8) :: RecvTime, SendTime, buffer(2), wtime(2)
244 real(r8),
dimension(LBi:UBi,LBj:UBj) :: Awrk
246 real(r8),
pointer :: A(:)
248 character (len=3 ),
dimension(2) :: op_handle
249 character (len=40) :: code
258 istr=
bounds(ng)%Istr(tile)
259 iend=
bounds(ng)%Iend(tile)
260 jstr=
bounds(ng)%Jstr(tile)
261 jend=
bounds(ng)%Jend(tile)
263 IF (
domain(ng)%Western_Edge(tile))
THEN
264 istrr=
bounds(ng)%Istr(tile)-1
266 istrr=
bounds(ng)%Istr(tile)
268 IF (
domain(ng)%Eastern_Edge(tile))
THEN
269 iendr=
bounds(ng)%Iend(tile)+1
271 iendr=
bounds(ng)%Iend(tile)
273 IF (
domain(ng)%Southern_Edge(tile))
THEN
274 jstrr=
bounds(ng)%Jstr(tile)-1
276 jstrr=
bounds(ng)%Jstr(tile)
278 IF (
domain(ng)%Northern_Edge(tile))
THEN
279 jendr=
bounds(ng)%Jend(tile)+1
281 jendr=
bounds(ng)%Jend(tile)
289 allocate ( a(asize) )
306 buffer(1)=my_wtime(wtime)
307 CALL mct_recv (wav2ocn_av, romstoswan, myerror)
308 recvtime=recvtime+my_wtime(wtime)-buffer(1)
309 IF (myerror.ne.0)
THEN
311 WRITE (
stdout,10)
'wave model, MyError = ', myerror
322 code=adjustl(
fields(id)%code)
325 add_offset=
fields(id)%AddOffset
327 SELECT CASE (trim(code))
331 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
334 a(i)=max(0.0_r8,a(i))
338 CALL roms_import2d (ng, tile, &
339 & id, gtype, scale, add_offset, &
341 & istrr, iendr, jstrr, jendr, &
342 & lbi, ubi, lbj, ubj, &
349 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
352 a(i)=max(0.0_r8,a(i))
356 CALL roms_import2d (ng, tile, &
357 & id, gtype, scale, add_offset, &
359 & istrr, iendr, jstrr, jendr, &
360 & lbi, ubi, lbj, ubj, &
367 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
370 a(i)=max(0.0_r8,a(i))
377 CALL roms_import2d (ng, tile, &
378 & id, gtype, scale, add_offset, &
380 & istrr, iendr, jstrr, jendr, &
381 & lbi, ubi, lbj, ubj, &
388 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
391 a(i)=max(0.0_r8,a(i))
395 CALL roms_import2d (ng, tile, &
396 & id, gtype, scale, add_offset, &
398 & istrr, iendr, jstrr, jendr, &
399 & lbi, ubi, lbj, ubj, &
406 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
409 a(i)=max(0.0_r8,a(i))
413 CALL roms_import2d (ng, tile, &
414 & id, gtype, scale, add_offset, &
416 & istrr, iendr, jstrr, jendr, &
417 & lbi, ubi, lbj, ubj, &
424 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
427 a(i)=max(0.0_r8,a(i))
431 CALL roms_import2d (ng, tile, &
432 & id, gtype, scale, add_offset, &
434 & istrr, iendr, jstrr, jendr, &
435 & lbi, ubi, lbj, ubj, &
437 &
forces(ng)%Wave_dissip, &
442 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
445 a(i)=max(0.0_r8,a(i))
449 CALL roms_import2d (ng, tile, &
450 & id, gtype, scale, add_offset, &
452 & istrr, iendr, jstrr, jendr, &
453 & lbi, ubi, lbj, ubj, &
458#ifdef SVENDSEN_ROLLER
462 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
465 a(i)=max(0.0_r8,a(i))
469 CALL roms_import2d (ng, tile, &
470 & id, gtype, scale, add_offset, &
472 & istrr, iendr, jstrr, jendr, &
473 & lbi, ubi, lbj, ubj, &
475 &
forces(ng)%Wave_break, &
490 code=adjustl(
fields(id)%code)
493 add_offset=
fields(id)%AddOffset
495 SELECT CASE (trim(code))
499 CALL roms_export2d (ng, tile, &
500 & id, gtype, scale, add_offset, &
501 & lbi, ubi, lbj, ubj, &
506 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
511 CALL roms_export2d (ng, tile, &
512 & id, gtype, scale, add_offset, &
513 & lbi, ubi, lbj, ubj, &
514 &
ocean(ng)%zeta(:,:,kout), &
518 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
523 CALL roms_export2d (ng, tile, &
524 & id, gtype, scale, add_offset, &
525 & lbi, ubi, lbj, ubj, &
527 &
ocean(ng)%u(:,:,
n(ng),nout), &
529 &
ocean(ng)%ubar(:,:,kout), &
534 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
539 CALL roms_export2d (ng, tile, &
540 & id, gtype, scale, add_offset, &
541 & lbi, ubi, lbj, ubj, &
543 &
ocean(ng)%v(:,:,
n(ng),nout), &
545 &
ocean(ng)%vbar(:,:,kout), &
550 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
558 awrk(i,j)=max(0.0001_r8, &
561 awrk(i,j)=max(0.0001_r8,
rdrg2(ng))
565 CALL roms_export2d (ng, tile, &
566 & id, gtype, scale, add_offset, &
567 & lbi, ubi, lbj, ubj, &
572 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
580 IF (iexport.gt.0)
THEN
581 buffer(2)=my_wtime(wtime)
582 CALL mct_send (ocn2wav_av, romstoswan, myerror)
583 sendtime=sendtime+my_wtime(wtime)-buffer(2)
584 IF (myerror.ne.0)
THEN
586 WRITE (
stdout,20)
'wave model, MyError = ', myerror
606 IF (
master.and.((iimport.gt.0).or.(iexport.gt.0)))
THEN
612 WRITE (
stdout,40)
'ROMS Import: ',trim(
fields(id)%name), &
617 WRITE (
stdout,40)
'ROMS Export: ',trim(
fields(id)%name), &
627 10
FORMAT (
' OCN2WAV_COUPLING - error while receiving fields from ', &
629 20
FORMAT (
' OCN2WAV_COUPLING - error while sending fields to: ', &
631 30
FORMAT (6x,
'OCN2WAV - (', i2.2,
') imported and (', i2.2, &
632 &
') exported fields,', t62,
't = ', a,/, 16x, &
633 &
'- ROMS coupling exchanges wait clock (s):',/, 19x, &
634 &
'(Recv= ', 1p,e14.8,0p,
' Send= ', 1p,e14.8,0p,
')')
635 40
FORMAT (16x,
'- ',a,a, &
636 & /,19x,
'(Min= ',1p,e15.8,0p,
' Max= ',1p,e15.8,0p,
')')
639 END SUBROUTINE ocn2wav_coupling_tile
641 SUBROUTINE finalize_ocn2wav_coupling
657 CALL router_clean (romstoswan, myerror)
658 CALL attrvect_clean (ocn2wav_av, myerror)
659 CALL globalsegmap_clean (gsmaproms, myerror)
663 END SUBROUTINE finalize_ocn2wav_coupling
type(t_integer), dimension(:), allocatable exportid
integer, dimension(:), allocatable nexport
type(t_field), dimension(maxnumberfields) fields
character(len=240), dimension(:), allocatable exportlist
integer, dimension(:), allocatable nthreads
integer, dimension(:), allocatable nimport
type(t_integer), dimension(:), allocatable importid
type(t_forces), dimension(:), allocatable forces
type(t_grid), dimension(:), allocatable grid
type(t_ocean), dimension(:), allocatable ocean
integer, dimension(:), allocatable n
type(t_bounds), dimension(:), allocatable bounds
type(t_domain), dimension(:), allocatable domain
integer, dimension(:), allocatable lm
real(dp), parameter deg2rad
character(len=22), dimension(:), allocatable time_code
real(r8), dimension(:), allocatable rdrg2
type(t_sedbed), dimension(:), allocatable sedbed
recursive subroutine wclock_off(ng, model, region, line, routine)
recursive subroutine wclock_on(ng, model, region, line, routine)