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: f90tst_parallel3.f90 2512 2016-09-29 01:29:37Z wkliao $ 8 9! This program tests PnetCDF parallel I/O from 10! fortran. It creates a file like this: 11 12! netcdf f90tst_parallel3 { 13! dimensions: 14! x = 16 ; 15! y = 16 ; 16! variables: 17! byte byte__(x, y) ; 18! short short_(x, y) ; 19! int int__(x, y) ; 20! float float_(x, y) ; 21! double double(x, y) ; 22! int64 int64_(x, y) ; 23 24 25program f90tst_parallel3 26 use mpi 27 use pnetcdf 28 implicit none 29 30 integer, parameter :: OneByteInt = selected_int_kind(2), & 31 TwoByteInt = selected_int_kind(4), & 32 FourByteInt = selected_int_kind(9), & 33 EightByteInt = selected_int_kind(18) 34 35 ! This is the name of the data file we will create. 36 character (len = *), parameter :: FILE_NAME = "f90tst_parallel3.nc" 37 integer, parameter :: MAX_DIMS = 2 38 integer, parameter :: NX = 16, NY = 16 39 integer, parameter :: HALF_NX = NX/2, HALF_NY = NY/2 40 integer, parameter :: NUM_PROC = 4 41 integer, parameter :: NUM_VARS = 6 42 integer, parameter :: CACHE_SIZE = 4194304, CACHE_NELEMS = 1013 43 integer, parameter :: CACHE_PREEMPTION = 79 44 character (len = *), parameter :: var_name(NUM_VARS) = & 45 (/ 'byte__', 'short_', 'int___', 'float_', 'double', 'int64_' /) 46 integer :: ncid, varid(NUM_VARS), dimids(MAX_DIMS) 47 integer :: var_type(NUM_VARS) = (/ nf90_byte, nf90_short, nf90_int, & 48 nf90_float, nf90_double, nf90_int64 /) 49 integer :: x_dimid, y_dimid 50 integer(kind=OneByteInt) :: byte_out(HALF_NY, HALF_NX), byte_in(HALF_NY, HALF_NX) 51 integer(kind=TwoByteInt) :: short_out(HALF_NY, HALF_NX), short_in(HALF_NY, HALF_NX) 52 integer :: int_out(HALF_NY, HALF_NX), int_in(HALF_NY, HALF_NX) 53 real :: areal_out(HALF_NY, HALF_NX), areal_in(HALF_NY, HALF_NX) 54 double precision :: double_out(HALF_NY, HALF_NX), double_in(HALF_NY, HALF_NX) 55 integer (kind=EightByteInt) :: int64_out(HALF_NY, HALF_NX), int64_in(HALF_NY, HALF_NX) 56 integer :: nvars, ngatts, ndims, unlimdimid, file_format 57 integer :: x, y, v 58 integer :: p, my_rank, err, ierr, get_args 59 integer(KIND=MPI_OFFSET_KIND) :: start(MAX_DIMS), count(MAX_DIMS) 60 integer :: cmode 61 integer(KIND=MPI_OFFSET_KIND) :: nx_ll, ny_ll 62 character(LEN=256) filename, cmd, msg 63 64 call MPI_Init(ierr) 65 call MPI_Comm_rank(MPI_COMM_WORLD, my_rank, ierr) 66 call MPI_Comm_size(MPI_COMM_WORLD, p, ierr) 67 68 ! take filename from command-line argument if there is any 69 if (my_rank .EQ. 0) then 70 filename = FILE_NAME 71 err = get_args(cmd, filename) 72 endif 73 call MPI_Bcast(err, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr) 74 if (err .EQ. 0) goto 999 75 76 call MPI_Bcast(filename, 256, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr) 77 78! if (p .ne. 4 .AND. my_rank .eq. 0) then 79! print *, 'Warning: ',trim(cmd),' is design to run on 4 processes.' 80! endif 81 82 ! Create some pretend data. 83 do x = 1, HALF_NX 84 do y = 1, HALF_NY 85 byte_out(y, x) = INT(my_rank,1) * (-1_1) 86 short_out(y, x) = INT(my_rank, TwoByteInt) * (-2_2) 87 int_out(y, x) = my_rank * (-4) 88 areal_out(y, x) = my_rank * 2.5 89 double_out(y, x) = my_rank * (-4.5) 90 int64_out(y, x) = my_rank * 4 91 end do 92 end do 93 94 ! Create the netCDF file. 95 cmode = IOR(NF90_CLOBBER, NF90_64BIT_DATA) 96 call check(nf90mpi_create(MPI_COMM_WORLD, filename, cmode, MPI_INFO_NULL, ncid)) 97 98 ! Define the dimensions. 99 nx_ll = NX 100 ny_ll = NY 101 call check(nf90mpi_def_dim(ncid, "x", nx_ll, x_dimid)) 102 call check(nf90mpi_def_dim(ncid, "y", ny_ll, y_dimid)) 103 dimids = (/ y_dimid, x_dimid /) 104 105 ! Define the variables. 106 do v = 1, NUM_VARS 107 call check(nf90mpi_def_var(ncid, var_name(v), var_type(v), dimids, varid(v))) 108 end do 109 110 ! This will be the last collective operation. 111 call check(nf90mpi_enddef(ncid)) 112 113 ! Determine what part of the variable will be written/read for this 114 ! processor. It's a checkerboard decomposition. 115 count = (/ HALF_NX, HALF_NY /) 116 if (my_rank .eq. 0) then 117 start = (/ 1, 1 /) 118 else if (my_rank .eq. 1) then 119 start = (/ HALF_NX + 1, 1 /) 120 else if (my_rank .eq. 2) then 121 start = (/ 1, HALF_NY + 1 /) 122 else if (my_rank .eq. 3) then 123 start = (/ HALF_NX + 1, HALF_NY + 1 /) 124 else 125 start = (/ 1, 1 /) 126 count = 0 127 endif 128 129 ! Write this processor's data, except for processor zero. 130 if (my_rank .EQ. 0) count = (/ 0, 0 /) 131 call check(nf90mpi_put_var_all(ncid, varid(1), byte_out, start = start, count = count)) 132 call check(nf90mpi_put_var_all(ncid, varid(2), short_out, start = start, count = count)) 133 call check(nf90mpi_put_var_all(ncid, varid(3), int_out, start = start, count = count)) 134 call check(nf90mpi_put_var_all(ncid, varid(4), areal_out, start = start, count = count)) 135 call check(nf90mpi_put_var_all(ncid, varid(5), double_out, start = start, count = count)) 136 call check(nf90mpi_put_var_all(ncid, varid(6), int64_out, start = start, count = count)) 137 138 ! Close the file. 139 call check(nf90mpi_close(ncid)) 140 141 ! Reopen the file. 142 call check(nf90mpi_open(MPI_COMM_WORLD, filename, nf90_nowrite, MPI_INFO_NULL, ncid)) 143 144 ! Check some stuff out. 145 call check(nf90mpi_inquire(ncid, ndims, nvars, ngatts, unlimdimid, file_format)) 146 if (ndims /= 2 .or. nvars /= NUM_VARS .or. ngatts /= 0 .or. unlimdimid /= -1 .or. & 147 file_format /= nf90_format_cdf5) stop 2 148 149 ! Read this processor's data. 150 if (my_rank .EQ. 0) count = (/ HALF_NX, HALF_NY /) 151 call check(nf90mpi_get_var_all(ncid, varid(1), byte_in, start = start, count = count)) 152 call check(nf90mpi_get_var_all(ncid, varid(2), short_in, start = start, count = count)) 153 call check(nf90mpi_get_var_all(ncid, varid(3), int_in, start = start, count = count)) 154 call check(nf90mpi_get_var_all(ncid, varid(4), areal_in, start = start, count = count)) 155 call check(nf90mpi_get_var_all(ncid, varid(5), double_in, start = start, count = count)) 156 call check(nf90mpi_get_var_all(ncid, varid(6), int64_in, start = start, count = count)) 157 158 ! Check the data. All the data from the processor zero are fill 159 ! value. 160 if (my_rank .LT. 4) then 161 do x = 1, HALF_NX 162 do y = 1, HALF_NY 163 if (my_rank .NE. 0) then 164 if (byte_in(y, x) .ne. (my_rank * (-1))) stop 13 165 if (short_in(y, x) .ne. (my_rank * (-2))) stop 14 166 if (int_in(y, x) .ne. (my_rank * (-4))) stop 15 167 if (areal_in(y, x) .ne. (my_rank * (2.5))) stop 16 168 if (double_in(y, x) .ne. (my_rank * (-4.5))) stop 17 169 if (int64_in(y, x) .ne. (my_rank * (4))) stop 20 170 endif 171 end do 172 end do 173 endif 174 175 ! Close the file. 176 call check(nf90mpi_close(ncid)) 177 178 if (my_rank .eq. 0) then 179 msg = '*** TESTING F90 '//trim(cmd) 180 call pass_fail(0, msg) 181 endif 182 183 999 call MPI_Finalize(ierr) 184 185contains 186! This subroutine handles errors by printing an error message and 187! exiting with a non-zero status. 188 subroutine check(errcode) 189 implicit none 190 integer, intent(in) :: errcode 191 192 if(errcode /= nf90_noerr) then 193 print *, 'Error: ', trim(nf90mpi_strerror(errcode)) 194 stop 99 195 endif 196 end subroutine check 197end program f90tst_parallel3 198 199