45 character(len=*),
intent(in) :: filename
46 character(len=*),
intent(in) :: fldname
50 integer(ip_i4_p) :: ncid
51 integer(ip_i4_p) :: varid
52 integer(ip_i4_p) :: status
53 character(len=*),
parameter :: subname =
'(oasis_io_varexists)' 59 inquire(file=trim(filename),exist=exists)
61 status = nf90_open(filename,nf90_nowrite,ncid)
62 status = nf90_inq_varid(ncid,trim(fldname),varid)
63 if (status == nf90_noerr)
then 66 status = nf90_close(ncid)
81 character(len=*),
intent(in) :: filename
82 type(mct_avect) ,
intent(inout) :: av
83 type(mct_gsmap) ,
intent(in) :: gsmap
84 integer(ip_i4_p),
intent(in) :: mpicom
85 character(len=*),
intent(in) :: avfld
86 character(len=*),
intent(in) :: filefld
87 character(len=*),
intent(in),
optional :: fldtype
90 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
91 integer(ip_i4_p) :: nx
92 integer(ip_i4_p) :: ny
93 type(mct_avect) :: av_g
94 integer(ip_i4_p) :: master_task,iam,ierr
95 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
96 integer(ip_i4_p) :: dlen
97 integer(ip_i4_p) :: status
99 real(ip_double_p),
allocatable :: array2(:,:)
100 integer(ip_i4_p) ,
allocatable :: array2i(:,:)
101 integer(ip_i4_p) :: ifldtype
103 character(len=*),
parameter :: subname =
'(oasis_io_read_avfld)' 109 IF (mpicom == mpi_comm_null)
return 114 if (len_trim(filename) < 1)
then 120 call mpi_comm_rank(mpicom,iam,ierr)
123 if (
present(fldtype))
then 125 if (trim(fldtype) ==
'int') ifldtype = 1
126 if (trim(fldtype) ==
'real') ifldtype = 2
127 if (ifldtype == 0)
then 128 WRITE(
nulprt,*) subname,
estr,
'in fldtype argument' 133 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
135 if (iam == master_task)
then 137 inquire(file=trim(filename),exist=exists)
139 status = nf90_open(trim(filename),nf90_nowrite,ncid)
140 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
143 write(
nulprt,*) subname,
estr,
'file missing ',trim(filename)
147 status = nf90_inq_varid(ncid,trim(filefld),varid)
148 if (status /= nf90_noerr)
then 149 write(
nulprt,*) subname,
':',trim(nf90_strerror(status))
150 WRITE(
nulprt,*) subname,
estr,
'filefld variable not found '//trim(filefld)
153 status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid2)
154 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
157 write(
nulprt,*) subname,
estr,
'variable ndims gt 2 ',trim(filefld),dlen
160 status = nf90_inquire_dimension(ncid,dimid2(1),len=nx)
161 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
165 status = nf90_inquire_dimension(ncid,dimid2(2),len=ny)
166 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
170 if (
size(av_g%rAttr,dim=2) /= nx*ny)
then 171 WRITE(
nulprt,*) subname,
estr,
'av gsize nx ny mismatch in file :',&
172 trim(filename),
SIZE(av_g%rAttr,dim=2),nx,ny
176 if (ifldtype == 1)
then 177 allocate(array2i(nx,ny))
178 status = nf90_get_var(ncid,varid,array2i)
179 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
182 n = mct_avect_indexia(av_g,trim(avfld))
187 av_g%iAttr(n,n1) = array2i(i,j)
192 allocate(array2(nx,ny))
193 status = nf90_get_var(ncid,varid,array2)
194 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
197 n = mct_avect_indexra(av_g,trim(avfld))
202 av_g%rAttr(n,n1) = array2(i,j)
208 status = nf90_close(ncid)
209 IF (status /= nf90_noerr)
THEN 211 trim(nf90_strerror(status))
216 call mct_avect_scatter(av_g,av,gsmap,master_task,mpicom)
217 if (iam == master_task)
then 218 call mct_avect_clean(av_g)
233 character(len=*),
intent(in) :: rstfile
234 type(mct_avect) ,
intent(in) :: av
235 type(mct_gsmap) ,
intent(in) :: gsmap
236 integer(ip_i4_p),
intent(in) :: mpicom
237 integer(ip_i4_p),
intent(in) :: nx
238 integer(ip_i4_p),
intent(in) :: ny
239 character(len=*),
intent(in),
optional :: nampre
242 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
243 integer(ip_i4_p) :: nxf,nyf
244 type(mct_avect) :: av_g
245 type(mct_string) :: mstring
246 character(ic_med) :: itemc
247 character(ic_med) :: lnampre
248 character(ic_med) :: lstring
249 integer(ip_i4_p) :: master_task,iam,ierr
250 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
251 integer(ip_i4_p) :: dlen
252 integer(ip_i4_p) :: status
254 real(ip_double_p),
allocatable :: array2(:,:)
256 character(len=*),
parameter :: subname =
'(oasis_io_write_avfile)' 262 IF (mpicom == mpi_comm_null)
return 268 if (len_trim(rstfile) < 1)
then 274 if (
present(nampre))
then 275 lnampre = trim(nampre)
279 call mpi_comm_rank(mpicom,iam,ierr)
281 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
283 if (iam == master_task)
then 284 if (
size(av_g%rAttr,dim=2) /= nx*ny)
then 285 write(
nulprt,*) subname,
estr,
'av gsize nx ny mismatch in file :',&
286 trim(lstring),
SIZE(av_g%rAttr,dim=2),nx,ny
290 inquire(file=trim(rstfile),exist=exists)
292 status = nf90_open(trim(rstfile),nf90_write,ncid)
293 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
295 status = nf90_redef(ncid)
297 status = nf90_create(trim(rstfile),nf90_clobber,ncid)
298 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
302 do n = 1,mct_avect_nrattr(av_g)
303 call mct_avect_getrlist(mstring,n,av_g)
304 itemc = mct_string_tochar(mstring)
305 itemc = trim(lnampre)//trim(itemc)
306 call mct_string_clean(mstring)
308 status = nf90_inq_dimid(ncid,trim(itemc)//
'_nx',dimid2(1))
309 if (status /= nf90_noerr)
then 310 status = nf90_def_dim(ncid,trim(itemc)//
'_nx',nx,dimid2(1))
313 status = nf90_inq_dimid(ncid,trim(itemc)//
'_ny',dimid2(2))
314 if (status /= nf90_noerr)
then 315 status = nf90_def_dim(ncid,trim(itemc)//
'_ny',ny,dimid2(2))
318 status = nf90_inquire_dimension(ncid,dimid2(1),len=dlen)
319 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
322 write(
nulprt,*) subname,
wstr,
'dlen ne nx ',dlen,nx
327 status = nf90_inquire_dimension(ncid,dimid2(2),len=dlen)
328 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
331 write(
nulprt,*) subname,
wstr,
'dlen ne ny ',dlen,ny
336 status = nf90_inq_varid(ncid,trim(itemc),varid)
337 if (status /= nf90_noerr)
then 338 status = nf90_def_var(ncid,trim(itemc),nf90_double,dimid2,varid)
339 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
345 status = nf90_enddef(ncid)
346 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
351 do n = 1,mct_avect_nrattr(av_g)
352 call mct_avect_getrlist(mstring,n,av_g)
353 itemc = mct_string_tochar(mstring)
354 itemc = trim(lnampre)//trim(itemc)
355 call mct_string_clean(mstring)
356 status = nf90_inq_varid(ncid,trim(itemc),varid)
357 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
360 status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid2)
361 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
363 status = nf90_inquire_dimension(ncid,dimid2(1),len=nxf)
364 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
366 status = nf90_inquire_dimension(ncid,dimid2(2),len=nyf)
367 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
369 if (dlen /= 2 .or. nx*ny /= nxf*nyf)
then 370 WRITE(
nulprt,*) subname,
estr,
'ndims and size does not match on file' 373 allocate(array2(nxf,nyf))
380 array2(i,j) = av_g%rAttr(n,n1)
384 status = nf90_put_var(ncid,varid,array2)
385 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
389 call mct_avect_clean(av_g)
391 status = nf90_close(ncid)
392 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
409 character(len=*),
intent(in) :: rstfile
410 type(mct_avect) ,
intent(inout) :: av
411 type(mct_gsmap) ,
intent(in) :: gsmap
412 integer(ip_i4_p),
intent(in) :: mpicom
413 logical ,
intent(in) ,
optional :: abort
414 character(len=*),
intent(in) ,
optional :: nampre
415 logical ,
intent(out),
optional :: didread
418 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
419 integer(ip_i4_p) :: nx
420 integer(ip_i4_p) :: ny
421 type(mct_avect) :: av_g
422 type(mct_string) :: mstring
423 character(ic_med) :: itemc
424 character(ic_med) :: lnampre
425 character(ic_med) :: lstring
426 integer(ip_i4_p) :: master_task,iam,ierr
427 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
428 integer(ip_i4_p) :: dlen
429 integer(ip_i4_p) :: status
432 real(ip_double_p),
allocatable :: array2(:,:)
434 character(len=*),
parameter :: subname =
'(oasis_io_read_avfile)' 440 IF (mpicom == mpi_comm_null)
return 444 if (
present(didread)) didread = .false.
448 if (len_trim(rstfile) < 1)
then 454 if (
present(abort))
then 459 if (
present(nampre))
then 460 lnampre = trim(nampre)
464 call mpi_comm_rank(mpicom,iam,ierr)
466 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
468 if (iam == master_task)
then 470 inquire(file=trim(rstfile),exist=exists)
471 if (.not.exists)
then 473 write(
nulprt,*) subname,
estr,
'file missing ',trim(rstfile)
476 write(
nulprt,*) subname,
wstr,
'file missing ',trim(rstfile)
480 status = nf90_open(trim(rstfile),nf90_nowrite,ncid)
481 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
484 do n = 1,mct_avect_nrattr(av_g)
485 call mct_avect_getrlist(mstring,n,av_g)
486 itemc = mct_string_tochar(mstring)
487 itemc = trim(lnampre)//trim(itemc)
488 call mct_string_clean(mstring)
490 status = nf90_inq_varid(ncid,trim(itemc),varid)
492 if (status /= nf90_noerr)
then 494 write(
nulprt,*) subname,
estr,
'var missing on file = ',trim(itemc),
':',trim(nf90_strerror(status))
502 status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid2)
503 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
506 write(
nulprt,*) subname,
estr,
'variable ndims gt 2 on file ',trim(itemc),dlen
509 status = nf90_inquire_dimension(ncid,dimid2(1),len=nx)
510 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
514 status = nf90_inquire_dimension(ncid,dimid2(2),len=ny)
515 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
519 if (
size(av_g%rAttr,dim=2) /= nx*ny)
then 520 WRITE(
nulprt,*) subname,
estr,
'av gsize nx ny mismatch in file = ',&
521 trim(rstfile),
SIZE(av_g%rAttr,dim=2),nx,ny
525 allocate(array2(nx,ny))
527 status = nf90_get_var(ncid,varid,array2)
528 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
535 av_g%rAttr(n,n1) = array2(i,j)
538 if (
present(didread)) didread = .true.
544 status = nf90_close(ncid)
545 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
551 call mct_avect_scatter(av_g,av,gsmap,master_task,mpicom)
552 if (iam == master_task)
then 553 call mct_avect_clean(av_g)
568 character(len=*),
intent(in) :: rstfile
569 integer(ip_i4_p),
intent(in) :: mpicom
570 integer(ip_i4_p),
intent(inout),
optional :: iarray(:)
571 character(len=*),
intent(in),
optional :: ivarname
572 real(ip_double_p),
intent(inout),
optional :: rarray(:)
573 character(len=*),
intent(in),
optional :: rvarname
574 logical ,
intent(in),
optional :: abort
577 integer(ip_i4_p) :: ncnt
578 integer(ip_i4_p) :: master_task,iam,ierr
579 integer(ip_i4_p) :: ncid,dimid,dimid1(1),varid
580 integer(ip_i4_p) :: dlen
581 integer(ip_i4_p) :: status
585 character(len=*),
parameter :: subname =
'(oasis_io_read_array)' 591 if (mpicom == mpi_comm_null)
return 597 if (len_trim(rstfile) < 1)
then 603 if (
present(abort))
then 608 call mpi_comm_rank(mpicom,iam,ierr)
610 if (iam == master_task)
then 612 inquire(file=trim(rstfile),exist=exists)
613 if (.not.exists)
then 615 write(
nulprt,*) subname,
estr,
'file missing ',trim(rstfile)
618 write(
nulprt,*) subname,
wstr,
'file missing ',trim(rstfile)
622 status = nf90_open(trim(rstfile),nf90_nowrite,ncid)
623 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
626 if (
present(iarray))
then 627 if (.not.
present(ivarname))
then 628 write(
nulprt,*) subname,
estr,
'iarray must have ivarname set' 634 status = nf90_inq_varid(ncid,trim(ivarname),varid)
635 if (status /= nf90_noerr)
then 637 write(
nulprt,*) subname,
estr,
'var missing on file = ',trim(ivarname),
':',trim(nf90_strerror(status))
644 status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid1)
645 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
648 write(
nulprt,*) subname,
estr,
'variable ndims ne 1 ',trim(ivarname),dlen
651 status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
652 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
655 if (ncnt /= dlen)
then 656 write(
nulprt,*) subname,
estr,
'iarray ncnt dlen mismatch ',ncnt,dlen
660 status = nf90_get_var(ncid,varid,iarray)
661 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
666 if (
present(rarray))
then 667 if (.not.
present(rvarname))
then 668 write(
nulprt,*) subname,
estr,
'rarray must have rvarname set' 674 status = nf90_inq_varid(ncid,trim(rvarname),varid)
675 if (status /= nf90_noerr)
then 677 write(
nulprt,*) subname,
estr,
'var missing on file = ',trim(rvarname),
':',trim(nf90_strerror(status))
684 status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid1)
685 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
688 write(
nulprt,*) subname,
estr,
'variable ndims ne 1 ',trim(rvarname),dlen
691 status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
692 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
695 if (ncnt /= dlen)
then 696 write(
nulprt,*) subname,
estr,
'rarray ncnt dlen mismatch ',ncnt,dlen
700 status = nf90_get_var(ncid,varid,rarray)
701 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
706 status = nf90_close(ncid)
707 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
713 if (
present(iarray))
then 714 call oasis_mpi_bcast(iarray,mpicom,trim(subname)//
':iarray',master_task)
717 if (
present(rarray))
then 718 call oasis_mpi_bcast(rarray,mpicom,trim(subname)//
':rarray',master_task)
733 character(len=*),
intent(in) :: rstfile
734 integer(ip_i4_p),
intent(in) :: mpicom
735 integer(ip_i4_p),
intent(in),
optional :: iarray(:)
736 character(len=*),
intent(in),
optional :: ivarname
737 real(ip_double_p),
intent(in),
optional :: rarray(:)
738 character(len=*),
intent(in),
optional :: rvarname
741 integer(ip_i4_p) :: ncnt
742 integer(ip_i4_p) :: master_task,iam,ierr
743 integer(ip_i4_p) :: ncid,dimid,dimid1(1),varid
744 integer(ip_i4_p) :: dlen
745 integer(ip_i4_p) :: status
748 character(len=*),
parameter :: subname =
'(oasis_io_write_array)' 754 IF (mpicom == mpi_comm_null)
return 760 if (len_trim(rstfile) < 1)
then 766 call mpi_comm_rank(mpicom,iam,ierr)
768 if (iam == master_task)
then 770 inquire(file=trim(rstfile),exist=exists)
772 status = nf90_open(trim(rstfile),nf90_write,ncid)
773 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
775 status = nf90_redef(ncid)
777 status = nf90_create(trim(rstfile),nf90_clobber,ncid)
778 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
782 if (
present(iarray))
then 783 if (.not.
present(ivarname))
then 784 write(
nulprt,*) subname,
estr,
'iarray must have ivarname set' 790 status = nf90_inq_dimid(ncid,trim(ivarname)//
'_ncnt',dimid1(1))
791 if (status /= nf90_noerr)
then 792 status = nf90_def_dim(ncid,trim(ivarname)//
'_ncnt',ncnt,dimid1(1))
795 status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
796 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
798 if (dlen /= ncnt)
then 799 write(
nulprt,*) subname,
estr,
'iarray dlen ne ncnt ',dlen,ncnt
803 status = nf90_inq_varid(ncid,trim(ivarname),varid)
804 if (status /= nf90_noerr)
then 805 status = nf90_def_var(ncid,trim(ivarname),nf90_int,dimid1,varid)
806 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
811 if (
present(rarray))
then 812 if (.not.
present(rvarname))
then 813 write(
nulprt,*) subname,
estr,
'rarray must have rvarname set' 819 status = nf90_inq_dimid(ncid,trim(rvarname)//
'_ncnt',dimid1(1))
820 if (status /= nf90_noerr)
then 821 status = nf90_def_dim(ncid,trim(rvarname)//
'_ncnt',ncnt,dimid1(1))
824 status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
825 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
827 if (dlen /= ncnt)
then 828 write(
nulprt,*) subname,
estr,
'rarray dlen ne ncnt ',dlen,ncnt
832 status = nf90_inq_varid(ncid,trim(rvarname),varid)
833 if (status /= nf90_noerr)
then 834 status = nf90_def_var(ncid,trim(rvarname),nf90_double,dimid1,varid)
835 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
840 status = nf90_enddef(ncid)
841 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
844 if (
present(iarray))
then 845 status = nf90_inq_varid(ncid,trim(ivarname),varid)
846 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
848 status = nf90_put_var(ncid,varid,iarray)
849 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
853 if (
present(rarray))
then 854 status = nf90_inq_varid(ncid,trim(rvarname),varid)
855 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
857 status = nf90_put_var(ncid,varid,rarray)
858 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
862 status = nf90_close(ncid)
863 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
885 type(mct_avect) ,
intent(in) :: av
886 type(mct_gsmap) ,
intent(in) :: gsmap
887 integer(ip_i4_p),
intent(in) :: mpicom
888 integer(ip_i4_p),
intent(in) :: nx
889 integer(ip_i4_p),
intent(in) :: ny
890 integer(ip_i4_p),
intent(in),
optional :: msec
891 character(len=*),
intent(in),
optional :: f_string
892 character(len=*),
intent(in),
optional :: filename
895 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
896 type(mct_avect) :: av_g
897 type(mct_string) :: mstring
898 character(ic_med) :: itemc
899 character(ic_med) :: lfn
900 character(ic_med) :: lstring
901 integer(ip_i4_p) :: master_task,iam,ierr
902 integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid
903 integer(ip_i4_p) :: start3(3),count3(3)
904 integer(ip_i4_p) :: start1(1),count1(1)
905 integer(ip_i4_p) :: lmsec(1)
906 integer(ip_i4_p) :: dlen
907 integer(ip_i4_p) :: status
910 logical :: whead,wdata
911 integer(ip_i4_p) ,
allocatable :: time(:)
912 real(ip_double_p),
allocatable :: array3(:,:,:)
913 real(ip_double_p) :: tbnds(2)
915 character(len=*),
parameter :: subname =
'(oasis_io_write_avfbf)' 921 IF (mpicom == mpi_comm_null)
return 926 if (
present(msec))
then 931 if (
present(f_string))
then 932 lstring = trim(f_string)
936 call mpi_comm_rank(mpicom,iam,ierr)
940 call oasis_ioshr_wopen(lfn,clobber=.true.,cdf64=.true.)
946 elseif (fk == 2)
then 954 call oasis_ioshr_write(lfn,&
955 time_units=
'seconds',time_cal=
'seconds',time_val=
real(lmsec,ip_double_p),&
956 nt=1,whead=whead,wdata=wdata)
958 call oasis_ioshr_write(lfn, gsmap, av,
'pout', &
959 whead=whead,wdata=wdata,nx=nx,ny=ny,nt=1, &
962 if (fk == 1)
call oasis_ioshr_enddef(lfn)
965 call oasis_ioshr_close(lfn)
968 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
969 if (iam == master_task)
then 970 if (
size(av_g%rAttr,dim=2) /= nx*ny)
then 971 WRITE(
nulprt,*) subname,
estr,
'av gsize nx ny mismatch in file :',&
972 trim(filename),
SIZE(av_g%rAttr,dim=2),nx,ny
976 allocate(array3(nx,ny,1))
977 do n = 1,mct_avect_nrattr(av_g)
978 call mct_avect_getrlist(mstring,n,av_g)
979 itemc = mct_string_tochar(mstring)
980 call mct_string_clean(mstring)
981 if (
present(filename))
then 984 lfn = trim(itemc)//trim(lstring)//
'.nc' 990 array3(i,j,1) = av_g%rAttr(n,n1)
1002 inquire(file=trim(lfn),exist=exists)
1005 status = nf90_open(lfn,nf90_write,ncid)
1006 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1008 status = nf90_inq_dimid(ncid,
'time',dimid)
1009 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1011 status = nf90_inquire_dimension(ncid,dimid,len=dlen)
1012 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1014 allocate(time(dlen))
1015 status = nf90_inq_varid(ncid,
'time',varid)
1016 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1018 status = nf90_get_var(ncid,varid,time)
1019 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1025 if (time(i) >= lmsec(1)) newfile=.true.
1028 if (.not.newfile)
then 1029 start1(1) = dlen + 1
1030 start3(3) = start1(1)
1035 status = nf90_create(lfn,nf90_clobber,ncid)
1036 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1038 status = nf90_def_dim(ncid,
'nx',nx,dimid3(1))
1039 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1041 status = nf90_def_dim(ncid,
'ny',ny,dimid3(2))
1042 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1044 status = nf90_def_dim(ncid,
'time',nf90_unlimited,dimid)
1045 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1048 status = nf90_def_var(ncid,
'time',nf90_int,dimid,varid)
1049 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1051 status = nf90_def_var(ncid,trim(itemc),nf90_double,dimid3,varid)
1052 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1054 status = nf90_enddef(ncid)
1055 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1059 status = nf90_inq_varid(ncid,
'time',varid)
1060 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1062 status = nf90_put_var(ncid,varid,lmsec,start1,count1)
1063 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1065 status = nf90_inq_varid(ncid,trim(itemc),varid)
1066 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1068 status = nf90_put_var(ncid,varid,array3,start3,count3)
1069 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1071 status = nf90_close(ncid)
1072 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1076 call mct_avect_clean(av_g)
1099 type(mct_avect) ,
intent(inout) :: av
1100 type(mct_gsmap) ,
intent(in) :: gsmap
1101 integer(ip_i4_p),
intent(in) :: mpicom
1102 integer(ip_i4_p),
intent(in),
optional :: msec
1103 character(len=*),
intent(in),
optional :: f_string
1104 character(len=*),
intent(in),
optional :: filename
1107 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
1108 integer(ip_i4_p) :: nx,ny
1109 type(mct_avect) :: av_g
1110 type(mct_string) :: mstring
1111 character(ic_med) :: itemc
1112 character(ic_med) :: lfn
1113 character(ic_med) :: lstring
1114 integer(ip_i4_p) :: master_task,iam,ierr
1115 integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid
1116 integer(ip_i4_p) :: start3(3),count3(3)
1117 integer(ip_i4_p) :: lmsec(1)
1118 integer(ip_i4_p) :: dlen
1119 integer(ip_i4_p) :: status
1121 logical :: whead,wdata
1122 real(ip_double_p),
allocatable :: array3(:,:,:)
1123 integer(ip_i4_p) ,
allocatable :: time(:)
1124 real(ip_double_p) :: tbnds(2)
1126 character(len=*),
parameter :: subname =
'(oasis_io_read_avfbf)' 1132 IF (mpicom == mpi_comm_null)
return 1137 if (
present(msec))
then 1142 if (
present(f_string))
then 1143 lstring = trim(f_string)
1147 call mpi_comm_rank(mpicom,iam,ierr)
1149 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
1150 if (iam == master_task)
then 1151 do n = 1,mct_avect_nrattr(av_g)
1152 call mct_avect_getrlist(mstring,n,av_g)
1153 itemc = mct_string_tochar(mstring)
1154 call mct_string_clean(mstring)
1155 if (
present(filename))
then 1156 lfn = trim(filename)
1158 lfn = trim(itemc)//trim(lstring)//
'.nc' 1161 inquire(file=trim(lfn),exist=exists)
1162 if (.not.exists)
then 1163 write(
nulprt,*) subname,
estr,
'file not found ',trim(lfn)
1167 status = nf90_open(lfn,nf90_nowrite,ncid)
1168 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1171 status = nf90_inq_dimid(ncid,
'time',dimid)
1172 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1174 status = nf90_inquire_dimension(ncid,dimid,len=dlen)
1175 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1177 allocate(time(dlen))
1178 status = nf90_inq_varid(ncid,
'time',varid)
1179 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1181 status = nf90_get_var(ncid,varid,time)
1182 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1186 if (time(j) == lmsec(1)) n1 = j
1190 write(
nulprt,*) subname,
estr,
'time not found on file ',trim(lfn),lmsec
1194 status = nf90_inq_varid(ncid,trim(itemc),varid)
1195 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1197 status = nf90_inquire_variable(ncid,varid,dimids=dimid3)
1198 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1200 status = nf90_inquire_dimension(ncid,dimid3(1),len=nx)
1201 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1203 status = nf90_inquire_dimension(ncid,dimid3(2),len=ny)
1204 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1207 if (
size(av_g%rAttr,dim=2) /= nx*ny)
then 1208 WRITE(
nulprt,*) subname,
estr,
'av gsize nx ny mismatch in file :',&
1209 trim(filename),
SIZE(av_g%rAttr,dim=2),nx,ny
1218 allocate(array3(nx,ny,1))
1220 status = nf90_get_var(ncid,varid,array3,start3,count3)
1221 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1223 status = nf90_close(ncid)
1224 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1231 av_g%rAttr(n,n1) = array3(i,j,1)
1240 call mct_avect_scatter(av_g,av,gsmap,master_task,mpicom)
1241 if (iam == master_task)
then 1242 call mct_avect_clean(av_g)
1261 character(len=*) ,
intent(in) :: filename
1262 character(len=*) ,
intent(in) :: fldname
1263 integer(ip_i4_p) ,
intent(inout),
optional :: ifld2(:,:)
1264 real(ip_realwp_p),
intent(inout),
optional :: fld2(:,:)
1265 real(ip_realwp_p),
intent(inout),
optional :: fld3(:,:,:)
1266 integer(ip_i4_p) ,
intent(inout),
optional :: nx
1267 integer(ip_i4_p) ,
intent(inout),
optional :: ny
1268 integer(ip_i4_p) ,
intent(inout),
optional :: nz
1271 integer(ip_i4_p) :: ncid,varid
1272 integer(ip_i4_p) :: n,ndims,xtype
1273 integer(ip_i4_p),
allocatable :: dimid(:),nd(:)
1274 integer(ip_i4_p) :: status
1275 integer(ip_i4_p) :: ind
1277 character(len=ic_med) :: gridname
1279 character(len=*),
parameter :: subname =
'(oasis_io_read_field_fromroot)' 1294 write(
nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname)
1297 inquire(file=trim(filename),exist=exists)
1299 status = nf90_open(filename,nf90_nowrite,ncid)
1300 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1303 write(
nulprt,*) subname,
estr,
'in filename ',trim(filename)
1307 status = nf90_inq_varid(ncid,trim(fldname),varid)
1308 if (status /= nf90_noerr)
then 1309 write(
nulprt,*) subname,
estr,
'in variable name ',trim(fldname)
1313 status = nf90_inquire_variable(ncid,varid,ndims=ndims,xtype=xtype)
1314 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1317 allocate(dimid(ndims),nd(ndims))
1319 status = nf90_inquire_variable(ncid,varid,dimids=dimid)
1320 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1323 status = nf90_inquire_dimension(ncid,dimid(n),len=nd(n))
1324 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1328 if (
present(ifld2) .or.
present(fld2) .or.
present(fld3))
then 1329 if (xtype == nf90_int .and. ndims == 2 .and.
present(ifld2))
then 1330 status = nf90_get_var(ncid,varid,ifld2)
1331 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1333 elseif (xtype /= nf90_int .and. ndims == 2 .and.
present(fld2))
then 1334 status = nf90_get_var(ncid,varid,fld2)
1335 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1337 elseif (xtype /= nf90_int .and. ndims == 3 .and.
present(fld3))
then 1338 status = nf90_get_var(ncid,varid,fld3)
1339 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1342 write(
nulprt,*) subname,
estr,
'mismatch in field and data' 1347 status = nf90_close(ncid)
1348 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1351 if (
present(nx))
then 1355 if (
present(ny))
then 1359 if (
present(nz))
then 1363 deallocate(dimid,nd)
1372 if (
present(ifld2) .or.
present(fld2) .or.
present(fld3))
then 1373 if (xtype == nf90_int .and. ndims == 2 .and.
present(ifld2))
then 1375 elseif (xtype /= nf90_int .and. ndims == 2 .and.
present(fld2))
then 1377 elseif (xtype /= nf90_int .and. ndims == 3 .and.
present(fld3))
then 1404 character(len=*),
intent(in) :: filename
1405 character(len=*),
intent(in) :: fldname
1406 real(ip_realwp_p),
intent(in) :: fld(:,:)
1407 integer(ip_i4_p),
intent(in) :: nx
1408 integer(ip_i4_p),
intent(in) :: ny
1411 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
1412 integer(ip_i4_p) :: status
1413 integer(ip_i4_p) :: ind
1415 character(len=ic_med) :: gridname
1417 character(len=*),
parameter :: subname =
'(oasis_io_write_2dgridfld_fromroot)' 1429 write(
nulprt,*) subname,
' write ',trim(filename),
' ',trim(fldname),nx,ny
1432 ind = index(trim(fldname),
'.')
1434 write(
nulprt,*) subname,
estr,
'in fldname ',trim(fldname)
1437 gridname = fldname(1:ind-1)
1439 inquire(file=trim(filename),exist=exists)
1441 status = nf90_open(filename,nf90_write,ncid)
1442 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1444 status = nf90_redef(ncid)
1446 status = nf90_create(filename,nf90_clobber,ncid)
1447 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1452 status = nf90_inq_dimid(ncid,
'x_'//trim(gridname),dimid2(1))
1453 if (status /= nf90_noerr)
then 1454 status = nf90_def_dim(ncid,
'x_'//trim(gridname),nx,dimid2(1))
1455 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1460 status = nf90_inq_dimid(ncid,
'y_'//trim(gridname),dimid2(2))
1461 if (status /= nf90_noerr)
then 1462 status = nf90_def_dim(ncid,
'y_'//trim(gridname),ny,dimid2(2))
1463 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1468 status = nf90_inq_varid(ncid,trim(fldname),varid)
1469 if (status /= nf90_noerr)
then 1470 status = nf90_def_var(ncid,trim(fldname),nf90_double,dimid2,varid)
1471 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1473 status = nf90_enddef(ncid)
1474 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1476 status = nf90_put_var(ncid,varid,fld)
1477 if (status /= nf90_noerr)
write(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1480 status = nf90_enddef(ncid)
1481 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1485 status = nf90_close(ncid)
1486 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1507 character(len=*),
intent(in) :: filename
1508 character(len=*),
intent(in) :: fldname
1509 integer(ip_i4_p),
intent(in) :: fld(:,:)
1510 integer(ip_i4_p),
intent(in) :: nx
1511 integer(ip_i4_p),
intent(in) :: ny
1514 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
1515 integer(ip_i4_p) :: status
1516 integer(ip_i4_p) :: ind
1518 character(len=ic_med) :: gridname
1520 character(len=*),
parameter :: subname =
'(oasis_io_write_2dgridint_fromroot)' 1532 write(
nulprt,*) subname,
' write ',trim(filename),
' ',trim(fldname),nx,ny
1535 ind = index(trim(fldname),
'.')
1537 write(
nulprt,*) subname,
estr,
'in fldname ',trim(fldname)
1540 gridname = fldname(1:ind-1)
1542 inquire(file=trim(filename),exist=exists)
1544 status = nf90_open(filename,nf90_write,ncid)
1545 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1547 status = nf90_redef(ncid)
1549 status = nf90_create(filename,nf90_clobber,ncid)
1550 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1555 status = nf90_inq_dimid(ncid,
'x_'//trim(gridname),dimid2(1))
1556 if (status /= nf90_noerr)
then 1557 status = nf90_def_dim(ncid,
'x_'//trim(gridname),nx,dimid2(1))
1558 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1563 status = nf90_inq_dimid(ncid,
'y_'//trim(gridname),dimid2(2))
1564 if (status /= nf90_noerr)
then 1565 status = nf90_def_dim(ncid,
'y_'//trim(gridname),ny,dimid2(2))
1566 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1571 status = nf90_inq_varid(ncid,trim(fldname),varid)
1572 if (status /= nf90_noerr)
then 1573 status = nf90_def_var(ncid,trim(fldname),nf90_int,dimid2,varid)
1574 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1576 status = nf90_enddef(ncid)
1577 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1579 status = nf90_put_var(ncid,varid,fld)
1580 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1583 status = nf90_enddef(ncid)
1584 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1588 status = nf90_close(ncid)
1589 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1610 character(len=*),
intent(in) :: filename
1611 character(len=*),
intent(in) :: fldname
1612 real(ip_realwp_p),
intent(in) :: fld(:,:,:)
1613 integer(ip_i4_p),
intent(in) :: nx
1614 integer(ip_i4_p),
intent(in) :: ny
1615 integer(ip_i4_p),
intent(in) :: nc
1618 integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid
1619 integer(ip_i4_p) :: status
1620 integer(ip_i4_p) :: ind
1622 character(len=ic_med) :: gridname
1624 character(len=*),
parameter :: subname =
'(oasis_io_write_3dgridfld_fromroot)' 1636 write(
nulprt,*) subname,
' write ',trim(filename),
' ',trim(fldname),nx,ny,nc
1639 ind = index(trim(fldname),
'.')
1641 write(
nulprt,*) subname,
estr,
'in fldname ',trim(fldname)
1644 gridname = fldname(1:ind-1)
1646 inquire(file=trim(filename),exist=exists)
1648 status = nf90_open(filename,nf90_write,ncid)
1649 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1651 status = nf90_redef(ncid)
1653 status = nf90_create(filename,nf90_clobber,ncid)
1654 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1659 status = nf90_inq_dimid(ncid,
'x_'//trim(gridname),dimid3(1))
1660 if (status /= nf90_noerr)
then 1661 status = nf90_def_dim(ncid,
'x_'//trim(gridname),nx,dimid3(1))
1662 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1667 status = nf90_inq_dimid(ncid,
'y_'//trim(gridname),dimid3(2))
1668 if (status /= nf90_noerr)
then 1669 status = nf90_def_dim(ncid,
'y_'//trim(gridname),ny,dimid3(2))
1670 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1675 status = nf90_inq_dimid(ncid,
'crn_'//trim(gridname),dimid3(3))
1676 if (status /= nf90_noerr)
then 1677 status = nf90_def_dim(ncid,
'crn_'//trim(gridname),nc,dimid3(3))
1678 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1683 status = nf90_inq_varid(ncid,trim(fldname),varid)
1684 if (status /= nf90_noerr)
then 1685 status = nf90_def_var(ncid,trim(fldname),nf90_double,dimid3,varid)
1686 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1688 status = nf90_enddef(ncid)
1689 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1691 status = nf90_put_var(ncid,varid,fld)
1692 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1695 status = nf90_enddef(ncid)
1696 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
1700 status = nf90_close(ncid)
1701 IF (status /= nf90_noerr)
WRITE(
nulprt,*) subname,
' model :',
compid,
' proc :',&
subroutine, public oasis_io_read_array(rstfile, mpicom, iarray, ivarname, rarray, rvarname, abort)
Reads an integer or real field from a file into an array.
Provides a common location for several OASIS variables.
Provides reusable IO routines for OASIS.
integer(kind=ip_i4_p) mpi_comm_map
subroutine, public oasis_io_write_avfile(rstfile, av, gsmap, mpicom, nx, ny, nampre)
Writes all fields from an attribute vector to a file.
integer(kind=ip_intwp_p) nulprt
Generic overloaded interface into MPI broadcast.
subroutine, public oasis_abort(id_compid, cd_routine, cd_message, file, line, rcode)
OASIS abort method, publically available to users.
subroutine, public oasis_io_write_2dgridint_fromroot(filename, fldname, fld, nx, ny)
Write an integer array named field from the root task to a file.
subroutine, public oasis_io_write_2dgridfld_fromroot(filename, fldname, fld, nx, ny)
Write a real array named field from the root task to a file.
integer(kind=ip_i4_p) compid
integer(kind=ip_i4_p) mpi_rank_local
Provides a generic and simpler interface into MPI calls for OASIS.
subroutine, public oasis_debug_enter(string)
Used when a subroutine is entered, write info to log file at some debug level.
subroutine, public oasis_io_write_3dgridfld_fromroot(filename, fldname, fld, nx, ny, nc)
Write a 3d real array named field from the root task to a file.
subroutine, public oasis_io_read_field_fromroot(filename, fldname, ifld2, fld2, fld3, nx, ny, nz)
Read a field on the root task from a file into an array.
Defines parameters for OASIS.
subroutine, public oasis_io_read_avfld(filename, av, gsmap, mpicom, avfld, filefld, fldtype)
Reads single field from a file into an attribute Vector.
subroutine, public oasis_io_write_avfbf(av, gsmap, mpicom, nx, ny, msec, f_string, filename)
Write each field in an attribute vector to an individual files.
integer(kind=ip_i4_p) oasis_debug
character(len= *), parameter, public estr
subroutine, public oasis_io_read_avfbf(av, gsmap, mpicom, msec, f_string, filename)
Read each field in an attribute vector from individual files.
logical function, public oasis_io_varexists(filename, fldname)
Checks whether the var fldname is in the file.
subroutine, public oasis_flush(nu)
Flushes output to file.
IO interfaces based on pio (not supported yet)
subroutine, public oasis_io_write_array(rstfile, mpicom, iarray, ivarname, rarray, rvarname)
Writes a real or integer array to a file.
subroutine, public oasis_debug_exit(string)
Used when a subroutine is exited, write info to log file at some debug level.
character(len= *), parameter, public wstr
subroutine, public oasis_io_read_avfile(rstfile, av, gsmap, mpicom, abort, nampre, didread)
Reads all fields for an attribute vector from a file.