1! 2! Copyright (C) 2013, Northwestern University and Argonne National Laboratory 3! See COPYRIGHT notice in top-level directory. 4! 5! $Id: test_intent.f90 2231 2015-12-16 21:08:02Z wkliao $ 6 7! 8! This program tests Fortran 90 INTENT modifier used in PnetCDF. 9! It also tests put_att on Little Endian when using parameter 10! buffer (read-only memory) 11! 12 subroutine check(err, message) 13 use mpi 14 use pnetcdf 15 implicit none 16 integer err 17 character(len=*) message 18 character(len=256) msg 19 20 ! It is a good idea to check returned value for possible error 21 if (err .NE. NF90_NOERR) then 22 write(6,*) trim(message), trim(nf90mpi_strerror(err)) 23 msg = '*** TESTING F90 test_intent.f90 ' 24 call pass_fail(1, msg) 25 call MPI_Abort(MPI_COMM_WORLD, -1, err) 26 end if 27 end subroutine check 28 29 program main 30 use mpi 31 use pnetcdf 32 implicit none 33 34 integer, parameter :: OneByteInt = selected_int_kind(2), & 35 TwoByteInt = selected_int_kind(4), & 36 FourByteInt = selected_int_kind(9), & 37 EightByteInt = selected_int_kind(18) 38 39 character(LEN=256) filename, cmd, msg 40 integer err, ierr, rank, get_args 41 integer cmode, ncid, varid, dimid(1), req(1), status(1) 42 integer(kind=MPI_OFFSET_KIND) start(1) 43 integer(kind=MPI_OFFSET_KIND) count(1) 44 integer(kind=MPI_OFFSET_KIND) bufsize 45 46 character(LEN=3) cbuf 47 integer(KIND=OneByteInt) i1buf(3) 48 integer(KIND=TwoByteInt) sbuf(3) 49 integer ibuf(3) 50 real rbuf(3) 51 double precision dbuf(3) 52 integer(KIND=EightByteInt) i8buf(3) 53 54 PARAMETER( cbuf="123") 55 PARAMETER(i1buf=(/1_OneByteInt,2_OneByteInt,3_OneByteInt/)) 56 PARAMETER( sbuf=(/1_TwoByteInt,2_TwoByteInt,3_TwoByteInt/)) 57 PARAMETER( ibuf=(/1,2,3/)) 58 PARAMETER( rbuf=(/1.0,2.0,3.0/)) 59 PARAMETER( dbuf=(/1.0,2.0,3.0/)) 60 PARAMETER(i8buf=(/1_EightByteInt,2_EightByteInt,3_EightByteInt/)) 61 62 call MPI_Init(ierr) 63 call MPI_Comm_rank(MPI_COMM_WORLD, rank, ierr) 64 65 ! take filename from command-line argument if there is any 66 if (rank .EQ. 0) then 67 filename = 'testfile.nc' 68 err = get_args(cmd, filename) 69 endif 70 call MPI_Bcast(err, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr) 71 if (err .EQ. 0) goto 999 72 73 call MPI_Bcast(filename, 256, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr) 74 75 ! create file, truncate it if exists 76 cmode = IOR(NF90_CLOBBER, NF90_64BIT_DATA) 77 err = nf90mpi_create(MPI_COMM_WORLD, filename, cmode, & 78 MPI_INFO_NULL, ncid) 79 call check(err, 'In nf90mpi_create: ') 80 81 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_text', cbuf) 82 call check(err, 'In nf90mpi_put_att: ') 83 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_int1', i1buf) 84 call check(err, 'In nf90mpi_put_att: ') 85 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_int2', sbuf) 86 call check(err, 'In nf90mpi_put_att: ') 87 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_int', ibuf) 88 call check(err, 'In nf90mpi_put_att: ') 89 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_real', rbuf) 90 call check(err, 'In nf90mpi_put_att: ') 91 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_double', dbuf) 92 call check(err, 'In nf90mpi_put_att: ') 93 err = nf90mpi_put_att(ncid, NF90_GLOBAL, 'nf90_attr_int8', i8buf) 94 call check(err, 'In nf90mpi_put_att: ') 95 96 err = nfmpi_put_att_text(ncid, NF90_GLOBAL, 'nf_attr_text', 3_MPI_OFFSET_KIND, cbuf) 97 call check(err, 'In nfmpi_put_att_text: ') 98 err = nfmpi_put_att_int1(ncid, NF90_GLOBAL, 'nf_attr_int1', NF90_INT1, 3_MPI_OFFSET_KIND, i1buf) 99 call check(err, 'In nfmpi_put_att_int1: ') 100 err = nfmpi_put_att_int2(ncid, NF90_GLOBAL, 'nf_attr_int2', NF90_INT2, 3_MPI_OFFSET_KIND, sbuf) 101 call check(err, 'In nfmpi_put_att_int2: ') 102 err = nfmpi_put_att_int(ncid, NF90_GLOBAL, 'nf_attr_int', NF90_INT, 3_MPI_OFFSET_KIND, ibuf) 103 call check(err, 'In nfmpi_put_att_int: ') 104 err = nfmpi_put_att_real(ncid, NF90_GLOBAL, 'nf_attr_real', NF90_FLOAT, 3_MPI_OFFSET_KIND, rbuf) 105 call check(err, 'In nfmpi_put_att_real: ') 106 err = nfmpi_put_att_double(ncid, NF90_GLOBAL, 'nf_attr_double', NF90_DOUBLE, 3_MPI_OFFSET_KIND, dbuf) 107 call check(err, 'In nfmpi_put_att_double: ') 108 err = nfmpi_put_att_int8(ncid, NF90_GLOBAL, 'nf_attr_int8', NF90_INT64, 3_MPI_OFFSET_KIND, i8buf) 109 call check(err, 'In nfmpi_put_att_int8: ') 110 111 ! define a variable of an integer array of size 3 in the nc file 112 err = nfmpi_def_dim(ncid, 'X', 3_MPI_OFFSET_KIND, dimid(1)) 113 call check(err, 'In nfmpi_def_dim: ') 114 115 err = nfmpi_def_var(ncid, 'var', NF90_INT, 1, dimid, varid) 116 call check(err, 'In nfmpi_def_var: ') 117 118 err = nfmpi_enddef(ncid) 119 call check(err, 'In nfmpi_enddef: ') 120 121 ! bufsize must be max of data type converted before and after 122 bufsize = 3*4 123 err = nfmpi_buffer_attach(ncid, bufsize) 124 call check(err, 'In nfmpi_buffer_attach: ') 125 126 start(1) = 1 127 count(1) = 3 128 err = nfmpi_bput_vara_int(ncid, varid, start, count, ibuf, req(1)) 129 call check(err, 'In nfmpi_bput_vara_int: ') 130 131 err = nfmpi_wait_all(ncid, 1, req, status) 132 call check(err, 'In nfmpi_wait_all: ') 133 134 if (status(1) .ne. NF90_NOERR) then 135 print*,'Error at bput status ', nfmpi_strerror(status(1)) 136 endif 137 138 err = nfmpi_buffer_detach(ncid) 139 call check(err, 'In nfmpi_buffer_detach: ') 140 141 ! close the file 142 err = nf90mpi_close(ncid) 143 call check(err, 'In nf90mpi_close: ') 144 145 msg = '*** TESTING F90 '//trim(cmd)//' for INTENT modifier ' 146 if (rank .eq. 0) call pass_fail(0, msg) 147 148 999 call MPI_Finalize(ierr) 149 end program main 150 151