ROMS
Loading...
Searching...
No Matches
mct_roms_wrf.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 wrf atmosphere model **
10** using the model coupling toolkit(mct). **
11** **
12************************************************************************
13*/
14
15 SUBROUTINE initialize_ocn2atm_coupling (ng, tile)
16!
17!=======================================================================
18! !
19! Initialize ocean and atmosphere models coupling stream. This is !
20! the training phase used to constuct MCT parallel interpolators !
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 atmosphere 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 (atm2ocn_av, rlist=trim(exportlist(iatmos)), &
117 & lsize=asize)
118!
119! Initialize attribute vector holding the export data code string of
120! the ocean model.
121!
122 CALL attrvect_init (ocn2atm_av, rlist=trim(exportlist(iocean)), &
123 & lsize=asize)
124 CALL attrvect_zero (ocn2atm_av)
125!
126! Initialize a router to the atmosphere model component.
127!
128 CALL router_init (atmid, gsmaproms, ocn_comm_world, romstowrf)
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_ocn2atm_coupling
141
142 SUBROUTINE ocn2atm_coupling (ng, tile)
143!
144!=======================================================================
145! !
146! This subroutine acquires the coupling data streams between ocean !
147! and atmosphere models. Currently, the following data streams are !
148! coded: !
149! !
150! (...) WRF units !
151! [...] ROMS units !
152! !
153! Fields imported WRF model: !
154! !
155! * Surface atmospheric pressure (Pa), [mb] !
156! * Surface air relative humidity (percent), [fraction] !
157! * Surface (2 m) air temperature (Celsius), [Celsius] !
158! * Surface (10 m) U-wind speed (m/s), [m/s] !
159! * Surface (10 m) V-wind speed (m/s), [m/s] !
160! * Cloud fraction (percent/100), [percent/100] !
161! * Precipitation (m/s), [kg/m2/s] !
162! * Shortwave radiation (Watts/m2), [Celsius m/s] !
163! * Long wave raditaion (Watts/m2), [Celsius m/s] !
164! * Latent heat flux (Watts/m2), [Celsius m/s] !
165! * Sensible heat flux (Watts/m2), [Celsius m/s] !
166! * Net surface heat flux (Watts/2), [Celsius m/s] !
167! * Surface U-wind stress (Pa), [m2/s2] !
168! * Surface V-wind stress (Pa), [m2/s2] !
169! !
170! Fields exported to WRF model: !
171! !
172! * Sea surface potential temperature (Celsius), [Celsius] !
173! !
174!=======================================================================
175!
176 USE mod_param
177!
178 implicit none
179!
180! Imported variable declarations.
181!
182 integer, intent(in) :: ng, tile
183!
184! Local variable declarations.
185!
186 character (len=*), parameter :: MyFile = &
187 & __FILE__
188
189#include "tile.h"
190!
191#ifdef PROFILE
192 CALL wclock_on (ng, inlm, 39, __line__, myfile)
193#endif
194 CALL ocn2atm_coupling_tile (ng, tile, &
195 & lbi, ubi, lbj, ubj)
196#ifdef PROFILE
197 CALL wclock_off (ng, inlm, 39, __line__, myfile)
198#endif
199
200 RETURN
201 END SUBROUTINE ocn2atm_coupling
202!
203!***********************************************************************
204 SUBROUTINE ocn2atm_coupling_tile (ng, tile, &
205 & LBi, UBi, LBj, UBj)
206!***********************************************************************
207!
208 USE mod_param
209 USE mod_parallel
210 USE mod_coupler
211 USE mod_forces
212 USE mod_ocean
213 USE mod_scalars
214 USE mod_stepping
215 USE mod_iounits
216!
217 USE distribute_mod, ONLY : mp_reduce
218 USE roms_import_mod, ONLY : roms_import2d
219 USE roms_export_mod, ONLY : roms_export2d
220!
221 implicit none
222!
223! Imported variable declarations.
224!
225 integer, intent(in) :: ng, tile
226 integer, intent(in) :: LBi, UBi, LBj, UBj
227!
228! Local variable declarations.
229!
230 integer :: Istr, Iend, Jstr, Jend
231 integer :: IstrR, IendR, JstrR, JendR, IstrU, JstrV
232 integer :: Asize, Iimport, Iexport, MyError
233 integer :: gtype, i, id, ifield, j, status
234
235 real(r8) :: add_offset, scale
236 real(r8) :: RecvTime, SendTime, buffer(2), wtime(2)
237
238 real(r8) :: my_wtime
239
240 real(r8), pointer :: A(:)
241
242 character (len=3 ), dimension(2) :: op_handle
243 character (len=40) :: code
244!
245!-----------------------------------------------------------------------
246! Compute lower and upper bounds over a particular domain partition or
247! tile for RHO-, U-, and V-variables. Notice that "set_bounds.h" is
248! not used here because of implementation of periodicity in other
249! models.
250!-----------------------------------------------------------------------
251!
252 istr=bounds(ng)%Istr(tile)
253 iend=bounds(ng)%Iend(tile)
254 jstr=bounds(ng)%Jstr(tile)
255 jend=bounds(ng)%Jend(tile)
256!
257 IF (domain(ng)%Western_Edge(tile)) THEN
258 istrr=bounds(ng)%Istr(tile)-1
259 ELSE
260 istrr=bounds(ng)%Istr(tile)
261 END IF
262 IF (domain(ng)%Eastern_Edge(tile)) THEN
263 iendr=bounds(ng)%Iend(tile)+1
264 ELSE
265 iendr=bounds(ng)%Iend(tile)
266 END IF
267 IF (domain(ng)%Southern_Edge(tile)) THEN
268 jstrr=bounds(ng)%Jstr(tile)-1
269 ELSE
270 jstrr=bounds(ng)%Jstr(tile)
271 END IF
272 IF (domain(ng)%Northern_Edge(tile)) THEN
273 jendr=bounds(ng)%Jend(tile)+1
274 ELSE
275 jendr=bounds(ng)%Jend(tile)
276 END IF
277!
278!-----------------------------------------------------------------------
279! Allocate communications array.
280!-----------------------------------------------------------------------
281!
282 asize=globalsegmap_lsize(gsmaproms, ocn_comm_world)
283 allocate ( a(asize) )
284 a=0.0_r8
285!
286! Initialize coupling wait time clocks.
287!
288 recvtime=0.0_r8
289 sendtime=0.0_r8
290!
291!-----------------------------------------------------------------------
292! Import fields from atmosphere model (WRF) to ocean model (ROMS).
293! Currently, both atmosphere and ocean model grids are the same.
294! We need to revisit this logic to allow interpolation.
295!-----------------------------------------------------------------------
296!
297! Schedule receiving fields from atmosphere model.
298!
299 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
300 buffer(1)=my_wtime(wtime)
301 CALL mct_recv (atm2ocn_av, romstowrf, myerror)
302 recvtime=recvtime+my_wtime(wtime)-buffer(1)
303 IF (myerror.ne.0) THEN
304 IF (master) THEN
305 WRITE (stdout,10) 'atmosphere model, MyError = ', myerror
306 END IF
307 exit_flag=2
308 RETURN
309 END IF
310!
311! Receive fields from atmosphere model.
312!
313 iimport=0
314 DO ifield=1,nimport(iocean)
315 id=importid(iocean)%val(ifield)
316 code=adjustl(fields(id)%code)
317 gtype=fields(id)%GridType
318 scale=fields(id)%scale
319 add_offset=fields(id)%AddOffset
320
321 SELECT CASE (trim(code))
322
323#if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
324
325 CASE ('Pair') ! surface air pressure
326
327 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
328 iimport=iimport+1
329 scale=0.01_r8 ! Pa to mb
330 add_offset=0.0_r8
331 CALL roms_import2d (ng, tile, &
332 & id, gtype, scale, add_offset, &
333 & asize, a, &
334 & istrr, iendr, jstrr, jendr, &
335 & lbi, ubi, lbj, ubj, &
336 & fields(id)%ImpMin, fields(id)%ImpMax, &
337 & forces(ng)%Pair, &
338 & status)
339#endif
340#if defined BULK_FLUXES || defined ECOSIM || \
341 (defined shortwave && defined ana_srflux)
342
343 CASE ('Hair') ! surface air humidity
344
345 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
346 iimport=iimport+1
347 scale=0.01_r8 ! percent to fraction
348 add_offset=0.0_r8
349 CALL roms_import2d (ng, tile, &
350 & id, gtype, scale, add_offset, &
351 & asize, a, &
352 & istrr, iendr, jstrr, jendr, &
353 & lbi, ubi, lbj, ubj, &
354 & fields(id)%ImpMin, fields(id)%ImpMax, &
355 & forces(ng)%Hair, &
356 & status)
357
358 CASE ('Tair') ! surface (2m) air temperature
359
360 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
361 iimport=iimport+1
362 scale=1.0_r8 ! Celsius
363 add_offset=0.0_r8
364 CALL roms_import2d (ng, tile, &
365 & id, gtype, scale, add_offset, &
366 & asize, a, &
367 & istrr, iendr, jstrr, jendr, &
368 & lbi, ubi, lbj, ubj, &
369 & fields(id)%ImpMin, fields(id)%ImpMax, &
370 & forces(ng)%Tair, &
371 & status)
372#endif
373#if defined BULK_FLUXES || defined ECOSIM
374
375 CASE ('Uwind') ! U-wind (10m) component
376
377 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
378 iimport=iimport+1
379 scale=1.0_r8 ! m/s
380 add_offset=0.0_r8
381 CALL roms_import2d (ng, tile, &
382 & id, gtype, scale, add_offset, &
383 & asize, a, &
384 & istrr, iendr, jstrr, jendr, &
385 & lbi, ubi, lbj, ubj, &
386 & fields(id)%ImpMin, fields(id)%ImpMax, &
387 & forces(ng)%Uwind, &
388 & status)
389
390 CASE ('Vwind') ! V-wind (10m) component
391
392 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
393 iimport=iimport+1
394 scale=1.0_r8 ! m/s
395 add_offset=0.0_r8
396 CALL roms_import2d (ng, tile, &
397 & id, gtype, scale, add_offset, &
398 & asize, a, &
399 & istrr, iendr, jstrr, jendr, &
400 & lbi, ubi, lbj, ubj, &
401 & fields(id)%ImpMin, fields(id)%ImpMax, &
402 & forces(ng)%Vwind, &
403 & status)
404#endif
405#ifdef CLOUDS
406
407 CASE ('cloud') ! cloud fraction
408
409 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
410 iimport=iimport+1
411 scale=1.0_r8 ! percent/100, so 0 to 1
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)%cloud, &
420 & status)
421#endif
422#ifdef BULK_FLUXES
423
424 CASE ('rain') ! precipitation
425
426 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
427 iimport=iimport+1
428 scale=rho0 ! kg/m2/s
429 add_offset=0.0_r8
430 CALL roms_import2d (ng, tile, &
431 & id, gtype, scale, add_offset, &
432 & asize, a, &
433 & istrr, iendr, jstrr, jendr, &
434 & lbi, ubi, lbj, ubj, &
435 & fields(id)%ImpMin, fields(id)%ImpMax, &
436 & forces(ng)%rain, &
437 & status)
438
439 CASE ('LWrad') ! longwave radiation
440
441 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
442 iimport=iimport+1
443 scale=-1.0_r8/(rho0*cp) ! Watts/m2 to Celsius m/s
444 add_offset=0.0_r8
445 CALL roms_import2d (ng, tile, &
446 & id, gtype, scale, add_offset, &
447 & asize, a, &
448 & istrr, iendr, jstrr, jendr, &
449 & lbi, ubi, lbj, ubj, &
450 & fields(id)%ImpMin, fields(id)%ImpMax, &
451 & forces(ng)%lrflx, &
452 & status)
453
454 CASE ('Lheat') ! latent heat flux
455
456 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
457 iimport=iimport+1
458 scale=1.0_r8 ! Watts/m2
459 add_offset=0.0_r8
460 CALL roms_import2d (ng, tile, &
461 & id, gtype, scale, add_offset, &
462 & asize, a, &
463 & istrr, iendr, jstrr, jendr, &
464 & lbi, ubi, lbj, ubj, &
465 & fields(id)%ImpMin, fields(id)%ImpMax, &
466 & forces(ng)%lhflx, &
467 & status)
468
469 CASE ('Sheat') ! sensible heat flux
470
471 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
472 iimport=iimport+1
473 scale=1.0_r8 ! Watts/m2
474 add_offset=0.0_r8
475 CALL roms_import2d (ng, tile, &
476 & id, gtype, scale, add_offset, &
477 & asize, a, &
478 & istrr, iendr, jstrr, jendr, &
479 & lbi, ubi, lbj, ubj, &
480 & fields(id)%ImpMin, fields(id)%ImpMax, &
481 & forces(ng)%shflx, &
482 & status)
483#endif
484#ifdef SHORTWAVE
485
486 CASE ('SWrad') ! shortwave radiation
487
488 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
489 iimport=iimport+1
490 scale=-1.0_r8/(rho0*cp) ! Watts/m2 to Celsius m/s
491 add_offset=0.0_r8
492 CALL roms_import2d (ng, tile, &
493 & id, gtype, scale, add_offset, &
494 & asize, a, &
495 & istrr, iendr, jstrr, jendr, &
496 & lbi, ubi, lbj, ubj, &
497 & fields(id)%ImpMin, fields(id)%ImpMax, &
498 & forces(ng)%srflx, &
499 & status)
500#endif
501#ifndef BULK_FLUXES
502
503 CASE ('heat') ! surface net heat flux
504
505 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
506 iimport=iimport+1
507 scale=-1.0_r8/(rho0*cp) ! Watts/m2 to Celsius m/s
508 add_offset=0.0_r8
509 CALL roms_import2d (ng, myrank, &
510 & id, gtype, scale, add_offset, &
511 & asize, a &
512 & istrr, iendr, jstrr, jendr, &
513 & lbi, ubi, lbj, ubj, &
514 & fields(id)%ImpMin, fields(id)%ImpMax, &
515 & forces(ng)%stflx(:,:,itemp), &
516 & status)
517
518 CASE ('Ustr') ! surface U-wind stress
519
520 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
521 iimport=iimport+1
522 scale=1.0_r8/rho0 ! Pa to m2/s2
523 add_offset=0.0_r8
524 CALL roms_import2d (ng, tile, &
525 & id, gtype, scale, add_offset, &
526 & asize, a, &
527 & istrr, iendr, jstrr, jendr, &
528 & lbi, ubi, lbj, ubj, &
529 & fields(id)%ImpMin, fields(id)%ImpMax, &
530 & forces(ng)%sustr, &
531 & status)
532
533 CASE ('Vstr') ! surface V-wind stress
534
535 CALL attrvect_exportrattr (atm2ocn_av, trim(code), a, asize)
536 iimport=iimport+1
537 scale=1.0_r8/rho0 ! Pa to m2/s2
538 add_offset=0.0_r8
539 CALL roms_import2d (ng, tile, &
540 & id, gtype, scale, add_offset, &
541 & asize, a, &
542 & istrr, iendr, jstrr, jendr, &
543 & lbi, ubi, lbj, ubj, &
544 & fields(id)%ImpMin, fields(id)%ImpMax, &
545 & forces(ng)%svstr, &
546 & status)
547#endif
548 END SELECT
549 END DO
550!
551!-----------------------------------------------------------------------
552! Export fields from ocean (ROMS) to atmosphere (WRF) model.
553!-----------------------------------------------------------------------
554!
555! Schedule sending fields to the atmosphere model.
556!
557 iexport=0
558 DO ifield=1,nexport(iocean)
559 id=exportid(iocean)%val(ifield)
560 code=adjustl(fields(id)%code)
561 gtype=fields(id)%GridType
562 scale=fields(id)%scale
563 add_offset=fields(id)%AddOffset
564
565 SELECT CASE (trim(code))
566
567 CASE ('SST') ! sea surface temperature
568
569 CALL roms_export2d (ng, tile, &
570 & id, gtype, scale, add_offset, &
571 & lbi, ubi, lbj, ubj, &
572 & ocean(ng)%t(:,:,n(ng),nout,itemp), &
573 & fields(id)%ExpMin, fields(id)%ExpMax, &
574 & asize, a, &
575 & status)
576 CALL attrvect_importrattr (ocn2atm_av, trim(code), a, asize)
577 iexport=iexport+1
578
579 END SELECT
580 END DO
581!
582! Send ocean fields to atmosphere model.
583!
584 IF (iexport.gt.0) THEN
585 buffer(2)=my_wtime(wtime)
586 CALL mct_send (ocn2atm_av, romstowrf, myerror)
587 sendtime=sendtime+my_wtime(wtime)-buffer(2)
588 IF (myerror.ne.0) THEN
589 IF (master) THEN
590 WRITE (stdout,20) 'atmosphere model, MyError = ', myerror
591 END IF
592 exit_flag=2
593 RETURN
594 END IF
595 END IF
596!
597!-----------------------------------------------------------------------
598! Report.
599!-----------------------------------------------------------------------
600!
601 IF (nthreads(iocean).gt.1) THEN
602 buffer(1)=recvtime
603 buffer(2)=sendtime
604 op_handle(1)='SUM'
605 op_handle(2)='SUM'
606 CALL mp_reduce (ng, inlm, 2, buffer, op_handle)
607 recvtime=buffer(1)
608 sendtime=buffer(2)
609 END IF
610 IF (master.and.((iimport.gt.0).or.(iexport.gt.0))) THEN
611 WRITE (stdout,30) iimport, iexport, time_code(ng), &
612 & recvtime, sendtime
613 IF (lreport) THEN
614 DO ifield=1,nimport(iocean)
615 id=importid(iocean)%val(ifield)
616 WRITE (stdout,40) 'ROMS Import: ',trim(fields(id)%name), &
617 & fields(id)%ImpMin, fields(id)%ImpMax
618 END DO
619 DO ifield=1,nexport(iocean)
620 id=exportid(iocean)%val(ifield)
621 WRITE (stdout,40) 'ROMS Export: ',trim(fields(id)%name), &
622 & fields(id)%ExpMin, fields(id)%ExpMax
623 END DO
624 END IF
625 END IF
626!
627! Deallocate communication arrays.
628!
629 deallocate (a)
630!
631 10 FORMAT (' OCN2ATM_COUPLING - error while receiving fields from ', &
632 & a, i4)
633 20 FORMAT (' OCN2ATM_COUPLING - error while sending fields to ', &
634 & a, i4)
635 30 FORMAT (6x,'OCN2ATM - (', i2.2, ') imported and (', i2.2, &
636 & ') exported fields,', t62, 't = ', a,/, 16x, &
637 & '- ROMS coupling exchages wait clock (s):',/, 19x, &
638 & '(Recv= ', 1p,e14.8,0p, ' Send= ', 1p,e14.8,0p,')')
639 40 FORMAT (16x,'- ',a,a, &
640 & /,19x,'(Min= ',1p,e15.8,0p,' Max= ',1p,e15.8,0p,')')
641
642 RETURN
643 END SUBROUTINE ocn2atm_coupling_tile
644
645 SUBROUTINE finalize_ocn2atm_coupling
646!
647!========================================================================
648! !
649! This routine finalizes ocean and atmosphere models coupling data !
650! streams. !
651! !
652!========================================================================
653!
654! Local variable declarations.
655!
656 integer :: MyError
657!
658!-----------------------------------------------------------------------
659! Deallocate MCT environment.
660!-----------------------------------------------------------------------
661!
662 CALL router_clean (romstowrf, myerror)
663 CALL attrvect_clean (ocn2atm_av, myerror)
664 CALL globalsegmap_clean (gsmaproms, myerror)
665
666 RETURN
667
668 END SUBROUTINE finalize_ocn2atm_coupling
integer iocean
integer iatmos
type(t_integer), dimension(:), allocatable exportid
logical lreport
Definition mod_coupler.F:94
integer, dimension(:), allocatable nexport
type(t_field), dimension(maxnumberfields) fields
integer atmid
Definition mod_coupler.F:88
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
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
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) cp
character(len=22), dimension(:), allocatable time_code
integer exit_flag
integer itemp
real(dp) rho0
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