ROMS
Loading...
Searching...
No Matches
mct_roms_swan.h
Go to the documentation of this file.
1/*
2** git $id$
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************************************************************************
8** **
9** these routines are use couple roms to swan wave model using **
10** the model coupling toolkit(mct). **
11** **
12************************************************************************
13*/
14
15 SUBROUTINE initialize_ocn2wav_coupling (ng, tile)
16!
17!=======================================================================
18! !
19! Initialize ocean and wave models coupling stream. This is the !
20! training phase used to constuct MCT parallel interpolators and !
21! and stablish communication patterns. !
22! !
23!=======================================================================
24!
25 USE mod_param
26 USE mod_parallel
27 USE mod_coupler
28 USE mod_forces
29 USE mod_kinds
30 USE mod_scalars
31!
32! Imported variable definitions.
33!
34 integer, intent(in) :: ng, tile
35!
36! Local variable declarations.
37!
38 integer :: Istr, Iend, Jstr, Jend
39 integer :: IstrR, IendR, JstrR, JendR, IstrU, JstrV
40 integer :: Asize, Jsize, MyError
41 integer :: j, jc, nprocs
42
43 integer, allocatable :: length(:)
44 integer, allocatable :: start(:)
45!
46!-----------------------------------------------------------------------
47! Compute lower and upper bounds over a particular domain partition or
48! tile for RHO-, U-, and V-variables. Notice that "set_bounds.h" is
49! not used here because of implementation of periodicity in other
50! models.
51!-----------------------------------------------------------------------
52!
53 istr=bounds(ng)%Istr(tile)
54 iend=bounds(ng)%Iend(tile)
55 jstr=bounds(ng)%Jstr(tile)
56 jend=bounds(ng)%Jend(tile)
57!
58 IF (domain(ng)%Western_Edge(tile)) THEN
59 istrr=bounds(ng)%Istr(tile)-1
60 ELSE
61 istrr=bounds(ng)%Istr(tile)
62 END IF
63 IF (domain(ng)%Eastern_Edge(tile)) THEN
64 iendr=bounds(ng)%Iend(tile)+1
65 ELSE
66 iendr=bounds(ng)%Iend(tile)
67 END IF
68 IF (domain(ng)%Southern_Edge(tile)) THEN
69 jstrr=bounds(ng)%Jstr(tile)-1
70 ELSE
71 jstrr=bounds(ng)%Jstr(tile)
72 END IF
73 IF (domain(ng)%Northern_Edge(tile)) THEN
74 jendr=bounds(ng)%Jend(tile)+1
75 ELSE
76 jendr=bounds(ng)%Jend(tile)
77 END IF
78!
79!-----------------------------------------------------------------------
80! Begin initialization phase.
81!-----------------------------------------------------------------------
82!
83! Get communicator local rank and size.
84!
85 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
86 CALL mpi_comm_size (ocn_comm_world, nprocs, myerror)
87!
88! Initialize MCT coupled model registry.
89!
90 CALL mctworld_init (nmodels, mpi_comm_world, ocn_comm_world, &
91 & ocnid)
92!
93! Determine start and lengths for domain decomposition.
94!
95 jsize=jendr-jstrr+1
96 IF (.not.allocated(start)) THEN
97 allocate ( start(jsize) )
98 END IF
99 IF (.not.allocated(length)) THEN
100 allocate ( length(jsize) )
101 END IF
102 jc=0
103 DO j=jstrr,jendr
104 jc=jc+1
105 start(jc)=j*(lm(ng)+2)+istrr+1
106 length(jc)=(iendr-istrr+1)
107 END DO
108 CALL globalsegmap_init (gsmaproms, start, length, 0, &
110!
111! Initialize attribute vector holding the export data code strings of
112! the wave model. The Asize is the number of grid point on this
113! processor.
114!
115 asize=globalsegmap_lsize(gsmaproms, ocn_comm_world)
116 CALL attrvect_init (wav2ocn_av, rlist=trim(exportlist(iwaves)), &
117 & lsize=asize)
118!
119! Initialize attribute vector holding the export data code string of
120! the ocean model.
121!
122 CALL attrvect_init (ocn2wav_av, rlist=trim(exportlist(iocean)), &
123 & lsize=asize)
124 CALL attrvect_zero (ocn2wav_av)
125!
126! Initialize a router to the wave model component.
127!
128 CALL router_init (wavid, gsmaproms, ocn_comm_world, romstoswan)
129!
130! Deallocate working arrays.
131!
132 IF (allocated(start)) THEN
133 deallocate (start)
134 END IF
135 IF (allocated(length)) THEN
136 deallocate (length)
137 END IF
138
139 RETURN
140 END SUBROUTINE initialize_ocn2wav_coupling
141
142 SUBROUTINE ocn2wav_coupling (ng, tile)
143!
144!=======================================================================
145! !
146! This routine acquires the coupling data streams between waves !
147! and ocean models. Currently, the following data streams are !
148! coded: !
149! !
150! (...) SWAN units !
151! [...] ROMS units !
152! !
153! Fields imported from SWAN model: !
154! !
155! * Wave direction (degrees), [radians] !
156! * Significant wave height (m), [m] !
157! * Average wave length (m), [m] !
158! * Surface wave relative peak period (s), [s] !
159! * Bottom wave period (s), [s] !
160! * Percent of breakig waves (nondimensional), [nondimensional] !
161! * Wave energy dissipation (W/m2), [m3/s3] !
162! * Wave bottom orbital velocity (m/s), [m/s] !
163! !
164! Fields exported to SWAN model: !
165! !
166! * Bathymetry, bottom elevation (m), [m] !
167! * Free-surface, water surface elevation (m), [m] !
168! * Depth integrated u-momentum (m/s), [m/s] !
169! * Depth integrated v-momentum (m/s), [m/s] !
170! !
171!=======================================================================
172!
173 USE mod_param
174!
175 implicit none
176!
177! Imported variable declarations.
178!
179 integer, intent(in) :: ng, tile
180!
181! Local variable declarations.
182!
183 character (len=*), parameter :: MyFile = &
184 & __FILE__
185
186#include "tile.h"
187!
188#ifdef PROFILE
189 CALL wclock_on (ng, inlm, 41, __line__, myfile)
190#endif
191 CALL ocn2wav_coupling_tile (ng, tile, &
192 & lbi, ubi, lbj, ubj)
193#ifdef PROFILE
194 CALL wclock_off (ng, inlm, 41, __line__, myfile)
195#endif
196
197 RETURN
198 END SUBROUTINE ocn2wav_coupling
199!
200!***********************************************************************
201 SUBROUTINE ocn2wav_coupling_tile (ng, tile, &
202 & LBi, UBi, LBj, UBj)
203!***********************************************************************
204!
205 USE mod_param
206 USE mod_parallel
207 USE mod_coupler
208 USE mod_forces
209 USE mod_grid
210 USE mod_ocean
211 USE mod_scalars
212 USE mod_stepping
213 USE mod_iounits
214#ifdef BBL_MODEL
215 USE mod_sedbed
216#endif
217 USE mod_sediment
218!
219 USE distribute_mod, ONLY : mp_reduce
220 USE roms_import_mod, ONLY : roms_import2d
221 USE roms_export_mod, ONLY : roms_export2d
222!
223 implicit none
224!
225! Imported variable declarations.
226!
227 integer, intent(in) :: ng, tile
228 integer, intent(in) :: LBi, UBi, LBj, UBj
229!
230! Local variable declarations.
231!
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
236
237 real(r8), parameter :: Lwave_max=500.0_r8
238
239 real(r8) :: add_offset, scale
240 real(r8) :: RecvTime, SendTime, buffer(2), wtime(2)
241
242 real(r8) :: my_wtime
243
244 real(r8), dimension(LBi:UBi,LBj:UBj) :: Awrk
245
246 real(r8), pointer :: A(:)
247
248 character (len=3 ), dimension(2) :: op_handle
249 character (len=40) :: code
250!
251!-----------------------------------------------------------------------
252! Compute lower and upper bounds over a particular domain partition or
253! tile for RHO-, U-, and V-variables. Notice that "set_bounds.h" is
254! not used here because of implementation of periodicity in other
255! models.
256!-----------------------------------------------------------------------
257!
258 istr=bounds(ng)%Istr(tile)
259 iend=bounds(ng)%Iend(tile)
260 jstr=bounds(ng)%Jstr(tile)
261 jend=bounds(ng)%Jend(tile)
262!
263 IF (domain(ng)%Western_Edge(tile)) THEN
264 istrr=bounds(ng)%Istr(tile)-1
265 ELSE
266 istrr=bounds(ng)%Istr(tile)
267 END IF
268 IF (domain(ng)%Eastern_Edge(tile)) THEN
269 iendr=bounds(ng)%Iend(tile)+1
270 ELSE
271 iendr=bounds(ng)%Iend(tile)
272 END IF
273 IF (domain(ng)%Southern_Edge(tile)) THEN
274 jstrr=bounds(ng)%Jstr(tile)-1
275 ELSE
276 jstrr=bounds(ng)%Jstr(tile)
277 END IF
278 IF (domain(ng)%Northern_Edge(tile)) THEN
279 jendr=bounds(ng)%Jend(tile)+1
280 ELSE
281 jendr=bounds(ng)%Jend(tile)
282 END IF
283!
284!-----------------------------------------------------------------------
285! Allocate communications array.
286!-----------------------------------------------------------------------
287!
288 asize=globalsegmap_lsize(gsmaproms, ocn_comm_world)
289 allocate ( a(asize) )
290 a=0.0_r8
291!
292! Initialize coupling wait time clocks.
293!
294 recvtime=0.0_r8
295 sendtime=0.0_r8
296!
297!-----------------------------------------------------------------------
298! Import fields from wave model (SWAN) to ocean model (ROMS).
299! Currently, both waves and ocean model grids are the same.
300! We need to revisit this logic to allow interpolation.
301!-----------------------------------------------------------------------
302!
303! Schedule receiving fields from wave model.
304!
305 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
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
310 IF (master) THEN
311 WRITE (stdout,10) 'wave model, MyError = ', myerror
312 END IF
313 exit_flag=2
314 RETURN
315 END IF
316!
317! Receive fields from wave model.
318!
319 iimport=0
320 DO ifield=1,nimport(iocean)
321 id=importid(iocean)%val(ifield)
322 code=adjustl(fields(id)%code)
323 gtype=fields(id)%GridType
324 scale=fields(id)%scale
325 add_offset=fields(id)%AddOffset
326
327 SELECT CASE (trim(code))
328
329 CASE ('Wdir') ! wave direction
330
331 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
332 iimport=iimport+1
333 DO i=1,asize
334 a(i)=max(0.0_r8,a(i))
335 END DO
336 scale=deg2rad ! degress to radians
337 add_offset=0.0_r8
338 CALL roms_import2d (ng, tile, &
339 & id, gtype, scale, add_offset, &
340 & asize, a, &
341 & istrr, iendr, jstrr, jendr, &
342 & lbi, ubi, lbj, ubj, &
343 & fields(id)%ImpMin, fields(id)%ImpMax, &
344 & forces(ng)%Dwave, &
345 & status)
346
347 CASE ('Wamp') ! significant wave hight
348
349 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
350 iimport=iimport+1
351 DO i=1,asize
352 a(i)=max(0.0_r8,a(i))
353 END DO
354 scale=1.0_r8 ! m
355 add_offset=0.0_r8
356 CALL roms_import2d (ng, tile, &
357 & id, gtype, scale, add_offset, &
358 & asize, a, &
359 & istrr, iendr, jstrr, jendr, &
360 & lbi, ubi, lbj, ubj, &
361 & fields(id)%ImpMin, fields(id)%ImpMax, &
362 & forces(ng)%Hwave, &
363 & status)
364
365 CASE ('Wlen') ! wave length
366
367 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
368 iimport=iimport+1
369 DO i=1,asize
370 a(i)=max(0.0_r8,a(i))
371 IF (a(i).eq.infinity) THEN
372 a(i)=lwave_max
373 END IF
374 END DO
375 scale=1.0_r8 ! m
376 add_offset=0.0_r8
377 CALL roms_import2d (ng, tile, &
378 & id, gtype, scale, add_offset, &
379 & asize, a, &
380 & istrr, iendr, jstrr, jendr, &
381 & lbi, ubi, lbj, ubj, &
382 & fields(id)%ImpMin, fields(id)%ImpMax, &
383 & forces(ng)%Lwave, &
384 & status)
385
386 CASE ('Wptop') ! peak surface wave period
387
388 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
389 iimport=iimport+1
390 DO i=1,asize
391 a(i)=max(0.0_r8,a(i))
392 END DO
393 scale=1.0_r8
394 add_offset=0.0_r8
395 CALL roms_import2d (ng, tile, &
396 & id, gtype, scale, add_offset, &
397 & asize, a, &
398 & istrr, iendr, jstrr, jendr, &
399 & lbi, ubi, lbj, ubj, &
400 & fields(id)%ImpMin, fields(id)%ImpMax, &
401 & forces(ng)%Pwave_top, &
402 & status)
403
404 CASE ('Wpbot') ! mean bottom wave period
405
406 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
407 iimport=iimport+1
408 DO i=1,asize
409 a(i)=max(0.0_r8,a(i))
410 END DO
411 scale=1.0_r8
412 add_offset=0.0_r8
413 CALL roms_import2d (ng, tile, &
414 & id, gtype, scale, add_offset, &
415 & asize, a, &
416 & istrr, iendr, jstrr, jendr, &
417 & lbi, ubi, lbj, ubj, &
418 & fields(id)%ImpMin, fields(id)%ImpMax, &
419 & forces(ng)%Pwave_bot, &
420 & status)
421
422 CASE ('Wdiss') ! wave dissipation
423
424 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
425 iimport=iimport+1
426 DO i=1,asize
427 a(i)=max(0.0_r8,a(i))
428 END DO
429 scale=1.0_r8/rho0
430 add_offset=0.0_r8
431 CALL roms_import2d (ng, tile, &
432 & id, gtype, scale, add_offset, &
433 & asize, a, &
434 & istrr, iendr, jstrr, jendr, &
435 & lbi, ubi, lbj, ubj, &
436 & fields(id)%ImpMin, fields(id)%ImpMax, &
437 & forces(ng)%Wave_dissip, &
438 & status)
439
440 CASE ('Wubot') ! bottom orbital velocity
441
442 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
443 iimport=iimport+1
444 DO i=1,asize
445 a(i)=max(0.0_r8,a(i))
446 END DO
447 scale=1.0_r8 ! m/s
448 add_offset=0.0_r8
449 CALL roms_import2d (ng, tile, &
450 & id, gtype, scale, add_offset, &
451 & asize, a, &
452 & istrr, iendr, jstrr, jendr, &
453 & lbi, ubi, lbj, ubj, &
454 & fields(id)%ImpMin, fields(id)%ImpMax, &
455 & forces(ng)%Ub_swan, &
456 & status)
457
458#ifdef SVENDSEN_ROLLER
459
460 CASE ('Wbrk') ! percent wave breaking
461
462 CALL attrvect_exportrattr (wav2ocn_av, trim(code), a, asize)
463 iimport=iimport+1
464 DO i=1,asize
465 a(i)=max(0.0_r8,a(i))
466 END DO
467 scale=1.0_r8
468 add_offset=0.0_r8
469 CALL roms_import2d (ng, tile, &
470 & id, gtype, scale, add_offset, &
471 & asize, a, &
472 & istrr, iendr, jstrr, jendr, &
473 & lbi, ubi, lbj, ubj, &
474 & fields(id)%ImpMin, fields(id)%ImpMax, &
475 & forces(ng)%Wave_break, &
476 & status)
477#endif
478 END SELECT
479 END DO
480!
481!-----------------------------------------------------------------------
482! Export fields from ocean (ROMS) to wave (SWAN) model.
483!-----------------------------------------------------------------------
484!
485! Schedule sending fields to the wave model.
486!
487 iexport=0
488 DO ifield=1,nexport(iocean)
489 id=exportid(iocean)%val(ifield)
490 code=adjustl(fields(id)%code)
491 gtype=fields(id)%GridType
492 scale=fields(id)%scale
493 add_offset=fields(id)%AddOffset
494
495 SELECT CASE (trim(code))
496
497 CASE ('bath') ! bathymetry (depth)
498
499 CALL roms_export2d (ng, tile, &
500 & id, gtype, scale, add_offset, &
501 & lbi, ubi, lbj, ubj, &
502 & grid(ng)%h, &
503 & fields(id)%ExpMin, fields(id)%ExpMax, &
504 & asize, a, &
505 & status)
506 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
507 iexport=iexport+1
508
509 CASE ('SSH') ! free-surface (water level)
510
511 CALL roms_export2d (ng, tile, &
512 & id, gtype, scale, add_offset, &
513 & lbi, ubi, lbj, ubj, &
514 & ocean(ng)%zeta(:,:,kout), &
515 & fields(id)%ExpMin, fields(id)%ExpMax, &
516 & asize, a, &
517 & status)
518 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
519 iexport=iexport+1
520
521 CASE ('Ubar') ! 2D U-momentum
522
523 CALL roms_export2d (ng, tile, &
524 & id, gtype, scale, add_offset, &
525 & lbi, ubi, lbj, ubj, &
526#ifdef SOLVE3D
527 & ocean(ng)%u(:,:,n(ng),nout), &
528#else
529 & ocean(ng)%ubar(:,:,kout), &
530#endif
531 & fields(id)%ExpMin, fields(id)%ExpMax, &
532 & asize, a, &
533 & status)
534 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
535 iexport=iexport+1
536
537 CASE ('Vbar') ! 2D V-momentum
538
539 CALL roms_export2d (ng, tile, &
540 & id, gtype, scale, add_offset, &
541 & lbi, ubi, lbj, ubj, &
542#ifdef SOLVE3D
543 & ocean(ng)%v(:,:,n(ng),nout), &
544#else
545 & ocean(ng)%vbar(:,:,kout), &
546#endif
547 & fields(id)%ExpMin, fields(id)%ExpMax, &
548 & asize, a, &
549 & status)
550 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
551 iexport=iexport+1
552
553 CASE ('ZO') ! bottom roughness
554
555 DO j=jstrr,jendr
556 DO i=istrr,iendr
557#ifdef BBL_MODEL
558 awrk(i,j)=max(0.0001_r8, &
559 & sedbed(ng)%bottom(i,j,iznik)*30.0_r8)
560#else
561 awrk(i,j)=max(0.0001_r8,rdrg2(ng))
562#endif
563 END DO
564 END DO
565 CALL roms_export2d (ng, tile, &
566 & id, gtype, scale, add_offset, &
567 & lbi, ubi, lbj, ubj, &
568 & awrk, &
569 & fields(id)%ExpMin, fields(id)%ExpMax, &
570 & asize, a, &
571 & status)
572 CALL attrvect_importrattr (ocn2wav_av, trim(code), a, asize)
573 iexport=iexport+1
574
575 END SELECT
576 END DO
577!
578! Send ocean fields to wave model.
579!
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
585 IF (master) THEN
586 WRITE (stdout,20) 'wave model, MyError = ', myerror
587 END IF
588 exit_flag=2
589 RETURN
590 END IF
591 END IF
592!
593!-----------------------------------------------------------------------
594! Report.
595!-----------------------------------------------------------------------
596!
597 IF (nthreads(iocean).gt.1) THEN
598 buffer(1)=recvtime
599 buffer(2)=sendtime
600 op_handle(1)='SUM'
601 op_handle(2)='SUM'
602 CALL mp_reduce (ng, inlm, 2, buffer, op_handle)
603 recvtime=buffer(1)
604 sendtime=buffer(2)
605 END IF
606 IF (master.and.((iimport.gt.0).or.(iexport.gt.0))) THEN
607 WRITE (stdout,30) iimport, iexport, time_code(ng), &
608 & recvtime, sendtime
609 IF (lreport) THEN
610 DO ifield=1,nimport(iocean)
611 id=importid(iocean)%val(ifield)
612 WRITE (stdout,40) 'ROMS Import: ',trim(fields(id)%name), &
613 & fields(id)%ImpMin, fields(id)%ImpMax
614 END DO
615 DO ifield=1,nexport(iocean)
616 id=exportid(iocean)%val(ifield)
617 WRITE (stdout,40) 'ROMS Export: ',trim(fields(id)%name), &
618 & fields(id)%ExpMin, fields(id)%ExpMax
619 END DO
620 END IF
621 END IF
622!
623! Deallocate communication arrays.
624!
625 deallocate (a)
626!
627 10 FORMAT (' OCN2WAV_COUPLING - error while receiving fields from ', &
628 & a, i4)
629 20 FORMAT (' OCN2WAV_COUPLING - error while sending fields to: ', &
630 & a, i4)
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,')')
637
638 RETURN
639 END SUBROUTINE ocn2wav_coupling_tile
640
641 SUBROUTINE finalize_ocn2wav_coupling
642!
643!========================================================================
644! !
645! This routine finalizes ocean and wave models coupling data streams. !
646! !
647!========================================================================
648!
649! Local variable declarations.
650!
651 integer :: MyError
652!
653!-----------------------------------------------------------------------
654! Deallocate MCT environment.
655!-----------------------------------------------------------------------
656!
657 CALL router_clean (romstoswan, myerror)
658 CALL attrvect_clean (ocn2wav_av, myerror)
659 CALL globalsegmap_clean (gsmaproms, myerror)
660
661 RETURN
662
663 END SUBROUTINE finalize_ocn2wav_coupling
integer iocean
integer wavid
Definition mod_coupler.F:89
type(t_integer), dimension(:), allocatable exportid
logical lreport
Definition mod_coupler.F:94
integer, dimension(:), allocatable nexport
type(t_field), dimension(maxnumberfields) fields
character(len=240), dimension(:), allocatable exportlist
integer, dimension(:), allocatable nthreads
integer ocnid
Definition mod_coupler.F:90
integer, dimension(:), allocatable nimport
integer nmodels
Definition mod_coupler.F:84
type(t_integer), dimension(:), allocatable importid
integer iwaves
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer stdout
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer ocn_comm_world
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable lm
Definition mod_param.F:455
real(dp), parameter deg2rad
character(len=22), dimension(:), allocatable time_code
integer exit_flag
real(r8), dimension(:), allocatable rdrg2
real(dp) infinity
real(dp) rho0
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, parameter iznik
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