1! 2! Copyright (C) 2013, Northwestern University and Argonne National Laboratory 3! See COPYRIGHT notice in top-level directory. 4! 5! This is part of the PnetCDF package. 6! 7! $Id: tst_io.f90 2147 2015-10-09 23:50:26Z wkliao $ 8 9! This program tests large files (> 4 GB) in PnetCDF. 10 11program tst_io 12 use mpi 13 use pnetcdf 14 implicit none 15 integer, parameter :: EightByteInt = selected_int_kind(18) 16 integer(kind=EightByteInt), parameter :: prsz1 = 50, prsz2 = 50, & 17 prsz3 = 50, prsz4 = 50, repct = 10 18 integer :: i1, i2, i3, i4 19 real :: psr 20 integer :: err, ierr, get_args 21 integer :: start, now, ncint1, ncint2, size 22 ! integer :: wrint1, wrint2, wrint3, ncint3 23 real, dimension (prsz1, prsz2, prsz3, prsz4) :: s, t, u, v, w, x, y, z 24 character(len = *), parameter :: nclFilenm1 = 'tst_io1.nc', & 25 nclFilenm2 = 'tst_io2.nc', nclFilenm3 = 'tst_io3.nc', & 26 nclFilenm4 = 'tst_io4.nc', nclFilenm5 = 'tst_io5.nc', & 27 nclFilenm6 = 'tst_io6.nc', nclFilenm7 = 'tst_io7.nc', & 28 nclFilenm8 = 'tst_io8.nc', nclFilenm9 = 'tst_io9.nc', & 29 nclFilenm10 = 'tst_io10.nc', nclFilenm11 = 'tst_io11.nc' 30 ! needed for netcdf 31 integer :: ncid, x1id, x2id, x3id, x4id, vrid 32 ! integer :: vrids, vridt, vridu, vridv, vridw, vridx, vridy, vridz 33 character(LEN=256) dirpath, cmd, msg 34 integer my_rank, p 35 36 call MPI_Init(ierr) 37 call MPI_Comm_rank(MPI_COMM_WORLD, my_rank, ierr) 38 call MPI_Comm_size(MPI_COMM_WORLD, p, ierr) 39 40 ! take filename from command-line argument if there is any 41 if (my_rank .EQ. 0) then 42 dirpath = '.' 43 err = get_args(cmd, dirpath) 44 endif 45 call MPI_Bcast(err, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr) 46 if (err .EQ. 0) goto 999 47 48 call MPI_Bcast(dirpath, 256, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr) 49 50! if (p .ne. 1 .AND. my_rank .eq. 0) then 51! print *, 'Warning: ',trim(cmd),' is design to run on 1 process' 52! endif 53 54 psr = 1.7/real(prsz1) 55 ! print *, "Starting data initialization." 56 size = (prsz1 * prsz2 * prsz3 * prsz4 )/ 250000 57 do i1 = 1, prsz1 58 do i2 = 1, prsz2 59 do i3 = 1, prsz3 ! Jackson Pollock it is not 60 do i4 = 1, prsz4 61 x(i1, i2, i3, i4) = sin(i1*psr)*(0.5 + cos(i2*psr))+(psr/i3)+ i4/(10.0*prsz4) 62 s(i1, i2, i3, i4) = x(i1, i2, i3, i4) - 5.0 63 t(i1, i2, i3, i4) = x(i1, i2, i3, i4) - 4.0 64 u(i1, i2, i3, i4) = x(i1, i2, i3, i4) - 3.0 65 v(i1, i2, i3, i4) = x(i1, i2, i3, i4) - 2.0 66 w(i1, i2, i3, i4) = x(i1, i2, i3, i4) - 1.0 67 y(i1, i2, i3, i4) = x(i1, i2, i3, i4) + 1.0 68 z(i1, i2, i3, i4) = x(i1, i2, i3, i4) + 2.0 69 enddo 70 enddo 71 enddo 72 enddo 73 74 ! call setupNetCDF (trim(dirpath)//'/'//nclFilenm1, ncid, vrid, x, prsz1, prsz2, prsz3, prsz4, & 75 call setupNetCDF (trim(dirpath)//'/'//nclFilenm1, ncid, vrid, prsz1, prsz2, prsz3, prsz4, & 76 x1id, x2id, x3id, x4id, NF90_CLOBBER, 20) 77 call system_clock(start) 78 call check(nfmpi_begin_indep_data(ncid), 11) 79 call check (NF90MPI_PUT_VAR(ncid, vrid, x), 18) 80 call check(nfmpi_end_indep_data(ncid), 11) 81 call system_clock(now) 82 ncint1 = now - start 83! print 3, size, "MB"," netcdf write = ", ncint1 * clockRate, & 84! real(ncint1)/real (wrint1) 85! 3 format("Time for", i5, a, a25, i7, " msec. Spd ratio = ", f5.2) 86 87 call check (NF90MPI_CLOSE(ncid), 14) 88 89 call system_clock(start) 90 do i1 = 1, repct 91 ! call setupNetCDF (trim(dirpath)//'/'//nclFilenm1, ncid, vrid, x, prsz1, prsz2, prsz3, prsz4, & 92 call setupNetCDF (trim(dirpath)//'/'//nclFilenm1, ncid, vrid, prsz1, prsz2, prsz3, prsz4, & 93 x1id, x2id, x3id, x4id, NF90_CLOBBER, 130) 94 call check(nfmpi_begin_indep_data(ncid), 11) 95 call check (NF90MPI_PUT_VAR(ncid, vrid, x), 23 + i1) 96 call check(nfmpi_end_indep_data(ncid), 11) 97 call check (NF90MPI_CLOSE(ncid), 15) 98 enddo 99 call system_clock(now) 100 ncint2 = now - start 101! print 4, size, repct, " repeated netcdf writes = ", ncint2 * clockRate, & 102! real(ncint2)/real(wrint2); 103! 4 format("Time for", i5, "MB", i3, a22, i7, " msec. Spd ratio = ", f5.2) 104 105! call system_clock(start) 106! call setupNetCDF (trim(dirpath)//'/'//nclFilenm3, ncid, vrids, s, prsz1, prsz2, prsz3, prsz4, & 107! x1id, x2id, x3id, x4id, NF90_CLOBBER, 20) 108! call setupNetCDF (trim(dirpath)//'/'//nclFilenm4, ncid, vridt, t, prsz1, prsz2, prsz3, prsz4, & 109! x1id, x2id, x3id, x4id, NF90_CLOBBER, 30) 110! call setupNetCDF (trim(dirpath)//'/'//nclFilenm5, ncid, vridu, u, prsz1, prsz2, prsz3, prsz4, & 111! x1id, x2id, x3id, x4id, NF90_CLOBBER, 40) 112! call setupNetCDF (trim(dirpath)//'/'//nclFilenm6, ncid, vridv, v, prsz1, prsz2, prsz3, prsz4, & 113! x1id, x2id, x3id, x4id, NF90_CLOBBER, 50) 114! call setupNetCDF (trim(dirpath)//'/'//nclFilenm7, ncid, vridw, w, prsz1, prsz2, prsz3, prsz4, & 115! x1id, x2id, x3id, x4id, NF90_CLOBBER, 60) 116! call setupNetCDF (trim(dirpath)//'/'//nclFilenm8, ncid, vridx, x, prsz1, prsz2, prsz3, prsz4, & 117! x1id, x2id, x3id, x4id, NF90_CLOBBER, 70) 118! call setupNetCDF (trim(dirpath)//'/'//nclFilenm9, ncid, vridy, y, prsz1, prsz2, prsz3, prsz4, & 119! x1id, x2id, x3id, x4id, NF90_CLOBBER, 80) 120! call setupNetCDF (trim(dirpath)//'/'//nclFilenm10, ncid, vridz, z, prsz1, prsz2, prsz3, prsz4, & 121! x1id, x2id, x3id, x4id, NF90_CLOBBER, 90) 122! call check(nfmpi_begin_indep_data(ncid), 11) 123! call check (NF90MPI_PUT_VAR(ncid, vrids, s), 118) 124! call check (NF90MPI_PUT_VAR(ncid, vridt, t), 119) 125! call check (NF90MPI_PUT_VAR(ncid, vridu, u), 120) 126! call check (NF90MPI_PUT_VAR(ncid, vridv, v), 121) 127! call check (NF90MPI_PUT_VAR(ncid, vridw, w), 122) 128! call check (NF90MPI_PUT_VAR(ncid, vridx, x), 123) 129! call check (NF90MPI_PUT_VAR(ncid, vridy, y), 124) 130! call check (NF90MPI_PUT_VAR(ncid, vridz, z), 125) 131! call check(nfmpi_end_indep_data(ncid), 11) 132! call system_clock(now) 133! ncint3 = now - start 134! call check (NF90MPI_CLOSE(ncid), 16) 135! print 4, size, 8, " netcdf file writes = ", ncint3 * clockRate, & 136! real(ncint3)/real(wrint3); 137 138 msg = '*** TESTING F90 '//trim(cmd) 139 if (my_rank .eq. 0) call pass_fail(0, msg) 140 141 999 call MPI_Finalize(ierr) 142 143contains 144 subroutine check (st, n) ! checks the return error code 145 integer, intent (in) :: st, n 146 if ((n < 10.and.st /= 0).or.(n > 10.and.st /= NF90_noerr))then 147 print *, "I/O error at", n, " status = ", st 148 stop 2 149 endif 150 end subroutine check 151 152! subroutine setupNetCDF(fn, nc, vr, vrnam, d1, d2, d3, d4, do1, do2, & 153! real, dimension (:, :, :, :), intent (in) :: vrnam 154 subroutine setupNetCDF(fn, nc, vr, d1, d2, d3, d4, do1, do2, & 155 do3, do4, stat, deb) 156 157 character(len = *), intent(in) :: fn 158 integer(kind=EightByteInt), intent(in) :: d1, d2, d3, d4 159 integer, intent(in) :: stat, deb 160 integer, intent(out) :: do1, do2, do3, do4, vr 161 integer, intent(inout) :: nc 162 integer, dimension(4) :: dimids (4) 163 164 call check (NF90MPI_CREATE (MPI_COMM_WORLD, fn, stat, MPI_INFO_NULL, nc), deb + 1) 165 call check (NF90MPI_DEF_DIM(nc, "d1", d1, do1), deb + 2) 166 call check (NF90MPI_DEF_DIM(nc, "d2", d2, do2), deb + 3) 167 call check (NF90MPI_DEF_DIM(nc, "d3", d3, do3), deb + 4) 168 call check (NF90MPI_DEF_DIM(nc, "d4", d4, do4), deb + 5) 169 170 dimids = (/ do1, do2, do3, do4 /) 171 call check (NF90MPI_DEF_VAR(nc, "data", NF90_REAL, dimids, vr), deb + 6) 172 call check (NF90MPI_ENDDEF (nc), deb + 7) 173 174 end subroutine setupNetCDF 175 176end program tst_io 177