CGpop mini-Application (2-sided MPI 1D data structure version) 0.1
io_serial.F90
Go to the documentation of this file.
00001 !==============================================================================
00002 ! Copyright (C) 2010, University Corporation for Atmospheric Research,
00003 !                     Colorado State University,
00004 !                     Los Alamos National Security, LLC,
00005 !                     United States Department of Energy
00006 !
00007 ! All rights reserved.  See ../COPYING for copyright details
00008 !==============================================================================
00009 
00010 !>
00011 !!  This module provides several methods to read and write 
00012 !!  data needed by the miniapp
00013 !<
00014 module io_serial
00015     use kinds_mod, only: i4, r8
00016     use simple_blocks, only: AppMD_t
00017     use netcdf
00018     use exit_mod
00019 
00020     implicit none
00021     private
00022 
00023     public :: write_AppMD, read_AppMD  ! read and write application metadata
00024 
00025     public :: read_tiles
00026 
00027     !>
00028     !! Reads in variables decomposed into tiles
00029     !<
00030     interface read_tiles
00031         module procedure read_tiles_int
00032         module procedure read_tiles_dbl
00033     end interface
00034 
00035     public :: write_tiles
00036 
00037     !>
00038     !! Writes out variables decomposed into tiles
00039     !<
00040     interface write_tiles
00041         module procedure write_tiles_int
00042         module procedure write_tiles_dbl
00043     end interface
00044 
00045     contains 
00046 
00047     !>
00048     !! @details Reads in a double precision variable from netCDF tilefile. Note 
00049     !! that this subroutine reads in all tiles. 
00050     !! @param fname Filename of the tilefile from which to read the variable.
00051     !! @param vname Name of the variable to read.
00052     !! @param var   Double precision variable that is returned.
00053     !<
00054     subroutine read_tiles_dbl(fname,vname,var)
00055         character(len=80), intent(in) :: fname
00056         character(len=*), intent(in)  :: vname
00057         real(r8), dimension(:,:,:) :: var
00058 
00059         integer(i4) :: mode,fh,varid,ierr
00060         integer(i4) :: nx,ny
00061         integer(i4) :: nx_read,ny_read
00062         integer(i4) :: dimid_nx, dimid_ny
00063         character(len=80) :: error
00064 
00065 
00066         mode = NF90_nowrite
00067         nx = SIZE(var,dim=1)
00068         ny = SIZE(var,dim=2)
00069         ierr = sio_open(fname,mode,fh)
00070 
00071         ierr = nf90_inq_dimid(fh,'nx_block',dimid_nx)
00072         ierr = nf90_inquire_dimension(fh,dimid_nx,len=nx_read)
00073         ierr = nf90_inq_dimid(fh,'ny_block',dimid_ny)
00074         ierr = nf90_inquire_dimension(fh,dimid_ny,len=ny_read)
00075         if((nx_read .ne. nx) .or. (ny_read .ne. ny)) then 
00076         print *,'Error: mismatch between size of blocks'
00077         print *,'nx:= ',nx,'nx_read:= ',nx_read
00078         print *,'ny:= ',ny,'ny_read:= ',ny_read
00079         stop 'in read_tiles:dbl'
00080         endif
00081 
00082 
00083         ierr = nf90_inq_varid(fh,TRIM(vname),varid)
00084         if(ierr /= nf90_noerr) &
00085             call check_netcdf(fh,ierr,'nf90_inq_varid: '//vname,__LINE__)
00086         ierr = nf90_get_var(fh,varid, var)
00087         if(ierr /= nf90_noerr) &
00088         call check_netcdf(fh,ierr,'nf90_get_var: '//vname,__LINE__)
00089 
00090         ierr = sio_close(fh) 
00091     end subroutine read_tiles_dbl
00092 
00093     !>
00094     !!
00095     !! Reads in a integer variable from netCDF tile file. Note 
00096     !! that this subroutine reads in all tiles. 
00097     !! 
00098     !! @param fname Filename of the tilefile from which to read the variable.
00099     !! @param vname Name of the variable to read 
00100     !! @param var   Integer variable that is returned.
00101     !<
00102     subroutine read_tiles_int(fname,vname,var)
00103         character(len=80), intent(in) :: fname
00104         character(len=*), intent(in)  :: vname
00105         integer(i4), dimension(:,:,:) :: var
00106 
00107         integer(i4) :: mode,fh,varid,ierr
00108         integer(i4) :: nx,ny
00109         integer(i4) :: nx_read,ny_read
00110         integer(i4) :: dimid_nx, dimid_ny
00111 
00112         mode = NF90_nowrite
00113         nx = SIZE(var,dim=1)
00114         ny = SIZE(var,dim=2)
00115         ierr = sio_open(fname,mode,fh)
00116 
00117         ierr = nf90_inq_dimid(fh,'nx_block',dimid_nx)
00118         ierr = nf90_inquire_dimension(fh,dimid_nx,len=nx_read)
00119         ierr = nf90_inq_dimid(fh,'ny_block',dimid_ny)
00120         ierr = nf90_inquire_dimension(fh,dimid_ny,len=ny_read)
00121         if((nx_read .ne. nx) .or. (ny_read .ne. ny)) then 
00122             print *,'Error: mismatch between size of blocks'
00123             print *,'nx:= ',nx,'nx_read:= ',nx_read
00124             print *,'ny:= ',ny,'ny_read:= ',ny_read
00125             stop 'in read_tiles:dbl'
00126         endif
00127 
00128         ierr = nf90_inq_varid(fh,TRIM(vname),varid)
00129         if(ierr /= nf90_noerr) &
00130             call check_netcdf(fh,ierr,'nf90_inq_varid: '//vname,__LINE__)
00131 
00132         ierr = nf90_get_var(fh,varid, var)
00133         if(ierr /= nf90_noerr) &
00134             call check_netcdf(fh,ierr,'nf90_get_var: '//vname,__LINE__)
00135 
00136         ierr = sio_close(fh) 
00137     end subroutine read_tiles_int
00138 
00139     !>
00140     !!
00141     !! Writes out an integer variable to a netCDF tile file. Note 
00142     !! that this subroutine writes out all tiles. 
00143     !! 
00144     !! @param fname Filename of the tilefile to which to write the variable.
00145     !! @param vname Name of the variable to write 
00146     !! @param var   Integer variable to write.
00147     !<
00148     subroutine write_tiles_int(fname,vname,var)
00149         character(len=80), intent(in) :: fname
00150         character(len=*), intent(in) :: vname
00151         integer(i4) :: var(:,:,:)
00152 
00153         integer(i4) :: dim_nx,dim_ny,dim_nb
00154         integer(i4) :: ncid, fh, ierr, mode
00155         integer(i4) :: nx_block,ny_block,nblocks
00156 
00157         mode = NF90_write
00158         ierr = sio_open(fname,mode,fh)
00159 
00160         nx_block = SIZE(var,dim=1)
00161         ny_block = SIZE(var,dim=2)
00162         nblocks  = SIZE(var,dim=3) 
00163 
00164         ierr = nf90_redef(fh)
00165         ierr = nf90_inq_dimid(fh,'nblocks',dimid=dim_nb)
00166         if(ierr /= nf90_noerr) then 
00167             ierr = nf90_def_dim(fh,name='nblocks',len=nblocks,dimid=dim_nb)
00168             if(ierr /= nf90_noerr) &
00169                 call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00170         endif
00171 
00172         ierr = nf90_inq_dimid(fh,'nx_block',dimid=dim_nx)
00173         if(ierr /= nf90_noerr) then 
00174             ierr = nf90_def_dim(fh,name='nx_block',len=nx_block,dimid=dim_nx)
00175             if(ierr /= nf90_noerr) &
00176                 call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00177         endif
00178 
00179         ierr = nf90_inq_dimid(fh,'ny_block',dimid=dim_ny)
00180         if(ierr /= nf90_noerr) then 
00181             ierr = nf90_def_dim(fh,name='ny_block',len=ny_block,dimid=dim_ny)
00182             if(ierr /= nf90_noerr) &
00183                 call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00184         endif
00185 
00186         ierr = nf90_def_var( &
00187             fh,name=TRIM(vname),xtype=NF90_INT,dimids= &
00188                 (/dim_nx,dim_ny,dim_nb/),varid=ncid)
00189         if(ierr /= nf90_noerr) &
00190             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00191 
00192         !------------------------
00193         ! end netCDF define mode
00194         !------------------------
00195         ierr = nf90_enddef(fh)
00196         if(ierr /= nf90_noerr) &
00197             call check_netcdf(fh,ierr,'nf90_enddef',__LINE__)
00198 
00199         !--------------------------------
00200         ! write the variables to the file
00201         !---------------------------------
00202         ierr = nf90_put_var(fh,ncid,var)
00203         if(ierr /= nf90_noerr) &
00204             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00205         ierr = nf90_close(fh)
00206     end subroutine write_tiles_int
00207 
00208     !>
00209     !!
00210     !! Writes out an double precision variable to a netCDF tile file. Note 
00211     !! that this subroutine writes out all tiles. 
00212     !! 
00213     !! @param fname Filename of the tilefile to which to write the variable.
00214     !! @param vname Name of the variable to write 
00215     !! @param var   Double precision variable to write.
00216     !<
00217     subroutine write_tiles_dbl(fname,vname,var)
00218         character(len=80), intent(in) :: fname
00219         character(len=*), intent(in) :: vname
00220         real(r8) :: var(:,:,:)
00221 
00222         integer(i4) :: dim_nx,dim_ny,dim_nb
00223         integer(i4) :: ncid,fh,ierr, mode
00224         integer(i4) :: nx_block,ny_block,nblocks
00225 
00226         mode = NF90_write
00227         ierr = sio_open(fname,mode,fh)
00228 
00229         nx_block = SIZE(var,dim=1)
00230         ny_block = SIZE(var,dim=2)
00231         nblocks  = SIZE(var,dim=3)
00232 
00233         ierr = nf90_redef(fh)
00234         ierr = nf90_inq_dimid(fh,'nblocks',dimid=dim_nb)
00235         if(ierr /= nf90_noerr) then
00236             ierr = nf90_def_dim(fh,name='nblocks',len=nblocks,dimid=dim_nb)
00237             if(ierr /= nf90_noerr) &
00238                 call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00239         endif
00240 
00241         ierr = nf90_inq_dimid(fh,'nx_block',dimid=dim_nx)
00242         if(ierr /= nf90_noerr) then
00243             ierr = nf90_def_dim(fh,name='nx_block',len=nx_block,dimid=dim_nx)
00244             if(ierr /= nf90_noerr) &
00245                 call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00246         endif
00247 
00248         ierr = nf90_inq_dimid(fh,'ny_block',dimid=dim_ny)
00249         if(ierr /= nf90_noerr) then
00250             ierr = nf90_def_dim(fh,name='ny_block',len=ny_block,dimid=dim_ny)
00251             if(ierr /= nf90_noerr) &
00252                 call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00253         endif
00254 
00255         ierr = nf90_def_var( &
00256             fh,name=TRIM(vname),xtype=NF90_DOUBLE, &
00257             dimids=(/dim_nx,dim_ny,dim_nb/),varid=ncid)
00258         if(ierr /= nf90_noerr) &
00259             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00260 
00261         !------------------------
00262         ! end netCDF define mode
00263         !------------------------
00264         ierr = nf90_enddef(fh)
00265         if(ierr /= nf90_noerr) &
00266             call check_netcdf(fh,ierr,'nf90_enddef',__LINE__)
00267 
00268         !--------------------------------
00269         ! write the variables to the file
00270         !---------------------------------
00271         ierr = nf90_put_var(fh,ncid,var)
00272         if(ierr /= nf90_noerr) &
00273             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00274         ierr = nf90_close(fh)
00275     end subroutine write_tiles_dbl
00276 
00277     !>
00278     !! This subroutine writes application specific metadata to the tilefiles.
00279     !! Note that this metadata is common among all implementations of the
00280     !! miniapp
00281     !!
00282     !! @param fname   Filename which to write the application metadata.   
00283     !! @param blocks  The blocks data structure which contains the application
00284     !!                metadata
00285     !! @param reorder An interger array which reorder the cannonical block
00286     !!                layout  into a space-filling curve based layout.  
00287     !<
00288     subroutine write_AppMD(fname,blocks,reorder)
00289         character(len=80), intent(in) :: fname
00290         type (AppMD_t) :: blocks(:)
00291         integer (i4) :: reorder(:)
00292 
00293         integer(i4) :: mode, fh, ierr
00294         integer(i4) :: dim_nb,dim_nx,dim_ny
00295         integer(i4) :: ncid_block_id,ncid_ib,ncid_ie,ncid_jb,ncid_je, 
00296             ncid_iblock,ncid_jbloc,ncid_jblock,ncid_npoints,ncid_iglob, 
00297             ncid_jglob, ncid_reorder
00298         integer(i4), allocatable :: iglob(:,:), jglob(:,:)
00299         integer(i4) :: nblocks, iblock, jblock, i, nx_block, ny_block
00300 
00301         mode = NF90_CLOBBER
00302         ierr = sio_create(fname,mode,fh)
00303 
00304         nblocks = SIZE(blocks)
00305         nx_block = SIZE(blocks(1)%i_glob(:))
00306         ny_block = SIZE(blocks(1)%j_glob(:))
00307 
00308         allocate(iglob(nx_block,nblocks))
00309         allocate(jglob(ny_block,nblocks))
00310 
00311         !---------------------------------------------------------
00312         ! define the dimensions for the file
00313         !---------------------------------------------------------
00314         ierr = nf90_def_dim(fh,name='nblocks',len=nblocks,dimid=dim_nb)
00315         if(ierr /= nf90_noerr) &
00316             call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00317         ierr = nf90_def_dim(fh,name='nx_block',len=nx_block,dimid=dim_nx)
00318         if(ierr /= nf90_noerr) &
00319             call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00320         ierr = nf90_def_dim(fh,name='ny_block',len=ny_block,dimid=dim_ny)
00321         if(ierr /= nf90_noerr) &
00322             call check_netcdf(fh,ierr,'nf90_def_dim',__LINE__)
00323 
00324         !---------------------------------------------------------
00325         ! define the variables which will be contained in the file
00326         !---------------------------------------------------------
00327         ierr = nf90_def_var( &
00328             fh,name='block_id', xtype = NF90_INT, &
00329             dimids=(/dim_nb/),varid=ncid_block_id)
00330         if(ierr /= nf90_noerr) &
00331             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00332 
00333         ierr = nf90_def_var( &
00334             fh,name='ib', xtype = NF90_INT, dimids=(/dim_nb/),varid=ncid_ib)
00335         if(ierr /= nf90_noerr) &
00336             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00337 
00338         ierr = nf90_def_var( &
00339             fh,name='ie', xtype = NF90_INT, dimids=(/dim_nb/),varid=ncid_ie)
00340         if(ierr /= nf90_noerr) &
00341             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00342 
00343         ierr = nf90_def_var( &
00344             fh,name='jb', xtype = NF90_INT, dimids=(/dim_nb/),varid=ncid_jb)
00345         if(ierr /= nf90_noerr)  &
00346             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00347 
00348         ierr = nf90_def_var( &
00349             fh,name='je', xtype = NF90_INT, dimids=(/dim_nb/),varid=ncid_je)
00350         if(ierr /= nf90_noerr)  &
00351             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00352 
00353         ierr = nf90_def_var( &
00354             fh,name='iblock', xtype = NF90_INT, dimids=(/dim_nb/), &
00355             varid=ncid_iblock)
00356         if(ierr /= nf90_noerr)  &
00357             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00358 
00359         ierr = nf90_def_var( &
00360             fh,name='jblock', xtype = NF90_INT, dimids=(/dim_nb/), &
00361             varid=ncid_jblock)
00362         if(ierr /= nf90_noerr)  &
00363             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00364 
00365         ierr = nf90_def_var( &
00366             fh,name='npoints', xtype = NF90_INT, dimids=(/dim_nb/), &
00367             varid=ncid_npoints)
00368         if(ierr /= nf90_noerr)  &
00369             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00370 
00371         ierr = nf90_def_var( &
00372             fh,name='i_glob', xtype = NF90_INT, dimids=(/dim_nx,dim_nb/), &
00373             varid=ncid_iglob)
00374         if(ierr /= nf90_noerr)  &
00375             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00376 
00377         ierr = nf90_def_var( &
00378             fh,name='j_glob', xtype = NF90_INT, dimids=(/dim_ny,dim_nb/), &
00379             varid=ncid_jglob)
00380         if(ierr /= nf90_noerr) &
00381             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00382 
00383         ierr = nf90_def_var( &
00384             fh,name='reorder', xtype = NF90_INT, dimids=(/dim_nb/), &
00385             varid=ncid_reorder)
00386         if(ierr /= nf90_noerr)  &
00387             call check_netcdf(fh,ierr,'nf90_def_var',__LINE__)
00388 
00389         !------------------------
00390         ! end netCDF define mode
00391         !------------------------
00392         ierr = nf90_enddef(fh)
00393         if(ierr /= nf90_noerr)  &
00394             call check_netcdf(fh,ierr,'nf90_enddef',__LINE__)
00395 
00396         !--------------------------------
00397         ! write the variables to the file
00398         !---------------------------------
00399         ierr = nf90_put_var(fh,ncid_block_id,blocks(:)%block_id)
00400         if(ierr /= nf90_noerr)  &
00401             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00402         ierr = nf90_put_var(fh,ncid_ib,blocks(:)%ib)
00403         if(ierr /= nf90_noerr)  &
00404             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00405         ierr = nf90_put_var(fh,ncid_ie,blocks(:)%ie)
00406         if(ierr /= nf90_noerr)  &
00407             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00408         ierr = nf90_put_var(fh,ncid_jb,blocks(:)%jb)
00409         if(ierr /= nf90_noerr)  &
00410             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00411         ierr = nf90_put_var(fh,ncid_je,blocks(:)%je)
00412         if(ierr /= nf90_noerr)  &
00413             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00414         ierr = nf90_put_var(fh,ncid_iblock,blocks(:)%iblock)
00415         if(ierr /= nf90_noerr)  &
00416             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00417         ierr = nf90_put_var(fh,ncid_jblock,blocks(:)%jblock)
00418         if(ierr /= nf90_noerr)  &
00419             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00420         ierr = nf90_put_var(fh,ncid_npoints,blocks(:)%npoints)
00421         if(ierr /= nf90_noerr)  &
00422             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00423 
00424         do i=1,nblocks
00425             iglob(:,i) = blocks(i)%i_glob(:)
00426             jglob(:,i) = blocks(i)%j_glob(:)
00427         enddo
00428 
00429         ierr = nf90_put_var(fh,ncid_iglob,iglob)        
00430         if(ierr /= nf90_noerr) &
00431             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00432         ierr = nf90_put_var(fh,ncid_jglob,jglob)        
00433         if(ierr /= nf90_noerr) &
00434             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00435         ierr = nf90_put_var(fh,ncid_reorder,reorder)        
00436         if(ierr /= nf90_noerr) &
00437             call check_netcdf(fh,ierr,'nf90_put_var',__LINE__)
00438         ierr = nf90_close(fh)
00439     end subroutine write_AppMD
00440 
00441     !>
00442     !! This subroutine reads application specific metadata from the tilefiles.
00443     !! Note that this metadata is common among all implementations of the
00444     !! miniapp
00445     !!
00446     !! @param fname   Filename from which to read the application metadata.   
00447     !! @param blocks  The blocks data structure which contains the application
00448     !! metadata
00449     !! @param reorder An interger array which reorders the cannonical block
00450     !! layout  into a space-filling curve based layout.  
00451     !<
00452     subroutine read_AppMD(fname,blocks,reorder)
00453         character(len=80), intent(in) :: fname
00454         type (AppMD_t), pointer, intent(inout)  :: blocks(:)
00455         integer (i4),pointer :: reorder(:)
00456 
00457         integer (i4) :: mode, fh, ierr, i, j
00458         integer (i4) :: nblocks, nx_block,ny_block
00459         integer (i4) :: dim_nb, dim_nx, dim_ny, varid
00460 
00461         integer (i4), allocatable :: iglob(:,:), jglob(:,:)
00462 
00463         mode = NF90_nowrite
00464         ierr = sio_open(fname,mode,fh)
00465 
00466         ierr = nf90_inq_dimid(fh,'nblocks',dimid=dim_nb)
00467         if(ierr /= nf90_noerr) &
00468             call check_netcdf(fh,ierr,'nf90_inq_dimid',__LINE__)
00469 
00470         ierr = nf90_inq_dimid(fh,'nx_block',dimid=dim_nx)
00471         if(ierr /= nf90_noerr) &
00472             call check_netcdf(fh,ierr,'nf90_inq_dimid',__LINE__)
00473 
00474         ierr = nf90_inq_dimid(fh,'ny_block',dimid=dim_ny)
00475         if(ierr /= nf90_noerr) &
00476             call check_netcdf(fh,ierr,'nf90_inq_dimid',__LINE__)
00477 
00478         ierr = nf90_inquire_dimension(fh,dim_nb,len=nblocks)
00479         if(ierr /= nf90_noerr) &
00480             call check_netcdf(fh,ierr,'nf90_inquire_dimension',__LINE__)
00481 
00482         ierr = nf90_inquire_dimension(fh,dim_nx,len=nx_block)
00483         if(ierr /= nf90_noerr) &
00484             call check_netcdf(fh,ierr,'nf90_inquire_dimension',__LINE__)
00485 
00486         ierr = nf90_inquire_dimension(fh,dim_ny,len=ny_block)
00487         if(ierr /= nf90_noerr) &
00488             call check_netcdf(fh,ierr,'nf90_inquire_dimension',__LINE__)
00489 
00490         allocate(blocks(nblocks))
00491 
00492         !--------------------------------------------------------------------
00493         ! Read and load the different componets of the AppMD_t data structure
00494         !--------------------------------------------------------------------
00495         ierr = nf90_inq_varid(fh,'block_id',varid)
00496         if(ierr /= nf90_noerr) &
00497             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00498 
00499         ierr = nf90_get_var(fh,varid, blocks(:)%block_id)
00500         if(ierr /= nf90_noerr) &
00501             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00502 
00503         ierr = nf90_inq_varid(fh,'ib',varid)
00504         if(ierr /= nf90_noerr) &
00505             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00506 
00507         ierr = nf90_get_var(fh,varid, blocks(:)%ib)
00508         if(ierr /= nf90_noerr) &
00509             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00510 
00511         ierr = nf90_inq_varid(fh,'ie',varid)
00512         if(ierr /= nf90_noerr) &
00513             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00514 
00515         ierr = nf90_get_var(fh,varid, blocks(:)%ie)
00516         if(ierr /= nf90_noerr) &
00517             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00518 
00519         ierr = nf90_inq_varid(fh,'jb',varid)
00520         if(ierr /= nf90_noerr) &
00521             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00522 
00523         ierr = nf90_get_var(fh,varid, blocks(:)%jb)
00524         if(ierr /= nf90_noerr) &
00525             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00526 
00527         ierr = nf90_inq_varid(fh,'je',varid)
00528         if(ierr /= nf90_noerr) &
00529             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00530 
00531         ierr = nf90_get_var(fh,varid, blocks(:)%je)
00532         if(ierr /= nf90_noerr) &
00533             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00534 
00535         ierr = nf90_inq_varid(fh,'iblock',varid)
00536         if(ierr /= nf90_noerr) &
00537             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00538 
00539         ierr = nf90_get_var(fh,varid, blocks(:)%iblock)
00540         if(ierr /= nf90_noerr) &
00541             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00542 
00543         ierr = nf90_inq_varid(fh,'jblock',varid)
00544         if(ierr /= nf90_noerr) &
00545             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00546 
00547         ierr = nf90_get_var(fh,varid, blocks(:)%jblock)
00548         if(ierr /= nf90_noerr) &
00549             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00550 
00551         ierr = nf90_inq_varid(fh,'npoints',varid)
00552         if(ierr /= nf90_noerr) &
00553             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00554 
00555         ierr = nf90_get_var(fh,varid, blocks(:)%npoints)
00556         if(ierr /= nf90_noerr) &
00557             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00558 
00559         allocate(iglob(nx_block,nblocks))
00560         allocate(jglob(ny_block,nblocks))
00561 
00562         ierr = nf90_inq_varid(fh,'i_glob',varid)
00563         if(ierr /= nf90_noerr) &
00564             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00565         ierr = nf90_get_var(fh,varid, iglob)
00566         if(ierr /= nf90_noerr) &
00567             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00568         ierr = nf90_inq_varid(fh,'j_glob',varid)
00569         if(ierr /= nf90_noerr) &
00570             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00571         ierr = nf90_get_var(fh,varid, jglob)
00572         if(ierr /= nf90_noerr) &
00573             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00574         
00575         do i=1,nblocks
00576             allocate(blocks(i)%i_glob(nx_block))
00577             allocate(blocks(i)%j_glob(ny_block))
00578             blocks(i)%i_glob(:) = iglob(:,i)
00579             blocks(i)%j_glob(:) = jglob(:,i)
00580         enddo
00581         
00582         deallocate(iglob,jglob)
00583         allocate(reorder(nblocks))
00584         
00585         ierr = nf90_inq_varid(fh,'reorder',varid)
00586         if(ierr /= nf90_noerr) &
00587             call check_netcdf(fh,ierr,'nf90_inq_varid',__LINE__)
00588         
00589         ierr = nf90_get_var(fh,varid, reorder)
00590         if(ierr /= nf90_noerr) &
00591             call check_netcdf(fh,ierr,'nf90_get_var',__LINE__)
00592     end subroutine read_AppMD
00593 
00594     !>
00595     !! Creates a netCDF file
00596     !!
00597     !! @param fname  Filename of the file to be created.
00598     !! @param mode   The netCDF file mode  
00599     !! @param fh     An integer which identifies the created netCDF file.
00600     !<
00601     function sio_create(fname,mode,fh)
00602         character(len=*), intent(in) :: fname
00603         integer(i4), intent(in) :: mode
00604         integer(i4), intent(out) :: fh
00605         integer(i4) :: sio_create
00606         integer(i4) :: fmode, ierr
00607         
00608         ierr = nf90_create(fname,mode,fh)
00609         ierr = nf90_set_fill(fh,NF90_NOFILL,fmode)
00610         
00611         sio_create = ierr
00612     end function sio_create
00613 
00614     !>
00615     !! Opens an existing netCDF file
00616     !!
00617     !! @param fname  The name of the file to ope.
00618     !! @param mode   The netCDF file mode  
00619     !! @param fh     An integer which identifies the netCDF file.
00620     !<
00621     function sio_open(fname,mode,fh)
00622         character(len=*), intent(in) :: fname
00623         integer(i4), intent(in) :: mode
00624         integer(i4), intent(out) :: fh
00625         integer(i4) :: sio_open
00626         integer(i4) :: ierr
00627 
00628         ierr = nf90_open(fname,mode,fh)
00629         if(ierr /= nf90_noerr) then
00630             call check_netcdf(fh,ierr,'nf90_open(): '//fname,__LINE__)
00631             call exit_POP(1, "Issue encountered while opening netCDF file")
00632         end if
00633         sio_open= ierr
00634     end function sio_open
00635 
00636     !>
00637     !! Closes an open netCDF file
00638     !!
00639     !! @param ncid The netCDF file identifier. 
00640     !<
00641     function sio_close(ncid)
00642         integer(i4), intent(inout) :: ncid
00643         integer(i4) :: sio_close
00644         integer(i4) :: ierr
00645 
00646         ierr = nf90_sync(ncid)
00647         ierr = nf90_close(ncid)
00648 
00649         sio_close = ierr
00650     end function sio_close
00651 
00652     !>
00653     !! This subroutine will print out a netCDF error message
00654     !!
00655     !! @param fh       The netCDF file identifier. 
00656     !! @param status   The error return code.
00657     !! @param filestr  A user supplied string to identify the location of the
00658     !!                 netCDF error
00659     !! @param line     The file line number which generated the error  
00660     !<
00661     subroutine check_netcdf(fh,status,filestr,line)
00662         integer(i4) :: fh
00663         integer(i4) :: status
00664         character(len=*), intent(in) :: filestr
00665         integer, intent(in) :: line
00666 
00667         print *,'ERROR: ',TRIM(filestr),' ',TRIM(nf90_strerror(status)), &
00668             ' on line: ',line 
00669     end subroutine check_netcdf
00670 end module io_serial
 All Classes Namespaces Files Functions Variables