#define __PIO_FILE__ 'piodarray' module piodarray use pio_types, only : file_desc_t, io_desc_t, var_desc_t, pio_noerr, iosystem_desc_t, & pio_iotype_pbinary, pio_iotype_binary, pio_iotype_direct_pbinary, & pio_iotype_netcdf, pio_iotype_pnetcdf, pio_iotype_netcdf4p, pio_iotype_netcdf4c, & PIO_MAX_VAR_DIMS, pio_iotype_vdc2 use pio_kinds use pio_support use pionfwrite_mod, only : write_nf use pionfread_mod, only : read_nf use iompi_mod use alloc_mod use rearrange #ifdef _COMPRESSION use piovdc #endif #ifdef TIMING use perf_mod, only : t_startf, t_stopf !_EXTERNAL #endif #ifndef NO_MPIMOD use mpi !_EXTERNAL #endif implicit none #ifdef NO_MPIMOD include 'mpif.h' !_EXTERNAL #endif private public :: pio_read_darray, pio_write_darray, darray_write_complete, pio_set_buffer_size_limit !> !! @defgroup PIO_write_darray PIO_write_darray !! @brief The overloaded PIO_write_darray writes a distributed array to disk. !< interface PIO_write_darray ! TYPE real,int,double ! DIMS 1,2,3,4,5,6,7 module procedure write_darray_{DIMS}d_{TYPE} end interface !> !! @defgroup PIO_read_darray PIO_read_darray !! @brief The overloaded PIO_read_darray function reads a distributed array from disk. !< interface PIO_read_darray ! TYPE real,int,double ! DIMS 1,2,3,4,5,6,7 module procedure read_darray_{DIMS}d_{TYPE} end interface !> !! @private !< interface add_data_to_buffer ! TYPE real,int,double module procedure add_data_to_buffer_{TYPE} end interface #ifdef _COMPRESSION interface subroutine WriteVDC2Var(iobuf, start, kount, iocomm, ts, lod, reflevel, iotasks, name ) bind(C) use, intrinsic :: iso_c_binding type(c_ptr), intent(in), value :: iobuf integer(c_int), intent(in) :: start(3), kount(3) integer(c_int), intent(in), value :: iocomm, ts, lod, reflevel, iotasks type(c_ptr), intent(in), value :: name end subroutine WriteVDC2Var end interface interface subroutine ReadVDC2Var(iobuf, start, kount, iocomm, ts, lod, reflevel, iotasks, name) bind(C) use, intrinsic :: iso_c_binding type(c_ptr), intent(in), value :: iobuf integer(c_int), intent(in) :: start(3), kount(3) integer(c_int), intent(in), value :: iocomm, ts, lod, reflevel, iotasks type(c_ptr), intent(in), value :: name end subroutine ReadVDC2Var end interface #endif character(len=*), parameter, private :: modName='piodarray' integer :: total_buffsize=0 integer :: pio_buffer_size_limit= 100000000 ! 100MB default #ifdef MEMCHK integer :: msize, rss, mshare, mtext, mstack, lastrss=0 #endif contains subroutine pio_set_buffer_size_limit(limit) integer, intent(in) :: limit if(limit<0) then call piodie(__PIO_FILE__,__LINE__,' bad value to buffer_size_limit',limit) end if pio_buffer_size_limit=limit end subroutine pio_set_buffer_size_limit ! TYPE real,int,double !> !! @public !! @ingroup PIO_write_darray !! @brief Writes a 1D array of type {TYPE}. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The data to be written !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !! @param fillval : An optional fill value to fill holes in the data written !< subroutine write_darray_1d_{TYPE} (File,varDesc,ioDesc, array, iostat, fillval) use pio_msg_mod, only : pio_msg_writedarray ! !DESCRIPTION: ! Writes a 2-d slab of TYPE to a netcdf file. ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! variable descriptor {VTYPE}, dimension(:), target, intent(in) :: & array ! array to be written {VTYPE}, optional, intent(in) :: fillval ! rearrange receiver fill value type(iosystem_desc_t), pointer :: ios integer(i4), intent(out) :: iostat integer :: msg, ierr character(len=*), parameter :: subName=modName//'::write_darray_{TYPE}' ios => file%iosystem if(ios%async_interface .and. .not. ios%ioproc) then msg = PIO_MSG_WRITEDARRAY if(debugasync) print *,__PIO_FILE__,__LINE__, iodesc%async_id if(ios%comp_rank==0) call mpi_send(msg, 1, mpi_integer, ios%ioroot, 1, ios%union_comm, ierr) call mpi_bcast(file%fh, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(vardesc%varid, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(vardesc%rec, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(vardesc%ndims, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(iodesc%async_id, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast({MPITYPE}, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) if(debugasync) print *,__PIO_FILE__,__LINE__, {MPITYPE} if(present(fillval)) then call mpi_bcast(1, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(fillval, 1, {MPITYPE}, ios%compmaster, ios%intercomm, ierr) else call mpi_bcast(0, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) end if if(debugasync) print *,__PIO_FILE__,__LINE__ endif if(debugasync .and. ios%ioproc) print *,__FILE__,__LINE__,iodesc%async_id select case(File%iotype) case(pio_iotype_pbinary, pio_iotype_direct_pbinary) if (present(fillval)) then call write_darray_bin_{TYPE}(File,varDesc,iodesc, array, iostat, fillval) else call write_darray_bin_{TYPE}(File,varDesc,iodesc, array, iostat) endif case(pio_iotype_pnetcdf, pio_iotype_netcdf, pio_iotype_netcdf4c, pio_iotype_netcdf4p) if (present(fillval)) then call write_darray_nf_{TYPE}(File,varDesc,iodesc, array, iostat, fillval) else call write_darray_nf_{TYPE}(File,varDesc,iodesc, array, iostat) endif case(pio_iotype_binary) print *, subName,': IO type not supported' #ifdef _COMPRESSION case(pio_iotype_vdc2) #if ( {ITYPE} == TYPEREAL ) call write_vdc2_real(File, Vardesc, iodesc, array, iostat) #endif #endif end select end subroutine write_darray_1d_{TYPE} ! TYPE real,int,double ! DIMS 2,3,4,5,6,7 !> !! @public !! @ingroup PIO_write_darray !! @brief Writes a {DIMS}D array of type {TYPE}. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The data to be written !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !! @param fillval : An optional fill value to fill holes in the data written !< subroutine write_darray_{DIMS}d_{TYPE} (File,varDesc,ioDesc, array, iostat, fillval) ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! variable descriptor {VTYPE}, intent(in) :: & array{DIMSTR} ! array to be written {VTYPE}, optional, intent(in) :: fillval ! rearrange receiver fill value integer(i4), intent(out) :: iostat {VTYPE} :: transvar(1), dumbvar(0) ! cannot call transfer function with a 0 sized array if(size(array)==0) then call write_darray_1d_{TYPE} (File, varDesc, iodesc, dumbvar, iostat) else if(present(fillval)) then call write_darray_1d_{TYPE} (File, varDesc, iodesc, transfer(array,transvar), iostat, fillval) else call write_darray_1d_{TYPE} (File, varDesc, iodesc, transfer(array,transvar), iostat) end if end subroutine write_darray_{DIMS}d_{TYPE} ! TYPE real,int,double !> !! @public !! @ingroup PIO_read_darray !! @brief !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The read data !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !< subroutine read_darray_1d_{TYPE} (File,varDesc, ioDesc, array, iostat) use pio_msg_mod, only : pio_msg_readdarray ! !DESCRIPTION: ! Reads a 2-d slab of TYPE to a netcdf file. ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! iodecomp descriptor {VTYPE}, dimension(:), intent(out) :: & array ! array to be read integer(i4), intent(out) :: iostat character(len=*), parameter :: subName=modName//'::read_darray_{TYPE}' type(iosystem_desc_t), pointer :: ios integer :: ierr, msg array = 0 ios => File%iosystem if(ios%async_interface .and. .not. ios%ioproc) then msg = PIO_MSG_READDARRAY if(DebugAsync) print *,__PIO_FILE__,__LINE__ if(ios%comp_rank==0) call mpi_send(msg, 1, mpi_integer, ios%ioroot, 1, ios%union_comm, ierr) call mpi_bcast(file%fh, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(vardesc%varid, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(vardesc%rec, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast(iodesc%async_id, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) call mpi_bcast({MPITYPE}, 1, mpi_integer, ios%compmaster, ios%intercomm, ierr) if(DebugAsync) print *,__PIO_FILE__,__LINE__, {MPITYPE} endif select case(File%iotype) case(pio_iotype_pbinary, pio_iotype_direct_pbinary) call read_darray_bin_{TYPE} (File,varDesc,iodesc,array, iostat) case(pio_iotype_pnetcdf, pio_iotype_netcdf, pio_iotype_netcdf4c, pio_iotype_netcdf4p) call read_darray_nf_{TYPE} (File,varDesc,iodesc,array, iostat) #ifdef _COMPRESSION case(pio_iotype_vdc2) #if ( {ITYPE} == TYPEREAL ) call read_vdc2_real(File, Vardesc, iodesc, array, iostat) #endif #endif case(pio_iotype_binary) print *, subName,': IO type not supported' end select end subroutine read_darray_1d_{TYPE} ! TYPE real,int,double ! DIMS 2,3,4,5,6,7 !> !! @public !! @ingroup PIO_read_darray !! @brief Reads a {DIMS}D array of type {TYPE}. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The read data !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !< subroutine read_darray_{DIMS}d_{TYPE} (File,varDesc,ioDesc, array, iostat) ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! iodecomp descriptor {VTYPE}, intent(out) :: array{DIMSTR} ! array to be read integer(i4), intent(out) :: iostat {VTYPE}, pointer :: tmpvar(:) call alloc_check(tmpvar,size(array)) call read_darray_1d_{TYPE} (File, varDesc, iodesc, tmpvar, iostat) array = reshape(tmpvar,shape(array)) call dealloc_check(tmpvar) end subroutine read_darray_{DIMS}d_{TYPE} ! TYPE real,int,double !> !! @private !! @brief Write a 1D array of type {TYPE} defined by varDesc using the decomposition described in iodesc to the netcdf or pnetcdf file File. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The data to be written !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !! @param fillval : An optional fill value to fill holes in the data written !< subroutine write_darray_nf_{TYPE} (File,varDesc,ioDesc,array, iostat, fillval) ! !DESCRIPTION: ! Writes a 2-d slab of TYPE to a netcdf file. ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! io decomp descriptor {VTYPE}, target, intent(in) :: & array(:) ! array to be written {VTYPE}, optional, intent(in) :: fillval ! rearrange receiver fill value integer(pio_offset), pointer :: start(:), count(:) integer(i4), intent(out) :: iostat integer :: request !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character(len=*), parameter :: subName=modName//'::write_darray_nf_{TYPE}' {VTYPE}, dimension(:), pointer :: & IOBUF ! local IO buffer logical (log_kind) :: IOproc ! true if IO processor integer (i4) :: len, &! length of IO decomp segmap iotype, &! type of IO to perform ndims ! number of variable dimensions logical(log_kind) :: UseRearranger #if DEBUG_REARR {VTYPE}, dimension(:), pointer :: array2 integer i #endif {VTYPE} :: rsum integer(i4) :: ierr #ifdef TIMING call t_startf("pio_write_darray") #endif ! ----------------------------------------------------- ! pull information from file_desc_t data structure ! ----------------------------------------------------- IOproc = File%iosystem%IOproc iotype = File%iotype UseRearranger = File%iosystem%UseRearranger ! --------------------------------------------------------- ! pull information out of the decomposition data structure ! --------------------------------------------------------- len = iodesc%IOmap%length if(Debug) print *,__PIO_FILE__,__LINE__,' NAME : IAM: ', & File%iosystem%comp_rank,' UseRearranger: ',UseRearranger,iodesc%glen, iodesc%iomap%start, len #ifdef TIMING call t_startf("pio_rearrange_write") #endif if(UseRearranger) then if (IOproc) then if(Debug) print *, subName,': IAM: ',File%iosystem%comp_rank,'Before call to allocate(IOBUF): ',len, iodesc%write%n_elemtype call alloc_check(IOBUF,len,' TYPE :IOBUF') if (present(fillval)) then IOBUF=fillval else IOBUF= -1.0_r8 endif !------------------------------------------------ ! set the IO buffer to a particular test pattern !------------------------------------------------ !JMD IOBUF(:) = real(File%iosystem%io_rank,kind=r8) if(Debug) print *, subName,': {comp,io}_rank: ',File%iosystem%comp_rank,File%iosystem%io_rank, & 'offset: ',iodesc%iomap%start,'len: ',len !,' IOBUF: ',IOBUF else call alloc_check(IOBUF,0) IOBUF= -1.0_r8 endif !------------------------------------ ! Rearrange data from comp->io decomp !------------------------------------ ! "array" is comp data call rearrange_comp2io(File%iosystem,iodesc, array, iobuf) #if DEBUG_REARR call alloc_check(array2,size(array),'array2') call rearrange_io2comp(File%iosystem,iodesc,IOBUF,array2) do i=1,size(array) if (array(i) /= array2(i)) then print *, subName,': error: write ping-pong test failed on index',i end if end do print *, subName,': passed write ping-pong test' call dealloc_check(array2) !!!!!!! end debug #endif !-------------------------------------------- ! End data rearrange !-------------------------------------------- else if(file%iotype==pio_iotype_pnetcdf) then allocate(iobuf(size(array))) iobuf=array else iobuf=>array end if endif ! if(UseRearranger) #ifdef TIMING call t_stopf("pio_rearrange_write") #endif if (IOproc) then !---------------------------------------------- ! write the global 2-d slice from IO processors !---------------------------------------------- if (DebugIO.and.userearranger.and.len>1) then print *,__PIO_FILE__,__LINE__, & File%iosystem%comp_rank,': write IOBUF r8', & IOBUF(1:2),' ...',IOBUF(len-1:len), & iodesc%glen,len ! write an ascii version write(10+File%iosystem%comp_rank,*) IOBUF(1:len) close(10+File%iosystem%comp_rank) endif ! this is a time dependent multidimensional array if(vardesc%rec>=0 .and. iodesc%start(1)>=0) then ndims = size(iodesc%start)+1 call alloc_check(start,ndims) call alloc_check(count,ndims) start(1:ndims-1)=iodesc%start count(1:ndims-1)=iodesc%count start(ndims:ndims)=vardesc%rec count(ndims:ndims)=1 if(Debug) print *, __PIO_FILE__,__LINE__,'start:',start, & ' count:',count,' ndims:',ndims, minval(iobuf),maxval(iobuf) ! this is a timedependent single value else if(vardesc%rec>=0) then call alloc_check(start,1) call alloc_check(count,1) start(1) = int(vardesc%rec,kind=PIO_Offset) count(1) = 1_PIO_Offset if(Debug) print *, __PIO_FILE__,__LINE__,'start:',start,' count:',count ! this is a non-timedependent array else ndims = size(iodesc%start) call alloc_check(start,ndims) call alloc_check(count,ndims) start=iodesc%start count=iodesc%count if(Debug) print *, __PIO_FILE__,__LINE__,'start:',start,' count:',count,' ndims:',ndims end if else ! some compilers have problems passing ! unassociated pointers when they are intent in call alloc_check(start, 0) call alloc_check(count, 0) endif #ifdef TIMING call t_startf("pio_write_nf") #endif ierr = write_nf(File,IOBUF,varDesc,iodesc,start,count, request) #ifdef TIMING call t_stopf("pio_write_nf") #endif call dealloc_check(start) call dealloc_check(count) if(IOPROC) then if(file%iotype==pio_iotype_pnetcdf) then call add_data_to_buffer(File, IOBUF, request) else if(Userearranger) then deallocate(iobuf) end if end if ! call MPI_Barrier(File%iosystem%comp_comm,ierr) !-------------------------- ! set the error return code !-------------------------- iostat=ierr !----------------------------------------------------------------------- !EOC #ifdef TIMING call t_stopf("pio_write_darray") #endif end subroutine write_darray_nf_{TYPE} ! TYPE real,int,double !> !! @private !! @brief Write a 1D array of type {TYPE}. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The data to be written !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !! @param fillval : An optional fill value to fill holes in the data written !< subroutine write_darray_bin_{TYPE}(File,varDesc,ioDesc,array, iostat, fillval) ! !DESCRIPTION: ! Writes a 2-d slab of integers to a file ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! varable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! IO decomp descriptor {VTYPE}, dimension(:), intent(in), target :: & array ! array to be written {VTYPE}, optional, intent(in) :: fillval integer (i4), intent(out) :: iostat !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character(len=*), parameter :: subName=modName//'::write_darray_bin_{TYPE}' {VTYPE}, dimension(:), pointer :: & IOBUF ! local IO buffer #if DEBUG_REARR {VTYPE}, dimension(:), pointer :: array2 integer i #endif logical (log_kind) :: IOproc ! true if IO processor integer (i4) :: len, &! length of IO decomp segmap iotype, &! type of IO to perform varID ! variable ID integer (i4) :: ierr logical(log_kind) :: UseRearranger #ifdef TIMING call t_startf("pio_write_darray") #endif ! ----------------------------------------------------- ! pull information from file_desc_t data structure ! ----------------------------------------------------- IOproc = File%iosystem%IOproc iotype = File%iotype UseRearranger = File%iosystem%UseRearranger ! ------------------------------------------------- ! Pull information about the IO decomposition ! ------------------------------------------------- if(Debug) print *, __PIO_FILE__,__LINE__,iodesc%Write%elemTYPE, iodesc%Write%fileTYPE len = iodesc%IOmap%length ierr = pio_noerr if (IOproc) then if(userearranger) then call alloc_check(IOBUF,len,' TYPE :IOBUF') if (present(fillval)) then IOBUF=fillval else IOBUF=-1.0_r8 endif else iobuf=>array end if !------------------------------------------------ ! set the IO buffer to a particular test pattern !------------------------------------------------ !JMD IOBUF(:) = File%iosystem%io_rank if(Debug) print *, subName,': {comp,io}_rank: ',File%iosystem%comp_rank,File%iosystem%io_rank, & 'offset: ',iodesc%iomap%start,'len: ',len ! ,' IOBUF: ',IOBUF else if(userearranger) call alloc_check(IOBUF,0,'write_darray_int:IOBUF') endif !----------------------------------------- !NEED HELP: ! ! Need a call to a data rearranger here ! ! call ESMF_rearrange() ! ! or ! ! call MCT_rearrange() !----------------------------------------- #ifdef TIMING call t_startf("pio_rearrange_write") #endif if(UseRearranger) then !------------------------------------ ! Rearrange data from comp->io decomp !------------------------------------ ! "array" is comp data call rearrange_comp2io(File%iosystem,iodesc,array,IOBUF) #if DEBUG_REARR call alloc_check(array2,size(array),'array2') call rearrange_io2comp(File%iosystem,iodesc,IOBUF,array2) do i=1,size(array) if (array(i) /= array2(i)) then print *, subName,': error: int write ping-pong test failed on index',i end if end do print *, subName,': passed int write ping-pong test' call dealloc_check(array2) !!!!!!! end debug #endif !-------------------------------------------- ! End data rearrange !-------------------------------------------- endif #ifdef TIMING call t_stopf("pio_rearrange_write") #endif if(IOProc) then #ifdef TIMING call t_startf("pio_write_bin") #endif !---------------------------------------------- ! write the global 2-d slice from IO processors !---------------------------------------------- ierr = write_mpiio(File,IOBUF,varDesc,iodesc) #ifdef TIMING call t_stopf("pio_write_bin") #endif endif !-------------------------- ! deallocate the IO buffer !-------------------------- if(userearranger) call dealloc_check(IOBUF) ! call MPI_Barrier(File%iosystem%comp_comm,ierr) !-------------------------- ! set the error return code !-------------------------- iostat=ierr !----------------------------------------------------------------------- !EOC #ifdef TIMING call t_stopf("pio_write_darray") #endif end subroutine write_darray_bin_{TYPE} ! TYPE real,int,double !> !! @private !! @brief Read a 1D array of type {TYPE} defined by varDesc using the decomposition described in ioDesc to the netcdf or pnetcdf file File. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The read data !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !< subroutine read_darray_nf_{TYPE} (File,varDesc,ioDesc,array, iostat) ! ! !DESCRIPTION: ! Reads a 2-d horizontal slice of integers from a binary file ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! info about data file type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! io decomp descriptor ! !INPUT/OUTPUT PARAMETERS: {VTYPE}, dimension(:), intent(out), target :: & array ! array to be read integer (i4), intent(out) :: iostat !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character(len=*), parameter :: subName=modName//'::read_darray_nf_{TYPE}' {VTYPE}, dimension(:), pointer :: & IOBUF ! local IO buffer logical (log_kind) :: IOproc ! true if an IO processor integer (i4) :: len, &! local length of IO decomp iotype, &! type of IO to perform ndims integer(pio_offset), dimension(PIO_MAX_VAR_DIMS) :: start, count logical(log_kind), parameter :: Debug = .FALSE. logical(log_kind) :: UseRearranger integer(i4) :: ierr #if DEBUG_REARR {VTYPE}, dimension(:), pointer :: iobuf2 integer i #endif #ifdef TIMING call t_startf("pio_read_darray") #endif ! ----------------------------------------------------- ! pull information from file_desc_t data structure ! ----------------------------------------------------- IOproc = File%iosystem%IOproc iotype = File%iotype UseRearranger = File%iosystem%UseRearranger ierr = PIO_NOERR ! ----------------------------------------------------- ! Pull out information of DecompMap_t data structure ! ----------------------------------------------------- len = iodesc%IOmap%length start=0 count=0 if (IOproc) then !----------------------------- ! allocate temporary IO buffer !----------------------------- if(userearranger) then if(File%iotype == pio_iotype_netcdf .and. file%iosystem%io_rank==0) then call alloc_check(IOBUF,iodesc%maxiobuflen,' TYPE :IOBUF') else call alloc_check(IOBUF,len,' TYPE :IOBUF') end if else iobuf=>array end if !---------------------------------------------- ! read the global 2-d slice to IO processors !---------------------------------------------- ! this is a time dependent multidimensional array if(vardesc%rec>=0 .and. iodesc%start(1)>=0) then ndims = size(iodesc%start)+1 start(1:ndims-1)=iodesc%start count(1:ndims-1)=iodesc%count start(ndims:ndims)=vardesc%rec count(ndims:ndims)=1 ! this is a timedependent single value else if(vardesc%rec>=0) then ndims=1 start(1) = int(vardesc%rec,kind=PIO_Offset) count(1) = 1_PIO_Offset ! this is a non-timedependent array else ndims = size(iodesc%start) start(1:ndims)=iodesc%start count(1:ndims)=iodesc%count end if else ndims=1 if(userearranger) then call alloc_check(IOBUF,0,'IOBUF') end if endif #ifdef TIMING call t_startf("pio_read_nf") #endif ierr = read_nf(File,IOBUF,varDesc,iodesc,start(1:ndims),count(1:ndims)) #ifdef TIMING call t_stopf("pio_read_nf") #endif if(DebugIO) print *, subName,': {comp,io}_rank: ',File%iosystem%comp_rank,File%iosystem%io_rank, & 'offset: ',iodesc%iomap%start,'len: ',len !,' IOBUF: ',IOBUF #ifdef TIMING call t_startf("pio_rearrange_read") #endif if(UseRearranger) then !------------------------------------ ! Rearrange data from io->comp decomp !------------------------------------ ! "array" is comp data call rearrange_io2comp(File%iosystem,iodesc,IOBUF,array) #if DEBUG_REARR call alloc_check(iobuf2,size(IOBUF),'iobuf2') call rearrange_comp2io(File%iosystem,iodesc,array,iobuf2) do i=1,size(iobuf) if (iobuf(i) /= iobuf2(i)) then print *, subName,': error: int read ping-pong test failed on index',i end if end do print *, subName,': passed int read ping-pong test' call dealloc_check(iobuf2) !!!!!!! end debug #endif ! -------------------------- ! deallocate IO buffer ! -------------------------- call dealloc_check(IOBUF) endif #ifdef TIMING call t_stopf("pio_rearrange_read") #endif !---------------- ! set errror code !---------------- iostat = ierr !----------------------------------------------------------------------- !EOC #ifdef TIMING call t_stopf("pio_read_darray") #endif end subroutine read_darray_nf_{TYPE} ! TYPE real,int,double !> !! @private !! @brief Read an array of type {TYPE} defined by varDesc using the decomposition !! described in ioDesc to the netcdf or pnetcdf file File. !! @details !! @param File @copydoc file_desc_t !! @param varDesc @copydoc var_desc_t !! @param ioDesc @copydoc io_desc_t !! @param array : The read data !! @param iostat : The status returned from this routine (see \ref PIO_seterrorhandling for details) !< subroutine read_darray_bin_{TYPE} (File,varDesc,ioDesc,array, iostat) ! ! !DESCRIPTION: ! Reads a 2-d horizontal slice of integers from a binary file ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! info about data file type (var_desc_t), intent(inout) :: & varDesc ! variable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! io decomp descriptor ! !INPUT/OUTPUT PARAMETERS: {VTYPE}, dimension(:), intent(out), target :: & array ! array to be read integer (i4), intent(out) :: iostat !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character(len=*), parameter :: subName=modName//'::read_darray_bin_{TYPE}' {VTYPE}, dimension(:), pointer :: & IOBUF ! local IO buffer logical (log_kind) :: IOproc ! true if an IO processor integer (i4) :: len, &! local length of IO decomp iotype ! type of IO to perform logical(log_kind), parameter :: Debug = .FALSE. logical(log_kind) :: UseRearranger integer(i4) :: ierr #if DEBUG_REARR {VTYPE}, dimension(:), pointer :: iobuf2 integer i #endif #ifdef TIMING call t_startf("pio_read_darray") #endif ! ----------------------------------------------------- ! pull information from file_desc_t data structure ! ----------------------------------------------------- IOproc = File%iosystem%IOproc iotype = File%iotype UseRearranger = File%iosystem%UseRearranger ! ----------------------------------------------------- ! Pull out information of DecompMap_t data structure ! ----------------------------------------------------- ! len = iodesc%IOmap%length len = iodesc%Read%n_words ierr = pio_noerr if (IOproc) then !----------------------------- ! allocate temporary IO buffer !----------------------------- if(userearranger) then call alloc_check(IOBUF,len,'read_darray_ :IOBUF') else iobuf=>array end if !---------------------------------------------- ! read the global 2-d slice to IO processors !---------------------------------------------- if(iotype.eq.pio_iotype_binary) then print *, subName,': TYPE : IO type not supported' iostat =-1 return end if #ifdef TIMING call t_startf("pio_read_bin") #endif ierr = read_mpiio(File,IOBUF,varDesc,iodesc) #ifdef TIMING call t_stopf("pio_read_bin") #endif if(DebugIO) print *, subName,': TYPE: {comp,io}_rank: ',File%iosystem%comp_rank,File%iosystem%io_rank, & 'len: ',len,' IOBUF: ',IOBUF(1:4) else if(userearranger) then call alloc_check(IOBUF,0,'IOBUF') endif #ifdef TIMING call t_startf("pio_rearrange_read") #endif if(UseRearranger) then !------------------------------------ ! Rearrange data from io->comp decomp !------------------------------------ ! "array" is comp data call rearrange_io2comp(File%iosystem,iodesc,IOBUF,array) #if DEBUG_REARR call alloc_check(iobuf2,size(IOBUF),'iobuf2') call rearrange_comp2io(File%iosystem,iodesc,array,iobuf2) do i=1,size(iobuf) if (iobuf(i) /= iobuf2(i)) then print *, subName,': error: int read ping-pong test failed on index',i end if end do print *, subName,': passed int read ping-pong test' call dealloc_check(iobuf2) !!!!!!! end debug #endif ! -------------------------- ! deallocate IO buffer ! -------------------------- call dealloc_check(IOBUF) endif #ifdef TIMING call t_stopf("pio_rearrange_read") #endif !---------------- ! set errror code !---------------- iostat = ierr !----------------------------------------------------------------------- !EOC #ifdef TIMING call t_stopf("pio_read_darray") #endif end subroutine read_darray_bin_{TYPE} ! TYPE real,int,double subroutine add_data_to_buffer_{TYPE} (File, IOBUF, request) use pio_types, only : io_data_list type(file_desc_t) :: File {VTYPE}, pointer :: IOBUF(:) integer, intent(in) :: request integer :: cnt, mpierr, maxbuffsize, this_buffsize type(io_data_list), pointer :: ptr #ifdef MEMCHK call GPTLget_memusage(msize, rss, mshare, mtext, mstack) if(rss>lastrss) then lastrss=rss print *,__PIO_FILE__,__LINE__,'mem=',rss end if #endif if(.not. associated(File%data_list_top)) then allocate(file%data_list_top) ptr => file%data_list_top else ptr => file%data_list_top do while(associated(ptr%next)) ptr => ptr%next end do allocate(ptr%next) ptr=>ptr%next nullify(ptr%next) end if ptr%request = request File%request_cnt=File%request_cnt+1 ptr%data_{TYPE} => IOBUF this_buffsize = size(iobuf)*sizeof(iobuf(1)) file%buffsize=file%buffsize+this_buffsize total_buffsize = total_buffsize+this_buffsize call MPI_ALLREDUCE(total_buffsize,maxbuffsize,1,MPI_INTEGER,MPI_MAX,file%iosystem%io_comm, mpierr) if(maxbuffsize > pio_buffer_size_limit) then call darray_write_complete(File) endif if(debug) print *,__FILE__,__LINE__,'buffsize = ',file%buffsize,file%request_cnt #ifdef MEMCHK call GPTLget_memusage(msize, rss, mshare, mtext, mstack) if(rss>lastrss) then lastrss=rss print *,__PIO_FILE__,__LINE__,'mem=',rss end if #endif end subroutine add_data_to_buffer_{TYPE} subroutine darray_write_complete(File) use pio_types, only : io_data_list #ifdef _PNETCDF # include /* _EXTERNAL */ #endif type(file_desc_t) :: File type(io_data_list), pointer :: ptr, prevptr integer :: cnt, ierr integer, pointer :: array_of_requests(:), status(:) if(associated(file%data_list_top)) then if(debug) print *,__FILE__,__LINE__, file%request_cnt #ifdef MEMCHK call GPTLget_memusage(msize, rss, mshare, mtext, mstack) if(rss>lastrss) then lastrss=rss print *,__PIO_FILE__,__LINE__,'mem=',rss end if #endif allocate(array_of_requests(file%request_cnt), status(file%request_cnt)) ptr=>file%data_list_top cnt=1 do while(associated(ptr)) array_of_requests(cnt)=ptr%request cnt=cnt+1 ptr=>ptr%next end do if(DEBUG) print *,__FILE__,__LINE__, array_of_requests, cnt, file%request_cnt #ifdef _PNETCDF ierr = nfmpi_wait_all(file%fh, file%request_cnt, array_of_requests, status) #endif if(DEBUG) print *,__FILE__,__LINE__,status, ierr, total_buffsize ptr=>file%data_list_top do while(associated(ptr)) if(associated(ptr%data_double)) then deallocate(ptr%data_double) else if(associated(ptr%data_real)) then deallocate(ptr%data_real) else if(associated(ptr%data_int)) then deallocate(ptr%data_int) end if prevptr=>ptr ptr => ptr%next deallocate(prevptr) end do nullify(file%data_list_top) total_buffsize=total_buffsize-file%buffsize file%buffsize=0 file%request_cnt=0 deallocate(array_of_requests) deallocate(status) end if #ifdef MEMCHK call GPTLget_memusage(msize, rss, mshare, mtext, mstack) if(rss>lastrss) then lastrss=rss print *,__PIO_FILE__,__LINE__,'mem=',rss end if #endif end subroutine darray_write_complete #ifdef _COMPRESSION subroutine write_vdc2_real(File, Vardesc, iodesc, array, iostat) use pio_support, only : piodie use, intrinsic :: iso_c_binding use C_interface_mod, only : F_C_STRING_DUP ! !DESCRIPTION: ! Writes a VDC2 vapor data collection ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! varable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! IO decomp descriptor real(r4), dimension(:), intent(in), target :: & array ! array to be written integer (i4), intent(out) :: iostat !EOP !BOC !!!!!!!!!!! locals !!!!!!!!!!!!!!!!!!!!! real(r4), dimension(:), pointer :: & IOBUF ! local IO buffer logical(log_kind) :: UseRearranger logical (log_kind) :: IOproc ! true if IO processor integer (i4) :: len, &! length of IO decomp segmap iotype, &! type of IO to perform ndims ! number of variable dimensions integer(4), pointer :: start(:), count(:) double precision :: timer; integer :: vlen type(c_ptr) :: cptr ! pull information from file_desc_t data structure IOproc = File%iosystem%IOproc iotype = File%iotype UseRearranger = File%iosystem%UseRearranger #ifdef MEMCHK call GPTLget_memusage(msize, rss, mshare, mtext, mstack) if(rss>lastrss) then lastrss=rss print *,__PIO_FILE__,__LINE__,'mem=',rss end if #endif ! pull information out of the decomposition data structure len = iodesc%IOmap%length ndims = size(iodesc%start) call alloc_check(start,ndims) call alloc_check(count,ndims) start = iodesc%start count = iodesc%count if(UseRearranger) then if (IOproc) then if(Debug) print *, 'write_vdc2_real: IAM: ',File%iosystem%comp_rank,'Before call to allocate(IOBUF): ',len call alloc_check(IOBUF,len,' TYPE :IOBUF') IOBUF= -1.0_r8 else call alloc_check(IOBUF,0) IOBUF= -1.0_r8 endif !------------------------------------ ! Rearrange data from comp->io decomp !------------------------------------ ! "array" is comp data call rearrange_comp2io(File%iosystem, iodesc, array, iobuf) !-------------------------------------------- ! End data rearrange !-------------------------------------------- else iobuf=>array endif ! if(UseRearranger) if (IOproc) then vlen = len_trim(vardesc%name) cptr = c_loc(iobuf(1)) call WriteVDC2Var(cptr, start, count, File%iosystem%IO_comm, Vardesc%rec, -1, -1, File%iosystem%num_iotasks, F_C_STRING_DUP(varDesc%name(1:vlen))) endif if(UseRearranger) call dealloc_check(IOBUF) #ifdef MEMCHK call GPTLget_memusage(msize, rss, mshare, mtext, mstack) if(rss>lastrss) then lastrss=rss print *,__PIO_FILE__,__LINE__,'mem=',rss end if #endif end subroutine write_vdc2_real subroutine read_vdc2_real(File, Vardesc, iodesc, array, iostat) use pio_support, only : piodie use, intrinsic :: iso_c_binding use C_interface_mod, only : F_C_STRING_DUP ! !DESCRIPTION: ! Writes a VDC2 vapor data collection ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (File_desc_t), intent(inout) :: & File ! file information type (var_desc_t), intent(inout) :: & varDesc ! varable descriptor type (io_desc_t), intent(inout) :: & ioDesc ! IO decomp descriptor real(r4), dimension(:), intent(inout), target :: & array ! array to be written integer (i4), intent(out) :: iostat !EOP !BOC !!!!!!!!!!! locals !!!!!!!!!!!!!!!!!!!!! real(r4), dimension(:), pointer :: & IOBUF ! local IO buffer logical(log_kind) :: UseRearranger logical (log_kind) :: IOproc ! true if IO processor integer (i4) :: len, &! length of IO decomp segmap iotype, &! type of IO to perform ndims ! number of variable dimensions integer(4), pointer :: start(:), count(:) double precision :: timer; type(c_ptr) :: cptr ! pull information from file_desc_t data structure IOproc = File%iosystem%IOproc iotype = File%iotype UseRearranger = File%iosystem%UseRearranger ! pull information out of the decomposition data structure len = iodesc%IOmap%length ndims = size(iodesc%start) call alloc_check(start,ndims) call alloc_check(count,ndims) start = iodesc%start count = iodesc%count if(UseRearranger) then if (IOproc) then if(Debug) print *, 'read_vdc2_real: IAM: ',File%iosystem%comp_rank,'Before call to allocate(IOBUF): ',len call alloc_check(IOBUF,len,' TYPE :IOBUF') IOBUF= -1.0_r8 cptr = c_loc(iobuf(1)) call ReadVDC2Var(cptr, start, count, File%iosystem%IO_comm, vardesc%rec, -1, -1, File%iosystem%num_iotasks, F_C_String_dup(varDesc%name)) else call alloc_check(IOBUF,0) IOBUF= -1.0_r8 endif !------------------------------------ ! Rearrange data from comp->io decomp !------------------------------------ ! "array" is comp data call rearrange_io2comp(File%iosystem, iodesc, iobuf, array) !-------------------------------------------- ! End data rearrange !-------------------------------------------- else iobuf=>array endif ! if(UseRearranger) if(UseRearranger) call dealloc_check(IOBUF) end subroutine read_vdc2_real #endif end module piodarray