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_vars4.f90 2131 2015-09-25 22:33:12Z wkliao $ 8 9! This program tests PnetCDF variable functions from fortran 90. 10 11program f90tst_vars4 12 use mpi 13 use pnetcdf 14 implicit none 15 16 ! This is the name of the data file we will create. 17 character (len = *), parameter :: FILE_NAME = "f90tst_vars4.nc" 18 19 integer, parameter :: MAX_DIMS = 2 20 integer, parameter :: NX = 40, NY = 4096 21 integer :: data_out(NY, NX), data_in(NY, NX) 22 23 ! We need these ids and other gunk for netcdf. 24 integer :: ncid, varid, dimids(MAX_DIMS) 25 integer :: x_dimid, y_dimid 26 integer :: mode_flag 27 integer :: nvars, ngatts, ndims, unlimdimid, file_format 28 integer :: x, y 29 integer, parameter :: CACHE_SIZE = 1000000 30 integer :: xtype_in, natts_in, dimids_in(MAX_DIMS) 31 character (len = NF90_MAX_NAME) :: name_in 32 integer :: err, ierr, get_args 33 integer(KIND=MPI_OFFSET_KIND) :: nx_ll, ny_ll 34 character(LEN=256) filename, cmd, msg 35 integer my_rank, p 36 37 call MPI_Init(ierr) 38 call MPI_Comm_rank(MPI_COMM_WORLD, my_rank, ierr) 39 call MPI_Comm_size(MPI_COMM_WORLD, p, ierr) 40 41 ! take filename from command-line argument if there is any 42 if (my_rank .EQ. 0) then 43 filename = FILE_NAME 44 err = get_args(cmd, filename) 45 endif 46 call MPI_Bcast(err, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr) 47 if (err .EQ. 0) goto 999 48 49 call MPI_Bcast(filename, 256, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr) 50 51! if (p .ne. 1 .AND. my_rank .eq. 0) then 52! print *, 'Warning: ',trim(cmd),' is design to run on 1 process' 53! endif 54 55 ! Create some pretend data. 56 do x = 1, NX 57 do y = 1, NY 58 data_out(y, x) = (x - 1) * NY + (y - 1) 59 end do 60 end do 61 62 ! Create the netCDF file. 63 mode_flag = IOR(NF90_CLOBBER, NF90_64BIT_DATA) 64 call handle_err(nf90mpi_create(MPI_COMM_WORLD, filename, mode_flag, MPI_INFO_NULL, ncid)) 65 66 ! Define the dimensions. 67 nx_ll = NX 68 ny_ll = NY 69 call handle_err(nf90mpi_def_dim(ncid, "x", nx_ll, x_dimid)) 70 call handle_err(nf90mpi_def_dim(ncid, "y", ny_ll, y_dimid)) 71 dimids = (/ y_dimid, x_dimid /) 72 73 ! Define the variable. 74 call handle_err(nf90mpi_def_var(ncid, 'data', NF90_INT, dimids, varid)) 75 76 ! enddef must be called. 77 call handle_err(nf90mpi_enddef(ncid)) 78 79 call handle_err(nf90mpi_begin_indep_data(ncid)) 80 81 ! Write the pretend data to the file. 82 call handle_err(nf90mpi_put_var(ncid, varid, data_out)) 83 84 ! Close the file. 85 call handle_err(nf90mpi_close(ncid)) 86 87 ! Reopen the file. 88 call handle_err(nf90mpi_open(MPI_COMM_WORLD, filename, nf90_nowrite, MPI_INFO_NULL, ncid)) 89 90 ! Check some stuff out. 91 call handle_err(nf90mpi_inquire(ncid, ndims, nvars, ngatts, unlimdimid, file_format)) 92 if (ndims /= 2 .or. nvars /= 1 .or. ngatts /= 0 .or. unlimdimid /= -1 .or. & 93 file_format /= nf90_format_cdf5) stop 3 94 95 call handle_err(nf90mpi_inquire_variable(ncid, varid, name_in, xtype_in, ndims, dimids_in, & 96 natts_in)) 97 if (name_in .ne. 'data' .or. xtype_in .ne. NF90_INT .or. ndims .ne. MAX_DIMS .or. & 98 dimids_in(1) /= y_dimid .or. dimids_in(2) /= x_dimid .or. & 99 natts_in .ne. 0) & 100 stop 4 101 102 ! Check the data. 103 call handle_err(nf90mpi_get_var_all(ncid, varid, data_in)) 104 do x = 1, NX 105 do y = 1, NY 106 if (data_out(y, x) .ne. data_in(y, x)) stop 3 107 end do 108 end do 109 110 ! Close the file. 111 call handle_err(nf90mpi_close(ncid)) 112 113 msg = '*** TESTING F90 '//trim(cmd)//' for def_var API' 114 if (my_rank .eq. 0) call pass_fail(0, msg) 115 116 999 call MPI_Finalize(ierr) 117 118contains 119! This subroutine handles errors by printing an error message and 120! exiting with a non-zero status. 121 subroutine handle_err(errcode) 122 implicit none 123 integer, intent(in) :: errcode 124 125 if(errcode /= nf90_noerr) then 126 print *, 'Error: ', trim(nf90mpi_strerror(errcode)) 127 stop 2 128 endif 129 end subroutine handle_err 130end program f90tst_vars4 131 132