! SVN:$Id: ice_restart.F90 1304 2015-11-12 19:28:14Z afrobert@nps.edu $ !======================================================================= ! Read and write ice model restart files using netCDF or binary ! interfaces. ! authors David A Bailey, NCAR module ice_restart use ice_broadcast use ice_exit, only: abort_ice use ice_kinds_mod use netcdf use ice_restart_shared, only: & restart, restart_ext, restart_dir, restart_file, pointer_file, & runid, runtype, use_restart_time, restart_format, lcdf64, lenstr implicit none private public :: init_restart_write, init_restart_read, & read_restart_field, write_restart_field, final_restart integer (kind=int_kind) :: ncid !======================================================================= contains !======================================================================= ! Sets up restart file for reading. ! author David A Bailey, NCAR subroutine init_restart_read(ice_ic) use ice_calendar, only: sec, month, mday, nyr, istep0, istep1, & time, time_forc, year_init, npt use ice_communicate, only: my_task, master_task use ice_domain, only: nblocks use ice_fileunits, only: nu_diag, nu_rst_pointer character(len=char_len_long), intent(in), optional :: ice_ic ! local variables character(len=char_len_long) :: & filename, filename0 integer (kind=int_kind) :: status if (present(ice_ic)) then filename = trim(ice_ic) else if (my_task == master_task) then open(nu_rst_pointer,file=pointer_file) read(nu_rst_pointer,'(a)') filename0 filename = trim(filename0) close(nu_rst_pointer) write(nu_diag,*) 'Read ',pointer_file(1:lenstr(pointer_file)) endif call broadcast_scalar(filename, master_task) endif if (my_task == master_task) then write(nu_diag,*) 'Using restart dump=', trim(filename) status = nf90_open(trim(filename), nf90_nowrite, ncid) if (status /= nf90_noerr) call abort_ice( & 'ice: Error reading restart ncfile '//trim(filename)) if (use_restart_time) then status = nf90_get_att(ncid, nf90_global, 'istep1', istep0) status = nf90_get_att(ncid, nf90_global, 'time', time) status = nf90_get_att(ncid, nf90_global, 'time_forc', time_forc) status = nf90_get_att(ncid, nf90_global, 'nyr', nyr) if (status == nf90_noerr) then status = nf90_get_att(ncid, nf90_global, 'month', month) status = nf90_get_att(ncid, nf90_global, 'mday', mday) status = nf90_get_att(ncid, nf90_global, 'sec', sec) endif endif ! use namelist values if use_restart_time = F write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc endif call broadcast_scalar(istep0,master_task) call broadcast_scalar(time,master_task) call broadcast_scalar(time_forc,master_task) istep1 = istep0 ! if runid is bering then need to correct npt for istep0 if (trim(runid) == 'bering') then npt = npt - istep0 endif end subroutine init_restart_read !======================================================================= ! Sets up restart file for writing. ! author David A Bailey, NCAR subroutine init_restart_write(filename_spec) use ice_blocks, only: nghost use ice_calendar, only: sec, month, mday, nyr, istep1, & time, time_forc, year_init use ice_communicate, only: my_task, master_task use ice_domain_size, only: nx_global, ny_global, ncat, nilyr, nslyr, & n_aero use ice_dyn_shared, only: kdyn use ice_fileunits, only: nu_diag, nu_rst_pointer use ice_ocean, only: oceanmixed_ice use ice_state, only: tr_iage, tr_FY, tr_lvl, tr_aero, tr_pond_cesm, & tr_pond_topo, tr_pond_lvl, tr_brine use ice_zbgc_shared, only: tr_bgc_N_sk, tr_bgc_C_sk, tr_bgc_Nit_sk, & tr_bgc_Sil_sk, tr_bgc_DMSPp_sk, tr_bgc_DMS_sk, & tr_bgc_chl_sk, tr_bgc_DMSPd_sk, tr_bgc_Am_sk, & skl_bgc character(len=char_len_long), intent(in), optional :: filename_spec ! local variables integer (kind=int_kind) :: & k, & ! index nx, ny, & ! global array size iyear, imonth, iday ! year, month, day character(len=char_len_long) :: filename integer (kind=int_kind), allocatable :: dims(:) integer (kind=int_kind) :: & dimid_ni, & ! netCDF identifiers dimid_nj, & ! dimid_ncat, & ! iflag, & ! netCDF creation flag status ! status variable from netCDF routine character (len=3) :: nchar ! construct path/file if (present(filename_spec)) then filename = trim(filename_spec) else iyear = nyr + year_init - 1 imonth = month iday = mday write(filename,'(a,a,a,i4.4,a,i2.2,a,i2.2,a,i5.5)') & restart_dir(1:lenstr(restart_dir)), & restart_file(1:lenstr(restart_file)),'.', & iyear,'-',month,'-',mday,'-',sec end if ! write pointer (path/file) if (my_task == master_task) then filename = trim(filename) // '.nc' open(nu_rst_pointer,file=pointer_file) write(nu_rst_pointer,'(a)') filename close(nu_rst_pointer) iflag = 0 if (lcdf64) iflag = nf90_64bit_offset status = nf90_create(trim(filename), iflag, ncid) if (status /= nf90_noerr) call abort_ice( & 'ice: Error creating restart ncfile '//trim(filename)) status = nf90_put_att(ncid,nf90_global,'istep1',istep1) status = nf90_put_att(ncid,nf90_global,'time',time) status = nf90_put_att(ncid,nf90_global,'time_forc',time_forc) status = nf90_put_att(ncid,nf90_global,'nyr',nyr) status = nf90_put_att(ncid,nf90_global,'month',month) status = nf90_put_att(ncid,nf90_global,'mday',mday) status = nf90_put_att(ncid,nf90_global,'sec',sec) nx = nx_global ny = ny_global if (restart_ext) then nx = nx_global + 2*nghost ny = ny_global + 2*nghost endif status = nf90_def_dim(ncid,'ni',nx,dimid_ni) status = nf90_def_dim(ncid,'nj',ny,dimid_nj) status = nf90_def_dim(ncid,'ncat',ncat,dimid_ncat) !----------------------------------------------------------------- ! 2D restart fields !----------------------------------------------------------------- allocate(dims(2)) dims(1) = dimid_ni dims(2) = dimid_nj call define_rest_field(ncid,'uvel',dims) call define_rest_field(ncid,'vvel',dims) #ifdef CCSMCOUPLED call define_rest_field(ncid,'coszen',dims) #endif call define_rest_field(ncid,'scale_factor',dims) call define_rest_field(ncid,'swvdr',dims) call define_rest_field(ncid,'swvdf',dims) call define_rest_field(ncid,'swidr',dims) call define_rest_field(ncid,'swidf',dims) call define_rest_field(ncid,'strocnxT',dims) call define_rest_field(ncid,'strocnyT',dims) call define_rest_field(ncid,'stressp_1',dims) call define_rest_field(ncid,'stressp_2',dims) call define_rest_field(ncid,'stressp_3',dims) call define_rest_field(ncid,'stressp_4',dims) call define_rest_field(ncid,'stressm_1',dims) call define_rest_field(ncid,'stressm_2',dims) call define_rest_field(ncid,'stressm_3',dims) call define_rest_field(ncid,'stressm_4',dims) call define_rest_field(ncid,'stress12_1',dims) call define_rest_field(ncid,'stress12_2',dims) call define_rest_field(ncid,'stress12_3',dims) call define_rest_field(ncid,'stress12_4',dims) call define_rest_field(ncid,'iceumask',dims) if (oceanmixed_ice) then call define_rest_field(ncid,'sst',dims) call define_rest_field(ncid,'frzmlt',dims) endif if (tr_FY) then call define_rest_field(ncid,'frz_onset',dims) endif if (kdyn == 2) then call define_rest_field(ncid,'a11_1',dims) call define_rest_field(ncid,'a11_2',dims) call define_rest_field(ncid,'a11_3',dims) call define_rest_field(ncid,'a11_4',dims) call define_rest_field(ncid,'a12_1',dims) call define_rest_field(ncid,'a12_2',dims) call define_rest_field(ncid,'a12_3',dims) call define_rest_field(ncid,'a12_4',dims) endif if (tr_pond_lvl) then call define_rest_field(ncid,'fsnow',dims) endif if (skl_bgc) then call define_rest_field(ncid,'algalN',dims) call define_rest_field(ncid,'nit' ,dims) if (tr_bgc_Am_sk) & call define_rest_field(ncid,'amm' ,dims) if (tr_bgc_Sil_sk) & call define_rest_field(ncid,'sil' ,dims) if (tr_bgc_DMSPp_sk) & call define_rest_field(ncid,'dmsp' ,dims) if (tr_bgc_DMS_sk) & call define_rest_field(ncid,'dms' ,dims) endif deallocate(dims) !----------------------------------------------------------------- ! 3D restart fields (ncat) !----------------------------------------------------------------- allocate(dims(3)) dims(1) = dimid_ni dims(2) = dimid_nj dims(3) = dimid_ncat call define_rest_field(ncid,'aicen',dims) call define_rest_field(ncid,'vicen',dims) call define_rest_field(ncid,'vsnon',dims) call define_rest_field(ncid,'Tsfcn',dims) if (tr_iage) then call define_rest_field(ncid,'iage',dims) end if if (tr_FY) then call define_rest_field(ncid,'FY',dims) end if if (tr_lvl) then call define_rest_field(ncid,'alvl',dims) call define_rest_field(ncid,'vlvl',dims) end if if (tr_pond_cesm) then call define_rest_field(ncid,'apnd',dims) call define_rest_field(ncid,'hpnd',dims) end if if (tr_pond_topo) then call define_rest_field(ncid,'apnd',dims) call define_rest_field(ncid,'hpnd',dims) call define_rest_field(ncid,'ipnd',dims) end if if (tr_pond_lvl) then call define_rest_field(ncid,'apnd',dims) call define_rest_field(ncid,'hpnd',dims) call define_rest_field(ncid,'ipnd',dims) call define_rest_field(ncid,'dhs',dims) call define_rest_field(ncid,'ffrac',dims) end if if (tr_brine) then call define_rest_field(ncid,'fbrn',dims) call define_rest_field(ncid,'first_ice',dims) endif if (skl_bgc) then call define_rest_field(ncid,'bgc_N_sk' ,dims) call define_rest_field(ncid,'bgc_Nit_sk' ,dims) if (tr_bgc_C_sk) & call define_rest_field(ncid,'bgc_C_sk' ,dims) if (tr_bgc_chl_sk) & call define_rest_field(ncid,'bgc_chl_sk' ,dims) if (tr_bgc_Am_sk) & call define_rest_field(ncid,'bgc_Am_sk' ,dims) if (tr_bgc_Sil_sk) & call define_rest_field(ncid,'bgc_Sil_sk' ,dims) if (tr_bgc_DMSPp_sk) & call define_rest_field(ncid,'bgc_DMSPp_sk',dims) if (tr_bgc_DMSPd_sk) & call define_rest_field(ncid,'bgc_DMSPd_sk',dims) if (tr_bgc_DMS_sk) & call define_rest_field(ncid,'bgc_DMS_sk' ,dims) endif !----------------------------------------------------------------- ! 4D restart fields, written as layers of 3D !----------------------------------------------------------------- do k=1,nilyr write(nchar,'(i3.3)') k call define_rest_field(ncid,'sice'//trim(nchar),dims) call define_rest_field(ncid,'qice'//trim(nchar),dims) enddo do k=1,nslyr write(nchar,'(i3.3)') k call define_rest_field(ncid,'qsno'//trim(nchar),dims) enddo if (tr_aero) then do k=1,n_aero write(nchar,'(i3.3)') k call define_rest_field(ncid,'aerosnossl'//trim(nchar),dims) call define_rest_field(ncid,'aerosnoint'//trim(nchar),dims) call define_rest_field(ncid,'aeroicessl'//trim(nchar),dims) call define_rest_field(ncid,'aeroiceint'//trim(nchar),dims) enddo endif deallocate(dims) status = nf90_enddef(ncid) write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif ! master_task end subroutine init_restart_write !======================================================================= ! Reads a single restart field ! author David A Bailey, NCAR subroutine read_restart_field(nu,nrec,work,atype,vname,ndim3, & diag, field_loc, field_type) use ice_blocks, only: nx_block, ny_block use ice_domain_size, only: max_blocks, ncat use ice_fileunits, only: nu_diag use ice_read_write, only: ice_read, ice_read_nc integer (kind=int_kind), intent(in) :: & nu , & ! unit number (not used for netcdf) ndim3 , & ! third dimension nrec ! record number (0 for sequential access) real (kind=dbl_kind), dimension(nx_block,ny_block,ndim3,max_blocks), & intent(inout) :: & work ! input array (real, 8-byte) character (len=4), intent(in) :: & atype ! format for output array ! (real/integer, 4-byte/8-byte) logical (kind=log_kind), intent(in) :: & diag ! if true, write diagnostic output character (len=*), intent(in) :: vname integer (kind=int_kind), optional, intent(in) :: & field_loc, & ! location of field on staggered grid field_type ! type of field (scalar, vector, angle) ! local variables integer (kind=int_kind) :: & n, & ! number of dimensions for variable varid, & ! variable id status ! status variable from netCDF routine real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) :: & work2 ! input array (real, 8-byte) if (present(field_loc)) then if (ndim3 == ncat) then if (restart_ext) then call ice_read_nc(ncid,1,vname,work,diag, & field_loc=field_loc,field_type=field_type,restart_ext=restart_ext) else call ice_read_nc(ncid,1,vname,work,diag,field_loc,field_type) endif elseif (ndim3 == 1) then if (restart_ext) then call ice_read_nc(ncid,1,vname,work2,diag, & field_loc=field_loc,field_type=field_type,restart_ext=restart_ext) else call ice_read_nc(ncid,1,vname,work2,diag,field_loc,field_type) endif work(:,:,1,:) = work2(:,:,:) else write(nu_diag,*) 'ndim3 not supported ',ndim3 endif else if (ndim3 == ncat) then if (restart_ext) then call ice_read_nc(ncid, 1, vname, work, diag, restart_ext=restart_ext) else call ice_read_nc(ncid, 1, vname, work, diag) endif elseif (ndim3 == 1) then if (restart_ext) then call ice_read_nc(ncid, 1, vname, work2, diag, restart_ext=restart_ext) else call ice_read_nc(ncid, 1, vname, work2, diag) endif work(:,:,1,:) = work2(:,:,:) else write(nu_diag,*) 'ndim3 not supported ',ndim3 endif endif end subroutine read_restart_field !======================================================================= ! Writes a single restart field. ! author David A Bailey, NCAR subroutine write_restart_field(nu,nrec,work,atype,vname,ndim3,diag) use ice_blocks, only: nx_block, ny_block use ice_domain_size, only: max_blocks, ncat use ice_fileunits, only: nu_diag use ice_read_write, only: ice_write, ice_write_nc integer (kind=int_kind), intent(in) :: & nu , & ! unit number ndim3 , & ! third dimension nrec ! record number (0 for sequential access) real (kind=dbl_kind), dimension(nx_block,ny_block,ndim3,max_blocks), & intent(in) :: & work ! input array (real, 8-byte) character (len=4), intent(in) :: & atype ! format for output array ! (real/integer, 4-byte/8-byte) logical (kind=log_kind), intent(in) :: & diag ! if true, write diagnostic output character (len=*), intent(in) :: vname ! local variables integer (kind=int_kind) :: & n, & ! dimension counter varid, & ! variable id status ! status variable from netCDF routine real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) :: & work2 ! input array (real, 8-byte) status = nf90_inq_varid(ncid,trim(vname),varid) if (ndim3 == ncat) then if (restart_ext) then call ice_write_nc(ncid, 1, varid, work, diag, restart_ext) else call ice_write_nc(ncid, 1, varid, work, diag) endif elseif (ndim3 == 1) then work2(:,:,:) = work(:,:,1,:) if (restart_ext) then call ice_write_nc(ncid, 1, varid, work2, diag, restart_ext) else call ice_write_nc(ncid, 1, varid, work2, diag) endif else write(nu_diag,*) 'ndim3 not supported',ndim3 endif end subroutine write_restart_field !======================================================================= ! Finalize the restart file. ! author David A Bailey, NCAR subroutine final_restart() use ice_calendar, only: istep1, time, time_forc use ice_communicate, only: my_task, master_task use ice_fileunits, only: nu_diag integer (kind=int_kind) :: status status = nf90_close(ncid) if (my_task == master_task) & write(nu_diag,*) 'Restart read/written ',istep1,time,time_forc end subroutine final_restart !======================================================================= ! Defines a restart field ! author David A Bailey, NCAR subroutine define_rest_field(ncid, vname, dims) character (len=*) , intent(in) :: vname integer (kind=int_kind), intent(in) :: dims(:) integer (kind=int_kind), intent(in) :: ncid integer (kind=int_kind) :: varid integer (kind=int_kind) :: & status ! status variable from netCDF routine status = nf90_def_var(ncid,trim(vname),nf90_double,dims,varid) end subroutine define_rest_field !======================================================================= end module ice_restart !=======================================================================