1! 2! Copyright (C) 2012, Northwestern University and Argonne National Laboratory 3! See COPYRIGHT notice in top-level directory. 4! 5! $Id: put_varn_real.f90 2476 2016-09-06 01:05:33Z wkliao $ 6 7! 8! This example shows how to use a single call of nf90mpi_put_varn_all() 9! to write a sequence of one-element requests with arbitrary array indices. 10! 11! The compile and run commands are given below, together with an ncmpidump of 12! the output file. 13! 14! % mpif90 -O2 -o put_varn_real put_varn_real.f90 -lpnetcdf 15! % mpiexec -n 4 ./put_varn_real /pvfs2/wkliao/testfile.nc 16! % ncmpidump /pvfs2/wkliao/testfile.nc 17! netcdf testfile { 18! // file format: CDF-5 (big variables) 19! dimensions: 20! Y = 4 ; 21! X = 10 ; 22! variables: 23! int var(Y, X) ; 24! data: 25! 26! var = 27! 3, 3, 3, 1, 1, 0, 0, 2, 1, 1, 28! 0, 2, 2, 2, 3, 1, 1, 2, 2, 2, 29! 1, 1, 2, 3, 3, 3, 0, 0, 1, 1, 30! 0, 0, 0, 2, 1, 1, 1, 3, 3, 3 ; 31! } 32! 33 34 subroutine check(err, message) 35 use mpi 36 use pnetcdf 37 implicit none 38 integer err 39 character(len=*) message 40 41 ! It is a good idea to check returned value for possible error 42 if (err .NE. NF90_NOERR) then 43 write(6,*) trim(message), trim(nf90mpi_strerror(err)) 44 call MPI_Abort(MPI_COMM_WORLD, -1, err) 45 end if 46 end subroutine check 47 48 program main 49 use mpi 50 use pnetcdf 51 implicit none 52 53 integer NDIMS 54 PARAMETER(NDIMS=2) 55 56 character(LEN=256) filename, cmd 57 integer rank, nprocs, err, num_reqs, ierr, get_args, dummy 58 integer ncid, cmode, varid, dimid(2), y, x 59 real buffer(13) 60 integer(kind=MPI_OFFSET_KIND) NY, NX 61 integer(kind=MPI_OFFSET_KIND) starts(NDIMS, 13) 62 integer(kind=MPI_OFFSET_KIND) counts(NDIMS, 13) 63 integer(kind=MPI_OFFSET_KIND) malloc_size, sum_size 64 logical verbose 65 66 NY = 4 67 NX = 10 68 69 call MPI_Init(err) 70 call MPI_Comm_rank(MPI_COMM_WORLD, rank, err) 71 call MPI_Comm_size(MPI_COMM_WORLD, nprocs, err) 72 73 ! take filename from command-line argument if there is any 74 if (rank .EQ. 0) then 75 verbose = .TRUE. 76 filename = "testfile.nc" 77 ierr = get_args(2, cmd, filename, verbose, dummy) 78 endif 79 call MPI_Bcast(ierr, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, err) 80 if (ierr .EQ. 0) goto 999 81 82 call MPI_Bcast(verbose, 1, MPI_LOGICAL, 0, MPI_COMM_WORLD, err) 83 call MPI_Bcast(filename, 256, MPI_CHARACTER, 0, MPI_COMM_WORLD, err) 84 85 if (nprocs .NE. 4 .AND. rank .EQ. 0 .AND. verbose) & 86 print*,'Warning: ',trim(cmd),' is intended to run on ', & 87 '4 processes' 88 89 ! create file, truncate it if exists 90 cmode = IOR(NF90_CLOBBER, NF90_64BIT_DATA) 91 err = nf90mpi_create(MPI_COMM_WORLD, filename, cmode, & 92 MPI_INFO_NULL, ncid) 93 call check(err, 'In nf90mpi_create: ') 94 95 ! create a global array of size NY * NX */ 96 err = nf90mpi_def_dim(ncid, "Y", NY, dimid(2)) 97 call check(err, 'In nf90mpi_def_dim Y: ') 98 err = nf90mpi_def_dim(ncid, "X", NX, dimid(1)) 99 call check(err, 'In nf90mpi_def_dim X: ') 100 err = nf90mpi_def_var(ncid, "var", NF90_FLOAT, dimid, varid) 101 call check(err, 'In nf90mpi_def_var var: ') 102 err = nf90mpi_enddef(ncid) 103 call check(err, 'In nf90mpi_enddef: ') 104 105 ! pick arbitrary numbers of requests for 4 processes 106 num_reqs = 0 107 if (rank .EQ. 0) then 108 num_reqs = 8 109 elseif (rank .EQ. 1) then 110 num_reqs = 13 111 elseif (rank .EQ. 2) then 112 num_reqs = 9 113 elseif (rank .EQ. 3) then 114 num_reqs = 10 115 endif 116 117 ! assign arbitrary starts 118 y=2 119 x=1 120 if (rank .EQ. 0) then 121 starts(y, 1) = 1 122 starts(x, 1) = 6 123 starts(y, 2) = 2 124 starts(x, 2) = 1 125 starts(y, 3) = 3 126 starts(x, 3) = 7 127 starts(y, 4) = 4 128 starts(x, 4) = 1 129 starts(y, 5) = 1 130 starts(x, 5) = 7 131 starts(y, 6) = 3 132 starts(x, 6) = 8 133 starts(y, 7) = 4 134 starts(x, 7) = 2 135 starts(y, 8) = 4 136 starts(x, 8) = 3 137 ! rank 0 is writing the following locations: ("-" means skip) 138 ! - - - - - 0 0 - - - 139 ! 0 - - - - - - - - - 140 ! - - - - - - 0 0 - - 141 ! 0 0 0 - - - - - - - 142 elseif (rank .EQ. 1) then 143 starts(y, 1) = 1 144 starts(x, 1) = 4 145 starts(y, 2) = 1 146 starts(x, 2) = 9 147 starts(y, 3) = 2 148 starts(x, 3) = 6 149 starts(y, 4) = 3 150 starts(x, 4) = 1 151 starts(y, 5) = 3 152 starts(x, 5) = 9 153 starts(y, 6) = 4 154 starts(x, 6) = 5 155 starts(y, 7) = 1 156 starts(x, 7) = 5 157 starts(y, 8) = 1 158 starts(x, 8) = 10 159 starts(y, 9) = 2 160 starts(x, 9) = 7 161 starts(y, 10) = 3 162 starts(x, 10) = 2 163 starts(y, 11) = 3 164 starts(x, 11) = 10 165 starts(y, 12) = 4 166 starts(x, 12) = 6 167 starts(y, 13) = 4 168 starts(x, 13) = 7 169 ! rank 1 is writing the following locations: ("-" means skip) 170 ! - - - 1 1 - - - 1 1 171 ! - - - - - 1 1 - - - 172 ! 1 1 - - - - - - 1 1 173 ! - - - - 1 1 1 - - - 174 elseif (rank .EQ. 2) then 175 starts(y, 1) = 1 176 starts(x, 1) = 8 177 starts(y, 2) = 2 178 starts(x, 2) = 2 179 starts(y, 3) = 2 180 starts(x, 3) = 8 181 starts(y, 4) = 3 182 starts(x, 4) = 3 183 starts(y, 5) = 4 184 starts(x, 5) = 4 185 starts(y, 6) = 2 186 starts(x, 6) = 3 187 starts(y, 7) = 2 188 starts(x, 7) = 9 189 starts(y, 8) = 2 190 starts(x, 8) = 4 191 starts(y, 9) = 2 192 starts(x, 9) = 10 193 ! rank 2 is writing the following locations: ("-" means skip) 194 ! - - - - - - - 2 - - 195 ! - 2 2 2 - - - 2 2 2 196 ! - - 2 - - - - - - - 197 ! - - - 2 - - - - - - 198 elseif (rank .EQ. 3) then 199 starts(y, 1) = 1 200 starts(x, 1) = 1 201 starts(y, 2) = 2 202 starts(x, 2) = 5 203 starts(y, 3) = 3 204 starts(x, 3) = 4 205 starts(y, 4) = 4 206 starts(x, 4) = 8 207 starts(y, 5) = 1 208 starts(x, 5) = 2 209 starts(y, 6) = 3 210 starts(x, 6) = 5 211 starts(y, 7) = 4 212 starts(x, 7) = 9 213 starts(y, 8) = 1 214 starts(x, 8) = 3 215 starts(y, 9) = 3 216 starts(x, 9) = 6 217 starts(y, 10) = 4 218 starts(x, 10) = 10 219 ! rank 3 is writing the following locations: ("-" means skip) 220 ! 3 3 3 - - - - - - - 221 ! - - - - 3 - - - - - 222 ! - - - 3 3 3 - - - - 223 ! - - - - - - - 3 3 3 224 endif 225 counts = 1 226 227 ! allocate I/O buffer and initialize its contents 228 buffer = rank 229 230 ! set the buffer pointers to different offsets to the I/O buffer 231 err = nf90mpi_put_varn_all(ncid, varid, buffer, num_reqs, & 232 starts, counts) 233 call check(err, 'In nf90mpi_put_varn_all: ') 234 235 err = nf90mpi_close(ncid); 236 call check(err, 'In nf90mpi_close: ') 237 238 ! check if there is any PnetCDF internal malloc residue 239 998 format(A,I13,A) 240 err = nf90mpi_inq_malloc_size(malloc_size) 241 if (err == NF90_NOERR) then 242 call MPI_Reduce(malloc_size, sum_size, 1, MPI_INTEGER8, & 243 MPI_SUM, 0, MPI_COMM_WORLD, err) 244 if (rank .EQ. 0 .AND. sum_size .GT. 0_MPI_OFFSET_KIND) print 998, & 245 'heap memory allocated by PnetCDF internally has ', & 246 sum_size/1048576, ' MiB yet to be freed' 247 endif 248 249 999 call MPI_Finalize(err) 250 ! call EXIT(0) ! EXIT() is a GNU extension 251 252 end program 253 254