ROMS
Loading...
Searching...
No Matches
wrt_state.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined ADJOINT && defined TANGENT
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This module writes requested tangent linear or adjoint state !
13! fields using the standard NetCDF library or the Parallel-IO (PIO) !
14! library into an output file described by IO structure S. The NetCDF !
15! is opened and closed, and file ID or descriptor is not updated. !
16! !
17! On Input: !
18! !
19! ng Nested grid number (integer) !
20! model Calling model identifier (integer) !
21! label Identification label (string) !
22! OutRec NetCDF file time record (integer) !
23! i2d 2D state variables time level index (integer) !
24! i3d 3D state variables time level index (integer) !
25! stime Time variable value (real; seconds) !
26! S File derived type structure (TYPE T_IO) !
27! CloseFile Switch to close NetCDF file (logical, OPTIONAL) !
28! (defaut = .TRUE.) !
29! On Output: !
30! !
31! ad_arrays ADM state arrays if model = iADM !
32! or, !
33! tl_arrays TLM state arrays if model = iTLM or iRPM !
34! !
35!=======================================================================
36!
37 USE mod_param
38 USE mod_parallel
39# ifdef ADJUST_BOUNDARY
40 USE mod_boundary
41# endif
42# ifdef SOLVE3D
43 USE mod_coupling
44# endif
45# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
46 USE mod_forces
47# endif
48 USE mod_grid
49 USE mod_iounits
50# ifdef SOLVE3D
51 USE mod_mixing
52# endif
53 USE mod_ncparam
54 USE mod_ocean
55 USE mod_scalars
56 USE mod_stepping
57!
58 USE dateclock_mod, ONLY : time_string
60# ifdef ADJUST_BOUNDARY
62# endif
63# ifdef SOLVE3D
65# ifdef ADJUST_BOUNDARY
67# endif
68# endif
69 USE strings_mod, ONLY : founderror
70!
71 implicit none
72!
73 PUBLIC :: wrt_state
74 PRIVATE :: wrt_state_nf90
75# if defined PIO_LIB && defined DISTRIBUTE
76 PRIVATE :: wrt_state_pio
77# endif
78!
79 CONTAINS
80!
81!***********************************************************************
82 SUBROUTINE wrt_state (ng, tile, model, label, &
83 & OutRec, i2d, i3d, stime, S)
84!***********************************************************************
85!
86! Imported variable declarations.
87!
88 integer, intent(in) :: ng, tile, model, outrec, i2d, i3d
89!
90 real(dp), intent(in) :: stime
91!
92 character (len=*), intent(in) :: label
93!
94 TYPE(t_io), intent(inout) :: s(ngrids)
95!
96! Local variable declarations.
97!
98# ifdef ADJUST_BOUNDARY
99 integer :: lbij, ubij
100# endif
101 integer :: lbi, ubi, lbj, ubj
102!
103 character (len=*), parameter :: myfile = &
104 & __FILE__
105!
106!-----------------------------------------------------------------------
107! Write out history fields according to IO type.
108!-----------------------------------------------------------------------
109!
110# ifdef ADJUST_BOUNDARY
111 lbij=bounds(ng)%LBij
112 ubij=bounds(ng)%UBij
113# endif
114 lbi=bounds(ng)%LBi(tile)
115 ubi=bounds(ng)%UBi(tile)
116 lbj=bounds(ng)%LBj(tile)
117 ubj=bounds(ng)%UBj(tile)
118!
119 SELECT CASE (his(ng)%IOtype)
120 CASE (io_nf90)
121 CALL wrt_state_nf90 (ng, tile, model, label, &
122 & outrec, i2d, i3d, stime, &
123# ifdef ADJUST_BOUNDARY
124 & lbij, ubij, &
125 & lbi, ubi, lbj, ubj, &
126# endif
127 & s)
128
129# if defined PIO_LIB && defined DISTRIBUTE
130 CASE (io_pio)
131 CALL wrt_state_pio (ng, tile, model, label, &
132 & outrec, i2d, i3d, stime, &
133# ifdef ADJUST_BOUNDARY
134 & lbij, ubij, &
135 & lbi, ubi, lbj, ubj, &
136# endif
137 & s)
138# endif
139 CASE DEFAULT
140 IF (master) THEN
141 WRITE (stdout,10) s(ng)%IOtype
142 END IF
143 exit_flag=3
144 END SELECT
145 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
146!
147 10 FORMAT (' WRT_STATE - Illegal output file type, io_type = ',i0, &
148 & /,13x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
149!
150 RETURN
151 END SUBROUTINE wrt_state
152!
153!***********************************************************************
154 SUBROUTINE wrt_state_nf90 (ng, tile, model, label, &
155 & OutRec, i2d, i3d, stime, &
156# ifdef ADJUST_BOUNDARY
157 & LBij, UBij, &
158 & LBi, UBi, LBj, UBj, &
159# endif
160 & S)
161!***********************************************************************
162!
163 USE mod_netcdf
164!
165! Imported variable declarations.
166!
167 integer, intent(in) :: ng, tile, model, outrec, i2d, i3d
168# ifdef ADJUST_BOUNDARY
169 integer, intent(in) :: lbij, ubij
170# endif
171 integer :: lbi, ubi, lbj, ubj
172!
173 real(dp), intent(in) :: stime
174!
175 character (len=*), intent(in) :: label
176!
177 TYPE(t_io), intent(inout) :: s(ngrids)
178!
179! Local variable declarations.
180!
181 integer :: sstr, send
182 integer :: fcount, gfactor, gtype, ncid, omode, status
183# ifdef SOLVE3D
184 integer :: i, itrc, j, k
185# endif
186!
187 real(r8) :: fmin, fmax
188 real(dp) :: scale
189 real(dp) :: tval(1)
190!
191 character (len=15) :: tstring
192 character (len=22) :: t_code
193 character (len=50) :: string
194
195 character (len=*), parameter :: myfile = &
196 & __FILE__//", wrt_state_nf90"
197!
198 sourcefile=myfile
199!
200!-----------------------------------------------------------------------
201! Write out tangent linear or adjoint state fields.
202!-----------------------------------------------------------------------
203# ifdef PROFILE
204!
205! Turn on time wall clock.
206!
207 CALL wclock_on (ng, model, 81, __line__, myfile)
208# endif
209!
210! Open output NetCDF file.
211!
212 omode=1 ! read and write access
213 CALL netcdf_open (ng, model, s(ng)%name, omode, ncid)
214 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
215!
216! Set grid type factor to write full (gfactor=1) fields or water
217! points (gfactor=-1) fields only.
218!
219# if defined WRITE_WATER && defined MASKING
220 gfactor=-1
221# else
222 gfactor=1
223# endif
224!
225! Set time record index.
226!
227 s(ng)%Rindex=s(ng)%Rindex+1
228 fcount=s(ng)%Fcount
229 s(ng)%Nrec(fcount)=s(ng)%Nrec(fcount)+1
230!
231! Report information.
232!
233 IF (master) THEN
234 CALL time_string (stime, t_code)
235 IF (model.eq.inlm) THEN
236 string='writing NLM state fields,'
237 ELSE IF (model.eq.itlm) THEN
238 string='writing TLM state fields,'
239 ELSE IF (model.eq.iadm) THEN
240 string='writing ADM state fields,'
241 ELSE IF (model.eq.irpm) THEN
242 string='writing RPM state fields,'
243 END IF
244 sstr=scan(calledfrom,'/',back=.true.)+1
245 send=len_trim(calledfrom)
246 WRITE (tstring,'(f15.4)') stime*sec2day
247 WRITE (stdout,10) trim(label), trim(string), t_code, &
248 & ng, trim(adjustl(tstring)), trim(s(ng)%name), &
249 & i3d, outrec, calledfrom(sstr:send)
250 END IF
251!
252! Write out model time (s).
253!
254 IF (lwrtper(ng)) THEN
255 tval(1)=real(outrec,dp)*day2sec
256 ELSE
257 tval(1)=stime
258 END IF
259 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
260 & trim(vname(1,idtime)), tval, &
261 & (/outrec/), (/1/), &
262 & ncid = ncid, &
263 & varid = s(ng)%Vid(idtime))
264 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
265!
266! Write out free-surface.
267!
268 scale=1.0_dp
269 gtype=gfactor*r2dvar
270 IF (model.eq.itlm) THEN
271 status=nf_fwrite2d(ng, model, ncid, idfsur, &
272 & s(ng)%Vid(idfsur), &
273 & outrec, gtype, &
274 & lbi, ubi, lbj, ubj, scale, &
275# ifdef MASKING
276 & grid(ng) % rmask, &
277# endif
278# ifdef WET_DRY
279 & ocean(ng) % tl_zeta(:,:,i2d), &
280 & setfillval = .false., &
281 & minvalue = fmin, &
282 & maxvalue = fmax)
283# else
284 & ocean(ng) % tl_zeta(:,:,i2d), &
285 & minvalue = fmin, &
286 & maxvalue = fmax)
287# endif
288 ELSE IF (model.eq.iadm) THEN
289 status=nf_fwrite2d(ng, model, ncid, idfsur, &
290 & s(ng)%Vid(idfsur), &
291 & outrec, gtype, &
292 & lbi, ubi, lbj, ubj, scale, &
293# ifdef MASKING
294 & grid(ng) % rmask, &
295# endif
296# ifdef WET_DRY
297 & ocean(ng) % ad_zeta(:,:,i2d), &
298 & setfillval = .false., &
299 & minvalue = fmin, &
300 & maxvalue = fmax)
301# else
302 & ocean(ng) % ad_zeta(:,:,i2d), &
303 & minvalue = fmin, &
304 & maxvalue = fmax)
305# endif
306 END IF
307 IF (status.ne.nf90_noerr) THEN
308 IF (master) THEN
309 WRITE (stdout,20) trim(vname(1,idfsur)), trim(label), &
310 & trim(s(ng)%name), outrec
311 END IF
312 exit_flag=3
313 ioerror=status
314 RETURN
315 ELSE
316 IF (master) THEN
317 WRITE (stdout,30) trim(vname(2,idfsur)), fmin, fmax
318 END IF
319 END IF
320
321# ifdef ADJUST_BOUNDARY
322!
323! Write out free-surface open boundaries.
324!
325 IF (any(lobc(:,isfsur,ng))) THEN
326 scale=1.0_dp
327 IF (model.eq.itlm) THEN
328 status=nf_fwrite2d_bry(ng, model, s(ng)%name, ncid, &
329 & vname(1,idsbry(isfsur)), &
330 & s(ng)%Vid(idsbry(isfsur)), &
331 & outrec, r2dvar, &
332 & lbij, ubij, nbrec(ng), scale, &
333 & boundary(ng) % tl_zeta_obc(lbij:,:,:, &
334 & lbout(ng)), &
335 & minvalue = fmin, &
336 & maxvalue = fmax)
337 ELSE IF (model.eq.iadm) THEN
338 status=nf_fwrite2d_bry(ng, model, s(ng)%name, ncid, &
339 & vname(1,idsbry(isfsur)), &
340 & s(ng)%Vid(idsbry(isfsur)), &
341 & outrec, r2dvar, &
342 & lbij, ubij, nbrec(ng), scale, &
343 & boundary(ng) % ad_zeta_obc(lbij:,:,:, &
344 & lbout(ng)), &
345 & minvalue = fmin, &
346 & maxvalue = fmax)
347 END IF
348 IF (status.ne.nf90_noerr) THEN
349 IF (master) THEN
350 WRITE (stdout,20) trim(vname(1,idsbry(isfsur))), &
351 & trim(label), trim(s(ng)%name), outrec
352 END IF
353 exit_flag=3
354 ioerror=status
355 RETURN
356 ELSE
357 IF (master) THEN
358 WRITE (stdout,30) trim(vname(2,idsbry(isfsur))), &
359 & fmin, fmax
360 END IF
361 END IF
362 END IF
363# endif
364!
365! Write out 2D U-momentum component (m/s).
366!
367 scale=1.0_dp
368 gtype=gfactor*u2dvar
369 IF (model.eq.itlm) THEN
370 status=nf_fwrite2d(ng, model, ncid, idubar, &
371 & s(ng)%Vid(idubar), &
372 & outrec, gtype, &
373 & lbi, ubi, lbj, ubj, scale, &
374# ifdef MASKING
375 & grid(ng) % umask, &
376# endif
377 & ocean(ng) % tl_ubar(:,:,i2d), &
378 & minvalue = fmin, &
379 & maxvalue = fmax)
380
381 ELSE IF (model.eq.iadm) THEN
382 status=nf_fwrite2d(ng, model, ncid, idubar, &
383 & s(ng)%Vid(idubar), &
384 & outrec, gtype, &
385 & lbi, ubi, lbj, ubj, scale, &
386# ifdef MASKING
387 & grid(ng) % umask, &
388# endif
389 & ocean(ng) % ad_ubar(:,:,i2d), &
390 & minvalue = fmin, &
391 & maxvalue = fmax)
392 END IF
393 IF (status.ne.nf90_noerr) THEN
394 IF (master) THEN
395 WRITE (stdout,20) trim(vname(1,idubar)), trim(label), &
396 & trim(s(ng)%name), outrec
397 END IF
398 exit_flag=3
399 ioerror=status
400 RETURN
401 ELSE
402 IF (master) THEN
403 WRITE (stdout,30) trim(vname(2,idubar)), fmin, fmax
404 END IF
405 END IF
406
407# ifdef ADJUST_BOUNDARY
408!
409! Write out 2D U-momentum component open boundaries.
410!
411 IF (any(lobc(:,isubar,ng))) THEN
412 scale=1.0_dp
413 IF (model.eq.itlm) THEN
414 status=nf_fwrite2d_bry(ng, model, s(ng)%name, ncid, &
415 & vname(1,idsbry(isubar)), &
416 & s(ng)%Vid(idsbry(isubar)), &
417 & outrec, u2dvar, &
418 & lbij, ubij, nbrec(ng), scale, &
419 & boundary(ng) % tl_ubar_obc(lbij:,:,:, &
420 & lbout(ng)), &
421 & minvalue = fmin, &
422 & maxvalue = fmax)
423 ELSE IF (model.eq.iadm) THEN
424 status=nf_fwrite2d_bry(ng, model, s(ng)%name, ncid, &
425 & vname(1,idsbry(isubar)), &
426 & s(ng)%Vid(idsbry(isubar)), &
427 & outrec, u2dvar, &
428 & lbij, ubij, nbrec(ng), scale, &
429 & boundary(ng) % ad_ubar_obc(lbij:,:,:, &
430 & lbout(ng)), &
431 & minvalue = fmin, &
432 & maxvalue = fmax)
433 END IF
434 IF (status.ne.nf90_noerr) THEN
435 IF (master) THEN
436 WRITE (stdout,20) trim(vname(1,idsbry(isubar))), &
437 & trim(label), trim(s(ng)%name), outrec
438 END IF
439 exit_flag=3
440 ioerror=status
441 RETURN
442 ELSE
443 IF (master) THEN
444 WRITE (stdout,30) trim(vname(2,idsbry(isubar))), &
445 & fmin, fmax
446 END IF
447 END IF
448 END IF
449# endif
450!
451! Write out 2D V-momentum component.
452!
453 scale=1.0_dp
454 gtype=gfactor*v2dvar
455 IF (model.eq.itlm) THEN
456 status=nf_fwrite2d(ng, model, ncid, idvbar, &
457 & s(ng)%Vid(idvbar), &
458 & outrec, gtype, &
459 & lbi, ubi, lbj, ubj, scale, &
460# ifdef MASKING
461 & grid(ng) % vmask, &
462# endif
463 & ocean(ng) % tl_vbar(:,:,i2d), &
464 & minvalue = fmin, &
465 & maxvalue = fmax)
466 ELSE IF (model.eq.iadm) THEN
467 status=nf_fwrite2d(ng, model, ncid, idvbar, &
468 & s(ng)%Vid(idvbar), &
469 & outrec, gtype, &
470 & lbi, ubi, lbj, ubj, scale, &
471# ifdef MASKING
472 & grid(ng) % vmask, &
473# endif
474 & ocean(ng) % ad_vbar(:,:,i2d), &
475 & minvalue = fmin, &
476 & maxvalue = fmax)
477 END IF
478 IF (status.ne.nf90_noerr) THEN
479 IF (master) THEN
480 WRITE (stdout,20) trim(vname(1,idvbar)), trim(label), &
481 & trim(s(ng)%name), outrec
482 END IF
483 exit_flag=3
484 ioerror=status
485 RETURN
486 ELSE
487 IF (master) THEN
488 WRITE (stdout,30) trim(vname(2,idvbar)), fmin, fmax
489 END IF
490 END IF
491
492# ifdef ADJUST_BOUNDARY
493!
494! Write out 2D V-momentum component open boundaries.
495!
496 IF (any(lobc(:,isvbar,ng))) THEN
497 scale=1.0_dp
498 IF (model.eq.itlm) THEN
499 status=nf_fwrite2d_bry(ng, model, s(ng)%name, ncid, &
500 & vname(1,idsbry(isvbar)), &
501 & s(ng)%Vid(idsbry(isvbar)), &
502 & outrec, v2dvar, &
503 & lbij, ubij, nbrec(ng), scale, &
504 & boundary(ng) % tl_vbar_obc(lbij:,:,:, &
505 & lbout(ng)), &
506 & minvalue = fmin, &
507 & maxvalue = fmax)
508 ELSE IF (model.eq.iadm) THEN
509 status=nf_fwrite2d_bry(ng, model, s(ng)%name, ncid, &
510 & vname(1,idsbry(isvbar)), &
511 & s(ng)%Vid(idsbry(isvbar)), &
512 & outrec, v2dvar, &
513 & lbij, ubij, nbrec(ng), scale, &
514 & boundary(ng) % ad_vbar_obc(lbij:,:,:, &
515 & lbout(ng)), &
516 & minvalue = fmin, &
517 & maxvalue = fmax)
518 END IF
519 IF (status.ne.nf90_noerr) THEN
520 IF (master) THEN
521 WRITE (stdout,20) trim(vname(1,idsbry(isvbar))), &
522 & trim(label), trim(s(ng)%name), outrec
523 END IF
524 exit_flag=3
525 ioerror=status
526 RETURN
527 ELSE
528 IF (master) THEN
529 WRITE (stdout,30) trim(vname(2,idsbry(isvbar))), &
530 & fmin, fmax
531 END IF
532 END IF
533 END IF
534# endif
535
536# ifdef ADJUST_WSTRESS
537!
538! Write out surface U-momentum stress. Notice that the stress has its
539! own fixed time-dimension (of size Nfrec) to allow 4DVAR adjustments
540! at other times in addition to initialization time.
541!
542 scale=1.0_dp
543 gtype=gfactor*u3dvar
544 IF (model.eq.itlm) THEN
545 status=nf_fwrite3d(ng, model, ncid, idusms, &
546 & s(ng)%Vid(idusms), &
547 & outrec, gtype, &
548 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
549# ifdef MASKING
550 & grid(ng) % umask, &
551# endif
552 & forces(ng) % tl_ustr(:,:,:,lfout(ng)), &
553 & minvalue = fmin, &
554 & maxvalue = fmax)
555 ELSE IF (model.eq.iadm) THEN
556 status=nf_fwrite3d(ng, model, ncid, idusms, &
557 & s(ng)%Vid(idusms), &
558 & outrec, gtype, &
559 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
560# ifdef MASKING
561 & grid(ng) % umask, &
562# endif
563 & forces(ng) % ad_ustr(:,:,:,lfout(ng)), &
564 & minvalue = fmin, &
565 & maxvalue = fmax)
566 END IF
567
568 IF (status.ne.nf90_noerr) THEN
569 IF (master) THEN
570 WRITE (stdout,20) trim(vname(1,idusms)), trim(label), &
571 & trim(s(ng)%name), lfout(ng)
572 END IF
573 exit_flag=3
574 ioerror=status
575 RETURN
576 ELSE
577 IF (master) THEN
578 WRITE (stdout,30) trim(vname(2,idusms)), fmin, fmax
579 END IF
580 END IF
581!
582! Write out surface V-momentum stress.
583!
584 scale=1.0_dp
585 gtype=gfactor*v3dvar
586 IF (model.eq.itlm) THEN
587 status=nf_fwrite3d(ng, model, ncid, idvsms, &
588 & s(ng)%Vid(idvsms), &
589 & outrec, gtype, &
590 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
591# ifdef MASKING
592 & grid(ng) % vmask, &
593# endif
594 & forces(ng) % tl_vstr(:,:,:,lfout(ng)), &
595 & minvalue = fmin, &
596 & maxvalue = fmax)
597
598 ELSE IF (model.eq.iadm) THEN
599 status=nf_fwrite3d(ng, model, ncid, idvsms, &
600 & s(ng)%Vid(idvsms), &
601 & outrec, gtype, &
602 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
603# ifdef MASKING
604 & grid(ng) % vmask, &
605# endif
606 & forces(ng) % ad_vstr(:,:,:,lfout(ng)), &
607 & minvalue = fmin, &
608 & maxvalue = fmax)
609 END IF
610 IF (status.ne.nf90_noerr) THEN
611 IF (master) THEN
612 WRITE (stdout,20) trim(vname(1,idvsms)), trim(label), &
613 & trim(s(ng)%name), lfout(ng)
614 END IF
615 exit_flag=3
616 ioerror=status
617 RETURN
618 ELSE
619 IF (master) THEN
620 WRITE (stdout,30) trim(vname(2,idvsms)), fmin, fmax
621 END IF
622 END IF
623# endif
624
625# ifdef SOLVE3D
626!
627! Write out 3D U-momentum component (m/s).
628!
629 scale=1.0_dp
630 gtype=gfactor*u3dvar
631 IF (model.eq.itlm) THEN
632 status=nf_fwrite3d(ng, model, ncid, iduvel, &
633 & s(ng)%Vid(iduvel), &
634 & outrec, gtype, &
635 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
636# ifdef MASKING
637 & grid(ng) % umask, &
638# endif
639 & ocean(ng) % tl_u(:,:,:,i3d), &
640 & minvalue = fmin, &
641 & maxvalue = fmax)
642
643 ELSE IF (model.eq.iadm) THEN
644 status=nf_fwrite3d(ng, model, ncid, iduvel, &
645 & s(ng)%Vid(iduvel), &
646 & outrec, gtype, &
647 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
648# ifdef MASKING
649 & grid(ng) % umask, &
650# endif
651 & ocean(ng) % ad_u(:,:,:,i3d), &
652 & minvalue = fmin, &
653 & maxvalue = fmax)
654 END IF
655 IF (status.ne.nf90_noerr) THEN
656 IF (master) THEN
657 WRITE (stdout,20) trim(vname(1,iduvel)), trim(label), &
658 & trim(s(ng)%name), outrec
659 END IF
660 exit_flag=3
661 ioerror=status
662 RETURN
663 ELSE
664 IF (master) THEN
665 WRITE (stdout,30) trim(vname(2,iduvel)), fmin, fmax
666 END IF
667 END IF
668
669# ifdef ADJUST_BOUNDARY
670!
671! Write out 3D U-momentum component open boundaries.
672!
673 IF (any(lobc(:,isuvel,ng))) THEN
674 scale=1.0_dp
675 IF (model.eq.itlm) THEN
676 status=nf_fwrite3d_bry(ng, model, s(ng)%name, ncid, &
677 & vname(1,idsbry(isuvel)), &
678 & s(ng)%Vid(idsbry(isuvel)), &
679 & outrec, u3dvar, &
680 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
681 & boundary(ng) % tl_u_obc(lbij:,:,:,:, &
682 & lbout(ng)), &
683 & minvalue = fmin, &
684 & maxvalue = fmax)
685 ELSE IF (model.eq.iadm) THEN
686 status=nf_fwrite3d_bry(ng, model, s(ng)%name, ncid, &
687 & vname(1,idsbry(isuvel)), &
688 & s(ng)%Vid(idsbry(isuvel)), &
689 & outrec, u3dvar, &
690 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
691 & boundary(ng) % ad_u_obc(lbij:,:,:,:, &
692 & lbout(ng)), &
693 & minvalue = fmin, &
694 & maxvalue = fmax)
695 END IF
696 IF (status.ne.nf90_noerr) THEN
697 IF (master) THEN
698 WRITE (stdout,20) trim(vname(1,idsbry(isuvel))), &
699 & trim(label), trim(s(ng)%name), outrec
700 END IF
701 exit_flag=3
702 ioerror=status
703 RETURN
704 ELSE
705 IF (master) THEN
706 WRITE (stdout,30) trim(vname(2,idsbry(isuvel))), &
707 & fmin, fmax
708 END IF
709 END IF
710 END IF
711# endif
712!
713! Write out 3D V-momentum component (m/s).
714!
715 scale=1.0_dp
716 gtype=gfactor*v3dvar
717 IF (model.eq.itlm) THEN
718 status=nf_fwrite3d(ng, model, ncid, idvvel, &
719 & s(ng)%Vid(idvvel), &
720 & outrec, gtype, &
721 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
722# ifdef MASKING
723 & grid(ng) % vmask, &
724# endif
725 & ocean(ng) % tl_v(:,:,:,i3d), &
726 & minvalue = fmin, &
727 & maxvalue = fmax)
728 ELSE IF (model.eq.iadm) THEN
729 status=nf_fwrite3d(ng, model, ncid, idvvel, &
730 & s(ng)%Vid(idvvel), &
731 & outrec, gtype, &
732 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
733# ifdef MASKING
734 & grid(ng) % vmask, &
735# endif
736 & ocean(ng) % ad_v(:,:,:,i3d), &
737 & minvalue = fmin, &
738 & maxvalue = fmax)
739 END IF
740 IF (status.ne.nf90_noerr) THEN
741 IF (master) THEN
742 WRITE (stdout,20) trim(vname(1,idvvel)), trim(label), &
743 & trim(s(ng)%name), outrec
744 END IF
745 exit_flag=3
746 ioerror=status
747 RETURN
748 ELSE
749 IF (master) THEN
750 WRITE (stdout,30) trim(vname(2,idvvel)), fmin, fmax
751 END IF
752 END IF
753
754# ifdef ADJUST_BOUNDARY
755!
756! Write out 3D V-momentum component open boundaries.
757!
758 IF (any(lobc(:,isvvel,ng))) THEN
759 scale=1.0_dp
760 IF (model.eq.itlm) THEN
761 status=nf_fwrite3d_bry(ng, model, s(ng)%name, ncid, &
762 & vname(1,idsbry(isvvel)), &
763 & s(ng)%Vid(idsbry(isvvel)), &
764 & outrec , v3dvar, &
765 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
766 & boundary(ng) % ad_v_obc(lbij:,:,:,:, &
767 & lbout(ng)), &
768 & minvalue = fmin, &
769 & maxvalue = fmax)
770 ELSE IF (model.eq.iadm) THEN
771 status=nf_fwrite3d_bry(ng, model, s(ng)%name, ncid, &
772 & vname(1,idsbry(isvvel)), &
773 & s(ng)%Vid(idsbry(isvvel)), &
774 & outrec , v3dvar, &
775 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
776 & boundary(ng) % ad_v_obc(lbij:,:,:,:, &
777 & lbout(ng)), &
778 & minvalue = fmin, &
779 & maxvalue = fmax)
780 END IF
781 IF (status.ne.nf90_noerr) THEN
782 IF (master) THEN
783 WRITE (stdout,20) trim(vname(1,idsbry(isvvel))), &
784 & trim(label), trim(s(ng)%name), outrec
785 END IF
786 exit_flag=3
787 ioerror=status
788 RETURN
789 ELSE
790 IF (master) THEN
791 WRITE (stdout,30) trim(vname(2,idsbry(isvvel))), &
792 & fmin, fmax
793 END IF
794 END IF
795 END IF
796# endif
797!
798! Write out tracer type variables.
799!
800 DO itrc=1,nt(ng)
801 scale=1.0_dp
802 gtype=gfactor*r3dvar
803 IF (model.eq.itlm) THEN
804 status=nf_fwrite3d(ng, model, ncid, idtvar(itrc), &
805 & s(ng)%Tid(itrc), &
806 & outrec, gtype, &
807 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
808# ifdef MASKING
809 & grid(ng) % rmask, &
810# endif
811 & ocean(ng) % tl_t(:,:,:,i3d,itrc), &
812 & minvalue = fmin, &
813 & maxvalue = fmax)
814 ELSE IF (model.eq.iadm) THEN
815 status=nf_fwrite3d(ng, model, ncid, idtvar(itrc), &
816 & s(ng)%Tid(itrc), &
817 & outrec, gtype, &
818 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
819# ifdef MASKING
820 & grid(ng) % rmask, &
821# endif
822 & ocean(ng) % ad_t(:,:,:,i3d,itrc), &
823 & minvalue = fmin, &
824 & maxvalue = fmax)
825 END IF
826 IF (status.ne.nf90_noerr) THEN
827 IF (master) THEN
828 WRITE (stdout,20) trim(vname(1,idtvar(itrc))), &
829 & trim(label), trim(s(ng)%name), outrec
830 END IF
831 exit_flag=3
832 ioerror=status
833 RETURN
834 ELSE
835 IF (master) THEN
836 WRITE (stdout,30) trim(vname(2,idtvar(itrc))), &
837 & fmin, fmax
838 END IF
839 END IF
840 END DO
841
842# ifdef ADJUST_BOUNDARY
843!
844! Write out tracers open boundaries.
845!
846 DO itrc=1,nt(ng)
847 IF (any(lobc(:,istvar(itrc),ng))) THEN
848 scale=1.0_dp
849 IF (model.eq.itlm) THEN
850 status=nf_fwrite3d_bry(ng, model, s(ng)%name, ncid, &
851 & vname(1,idsbry(istvar(itrc))), &
852 & s(ng)%Vid(idsbry(istvar(itrc))), &
853 & outrec, r3dvar, &
854 & lbij, ubij, 1, n(ng), nbrec(ng), &
855 & scale, &
856 & boundary(ng) % tl_t_obc(lbij:,:,:,:, &
857 & lbout(ng),itrc), &
858 & minvalue = fmin, &
859 & maxvalue = fmax)
860 ELSE IF (model.eq.iadm) THEN
861 status=nf_fwrite3d_bry(ng, model, s(ng)%name, ncid, &
862 & vname(1,idsbry(istvar(itrc))), &
863 & s(ng)%Vid(idsbry(istvar(itrc))), &
864 & outrec, r3dvar, &
865 & lbij, ubij, 1, n(ng), nbrec(ng), &
866 & scale, &
867 & boundary(ng) % ad_t_obc(lbij:,:,:,:, &
868 & lbout(ng),itrc), &
869 & minvalue = fmin, &
870 & maxvalue = fmax)
871 END IF
872 IF (status.ne.nf90_noerr) THEN
873 IF (master) THEN
874 WRITE (stdout,20) trim(vname(1,idsbry(istvar(itrc)))), &
875 & trim(label), trim(s(ng)%name), outrec
876 END IF
877 exit_flag=3
878 ioerror=status
879 RETURN
880 ELSE
881 IF (master) THEN
882 WRITE (stdout,30) trim(vname(2,idsbry(istvar(itrc)))), &
883 & fmin, fmax
884 END IF
885 END IF
886 END IF
887 END DO
888# endif
889
890# ifdef ADJUST_STFLUX
891!
892! Write out surface net tracers fluxes. Notice that fluxes have their
893! own fixed time-dimension (of size Nfrec) to allow 4DVAR adjustments
894! at other times in addition to initialization time.
895!
896 DO itrc=1,nt(ng)
897 IF (lstflux(itrc,ng)) THEN
898 scale=1.0_dp ! kinematic flux units
899 gtype=gfactor*r3dvar
900 IF (model.eq.itlm) THEN
901 status=nf_fwrite3d(ng, itlm, ncid, idtsur(itrc), &
902 & s(ng)%Vid(idtsur(itrc)), &
903 & outrec, gtype, &
904 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
905# ifdef MASKING
906 & grid(ng) % rmask, &
907# endif
908 & forces(ng)% tl_tflux(:,:,:, &
909 & lfout(ng),itrc), &
910 & minvalue = fmin, &
911 & maxvalue = fmax)
912 ELSE IF (model.eq.iadm) THEN
913 status=nf_fwrite3d(ng, itlm, ncid, idtsur(itrc), &
914 & s(ng)%Vid(idtsur(itrc)), &
915 & outrec, gtype, &
916 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
917# ifdef MASKING
918 & grid(ng) % rmask, &
919# endif
920 & forces(ng)% ad_tflux(:,:,:, &
921 & lfout(ng),itrc), &
922 & minvalue = fmin, &
923 & maxvalue = fmax)
924 END IF
925 IF (status.ne.nf90_noerr) THEN
926 IF (master) THEN
927 WRITE (stdout,20) trim(vname(1,idtsur(itrc))), &
928 & trim(label), trim(s(ng)%name), lfout(ng)
929 END IF
930 exit_flag=3
931 ioerror=status
932 RETURN
933 ELSE
934 IF (master) THEN
935 WRITE (stdout,30) trim(vname(2,idtsur(itrc))), &
936 & fmin, fmax
937 END IF
938 END IF
939 END IF
940 END DO
941# endif
942# endif
943!
944!-----------------------------------------------------------------------
945! Synchronize tangent NetCDF file to disk to allow other processes
946! to access data immediately after it is written. Then, close NetCDF
947! file.
948!-----------------------------------------------------------------------
949!
950 CALL netcdf_sync (ng, itlm, s(ng)%name, ncid)
951 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
952!
953 CALL netcdf_close (ng, itlm, ncid, s(ng)%name, .false.)
954
955# ifdef PROFILE
956!
957! Turn off time wall clock.
958!
959 CALL wclock_off (ng, model, 81, __line__, myfile)
960# endif
961!
962 10 FORMAT (1x,a,': WRT_STATE_NF90 - ',a,t75,a, &
963 & /,24x,'(Grid ',i2.2,', t = ',a,', File: ',a, &
964 & ', Index=',i1,', Rec=',i0,')', &
965 & /,24x,'Called from ''',a,'''')
966 20 FORMAT (/,' WRT_STATE_NF90 - error while writing variable: ',a, &
967 & /,18x,'into ',a,' NetCDF file: ',a, &
968 & /,18x,'for time record: ',i0)
969 30 FORMAT (21x,'- ',a,/,24x,'(Min = ',1p,e15.8, &
970 & ' Max = ',1p,e15.8,')')
971!
972 RETURN
973 END SUBROUTINE wrt_state_nf90
974
975# if defined PIO_LIB && defined DISTRIBUTE
976!
977!***********************************************************************
978 SUBROUTINE wrt_state_pio (ng, tile, model, label, &
979 & OutRec, i2d, i3d, stime, &
980# ifdef ADJUST_BOUNDARY
981 & LBij, UBij, &
982 & LBi, UBi, LBj, UBj, &
983# endif
984 & S)
985!***********************************************************************
986!
988!
989! Imported variable declarations.
990!
991 integer, intent(in) :: ng, tile, model, outrec, i2d, i3d
992# ifdef ADJUST_BOUNDARY
993 integer, intent(in) :: lbij, ubij
994# endif
995 integer :: lbi, ubi, lbj, ubj
996!
997 real(dp), intent(in) :: stime
998!
999 character (len=*), intent(in) :: label
1000!
1001 TYPE(t_io), intent(inout) :: s(ngrids)
1002!
1003! Local variable declarations.
1004!
1005 integer :: sstr, send
1006 integer :: fcount, ifield, omode, status
1007# ifdef SOLVE3D
1008 integer :: i, itrc, j, k
1009# endif
1010!
1011 real(r8) :: fmin, fmax
1012 real(dp) :: scale
1013 real(dp) :: tval(1)
1014!
1015 character (len=15) :: tstring
1016 character (len=22) :: t_code
1017 character (len=50) :: string
1018
1019 character (len=*), parameter :: myfile = &
1020 & __FILE__//", wrt_state_pio"
1021!
1022 TYPE (io_desc_t), pointer :: iodesc
1023 TYPE (file_desc_t) :: piofile
1024!
1025 sourcefile=myfile
1026!
1027!-----------------------------------------------------------------------
1028! Write out tangent linear or adjoint state fields.
1029!-----------------------------------------------------------------------
1030# ifdef PROFILE
1031!
1032! Turn on time wall clock.
1033!
1034 CALL wclock_on (ng, model, 81, __line__, myfile)
1035# endif
1036!
1037! Open output NetCDF file.
1038!
1039 omode=1 ! read and write access
1040 CALL pio_netcdf_open (ng, model, s(ng)%name, omode, piofile)
1041 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1042!
1043! Set time record index.
1044!
1045 s(ng)%Rindex=s(ng)%Rindex+1
1046 fcount=s(ng)%Fcount
1047 s(ng)%Nrec(fcount)=s(ng)%Nrec(fcount)+1
1048!
1049! Report information.
1050!
1051 IF (master) THEN
1052 CALL time_string (stime, t_code)
1053 IF (model.eq.inlm) THEN
1054 string='writing NLM state fields,'
1055 ELSE IF (model.eq.itlm) THEN
1056 string='writing TLM state fields,'
1057 ELSE IF (model.eq.iadm) THEN
1058 string='writing ADM state fields,'
1059 ELSE IF (model.eq.irpm) THEN
1060 string='writing RPM state fields,'
1061 END IF
1062 sstr=scan(calledfrom,'/',back=.true.)+1
1063 send=len_trim(calledfrom)
1064 WRITE (tstring,'(f15.4)') stime*sec2day
1065 WRITE (stdout,10) trim(label), trim(string), t_code, &
1066 & ng, trim(adjustl(tstring)), trim(s(ng)%name), &
1067 & i3d, outrec, calledfrom(sstr:send)
1068 END IF
1069!
1070! Write out model time (s).
1071!
1072 IF (lwrtper(ng)) THEN
1073 tval(1)=real(outrec,dp)*day2sec
1074 ELSE
1075 tval(1)=stime
1076 END IF
1077 CALL pio_netcdf_put_fvar (ng, model, s(ng)%name, &
1078 & trim(vname(1,idtime)), tval, &
1079 & (/outrec/), (/1/), &
1080 & piofile = piofile, &
1081 & piovar = s(ng)%pioVar(idtime)%vd)
1082 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1083!
1084! Write out free-surface.
1085!
1086 scale=1.0_dp
1087 IF (s(ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
1088 iodesc => iodesc_dp_r2dvar(ng)
1089 ELSE
1090 iodesc => iodesc_sp_r2dvar(ng)
1091 END IF
1092!
1093 IF (model.eq.itlm) THEN
1094 status=nf_fwrite2d(ng, model, piofile, idfsur, &
1095 & s(ng)%pioVar(idfsur), outrec, &
1096 & iodesc, &
1097 & lbi, ubi, lbj, ubj, scale, &
1098# ifdef MASKING
1099 & grid(ng) % rmask, &
1100# endif
1101# ifdef WET_DRY
1102 & ocean(ng) % tl_zeta(:,:,i2d), &
1103 & setfillval = .false., &
1104 & minvalue = fmin, &
1105 & maxvalue = fmax)
1106# else
1107 & ocean(ng) % tl_zeta(:,:,i2d), &
1108 & minvalue = fmin, &
1109 & maxvalue = fmax)
1110# endif
1111 ELSE IF (model.eq.iadm) THEN
1112 status=nf_fwrite2d(ng, model, piofile, idfsur, &
1113 & s(ng)%pioVar(idfsur), outrec, &
1114 & iodesc, &
1115 & lbi, ubi, lbj, ubj, scale, &
1116# ifdef MASKING
1117 & grid(ng) % rmask, &
1118# endif
1119# ifdef WET_DRY
1120 & ocean(ng) % ad_zeta(:,:,i2d), &
1121 & setfillval = .false., &
1122 & minvalue = fmin, &
1123 & maxvalue = fmax)
1124# else
1125 & ocean(ng) % ad_zeta(:,:,i2d), &
1126 & minvalue = fmin, &
1127 & maxvalue = fmax)
1128# endif
1129 END IF
1130 IF (status.ne.pio_noerr) THEN
1131 IF (master) THEN
1132 WRITE (stdout,20) trim(vname(1,idfsur)), trim(label), &
1133 & trim(s(ng)%name), outrec
1134 END IF
1135 exit_flag=3
1136 ioerror=status
1137 RETURN
1138 ELSE
1139 IF (master) THEN
1140 WRITE (stdout,30) trim(vname(2,idfsur)), fmin, fmax
1141 END IF
1142 END IF
1143
1144# ifdef ADJUST_BOUNDARY
1145!
1146! Write out free-surface open boundaries.
1147!
1148 IF (any(lobc(:,isfsur,ng))) THEN
1149 scale=1.0_dp
1150 ifield=idsbry(isfsur)
1151 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1152 iodesc => iodesc_dp_r2dobc(ng)
1153 ELSE
1154 iodesc => iodesc_sp_r2dobc(ng)
1155 END IF
1156!
1157 IF (model.eq.itlm) THEN
1158 status=nf_fwrite2d_bry(ng, model, s(ng)%name, piofile, &
1159 & vname(1,ifield), &
1160 & s(ng)%pioVar(ifield), outrec, &
1161 & iodesc, &
1162 & lbij, ubij, nbrec(ng), scale, &
1163 & boundary(ng) % tl_zeta_obc(lbij:,:,:, &
1164 & lbout(ng)), &
1165 & minvalue = fmin, &
1166 & maxvalue = fmax)
1167 ELSE IF (model.eq.iadm) THEN
1168 status=nf_fwrite2d_bry(ng, model, s(ng)%name, piofile, &
1169 & vname(1,ifield), &
1170 & s(ng)%pioVar(ifield), outrec, &
1171 & iodesc, &
1172 & lbij, ubij, nbrec(ng), scale, &
1173 & boundary(ng) % ad_zeta_obc(lbij:,:,:, &
1174 & lbout(ng)), &
1175 & minvalue = fmin, &
1176 & maxvalue = fmax)
1177 END IF
1178 IF (status.ne.pio_noerr) THEN
1179 IF (master) THEN
1180 WRITE (stdout,20) trim(vname(1,ifield)), &
1181 & trim(label), trim(s(ng)%name), outrec
1182 END IF
1183 exit_flag=3
1184 ioerror=status
1185 RETURN
1186 ELSE
1187 IF (master) THEN
1188 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1189 END IF
1190 END IF
1191 END IF
1192# endif
1193!
1194! Write out 2D U-momentum component (m/s).
1195!
1196 scale=1.0_dp
1197 IF (s(ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
1198 iodesc => iodesc_dp_u2dvar(ng)
1199 ELSE
1200 iodesc => iodesc_sp_u2dvar(ng)
1201 END IF
1202!
1203 IF (model.eq.itlm) THEN
1204 status=nf_fwrite2d(ng, model, piofile, idubar, &
1205 & s(ng)%pioVar(idubar), outrec, &
1206 & iodesc, &
1207 & lbi, ubi, lbj, ubj, scale, &
1208# ifdef MASKING
1209 & grid(ng) % umask, &
1210# endif
1211 & ocean(ng) % tl_ubar(:,:,i2d), &
1212 & minvalue = fmin, &
1213 & maxvalue = fmax)
1214
1215 ELSE IF (model.eq.iadm) THEN
1216 status=nf_fwrite2d(ng, model, piofile, idubar, &
1217 & s(ng)%pioVar(idubar), outrec, &
1218 & iodesc, &
1219 & lbi, ubi, lbj, ubj, scale, &
1220# ifdef MASKING
1221 & grid(ng) % umask, &
1222# endif
1223 & ocean(ng) % ad_ubar(:,:,i2d), &
1224 & minvalue = fmin, &
1225 & maxvalue = fmax)
1226 END IF
1227 IF (status.ne.pio_noerr) THEN
1228 IF (master) THEN
1229 WRITE (stdout,20) trim(vname(1,idubar)), trim(label), &
1230 & trim(s(ng)%name), outrec
1231 END IF
1232 exit_flag=3
1233 ioerror=status
1234 RETURN
1235 ELSE
1236 IF (master) THEN
1237 WRITE (stdout,30) trim(vname(2,idubar)), fmin, fmax
1238 END IF
1239 END IF
1240
1241# ifdef ADJUST_BOUNDARY
1242!
1243! Write out 2D U-momentum component open boundaries.
1244!
1245 IF (any(lobc(:,isubar,ng))) THEN
1246 scale=1.0_dp
1247 ifield=idsbry(isubar)
1248 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1249 iodesc => iodesc_dp_u2dobc(ng)
1250 ELSE
1251 iodesc => iodesc_sp_u2dobc(ng)
1252 END IF
1253!
1254 IF (model.eq.itlm) THEN
1255 status=nf_fwrite2d_bry(ng, model, s(ng)%name, piofile, &
1256 & vname(1,ifield), &
1257 & s(ng)%pioVar(ifield), outrec, &
1258 & iodesc, &
1259 & lbij, ubij, nbrec(ng), scale, &
1260 & boundary(ng) % tl_ubar_obc(lbij:,:,:, &
1261 & lbout(ng)), &
1262 & minvalue = fmin, &
1263 & maxvalue = fmax)
1264 ELSE IF (model.eq.iadm) THEN
1265 status=nf_fwrite2d_bry(ng, model, s(ng)%name, piofile, &
1266 & vname(1,ifield), &
1267 & s(ng)%pioVar(ifield), outrec, &
1268 & iodesc, &
1269 & lbij, ubij, nbrec(ng), scale, &
1270 & boundary(ng) % ad_ubar_obc(lbij:,:,:, &
1271 & lbout(ng)), &
1272 & minvalue = fmin, &
1273 & maxvalue = fmax)
1274 END IF
1275 IF (status.ne.pio_noerr) THEN
1276 IF (master) THEN
1277 WRITE (stdout,20) trim(vname(1,ifield)), &
1278 & trim(label), trim(s(ng)%name), outrec
1279 END IF
1280 exit_flag=3
1281 ioerror=status
1282 RETURN
1283 ELSE
1284 IF (master) THEN
1285 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1286 END IF
1287 END IF
1288 END IF
1289# endif
1290!
1291! Write out 2D V-momentum component.
1292!
1293 scale=1.0_dp
1294 IF (s(ng)%pioVar(idvbar)%dkind.eq.pio_double) THEN
1295 iodesc => iodesc_dp_v2dvar(ng)
1296 ELSE
1297 iodesc => iodesc_sp_v2dvar(ng)
1298 END IF
1299!
1300 IF (model.eq.itlm) THEN
1301 status=nf_fwrite2d(ng, model, piofile, idvbar, &
1302 & s(ng)%pioVar(idvbar), outrec, &
1303 & iodesc, &
1304 & lbi, ubi, lbj, ubj, scale, &
1305# ifdef MASKING
1306 & grid(ng) % vmask, &
1307# endif
1308 & ocean(ng) % tl_vbar(:,:,i2d), &
1309 & minvalue = fmin, &
1310 & maxvalue = fmax)
1311 ELSE IF (model.eq.iadm) THEN
1312 status=nf_fwrite2d(ng, model, piofile, idvbar, &
1313 & s(ng)%pioVar(idvbar), outrec, &
1314 & iodesc, &
1315 & lbi, ubi, lbj, ubj, scale, &
1316# ifdef MASKING
1317 & grid(ng) % vmask, &
1318# endif
1319 & ocean(ng) % ad_vbar(:,:,i2d), &
1320 & minvalue = fmin, &
1321 & maxvalue = fmax)
1322 END IF
1323 IF (status.ne.pio_noerr) THEN
1324 IF (master) THEN
1325 WRITE (stdout,20) trim(vname(1,idvbar)), trim(label), &
1326 & trim(s(ng)%name), outrec
1327 END IF
1328 exit_flag=3
1329 ioerror=status
1330 RETURN
1331 ELSE
1332 IF (master) THEN
1333 WRITE (stdout,30) trim(vname(2,idvbar)), fmin, fmax
1334 END IF
1335 END IF
1336
1337# ifdef ADJUST_BOUNDARY
1338!
1339! Write out 2D V-momentum component open boundaries.
1340!
1341 IF (any(lobc(:,isvbar,ng))) THEN
1342 scale=1.0_dp
1343 ifield=idsbry(isvbar)
1344 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1345 iodesc => iodesc_dp_v2dobc(ng)
1346 ELSE
1347 iodesc => iodesc_sp_v2dobc(ng)
1348 END IF
1349!
1350 IF (model.eq.itlm) THEN
1351 status=nf_fwrite2d_bry(ng, model, s(ng)%name, piofile, &
1352 & vname(1,ifield), &
1353 & s(ng)%pioVar(ifield), outrec, &
1354 & iodesc, &
1355 & lbij, ubij, nbrec(ng), scale, &
1356 & boundary(ng) % tl_vbar_obc(lbij:,:,:, &
1357 & lbout(ng)), &
1358 & minvalue = fmin, &
1359 & maxvalue = fmax)
1360 ELSE IF (model.eq.iadm) THEN
1361 status=nf_fwrite2d_bry(ng, model, s(ng)%name, piofile, &
1362 & vname(1,ifield), &
1363 & s(ng)%pioVar(ifield), outrec, &
1364 & iodesc, &
1365 & lbij, ubij, nbrec(ng), scale, &
1366 & boundary(ng) % ad_vbar_obc(lbij:,:,:, &
1367 & lbout(ng)), &
1368 & minvalue = fmin, &
1369 & maxvalue = fmax)
1370 END IF
1371 IF (status.ne.pio_noerr) THEN
1372 IF (master) THEN
1373 WRITE (stdout,20) trim(vname(1,ifield)), &
1374 & trim(label), trim(s(ng)%name), outrec
1375 END IF
1376 exit_flag=3
1377 ioerror=status
1378 RETURN
1379 ELSE
1380 IF (master) THEN
1381 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1382 END IF
1383 END IF
1384 END IF
1385# endif
1386
1387# ifdef ADJUST_WSTRESS
1388!
1389! Write out surface U-momentum stress. Notice that the stress has its
1390! own fixed time-dimension (of size Nfrec) to allow 4DVAR adjustments
1391! at other times in addition to initialization time.
1392!
1393 scale=1.0_dp
1394 IF (s(ng)%pioVar(idusms)%dkind.eq.pio_double) THEN
1395 iodesc => iodesc_dp_u2dfrc(ng)
1396 ELSE
1397 iodesc => iodesc_sp_u2dfrc(ng)
1398 END IF
1399!
1400 IF (model.eq.itlm) THEN
1401 status=nf_fwrite3d(ng, model, piofile, idusms, &
1402 & s(ng)%pioVar(idusms), outrec, &
1403 & iodesc, &
1404 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1405# ifdef MASKING
1406 & grid(ng) % umask, &
1407# endif
1408 & forces(ng) % tl_ustr(:,:,:,lfout(ng)), &
1409 & minvalue = fmin, &
1410 & maxvalue = fmax)
1411 ELSE IF (model.eq.iadm) THEN
1412 status=nf_fwrite3d(ng, model, piofile, idusms, &
1413 & s(ng)%pioVar(idusms), outrec, &
1414 & iodesc, &
1415 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1416# ifdef MASKING
1417 & grid(ng) % umask, &
1418# endif
1419 & forces(ng) % ad_ustr(:,:,:,lfout(ng)), &
1420 & minvalue = fmin, &
1421 & maxvalue = fmax)
1422 END IF
1423 IF (status.ne.pio_noerr) THEN
1424 IF (master) THEN
1425 WRITE (stdout,20) trim(vname(1,idusms)), trim(label), &
1426 & trim(s(ng)%name), lfout(ng)
1427 END IF
1428 exit_flag=3
1429 ioerror=status
1430 RETURN
1431 ELSE
1432 IF (master) THEN
1433 WRITE (stdout,30) trim(vname(2,idusms)), fmin, fmax
1434 END IF
1435 END IF
1436!
1437! Write out surface V-momentum stress.
1438!
1439 scale=1.0_dp
1440 IF (s(ng)%pioVar(idvsms)%dkind.eq.pio_double) THEN
1441 iodesc => iodesc_dp_v2dfrc(ng)
1442 ELSE
1443 iodesc => iodesc_sp_v2dfrc(ng)
1444 END IF
1445!
1446 IF (model.eq.itlm) THEN
1447 status=nf_fwrite3d(ng, model, piofile, idvsms, &
1448 & s(ng)%pioVar(idvsms), outrec, &
1449 & iodesc, &
1450 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1451# ifdef MASKING
1452 & grid(ng) % vmask, &
1453# endif
1454 & forces(ng) % tl_vstr(:,:,:,lfout(ng)), &
1455 & minvalue = fmin, &
1456 & maxvalue = fmax)
1457
1458 ELSE IF (model.eq.iadm) THEN
1459 status=nf_fwrite3d(ng, model, piofile, idvsms, &
1460 & s(ng)%pioVar(idvsms), outrec, &
1461 & iodesc, &
1462 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1463# ifdef MASKING
1464 & grid(ng) % vmask, &
1465# endif
1466 & forces(ng) % ad_vstr(:,:,:,lfout(ng)), &
1467 & minvalue = fmin, &
1468 & maxvalue = fmax)
1469 END IF
1470 IF (status.ne.pio_noerr) THEN
1471 IF (master) THEN
1472 WRITE (stdout,20) trim(vname(1,idvsms)), trim(label), &
1473 & trim(s(ng)%name), lfout(ng)
1474 END IF
1475 exit_flag=3
1476 ioerror=status
1477 RETURN
1478 ELSE
1479 IF (master) THEN
1480 WRITE (stdout,30) trim(vname(2,idvsms)), fmin, fmax
1481 END IF
1482 END IF
1483# endif
1484
1485# ifdef SOLVE3D
1486!
1487! Write out 3D U-momentum component (m/s).
1488!
1489 scale=1.0_dp
1490 IF (s(ng)%pioVar(idfsur)%dkind.eq.pio_double) THEN
1491 iodesc => iodesc_dp_u3dvar(ng)
1492 ELSE
1493 iodesc => iodesc_sp_u3dvar(ng)
1494 END IF
1495!
1496 IF (model.eq.itlm) THEN
1497 status=nf_fwrite3d(ng, model, piofile, iduvel, &
1498 & s(ng)%pioVar(iduvel), outrec, &
1499 & iodesc, &
1500 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1501# ifdef MASKING
1502 & grid(ng) % umask, &
1503# endif
1504 & ocean(ng) % tl_u(:,:,:,i3d), &
1505 & minvalue = fmin, &
1506 & maxvalue = fmax)
1507
1508 ELSE IF (model.eq.iadm) THEN
1509 status=nf_fwrite3d(ng, model, piofile, iduvel, &
1510 & s(ng)%pioVar(iduvel), outrec, &
1511 & iodesc, &
1512 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1513# ifdef MASKING
1514 & grid(ng) % umask, &
1515# endif
1516 & ocean(ng) % ad_u(:,:,:,i3d), &
1517 & minvalue = fmin, &
1518 & maxvalue = fmax)
1519 END IF
1520 IF (status.ne.pio_noerr) THEN
1521 IF (master) THEN
1522 WRITE (stdout,20) trim(vname(1,iduvel)), trim(label), &
1523 & trim(s(ng)%name), outrec
1524 END IF
1525 exit_flag=3
1526 ioerror=status
1527 RETURN
1528 ELSE
1529 IF (master) THEN
1530 WRITE (stdout,30) trim(vname(2,iduvel)), fmin, fmax
1531 END IF
1532 END IF
1533
1534# ifdef ADJUST_BOUNDARY
1535!
1536! Write out 3D U-momentum component open boundaries.
1537!
1538 IF (any(lobc(:,isuvel,ng))) THEN
1539 scale=1.0_dp
1540 ifield=idsbry(isuvel)
1541 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1542 iodesc => iodesc_dp_u3dobc(ng)
1543 ELSE
1544 iodesc => iodesc_sp_u3dobc(ng)
1545 END IF
1546!
1547 IF (model.eq.itlm) THEN
1548 status=nf_fwrite3d_bry(ng, model, s(ng)%name, piofile, &
1549 & vname(1,ifield), &
1550 & s(ng)%pioVar(ifield), outrec, &
1551 & iodesc, &
1552 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
1553 & boundary(ng) % tl_u_obc(lbij:,:,:,:, &
1554 & lbout(ng)), &
1555 & minvalue = fmin, &
1556 & maxvalue = fmax)
1557 ELSE IF (model.eq.iadm) THEN
1558 status=nf_fwrite3d_bry(ng, model, s(ng)%name, piofile, &
1559 & vname(1,ifield), &
1560 & s(ng)%pioVar(ifield), outrec, &
1561 & iodesc, &
1562 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
1563 & boundary(ng) % ad_u_obc(lbij:,:,:,:, &
1564 & lbout(ng)), &
1565 & minvalue = fmin, &
1566 & maxvalue = fmax)
1567 END IF
1568 IF (status.ne.pio_noerr) THEN
1569 IF (master) THEN
1570 WRITE (stdout,20) trim(vname(1,ifield)), &
1571 & trim(label), trim(s(ng)%name), outrec
1572 END IF
1573 exit_flag=3
1574 ioerror=status
1575 RETURN
1576 ELSE
1577 IF (master) THEN
1578 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1579 END IF
1580 END IF
1581 END IF
1582# endif
1583!
1584! Write out 3D V-momentum component (m/s).
1585!
1586 scale=1.0_dp
1587 IF (s(ng)%pioVar(idvvel)%dkind.eq.pio_double) THEN
1588 iodesc => iodesc_dp_v3dvar(ng)
1589 ELSE
1590 iodesc => iodesc_sp_v3dvar(ng)
1591 END IF
1592!
1593 IF (model.eq.itlm) THEN
1594 status=nf_fwrite3d(ng, model, piofile, idvvel, &
1595 & s(ng)%pioVar(idvvel), outrec, &
1596 & iodesc, &
1597 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1598# ifdef MASKING
1599 & grid(ng) % vmask, &
1600# endif
1601 & ocean(ng) % tl_v(:,:,:,i3d), &
1602 & minvalue = fmin, &
1603 & maxvalue = fmax)
1604 ELSE IF (model.eq.iadm) THEN
1605 status=nf_fwrite3d(ng, model, piofile, idvvel, &
1606 & s(ng)%pioVar(idvvel), outrec, &
1607 & iodesc, &
1608 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1609# ifdef MASKING
1610 & grid(ng) % vmask, &
1611# endif
1612 & ocean(ng) % ad_v(:,:,:,i3d), &
1613 & minvalue = fmin, &
1614 & maxvalue = fmax)
1615 END IF
1616 IF (status.ne.pio_noerr) THEN
1617 IF (master) THEN
1618 WRITE (stdout,20) trim(vname(1,idvvel)), trim(label), &
1619 & trim(s(ng)%name), outrec
1620 END IF
1621 exit_flag=3
1622 ioerror=status
1623 RETURN
1624 ELSE
1625 IF (master) THEN
1626 WRITE (stdout,30) trim(vname(2,idvvel)), fmin, fmax
1627 END IF
1628 END IF
1629
1630# ifdef ADJUST_BOUNDARY
1631!
1632! Write out 3D V-momentum component open boundaries.
1633!
1634 IF (any(lobc(:,isvvel,ng))) THEN
1635 scale=1.0_dp
1636 ifield=idsbry(isvvel)
1637 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1638 iodesc => iodesc_dp_v3dobc(ng)
1639 ELSE
1640 iodesc => iodesc_sp_v3dobc(ng)
1641 END IF
1642!
1643 IF (model.eq.itlm) THEN
1644 status=nf_fwrite3d_bry(ng, model, s(ng)%name, piofile, &
1645 & vname(1,ifield), &
1646 & s(ng)%pioVar(ifield), outrec, &
1647 & iodesc, &
1648 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
1649 & boundary(ng) % ad_v_obc(lbij:,:,:,:, &
1650 & lbout(ng)), &
1651 & minvalue = fmin, &
1652 & maxvalue = fmax)
1653 ELSE IF (model.eq.iadm) THEN
1654 status=nf_fwrite3d_bry(ng, model, s(ng)%name, piofile, &
1655 & vname(1,ifield), &
1656 & s(ng)%pioVar(ifield), outrec, &
1657 & iodesc, &
1658 & lbij, ubij, 1, n(ng), nbrec(ng), scale,&
1659 & boundary(ng) % ad_v_obc(lbij:,:,:,:, &
1660 & lbout(ng)), &
1661 & minvalue = fmin, &
1662 & maxvalue = fmax)
1663 END IF
1664 IF (status.ne.pio_noerr) THEN
1665 IF (master) THEN
1666 WRITE (stdout,20) trim(vname(1,ifield)), &
1667 & trim(label), trim(s(ng)%name), outrec
1668 END IF
1669 exit_flag=3
1670 ioerror=status
1671 RETURN
1672 ELSE
1673 IF (master) THEN
1674 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1675 END IF
1676 END IF
1677 END IF
1678# endif
1679!
1680! Write out tracer type variables.
1681!
1682 DO itrc=1,nt(ng)
1683 scale=1.0_dp
1684 IF (s(ng)%pioTrc(itrc)%dkind.eq.pio_double) THEN
1685 iodesc => iodesc_dp_r3dvar(ng)
1686 ELSE
1687 iodesc => iodesc_sp_r3dvar(ng)
1688 END IF
1689!
1690 IF (model.eq.itlm) THEN
1691 status=nf_fwrite3d(ng, model, piofile, idtvar(itrc), &
1692 & s(ng)%pioTrc(itrc), outrec, &
1693 & iodesc, &
1694 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1695# ifdef MASKING
1696 & grid(ng) % rmask, &
1697# endif
1698 & ocean(ng) % tl_t(:,:,:,i3d,itrc), &
1699 & minvalue = fmin, &
1700 & maxvalue = fmax)
1701 ELSE IF (model.eq.iadm) THEN
1702 status=nf_fwrite3d(ng, model, piofile, idtvar(itrc), &
1703 & s(ng)%pioTrc(itrc), outrec, &
1704 & iodesc, &
1705 & lbi, ubi, lbj, ubj, 1, n(ng), scale, &
1706# ifdef MASKING
1707 & grid(ng) % rmask, &
1708# endif
1709 & ocean(ng) % ad_t(:,:,:,i3d,itrc), &
1710 & minvalue = fmin, &
1711 & maxvalue = fmax)
1712 END IF
1713 IF (status.ne.pio_noerr) THEN
1714 IF (master) THEN
1715 WRITE (stdout,20) trim(vname(1,idtvar(itrc))), &
1716 & trim(label), trim(s(ng)%name), outrec
1717 END IF
1718 exit_flag=3
1719 ioerror=status
1720 RETURN
1721 ELSE
1722 IF (master) THEN
1723 WRITE (stdout,30) trim(vname(2,idtvar(itrc))), &
1724 & fmin, fmax
1725 END IF
1726 END IF
1727 END DO
1728
1729# ifdef ADJUST_BOUNDARY
1730!
1731! Write out tracers open boundaries.
1732!
1733 DO itrc=1,nt(ng)
1734 IF (any(lobc(:,istvar(itrc),ng))) THEN
1735 scale=1.0_dp
1736 ifield=idsbry(istvar(itrc))
1737 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1738 iodesc => iodesc_dp_r3dobc(ng)
1739 ELSE
1740 iodesc => iodesc_sp_r3dobc(ng)
1741 END IF
1742!
1743 IF (model.eq.itlm) THEN
1744 status=nf_fwrite3d_bry(ng, model, s(ng)%name, piofile, &
1745 & vname(1,ifield), &
1746 & s(ng)%pioVar(ifield), &
1747 & outrec, iodesc, &
1748 & lbij, ubij, 1, n(ng), nbrec(ng), &
1749 & scale, &
1750 & boundary(ng) % tl_t_obc(lbij:,:,:,:,&
1751 & lbout(ng),itrc), &
1752 & minvalue = fmin, &
1753 & maxvalue = fmax)
1754 ELSE IF (model.eq.iadm) THEN
1755 status=nf_fwrite3d_bry(ng, model, s(ng)%name, piofile, &
1756 & vname(1,ifield), &
1757 & s(ng)%pioVar(ifield), &
1758 & outrec, iodesc, &
1759 & lbij, ubij, 1, n(ng), nbrec(ng), &
1760 & scale, &
1761 & boundary(ng) % ad_t_obc(lbij:,:,:,:,&
1762 & lbout(ng),itrc), &
1763 & minvalue = fmin, &
1764 & maxvalue = fmax)
1765 END IF
1766 IF (status.ne.pio_noerr) THEN
1767 IF (master) THEN
1768 WRITE (stdout,20) trim(vname(1,ifield)), &
1769 & trim(label), trim(s(ng)%name), outrec
1770 END IF
1771 exit_flag=3
1772 ioerror=status
1773 RETURN
1774 ELSE
1775 IF (master) THEN
1776 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1777 END IF
1778 END IF
1779 END IF
1780 END DO
1781# endif
1782
1783# ifdef ADJUST_STFLUX
1784!
1785! Write out surface net tracers fluxes. Notice that fluxes have their
1786! own fixed time-dimension (of size Nfrec) to allow 4DVAR adjustments
1787! at other times in addition to initialization time.
1788!
1789 DO itrc=1,nt(ng)
1790 IF (lstflux(itrc,ng)) THEN
1791 scale=1.0_dp ! kinematic flux units
1792 ifield=idtsur(itrc)
1793 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1794 iodesc => iodesc_dp_r2dfrc(ng)
1795 ELSE
1796 iodesc => iodesc_sp_r2dfrc(ng)
1797 END IF
1798!
1799 IF (model.eq.itlm) THEN
1800 status=nf_fwrite3d(ng, itlm, piofile, ifield, &
1801 & s(ng)%pioVar(ifield), &
1802 & outrec, iodesc, &
1803 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1804# ifdef MASKING
1805 & grid(ng) % rmask, &
1806# endif
1807 & forces(ng)% tl_tflux(:,:,:, &
1808 & lfout(ng),itrc), &
1809 & minvalue = fmin, &
1810 & maxvalue = fmax)
1811 ELSE IF (model.eq.iadm) THEN
1812 status=nf_fwrite3d(ng, itlm, piofile, ifield, &
1813 & s(ng)%pioVar(ifield), &
1814 & outrec, iodesc, &
1815 & lbi, ubi, lbj, ubj, 1, nfrec(ng), scale, &
1816# ifdef MASKING
1817 & grid(ng) % rmask, &
1818# endif
1819 & forces(ng)% ad_tflux(:,:,:, &
1820 & lfout(ng),itrc), &
1821 & minvalue = fmin, &
1822 & maxvalue = fmax)
1823 END IF
1824 IF (status.ne.pio_noerr) THEN
1825 IF (master) THEN
1826 WRITE (stdout,20) trim(vname(1,ifield)), &
1827 & trim(label), trim(s(ng)%name), lfout(ng)
1828 END IF
1829 exit_flag=3
1830 ioerror=status
1831 RETURN
1832 ELSE
1833 IF (master) THEN
1834 WRITE (stdout,30) trim(vname(2,ifield)), fmin, fmax
1835 END IF
1836 END IF
1837 END IF
1838 END DO
1839# endif
1840# endif
1841!
1842!-----------------------------------------------------------------------
1843! Synchronize tangent NetCDF file to disk to allow other processes
1844! to access data immediately after it is written. Then, close NetCDF
1845! file.
1846!-----------------------------------------------------------------------
1847!
1848 CALL pio_netcdf_sync (ng, itlm, s(ng)%name, piofile)
1849 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1850!
1851 CALL pio_netcdf_close (ng, itlm, piofile, s(ng)%name, .false.)
1852
1853# ifdef PROFILE
1854!
1855! Turn off time wall clock.
1856!
1857 CALL wclock_off (ng, model, 81, __line__, myfile)
1858# endif
1859!
1860 10 FORMAT (1x,a,': WRT_STATE_PIO - ',a,t75,a, &
1861 & /,23x,'(Grid ',i2.2,', t = ',a,', File: ',a, &
1862 & ', Index=',i1,', Rec=',i0,')', &
1863 & /,23x,'Called from ''',a,'''')
1864 20 FORMAT (/,' WRT_STATE_PIO - error while writing variable: ',a, &
1865 & /,17x,'into ',a,' NetCDF file: ',a, &
1866 & /,17x,'for time record: ',i0)
1867 30 FORMAT (20x,'- ',a,/,23x,'(Min = ',1p,e15.8, &
1868 & ' Max = ',1p,e15.8,')')
1869!
1870 RETURN
1871 END SUBROUTINE wrt_state_pio
1872# endif
1873#endif
1874 END MODULE wrt_state_mod
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
type(t_boundary), dimension(:), allocatable boundary
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
character(len=256) calledfrom
integer ioerror
type(t_io), dimension(:), allocatable his
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer idubar
integer idvvel
integer idvsms
integer isvvel
integer, dimension(:), allocatable idsbry
integer, parameter io_pio
Definition mod_ncparam.F:96
integer isvbar
integer, dimension(:), allocatable idtsur
integer idfsur
integer, dimension(:), allocatable idtvar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer iduvel
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer isubar
integer idusms
integer idvbar
subroutine, public netcdf_close(ng, model, ncid, ncname, lupdate)
subroutine, public netcdf_open(ng, model, ncname, omode, ncid)
subroutine, public netcdf_sync(ng, model, ncname, ncid)
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer, parameter inlm
Definition mod_param.F:662
integer, parameter irpm
Definition mod_param.F:664
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter r3dvar
Definition mod_param.F:721
integer, parameter iadm
Definition mod_param.F:665
integer, parameter u3dvar
Definition mod_param.F:722
integer, parameter u2dvar
Definition mod_param.F:718
integer ngrids
Definition mod_param.F:113
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
integer, parameter v3dvar
Definition mod_param.F:723
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dobc
subroutine, public pio_netcdf_sync(ng, model, ncname, piofile)
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dobc
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dobc
subroutine, public pio_netcdf_close(ng, model, piofile, ncname, lupdate)
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dobc
subroutine, public pio_netcdf_open(ng, model, ncname, omode, piofile)
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dfrc
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dobc
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dfrc
real(dp), parameter day2sec
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
real(dp), parameter sec2day
integer exit_flag
logical, dimension(:), allocatable lwrtper
integer, dimension(:), allocatable nbrec
integer noerror
integer, dimension(:), allocatable lbout
integer, dimension(:), allocatable lfout
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, private wrt_state_nf90(ng, tile, model, label, outrec, i2d, i3d, stime, lbij, ubij, lbi, ubi, lbj, ubj, s)
Definition wrt_state.F:161
subroutine, private wrt_state_pio(ng, tile, model, label, outrec, i2d, i3d, stime, ifdef adjust_boundary
Definition wrt_state.F:981
subroutine, public wrt_state(ng, tile, model, label, outrec, i2d, i3d, stime, s)
Definition wrt_state.F:84
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