diff --git a/cice.setup b/cice.setup index 3efe94827..8511f4cae 100755 --- a/cice.setup +++ b/cice.setup @@ -1094,7 +1094,7 @@ cd ${testname_base} source ./cice.settings if (\${dobuild} == true) then if (\${doreuse} == true) then - set ciceexe = "../ciceexe.\${ICE_ENVNAME}.\${ICE_COMMDIR}.\${ICE_BLDDEBUG}.\${ICE_THREADED}.\${ICE_IOTYPE}" + set ciceexe = "../ciceexe.\${ICE_TARGET}.\${ICE_ENVNAME}.\${ICE_COMMDIR}.\${ICE_BLDDEBUG}.\${ICE_THREADED}.\${ICE_IOTYPE}" ./cice.build --exe \${ciceexe} if !(-e \${ciceexe}) cp -p \${ICE_RUNDIR}/cice \${ciceexe} else diff --git a/cicecore/cicedynB/analysis/ice_history.F90 b/cicecore/cicedynB/analysis/ice_history.F90 index 1aa2515a4..f91562449 100644 --- a/cicecore/cicedynB/analysis/ice_history.F90 +++ b/cicecore/cicedynB/analysis/ice_history.F90 @@ -1713,7 +1713,7 @@ subroutine accum_hist (dt) use ice_domain_size, only: nfsd use ice_grid, only: tmask, lmask_n, lmask_s, dxu, dyu use ice_calendar, only: new_year, write_history, & - write_ic, time, histfreq, nstreams, month, & + write_ic, timesecs, histfreq, nstreams, mmonth, & new_month use ice_dyn_eap, only: a11, a12, e11, e12, e22, s11, s12, s22, & yieldstress11, yieldstress12, yieldstress22 @@ -1864,7 +1864,7 @@ subroutine accum_hist (dt) avgct(ns) = avgct(ns) + c1 ! if (avgct(ns) == c1) time_beg(ns) = (time-dt)/int(secday) if (avgct(ns) == c1) then - time_beg(ns) = (time-dt)/int(secday) + time_beg(ns) = (timesecs-dt)/int(secday) time_beg(ns) = real(time_beg(ns),kind=real_kind) endif endif @@ -3966,7 +3966,7 @@ subroutine accum_hist (dt) enddo ! iblk !$OMP END PARALLEL DO - time_end(ns) = time/int(secday) + time_end(ns) = timesecs/int(secday) time_end(ns) = real(time_end(ns),kind=real_kind) !--------------------------------------------------------------- @@ -4057,7 +4057,7 @@ subroutine accum_hist (dt) enddo endif ! new_year - if ( (month .eq. 7) .and. new_month ) then + if ( (mmonth .eq. 7) .and. new_month ) then do j=jlo,jhi do i=ilo,ihi ! reset SH Jul 1 diff --git a/cicecore/cicedynB/analysis/ice_history_shared.F90 b/cicecore/cicedynB/analysis/ice_history_shared.F90 index ce177ad1e..52d268990 100644 --- a/cicecore/cicedynB/analysis/ice_history_shared.F90 +++ b/cicecore/cicedynB/analysis/ice_history_shared.F90 @@ -653,9 +653,9 @@ module ice_history_shared subroutine construct_filename(ncfile,suffix,ns) - use ice_calendar, only: sec, nyr, month, daymo, & + use ice_calendar, only: msec, myear, mmonth, daymo, & mday, write_ic, histfreq, histfreq_n, & - year_init, new_year, new_month, new_day, & + new_year, new_month, new_day, & dt use ice_restart_shared, only: lenstr @@ -667,12 +667,12 @@ subroutine construct_filename(ncfile,suffix,ns) character (len=1) :: cstream character(len=*), parameter :: subname = '(construct_filename)' - iyear = nyr + year_init - 1 ! set year_init=1 in ice_in to get iyear=nyr - imonth = month + iyear = myear + imonth = mmonth iday = mday - isec = sec - dt + isec = msec - dt - if (write_ic) isec = sec + if (write_ic) isec = msec ! construct filename if (write_ic) then write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & @@ -688,7 +688,7 @@ subroutine construct_filename(ncfile,suffix,ns) imonth = 12 iday = daymo(imonth) elseif (new_month) then - imonth = month - 1 + imonth = mmonth - 1 iday = daymo(imonth) elseif (new_day) then iday = iday - 1 @@ -703,7 +703,7 @@ subroutine construct_filename(ncfile,suffix,ns) if (histfreq(ns) == '1') then ! instantaneous, write every dt write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & history_file(1:lenstr(history_file))//trim(cstream),'_inst.', & - iyear,'-',imonth,'-',iday,'-',sec,'.',suffix + iyear,'-',imonth,'-',iday,'-',msec,'.',suffix elseif (hist_avg) then ! write averaged data @@ -714,7 +714,7 @@ subroutine construct_filename(ncfile,suffix,ns) elseif (histfreq(ns) == 'h'.or.histfreq(ns) == 'H') then ! hourly write(ncfile,'(a,a,i2.2,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & history_file(1:lenstr(history_file))//trim(cstream),'_', & - histfreq_n(ns),'h.',iyear,'-',imonth,'-',iday,'-',sec,'.',suffix + histfreq_n(ns),'h.',iyear,'-',imonth,'-',iday,'-',msec,'.',suffix elseif (histfreq(ns) == 'm'.or.histfreq(ns) == 'M') then ! monthly write(ncfile,'(a,a,i4.4,a,i2.2,a,a)') & history_file(1:lenstr(history_file))//trim(cstream),'.', & @@ -728,7 +728,7 @@ subroutine construct_filename(ncfile,suffix,ns) else ! instantaneous with histfreq > dt write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & history_file(1:lenstr(history_file)),'_inst.', & - iyear,'-',imonth,'-',iday,'-',sec,'.',suffix + iyear,'-',imonth,'-',iday,'-',msec,'.',suffix endif endif diff --git a/cicecore/cicedynB/general/ice_forcing.F90 b/cicecore/cicedynB/general/ice_forcing.F90 index edbba8101..53695ada0 100644 --- a/cicecore/cicedynB/general/ice_forcing.F90 +++ b/cicecore/cicedynB/general/ice_forcing.F90 @@ -22,8 +22,8 @@ module ice_forcing use ice_blocks, only: nx_block, ny_block use ice_domain_size, only: ncat, max_blocks, nx_global, ny_global use ice_communicate, only: my_task, master_task - use ice_calendar, only: istep, istep1, time, time_forc, & - sec, mday, month, nyr, yday, daycal, dayyr, & + use ice_calendar, only: istep, istep1, & + msec, mday, mmonth, myear, yday, daycal, & daymo, days_per_year, hc_jday use ice_fileunits, only: nu_diag, nu_forcing use ice_exit, only: abort_ice @@ -31,7 +31,7 @@ module ice_forcing ice_get_ncvarsize, ice_read_vec_nc, & ice_open_nc, ice_read_nc, ice_close_nc use ice_timers, only: ice_timer_start, ice_timer_stop, timer_readwrite, & - timer_bound + timer_bound, timer_forcing use ice_arrays_column, only: oceanmixed_ice, restore_bgc use ice_constants, only: c0, c1, c2, c3, c4, c5, c8, c10, c12, c15, c20, & c180, c360, c365, c1000, c3600 @@ -53,10 +53,10 @@ module ice_forcing read_data_nc_point, interp_coeff integer (kind=int_kind), public :: & - ycycle , & ! number of years in forcing cycle - fyear_init , & ! first year of data in forcing cycle - fyear , & ! current year in forcing cycle - fyear_final ! last year in cycle + ycycle , & ! number of years in forcing cycle, set by namelist + fyear_init , & ! first year of data in forcing cycle, set by namelist + fyear , & ! current year in forcing cycle, varying during the run + fyear_final ! last year in cycle, computed at init character (char_len_long) :: & ! input data file names uwind_file, & @@ -80,8 +80,7 @@ module ice_forcing botmelt_file real (kind=dbl_kind), public :: & - c1intp, c2intp , & ! interpolation coefficients - ftime ! forcing time (for restart) + c1intp, c2intp ! interpolation coefficients integer (kind=int_kind) :: & oldrecnum = 0 , & ! old record number (save between steps) @@ -167,6 +166,12 @@ module ice_forcing integer (kind=int_kind), public :: & Njday_atm ! Number of atm forcing timesteps + + ! PRIVATE: + + logical (kind=log_kind), parameter :: & + forcing_debug = .false. ! local debug flag + !======================================================================= contains @@ -177,6 +182,9 @@ module ice_forcing ! subroutine alloc_forcing integer (int_kind) :: ierr + character(len=*), parameter :: subname = '(alloc_forcing)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' allocate ( & cldf(nx_block,ny_block, max_blocks), & ! cloud fraction @@ -221,14 +229,23 @@ subroutine init_forcing_atmo use ice_calendar, only: use_leap_years + integer (kind=int_kind) :: modadj ! adjustment for mod function character(len=*), parameter :: subname = '(init_forcing_atmo)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + ! Allocate forcing arrays call alloc_forcing() - fyear = fyear_init + mod(nyr-1,ycycle) ! current year + modadj = abs((min(0,myear-fyear_init)/ycycle+1)*ycycle) + fyear = fyear_init + mod(myear-fyear_init+modadj,ycycle) fyear_final = fyear_init + ycycle - 1 ! last year in forcing cycle + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg fyear = ',fyear,fyear_init,fyear_final + write(nu_diag,*) subname,'fdbg atm_data_type = ',trim(atm_data_type) + endif + if (trim(atm_data_type) /= 'default' .and. & my_task == master_task) then write (nu_diag,*) ' Initial forcing data year = ',fyear_init @@ -327,6 +344,8 @@ subroutine init_forcing_ocn(dt) character(len=*), parameter :: subname = '(init_forcing_ocn)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & @@ -415,7 +434,7 @@ subroutine init_forcing_ocn(dt) if (my_task == master_task) & call ice_open (nu_forcing, sst_file, nbits) - call ice_read (nu_forcing, month, sst, 'rda8', dbug, & + call ice_read (nu_forcing, mmonth, sst, 'rda8', dbug, & field_loc_center, field_type_scalar) if (my_task == master_task) close(nu_forcing) @@ -451,7 +470,7 @@ subroutine init_forcing_ocn(dt) endif fieldname='sst' - call ice_read_nc(fid,month,fieldname,sst,diag) + call ice_read_nc(fid,mmonth,fieldname,sst,diag) if (my_task == master_task) call ice_close_nc(fid) @@ -499,6 +518,8 @@ subroutine ocn_freezing_temperature character(len=*), parameter :: subname = '(ocn_freezing_temperature)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + !$OMP PARALLEL DO PRIVATE(iblk,i,j) do iblk = 1, nblocks do j = 1, ny_block @@ -533,7 +554,8 @@ subroutine get_forcing_atmo integer (kind=int_kind) :: & iblk, & ! block index ilo,ihi,jlo,jhi, & ! beginning and end of physical domain - fyear_old, & ! prior fyear value + modadj, & ! adjustment to make mod a postive number + fyear_old, & ! fyear setting on last timestep nt_Tsfc type (block) :: & @@ -541,12 +563,17 @@ subroutine get_forcing_atmo character(len=*), parameter :: subname = '(get_forcing_atmo)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + + call ice_timer_start(timer_forcing) + fyear_old = fyear - fyear = fyear_init + mod(nyr-1,ycycle) ! current year + modadj = abs((min(0,myear-fyear_init)/ycycle+1)*ycycle) + fyear = fyear_init + mod(myear-fyear_init+modadj,ycycle) if (trim(atm_data_type) /= 'default' .and. & (istep <= 1 .or. fyear /= fyear_old)) then if (my_task == master_task) then - write (nu_diag,*) ' Current forcing data year = ',fyear + write (nu_diag,*) ' Set current forcing data year = ',fyear endif endif @@ -555,23 +582,25 @@ subroutine get_forcing_atmo if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) - ftime = time ! forcing time - time_forc = ftime ! for restarting + !------------------------------------------------------------------- + ! Read and interpolate atmospheric data + !------------------------------------------------------------------- - !------------------------------------------------------------------- - ! Read and interpolate atmospheric data - !------------------------------------------------------------------- + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg fyear = ',fyear + write(nu_diag,*) subname,'fdbg atm_data_type = ',trim(atm_data_type) + endif if (trim(atm_data_type) == 'ncar') then call ncar_data elseif (trim(atm_data_type) == 'LYq') then call LY_data elseif (trim(atm_data_type) == 'JRA55_gx1') then - call JRA55_data(fyear) + call JRA55_data elseif (trim(atm_data_type) == 'JRA55_gx3') then - call JRA55_data(fyear) + call JRA55_data elseif (trim(atm_data_type) == 'JRA55_tx1') then - call JRA55_data(fyear) + call JRA55_data elseif (trim(atm_data_type) == 'hadgem') then call hadgem_data elseif (trim(atm_data_type) == 'monthly') then @@ -586,9 +615,9 @@ subroutine get_forcing_atmo return endif - !------------------------------------------------------------------- - ! Convert forcing data to fields needed by ice model - !------------------------------------------------------------------- + !------------------------------------------------------------------- + ! Convert forcing data to fields needed by ice model + !------------------------------------------------------------------- !$OMP PARALLEL DO PRIVATE(iblk,ilo,ihi,jlo,jhi,this_block) do iblk = 1, nblocks @@ -640,6 +669,8 @@ subroutine get_forcing_atmo field_loc_center, field_type_scalar) call ice_timer_stop(timer_bound) + call ice_timer_stop(timer_forcing) + end subroutine get_forcing_atmo !======================================================================= @@ -655,6 +686,15 @@ subroutine get_forcing_ocn (dt) character(len=*), parameter :: subname = '(get_forcing_ocn)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + + call ice_timer_start(timer_forcing) + + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg fyear = ',fyear + write(nu_diag,*) subname,'fdbg ocn_data_type = ',trim(ocn_data_type) + endif + if (trim(ocn_data_type) == 'clim') then call ocn_data_clim(dt) elseif (trim(ocn_data_type) == 'ncar' .or. & @@ -670,6 +710,8 @@ subroutine get_forcing_ocn (dt) !MHRI: NOT IMPLEMENTED YET endif + call ice_timer_stop(timer_forcing) + end subroutine get_forcing_ocn !======================================================================= @@ -726,6 +768,8 @@ subroutine read_data (flag, recd, yr, ixm, ixx, ixp, & character(len=*), parameter :: subname = '(read_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_readwrite) ! reading/writing nbits = 64 ! double precision data @@ -860,16 +904,14 @@ subroutine read_data_nc (flag, recd, yr, ixm, ixx, ixp, & fieldname ! field name in netCDF file integer (kind=int_kind), intent(in) :: & - field_loc, & ! location of field on staggered grid - field_type ! type of field (scalar, vector, angle) + field_loc, & ! location of field on staggered grid + field_type ! type of field (scalar, vector, angle) real (kind=dbl_kind), dimension(nx_block,ny_block,2,max_blocks), intent(out) :: & field_data ! 2 values needed for interpolation ! local variables - character(len=*), parameter :: subname = '(read_data_nc)' - integer (kind=int_kind) :: & nrec , & ! record number to read n2, n4 , & ! like ixm and ixp, but @@ -877,6 +919,10 @@ subroutine read_data_nc (flag, recd, yr, ixm, ixx, ixp, & arg , & ! value of time argument in field_data fid ! file id for netCDF routines + character(len=*), parameter :: subname = '(read_data_nc)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_readwrite) ! reading/writing if (istep1 > check_step) dbug = .true. !! debugging @@ -1011,6 +1057,10 @@ subroutine read_data_nc_hycom (flag, recd, & integer (kind=int_kind) :: & fid ! file id for netCDF routines + character(len=*), parameter :: subname = '(read_data_nc_hycom)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_readwrite) ! reading/writing if (istep1 > check_step) dbug = .true. !! debugging @@ -1079,6 +1129,8 @@ subroutine read_clim_data (readflag, recd, ixm, ixx, ixp, & character(len=*), parameter :: subname = '(read_clim_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_readwrite) ! reading/writing nbits = 64 ! double precision data @@ -1164,6 +1216,8 @@ subroutine read_clim_data_nc (readflag, recd, ixm, ixx, ixp, & character(len=*), parameter :: subname = '(read_clim_data_nc)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_readwrite) ! reading/writing if (istep1 > check_step) dbug = .true. !! debugging @@ -1222,14 +1276,16 @@ subroutine interp_coeff_monthly (recslot) real (kind=dbl_kind) :: & secday , & ! seconds in day - tt , & ! seconds elapsed in current year - t1, t2 ! seconds elapsed at month midpoint + tt , & ! days elapsed in current year + t1, t2 ! days elapsed at month midpoint real (kind=dbl_kind) :: & daymid(0:13) ! month mid-points character(len=*), parameter :: subname = '(interp_coeff_monthly)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & @@ -1238,21 +1294,27 @@ subroutine interp_coeff_monthly (recslot) daymid(1:13) = 14._dbl_kind ! time frame ends 0 sec into day 15 daymid(0) = 14._dbl_kind - daymo(12) ! Dec 15, 0 sec - ! make time cyclic - tt = mod(ftime/secday,dayyr) + ! compute days since Jan 1, 00h, yday is the day counter for the year + tt = real(yday-1,kind=dbl_kind) + real(msec,kind=dbl_kind)/secday ! Find neighboring times if (recslot==2) then ! first half of month - t2 = daycal(month) + daymid(month) ! midpoint, current month - if (month == 1) then + t2 = daycal(mmonth) + daymid(mmonth) ! midpoint, current month + if (mmonth == 1) then t1 = daymid(0) ! Dec 15 (0 sec) else - t1 = daycal(month-1) + daymid(month-1) ! midpoint, previous month + t1 = daycal(mmonth-1) + daymid(mmonth-1) ! midpoint, previous month endif else ! second half of month - t1 = daycal(month) + daymid(month) ! midpoint, current month - t2 = daycal(month+1) + daymid(month+1)! day 15 of next month (0 sec) + t1 = daycal(mmonth) + daymid(mmonth) ! midpoint, current month + t2 = daycal(mmonth+1) + daymid(mmonth+1)! day 15 of next month (0 sec) + endif + + if (tt < t1 .or. tt > t2) then + write(nu_diag,*) subname,' ERROR in tt',tt,t1,t2 + call abort_ice (error_message=subname//' ERROR in tt', & + file=__FILE__, line=__LINE__) endif ! Compute coefficients @@ -1282,8 +1344,7 @@ subroutine interp_coeff (recnum, recslot, secint, dataloc) ! local variables real (kind=dbl_kind) :: & - secday, & ! seconds in a day - secyr ! seconds in a year + secday ! seconds in a day real (kind=dbl_kind) :: & tt , & ! seconds elapsed in current year @@ -1292,13 +1353,15 @@ subroutine interp_coeff (recnum, recslot, secint, dataloc) character(len=*), parameter :: subname = '(interp_coeff)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) - secyr = dayyr * secday ! seconds in a year - tt = mod(ftime,secyr) + ! compute seconds since Jan 1, 00h, yday is the day counter for the year + tt = real(yday-1,kind=dbl_kind)*secday + real(msec,kind=dbl_kind) ! Find neighboring times rcnum = real(recnum,kind=dbl_kind) @@ -1322,6 +1385,12 @@ subroutine interp_coeff (recnum, recslot, secint, dataloc) c1intp = abs((t2 - tt) / (t2 - t1)) c2intp = c1 - c1intp + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg yday,sec = ',yday,msec + write(nu_diag,*) subname,'fdbg tt = ',tt + write(nu_diag,*) subname,'fdbg c12intp = ',c1intp,c2intp + endif + end subroutine interp_coeff !======================================================================= @@ -1335,6 +1404,9 @@ subroutine interp_coeff2 (tt, t1, t2) real (kind=dbl_kind), intent(in) :: & tt , & ! current decimal daynumber t1, t2 ! first+last decimal daynumber + character(len=*), parameter :: subname = '(interp_coeff2)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' ! Compute coefficients c1intp = abs((t2 - tt) / (t2 - t1)) @@ -1364,6 +1436,8 @@ subroutine interpolate_data (field_data, field) character(len=*), parameter :: subname = '(interpolate data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + !$OMP PARALLEL DO PRIVATE(iblk,i,j) do iblk = 1, nblocks do j = 1, ny_block @@ -1395,6 +1469,8 @@ subroutine file_year (data_file, yr) character(len=*), parameter :: subname = '(file_year)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + if (trim(atm_data_type) == 'hadgem') then ! netcdf i = index(data_file,'.nc') - 5 tmpname = data_file @@ -1481,6 +1557,8 @@ subroutine prepare_forcing (nx_block, ny_block, & character(len=*), parameter :: subname = '(prepare_forcing)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Tffresh_out=Tffresh, puny_out=puny) call icepack_query_parameters(secday_out=secday) call icepack_query_parameters(calc_strair_out=calc_strair) @@ -1579,7 +1657,7 @@ subroutine prepare_forcing (nx_block, ny_block, & ! convert precipitation units to kg/m^2 s if (trim(precip_units) == 'mm_per_month') then - precip_factor = c12/(secday*days_per_year) + precip_factor = c12/(secday*real(days_per_year,kind=dbl_kind)) elseif (trim(precip_units) == 'mm_per_day') then precip_factor = c1/secday elseif (trim(precip_units) == 'mm_per_sec' .or. & @@ -1699,6 +1777,8 @@ subroutine longwave_parkinson_washington(Tair, cldf, flw) character(len=*), parameter :: subname = '(longwave_parkinson_washington)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Tffresh_out=Tffresh, & stefan_boltzmann_out=stefan_boltzmann) call icepack_warnings_flush(nu_diag) @@ -1749,6 +1829,8 @@ subroutine longwave_rosati_miyakoda(cldf, Tsfc, & character(len=*), parameter :: subname = '(longwave_rosati_miyakoda)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Tffresh_out=Tffresh, & stefan_boltzmann_out=stefan_boltzmann, & emissivity_out=emissivity) @@ -1786,6 +1868,8 @@ subroutine ncar_files (yr) character(len=*), parameter :: subname = '(ncar_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + fsw_file = & trim(atm_data_dir)//'/MONTHLY/swdn.1996.dat' call file_year(fsw_file,yr) @@ -1857,6 +1941,8 @@ subroutine ncar_data character(len=*), parameter :: subname = '(ncar_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & @@ -1870,12 +1956,12 @@ subroutine ncar_data !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -1892,29 +1978,29 @@ subroutine ncar_data ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. if (trim(atm_data_format) == 'bin') then - call read_data (readm, 0, fyear, ixm, month, ixp, & + call read_data (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, fsw_file, fsw_data, & field_loc_center, field_type_scalar) - call read_data (readm, 0, fyear, ixm, month, ixp, & + call read_data (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, flw_file, cldf_data, & field_loc_center, field_type_scalar) - call read_data (readm, 0, fyear, ixm, month, ixp, & + call read_data (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, rain_file, fsnow_data, & field_loc_center, field_type_scalar) else call abort_ice (error_message=subname//'nonbinary atm_data_format unavailable', & file=__FILE__, line=__LINE__) ! The routine exists, for example: -! call read_data_nc (readm, 0, fyear, ixm, month, ixp, & +! call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & ! maxrec, fsw_file, 'fsw', fsw_data, & ! field_loc_center, field_type_scalar) -! call read_data_nc (readm, 0, fyear, ixm, month, ixp, & +! call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & ! maxrec, flw_file, 'cldf',cldf_data, & ! field_loc_center, field_type_scalar) -! call read_data_nc (readm, 0, fyear, ixm, month, ixp, & +! call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & ! maxrec, rain_file,'prec',fsnow_data, & ! field_loc_center, field_type_scalar) endif @@ -1937,7 +2023,7 @@ subroutine ncar_data maxrec = 1460 ! 365*4 ! current record number - recnum = 4*int(yday) - 3 + int(real(sec,kind=dbl_kind)/sec6hr) + recnum = 4*int(yday) - 3 + int(real(msec,kind=dbl_kind)/sec6hr) ! Compute record numbers for surrounding data @@ -2009,6 +2095,8 @@ subroutine LY_files (yr) character(len=*), parameter :: subname = '(LY_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + flw_file = & trim(atm_data_dir)//'/MONTHLY/cldf.omip.dat' @@ -2044,6 +2132,9 @@ subroutine LY_files (yr) endif ! master_task end subroutine LY_files + +!======================================================================= + subroutine JRA55_gx1_files(yr) ! integer (kind=int_kind), intent(in) :: & @@ -2051,6 +2142,8 @@ subroutine JRA55_gx1_files(yr) character(len=*), parameter :: subname = '(JRA55_gx1_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + uwind_file = & trim(atm_data_dir)//'/8XDAILY/JRA55_03hr_forcing_2005.nc' call file_year(uwind_file,yr) @@ -2060,6 +2153,9 @@ subroutine JRA55_gx1_files(yr) write (nu_diag,*) trim(uwind_file) endif end subroutine JRA55_gx1_files + +!======================================================================= + subroutine JRA55_tx1_files(yr) ! integer (kind=int_kind), intent(in) :: & @@ -2067,6 +2163,8 @@ subroutine JRA55_tx1_files(yr) character(len=*), parameter :: subname = '(JRA55_tx1_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + uwind_file = & trim(atm_data_dir)//'/8XDAILY/JRA55_03hr_forcing_tx1_2005.nc' call file_year(uwind_file,yr) @@ -2076,6 +2174,9 @@ subroutine JRA55_tx1_files(yr) write (nu_diag,*) trim(uwind_file) endif end subroutine JRA55_tx1_files + +!======================================================================= + subroutine JRA55_gx3_files(yr) ! integer (kind=int_kind), intent(in) :: & @@ -2083,6 +2184,8 @@ subroutine JRA55_gx3_files(yr) character(len=*), parameter :: subname = '(JRA55_gx3_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + uwind_file = & trim(atm_data_dir)//'/8XDAILY/JRA55_gx3_03hr_forcing_2005.nc' call file_year(uwind_file,yr) @@ -2092,6 +2195,7 @@ subroutine JRA55_gx3_files(yr) write (nu_diag,*) trim(uwind_file) endif end subroutine JRA55_gx3_files + !======================================================================= ! ! read Large and Yeager atmospheric data @@ -2131,6 +2235,8 @@ subroutine LY_data character(len=*), parameter :: subname = '(LY_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Tffresh_out=Tffresh) call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) @@ -2145,12 +2251,12 @@ subroutine LY_data !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -2167,11 +2273,11 @@ subroutine LY_data ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & flw_file, cldf_data, field_loc_center, field_type_scalar) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & rain_file, fsnow_data, field_loc_center, field_type_scalar) call interpolate_data (cldf_data, cldf) @@ -2190,7 +2296,7 @@ subroutine LY_data maxrec = 1460 ! 365*4 ! current record number - recnum = 4*int(yday) - 3 + int(real(sec,kind=dbl_kind)/sec6hr) + recnum = 4*int(yday) - 3 + int(real(msec,kind=dbl_kind)/sec6hr) ! Compute record numbers for surrounding data (2 on each side) @@ -2316,7 +2422,256 @@ end subroutine LY_data !======================================================================= - subroutine JRA55_data (yr) + subroutine JRA55_data + + use ice_blocks, only: block, get_block + use ice_global_reductions, only: global_minval, global_maxval + use ice_domain, only: nblocks, distrb_info, blocks_ice + use ice_flux, only: fsnow, Tair, uatm, vatm, Qa, fsw, flw + use ice_grid, only: hm, tlon, tlat, tmask, umask + use ice_state, only: aice + use ice_calendar, only: days_per_year + + integer (kind=int_kind) :: & + ncid , & ! netcdf file id + i, j, n1 , & + lfyear , & ! local year value + recnum , & ! record number + maxrec , & ! maximum record number + iblk ! block index + + integer (kind=int_kind), save :: & + frec_info(2,2) = -99 ! remember prior values to reduce reading + ! first dim is yr, recnum + ! second dim is data1 data2 + + real (kind=dbl_kind) :: & + sec3hr , & ! number of seconds in 3 hours + secday , & ! number of seconds in day + eps, tt , & ! for interpolation coefficients + Tffresh , & + vmin, vmax + + character(len=64) :: fieldname !netcdf field name + character (char_len_long) :: uwind_file_old + character(len=*), parameter :: subname = '(JRA55_data)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + + call icepack_query_parameters(Tffresh_out=Tffresh) + call icepack_query_parameters(secday_out=secday) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + sec3hr = secday/c8 ! seconds in 3 hours + maxrec = days_per_year * 8 + + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg dpy, maxrec = ',days_per_year,maxrec + endif + + !------------------------------------------------------------------- + ! 3-hourly data + ! states are instantaneous, 1st record is 00z Jan 1 + ! fluxes are 3 hour averages, 1st record is 00z-03z Jan 1 + ! interpolate states, do not interpolate fluxes + !------------------------------------------------------------------- + ! File is NETCDF with winds in NORTH and EAST direction + ! file variable names are: + ! glbrad (shortwave W/m^2), 3 hr average + ! dlwsfc (longwave W/m^2), 3 hr average + ! wndewd (eastward wind m/s), instantaneous + ! wndnwd (northward wind m/s), instantaneous + ! airtmp (air temperature K), instantaneous + ! spchmd (specific humidity kg/kg), instantaneous + ! ttlpcp (precipitation kg/m s-1), 3 hr average + !------------------------------------------------------------------- + + uwind_file_old = uwind_file + if (uwind_file /= uwind_file_old .and. my_task == master_task) then + write(nu_diag,*) subname,' reading forcing file = ',trim(uwind_file) + endif + + call ice_open_nc(uwind_file,ncid) + + do n1 = 1, 2 + + lfyear = fyear + call file_year(uwind_file,lfyear) + if (n1 == 1) then + recnum = 8*int(yday) - 7 + int(real(msec,kind=dbl_kind)/sec3hr) + if (my_task == master_task .and. (recnum <= 2 .or. recnum >= maxrec-1)) then + write(nu_diag,*) subname,' reading forcing file 1st ts = ',trim(uwind_file) + endif + elseif (n1 == 2) then + recnum = 8*int(yday) - 7 + int(real(msec,kind=dbl_kind)/sec3hr) + 1 + if (recnum > maxrec) then + lfyear = fyear + 1 ! next year + if (lfyear > fyear_final) lfyear = fyear_init + recnum = 1 + call file_year(uwind_file,lfyear) + if (my_task == master_task) then + write(nu_diag,*) subname,' reading forcing file 2nd ts = ',trim(uwind_file) + endif + call ice_close_nc(ncid) + call ice_open_nc(uwind_file,ncid) + endif + endif + + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg read recnum = ',recnum,n1 + endif + + ! to reduce reading, check whether it's the same data as last read + + if (lfyear /= frec_info(1,n1) .or. recnum /= frec_info(2,n1)) then + + ! check whether we can copy values from 2 to 1, should be faster than reading + ! can only do this from 2 to 1 or 1 to 2 without setting up a temporary + ! it's more likely that the values from data2 when time advances are needed in data1 + ! compare n1=1 year/record with data from last timestep at n1=2 + + if (n1 == 1 .and. lfyear == frec_info(1,2) .and. recnum == frec_info(2,2)) then + Tair_data(:,:,1,:) = Tair_data(:,:,2,:) + uatm_data(:,:,1,:) = uatm_data(:,:,2,:) + vatm_data(:,:,1,:) = vatm_data(:,:,2,:) + Qa_data(:,:,1,:) = Qa_data(:,:,2,:) + fsw_data(:,:,1,:) = fsw_data(:,:,2,:) + flw_data(:,:,1,:) = flw_data(:,:,2,:) + fsnow_data(:,:,1,:) = fsnow_data(:,:,2,:) + else + + fieldname = 'airtmp' + call ice_read_nc(ncid,recnum,fieldname,Tair_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + + fieldname = 'wndewd' + call ice_read_nc(ncid,recnum,fieldname,uatm_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + + fieldname = 'wndnwd' + call ice_read_nc(ncid,recnum,fieldname,vatm_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + + fieldname = 'spchmd' + call ice_read_nc(ncid,recnum,fieldname,Qa_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + + fieldname = 'glbrad' + call ice_read_nc(ncid,recnum,fieldname,fsw_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + + fieldname = 'dlwsfc' + call ice_read_nc(ncid,recnum,fieldname,flw_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + + fieldname = 'ttlpcp' + call ice_read_nc(ncid,recnum,fieldname,fsnow_data(:,:,n1,:),forcing_debug, & + field_loc=field_loc_center, & + field_type=field_type_scalar) + endif ! copy data from n1=2 from last timestep to n1=1 + endif ! input data is same as last timestep + + frec_info(1,n1) = lfyear + frec_info(2,n1) = recnum + + enddo ! n1 + + call ice_close_nc(ncid) + + ! reset uwind_file to original year + call file_year(uwind_file,fyear) + + ! Compute interpolation coefficients + eps = 1.0e-6 + tt = real(mod(msec,nint(sec3hr)),kind=dbl_kind) + c2intp = tt / sec3hr + if (c2intp < c0 .and. c2intp > c0-eps) c2intp = c0 + if (c2intp > c1 .and. c2intp < c1+eps) c2intp = c1 + c1intp = 1.0_dbl_kind - c2intp + if (c2intp < c0 .or. c2intp > c1) then + write(nu_diag,*) subname,' ERROR: c2intp = ',c2intp + call abort_ice (error_message=subname//' ERROR: c2intp out of range', & + file=__FILE__, line=__LINE__) + endif + if (forcing_debug .and. my_task == master_task) then + write(nu_diag,*) subname,'fdbg c12intp = ',c1intp,c2intp + endif + + ! Interpolate + call interpolate_data (Tair_data, Tair) + call interpolate_data (uatm_data, uatm) + call interpolate_data (vatm_data, vatm) + call interpolate_data (Qa_data, Qa) + ! use 3 hr average for heat flux and precip fields, no interpolation +! call interpolate_data (fsw_data, fsw) +! call interpolate_data (flw_data, flw) +! call interpolate_data (fsnow_data, fsnow) + fsw(:,:,:) = fsw_data(:,:,1,:) + flw(:,:,:) = flw_data(:,:,1,:) + fsnow(:,:,:) = fsnow_data(:,:,1,:) + + !$OMP PARALLEL DO PRIVATE(iblk,i,j) + do iblk = 1, nblocks + ! limit summer Tair values where ice is present + do j = 1, ny_block + do i = 1, nx_block + if (aice(i,j,iblk) > p1) Tair(i,j,iblk) = min(Tair(i,j,iblk), Tffresh+p1) + enddo + enddo + + do j = 1, ny_block + do i = 1, nx_block + Qa (i,j,iblk) = Qa (i,j,iblk) * hm(i,j,iblk) + Tair(i,j,iblk) = Tair(i,j,iblk) * hm(i,j,iblk) + uatm(i,j,iblk) = uatm(i,j,iblk) * hm(i,j,iblk) + vatm(i,j,iblk) = vatm(i,j,iblk) * hm(i,j,iblk) + fsw (i,j,iblk) = fsw (i,j,iblk) * hm(i,j,iblk) + flw (i,j,iblk) = flw (i,j,iblk) * hm(i,j,iblk) + fsnow(i,j,iblk) = fsnow (i,j,iblk) * hm(i,j,iblk) + enddo + enddo + + enddo ! iblk + !$OMP END PARALLEL DO + + if (dbug .or. forcing_debug) then + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg JRA55_bulk_data' + vmin = global_minval(fsw,distrb_info,tmask) + vmax = global_maxval(fsw,distrb_info,tmask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg fsw',vmin,vmax + vmin = global_minval(flw,distrb_info,tmask) + vmax = global_maxval(flw,distrb_info,tmask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg flw',vmin,vmax + vmin =global_minval(fsnow,distrb_info,tmask) + vmax =global_maxval(fsnow,distrb_info,tmask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg fsnow',vmin,vmax + vmin = global_minval(Tair,distrb_info,tmask) + vmax = global_maxval(Tair,distrb_info,tmask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg Tair',vmin,vmax + vmin = global_minval(uatm,distrb_info,umask) + vmax = global_maxval(uatm,distrb_info,umask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg uatm',vmin,vmax + vmin = global_minval(vatm,distrb_info,umask) + vmax = global_maxval(vatm,distrb_info,umask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg vatm',vmin,vmax + vmin = global_minval(Qa,distrb_info,tmask) + vmax = global_maxval(Qa,distrb_info,tmask) + if (my_task.eq.master_task) write (nu_diag,*) subname,'fdbg Qa',vmin,vmax + endif ! dbug + + end subroutine JRA55_data + +!======================================================================= + + subroutine Jra55_data_old (yr) use ice_blocks, only: block, get_block use ice_global_reductions, only: global_minval, global_maxval @@ -2350,7 +2705,9 @@ subroutine JRA55_data (yr) character (char_len_long) :: uwind_file_old character(len=64) :: fieldname !netcdf field name - character(len=*), parameter :: subname = '(JRA55_data)' + character(len=*), parameter :: subname = '(Jra55_data_old)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' call icepack_query_parameters(Tffresh_out=Tffresh) call icepack_query_parameters(secday_out=secday) @@ -2397,14 +2754,14 @@ subroutine JRA55_data (yr) do n1 = 1,2 if (n1 == 1) then - recnum = 8*int(yday) - 7 + int(real(sec,kind=dbl_kind)/sec3hr) + recnum = 8*int(yday) - 7 + int(real(msec,kind=dbl_kind)/sec3hr) if (my_task == master_task .and. (recnum <= 2 .or. recnum >= maxrec-1)) then write(nu_diag,*) subname,' reading forcing file 1st ts = ',trim(uwind_file) endif elseif (n1 == 2) then - recnum = 8*int(yday) - 7 + int(real(sec,kind=dbl_kind)/sec3hr) + 1 + recnum = 8*int(yday) - 7 + int(real(msec,kind=dbl_kind)/sec3hr) + 1 if (recnum > maxrec) then - yrp = fyear_init + mod(nyr,ycycle) ! next year + yrp = fyear_init + mod(myear,ycycle) ! next year recnum = 1 call file_year(uwind_file,yrp) if (my_task == master_task) then @@ -2466,7 +2823,7 @@ subroutine JRA55_data (yr) ! Compute interpolation coefficients eps = 1.0e-6 - tt = real(mod(sec,nint(sec3hr)),kind=dbl_kind) + tt = real(mod(msec,nint(sec3hr)),kind=dbl_kind) c2intp = tt / sec3hr if (c2intp < c0 .and. c2intp > c0-eps) c2intp = c0 if (c2intp > c1 .and. c2intp < c1+eps) c2intp = c1 @@ -2551,7 +2908,7 @@ subroutine JRA55_data (yr) endif ! dbug - end subroutine JRA55_data + end subroutine Jra55_data_old !======================================================================= ! @@ -2596,6 +2953,8 @@ subroutine compute_shortwave(nx_block, ny_block, & character(len=*), parameter :: subname = '(compute_shortwave)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(secday_out=secday, pi_out=pi) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & @@ -2604,7 +2963,7 @@ subroutine compute_shortwave(nx_block, ny_block, & do j=jlo,jhi do i=ilo,ihi deg2rad = pi/c180 -! solar_time = mod(real(sec,kind=dbl_kind),secday)/c3600 & +! solar_time = mod(real(msec,kind=dbl_kind),secday)/c3600 & ! + c12*sin(p5*TLON(i,j)) ! Convert longitude to range of -180 to 180 for LST calculation @@ -2613,7 +2972,7 @@ subroutine compute_shortwave(nx_block, ny_block, & if (lontmp .gt. c180) lontmp = lontmp - c360 if (lontmp .lt. -c180) lontmp = lontmp + c360 - solar_time = mod(real(sec,kind=dbl_kind),secday)/c3600 & + solar_time = mod(real(msec,kind=dbl_kind),secday)/c3600 & + lontmp/c15 if (solar_time .ge. 24._dbl_kind) solar_time = solar_time - 24._dbl_kind hour_angle = (c12 - solar_time)*pi/c12 @@ -2658,6 +3017,8 @@ subroutine Qa_fixLY(nx_block, ny_block, Tair, Qa) character(len=*), parameter :: subname = '(Qa_fixLY)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Tffresh_out=Tffresh, puny_out=puny) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & @@ -2700,6 +3061,8 @@ subroutine hadgem_files (yr) character(len=*), parameter :: subname = '(hadgem_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(calc_strair_out=calc_strair, & calc_Tsfc_out=calc_Tsfc) call icepack_warnings_flush(nu_diag) @@ -2898,6 +3261,8 @@ subroutine hadgem_data character(len=*), parameter :: subname = '(hadgem_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Lsub_out=Lsub) call icepack_query_parameters(calc_strair_out=calc_strair, & calc_Tsfc_out=calc_Tsfc) @@ -2913,12 +3278,12 @@ subroutine hadgem_data !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -2935,18 +3300,18 @@ subroutine hadgem_data ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. ! ----------------------------------------------------------- ! Rainfall and snowfall ! ----------------------------------------------------------- fieldname='rainfall' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, rain_file, fieldname, frain_data, & field_loc_center, field_type_scalar) fieldname='snowfall' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, snow_file, fieldname, fsnow_data, & field_loc_center, field_type_scalar) @@ -2961,11 +3326,11 @@ subroutine hadgem_data ! -------------------------------------------------------- fieldname='u_10' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, uwind_file, fieldname, uatm_data, & field_loc_center, field_type_vector) fieldname='v_10' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, vwind_file, fieldname, vatm_data, & field_loc_center, field_type_vector) @@ -2980,11 +3345,11 @@ subroutine hadgem_data ! -------------------------------------------------------- fieldname='taux' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, strax_file, fieldname, strax_data, & field_loc_center, field_type_vector) fieldname='tauy' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, stray_file, fieldname, stray_data, & field_loc_center, field_type_vector) @@ -2999,7 +3364,7 @@ subroutine hadgem_data ! -------------------------------------------------- fieldname='wind_10' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, wind_file, fieldname, wind_data, & field_loc_center, field_type_scalar) @@ -3022,23 +3387,23 @@ subroutine hadgem_data if (calc_Tsfc .or. oceanmixed_ice .or. calc_strair) then fieldname='SW_incoming' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, fsw_file, fieldname, fsw_data, & field_loc_center, field_type_scalar) fieldname='LW_incoming' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, flw_file, fieldname, flw_data, & field_loc_center, field_type_scalar) fieldname='t_10' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, tair_file, fieldname, Tair_data, & field_loc_center, field_type_scalar) fieldname='rho_10' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, rhoa_file, fieldname, rhoa_data, & field_loc_center, field_type_scalar) fieldname='q_10' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, humid_file, fieldname, Qa_data, & field_loc_center, field_type_scalar) @@ -3059,7 +3424,7 @@ subroutine hadgem_data ! ------------------------------------------------------ fieldname='sublim' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, sublim_file, fieldname, sublim_data, & field_loc_center, field_type_scalar) @@ -3068,12 +3433,12 @@ subroutine hadgem_data do n = 1, ncat write(fieldname, '(a,i1)') 'topmeltn',n - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, topmelt_file(n), fieldname, topmelt_data(:,:,:,:,n), & field_loc_center, field_type_scalar) write(fieldname, '(a,i1)') 'botmeltn',n - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, botmelt_file(n), fieldname, botmelt_data(:,:,:,:,n), & field_loc_center, field_type_scalar) @@ -3127,6 +3492,8 @@ subroutine monthly_files (yr) character(len=*), parameter :: subname = '(monthly_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + flw_file = & trim(atm_data_dir)//'/MONTHLY/cldf.omip.dat' @@ -3198,6 +3565,8 @@ subroutine monthly_data character(len=*), parameter :: subname = '(monthly_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + !------------------------------------------------------------------- ! monthly data ! @@ -3206,12 +3575,12 @@ subroutine monthly_data !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -3228,27 +3597,27 @@ subroutine monthly_data ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & flw_file, cldf_data, & field_loc_center, field_type_scalar) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & rain_file, fsnow_data, & field_loc_center, field_type_scalar) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & tair_file, Tair_data, & field_loc_center, field_type_scalar) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & humid_file, Qa_data, & field_loc_center, field_type_scalar) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & wind_file, wind_data, & field_loc_center, field_type_scalar) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & strax_file, strax_data, & field_loc_center, field_type_vector) - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & stray_file, stray_data, & field_loc_center, field_type_vector) @@ -3377,6 +3746,8 @@ subroutine oned_data character(len=*), parameter :: subname = '(oned_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + diag = .false. ! write diagnostic information if (trim(atm_data_format) == 'nc') then ! read nc file @@ -3452,6 +3823,8 @@ subroutine oned_files character(len=*), parameter :: subname = '(oned_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + fsw_file = & trim(atm_data_dir)//'/hourlysolar_brw1989_5yr.nc' @@ -3517,6 +3890,8 @@ subroutine ocn_data_clim (dt) character(len=*), parameter :: subname = '(ocn_data_clim)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + if (my_task == master_task .and. istep == 1) then if (trim(ocn_data_type)=='clim') then write (nu_diag,*) ' ' @@ -3540,12 +3915,12 @@ subroutine ocn_data_clim (dt) if (trim(ocn_data_type)=='clim') then midmonth = 15 ! data is given on 15th of every month -!!! midmonth = fix(p5 * real(daymo(month))) ! exact middle +!!! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -3561,14 +3936,14 @@ subroutine ocn_data_clim (dt) call interp_coeff_monthly (recslot) readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. !------------------------------------------------------------------- ! Read two monthly SSS values and interpolate. ! Note: SSS is restored instantaneously to data. !------------------------------------------------------------------- - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & sss_file, sss_data, & field_loc_center, field_type_scalar) call interpolate_data (sss_data, sss) @@ -3592,7 +3967,7 @@ subroutine ocn_data_clim (dt) !------------------------------------------------------------------- if (trim(ocn_data_type)=='clim') then - call read_clim_data (readm, 0, ixm, month, ixp, & + call read_clim_data (readm, 0, ixm, mmonth, ixp, & sst_file, sst_data, & field_loc_center, field_type_scalar) call interpolate_data (sst_data, sstdat) @@ -3673,6 +4048,8 @@ subroutine ocn_data_ncar_init character(len=*), parameter :: subname = '(ocn_data_ncar_init)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + if (my_task == master_task) then write (nu_diag,*) 'WARNING: evp_prep calculates surface tilt' @@ -3830,6 +4207,8 @@ subroutine ocn_data_ncar_init_3D character(len=*), parameter :: subname = '(ocn_data_ncar_init_3D)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + if (my_task == master_task) then write (nu_diag,*) 'WARNING: evp_prep calculates surface tilt' @@ -3967,6 +4346,8 @@ subroutine ocn_data_ncar(dt) character(len=*), parameter :: subname = '(ocn_data_ncar)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + !------------------------------------------------------------------- ! monthly data ! @@ -3975,12 +4356,12 @@ subroutine ocn_data_ncar(dt) !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month),kind=dbl_kind)) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth),kind=dbl_kind)) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -4001,9 +4382,9 @@ subroutine ocn_data_ncar(dt) ! use sst_data arrays as temporary work space until n=1 if (ixm /= -99) then ! first half of month sst_data(:,:,1,iblk) = ocn_frc_m(:,:,iblk,n,ixm) - sst_data(:,:,2,iblk) = ocn_frc_m(:,:,iblk,n,month) + sst_data(:,:,2,iblk) = ocn_frc_m(:,:,iblk,n,mmonth) else ! second half of month - sst_data(:,:,1,iblk) = ocn_frc_m(:,:,iblk,n,month) + sst_data(:,:,1,iblk) = ocn_frc_m(:,:,iblk,n,mmonth) sst_data(:,:,2,iblk) = ocn_frc_m(:,:,iblk,n,ixp) endif enddo @@ -4125,6 +4506,8 @@ subroutine ocn_data_oned character(len=*), parameter :: subname = '(ocn_data_oned)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + sss (:,:,:) = 34.0_dbl_kind ! sea surface salinity (ppt) call ocn_freezing_temperature @@ -4180,6 +4563,8 @@ subroutine ocn_data_hadgem(dt) character(len=*), parameter :: subname = '(ocn_data_hadgem)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + !------------------------------------------------------------------- ! monthly data ! @@ -4188,12 +4573,12 @@ subroutine ocn_data_hadgem(dt) !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -4210,7 +4595,7 @@ subroutine ocn_data_hadgem(dt) ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. if (my_task == master_task .and. istep == 1) then write (nu_diag,*) ' ' @@ -4231,7 +4616,7 @@ subroutine ocn_data_hadgem(dt) ! ----------------------------------------------------------- sst_file = trim(ocn_data_dir)//'/MONTHLY/sst.1997.nc' fieldname='sst' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, sst_file, fieldname, sst_data, & field_loc_center, field_type_scalar) @@ -4265,7 +4650,7 @@ subroutine ocn_data_hadgem(dt) filename = trim(ocn_data_dir)//'/MONTHLY/uocn.1997.nc' fieldname='uocn' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, filename, fieldname, uocn_data, & field_loc_center, field_type_vector) @@ -4274,7 +4659,7 @@ subroutine ocn_data_hadgem(dt) filename = trim(ocn_data_dir)//'/MONTHLY/vocn.1997.nc' fieldname='vocn' - call read_data_nc (readm, 0, fyear, ixm, month, ixp, & + call read_data_nc (readm, 0, fyear, ixm, mmonth, ixp, & maxrec, filename, fieldname, vocn_data, & field_loc_center, field_type_vector) @@ -4334,6 +4719,10 @@ subroutine ocn_data_hycom_init character (char_len) :: & fieldname ! field name in netcdf file + character(len=*), parameter :: subname = '(ocn_data_hycom_init)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + if (trim(ocn_data_type) == 'hycom') then sss_file = trim(ocn_data_dir)//'ice.restart.surf.nc' @@ -4387,6 +4776,9 @@ subroutine hycom_atm_files fid ! File id character (char_len) :: & varname ! variable name in netcdf file + character(len=*), parameter :: subname = '(hycom_atm_files)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' fsw_file = trim(atm_data_dir)//'/forcing.shwflx.nc' flw_file = trim(atm_data_dir)//'/forcing.radflx.nc' @@ -4430,7 +4822,6 @@ subroutine hycom_atm_data use ice_flux, only: fsw, fsnow, Tair, uatm, vatm, Qa, flw use ice_domain, only: nblocks - use ice_calendar, only: year_init integer (kind=int_kind) :: & recnum ! record number @@ -4450,11 +4841,13 @@ subroutine hycom_atm_data character(len=*), parameter :: subname = '(hycom_atm_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(Tffresh_out=Tffresh) call icepack_query_parameters(secday_out=secday) ! current time in HYCOM jday units - hcdate = hc_jday(nyr+year_init-1,0,0)+ yday+sec/secday + hcdate = hc_jday(myear,0,0)+ yday+msec/secday ! Init recnum try recnum=min(max(oldrecnum,1),Njday_atm-1) @@ -4477,7 +4870,7 @@ subroutine hycom_atm_data write (nu_diag,*) & 'ERROR: CICE: Atm forcing not available at hcdate =',hcdate write (nu_diag,*) & - 'ERROR: CICE: nyr, year_init, yday ,sec = ',nyr, year_init, yday, sec + 'ERROR: CICE: myear, yday ,msec = ',myear, yday, msec call abort_ice ('ERROR: CICE stopped') endif @@ -4605,8 +4998,6 @@ subroutine read_data_nc_point (flag, recd, yr, ixm, ixx, ixp, & real (kind=dbl_kind), dimension(2), intent(inout) :: & field_data ! 2 values needed for interpolation - character(len=*), parameter :: subname = '(read_data_nc_point)' - integer (kind=int_kind) :: & nrec , & ! record number to read n2, n4 , & ! like ixm and ixp, but @@ -4614,6 +5005,10 @@ subroutine read_data_nc_point (flag, recd, yr, ixm, ixx, ixp, & arg , & ! value of time argument in field_data fid ! file id for netCDF routines + character(len=*), parameter :: subname = '(read_data_nc_point)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_readwrite) ! reading/writing field_data = c0 ! to satisfy intent(out) attribute @@ -4726,6 +5121,8 @@ subroutine ISPOL_files character(len=*), parameter :: subname = '(ISPOL_files)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + fsw_file = & trim(atm_data_dir)//'/fsw_sfc_4Xdaily.nc' @@ -4817,6 +5214,8 @@ subroutine ISPOL_data character(len=*), parameter :: subname = '(ISPOL_data)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & @@ -4914,7 +5313,7 @@ subroutine ISPOL_data maxrec = 1460 ! 366*4 ! current record number - recnum4X = 4*int(yday) - 3 + int(real(sec,kind=dbl_kind)/sec1hr) + recnum4X = 4*int(yday) - 3 + int(real(msec,kind=dbl_kind)/sec1hr) ! Compute record numbers for surrounding data (2 on each side) ixm = mod(recnum4X+maxrec-2,maxrec) + 1 @@ -5015,6 +5414,8 @@ subroutine ocn_data_ispol_init character(len=*), parameter :: subname = '(ocn_data_ispol_init)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + if (my_task == master_task) then if (restore_ocn) write (nu_diag,*) & @@ -5074,6 +5475,7 @@ subroutine box2001_data use ice_domain, only: nblocks use ice_domain_size, only: max_blocks + use ice_calendar, only: timesecs use ice_blocks, only: nx_block, ny_block, nghost use ice_flux, only: uocn, vocn, uatm, vatm, wind, rhoa, strax, stray use ice_grid, only: uvm, to_ugrid @@ -5089,6 +5491,11 @@ subroutine box2001_data real (kind=dbl_kind) :: & secday, pi , puny, period, pi2, tau + + character(len=*), parameter :: subname = '(box2001_data)' + + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call icepack_query_parameters(pi_out=pi, pi2_out=pi2, puny_out=puny) call icepack_query_parameters(secday_out=secday) @@ -5111,12 +5518,12 @@ subroutine box2001_data vocn(i,j,iblk) = vocn(i,j,iblk) * uvm(i,j,iblk) ! wind components - uatm(i,j,iblk) = c5 + (sin(pi2*time/period)-c3) & + uatm(i,j,iblk) = c5 + (sin(pi2*timesecs/period)-c3) & * sin(pi2*real(i-nghost, kind=dbl_kind) & /real(nx_global,kind=dbl_kind)) & * sin(pi *real(j-nghost, kind=dbl_kind) & /real(ny_global,kind=dbl_kind)) - vatm(i,j,iblk) = c5 + (sin(pi2*time/period)-c3) & + vatm(i,j,iblk) = c5 + (sin(pi2*timesecs/period)-c3) & * sin(pi *real(i-nghost, kind=dbl_kind) & /real(nx_global,kind=dbl_kind)) & * sin(pi2*real(j-nghost, kind=dbl_kind) & @@ -5180,6 +5587,8 @@ subroutine get_wave_spec logical (kind=log_kind) :: wave_spec character(len=*), parameter :: subname = '(get_wave_spec)' + if (forcing_debug .and. my_task == master_task) write(nu_diag,*) subname,'fdbg start' + call ice_timer_start(timer_fsd) call icepack_query_parameters(wave_spec_out=wave_spec, & diff --git a/cicecore/cicedynB/general/ice_forcing_bgc.F90 b/cicecore/cicedynB/general/ice_forcing_bgc.F90 index e5ef851fa..d9408c304 100644 --- a/cicecore/cicedynB/general/ice_forcing_bgc.F90 +++ b/cicecore/cicedynB/general/ice_forcing_bgc.F90 @@ -14,7 +14,7 @@ module ice_forcing_bgc use ice_blocks, only: nx_block, ny_block use ice_domain_size, only: max_blocks use ice_communicate, only: my_task, master_task - use ice_calendar, only: dt, istep, sec, mday, month + use ice_calendar, only: dt, istep, msec, mday, mmonth use ice_fileunits, only: nu_diag use ice_arrays_column, only: restore_bgc, & bgc_data_dir, fe_data_type @@ -163,12 +163,12 @@ subroutine get_forcing_bgc !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -!!! midmonth = fix(p5 * real(daymo(month))) ! exact middle +!!! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -184,7 +184,7 @@ subroutine get_forcing_bgc call interp_coeff_monthly (recslot) readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. endif ! 'clim prep' @@ -194,11 +194,11 @@ subroutine get_forcing_bgc !------------------------------------------------------------------- if (trim(bgc_data_type)=='clim' .AND. tr_bgc_Sil) then - ! call read_clim_data (readm, 0, ixm, month, ixp, & + ! call read_clim_data (readm, 0, ixm, mmonth, ixp, & ! sil_file, sil_data, & ! field_loc_center, field_type_scalar) fieldname = 'silicate' - call read_clim_data_nc (readm, 0, ixm, month, ixp, & + call read_clim_data_nc (readm, 0, ixm, mmonth, ixp, & sil_file, fieldname, sil_data, & field_loc_center, field_type_scalar) call interpolate_data (sil_data, sildat) @@ -276,11 +276,11 @@ subroutine get_forcing_bgc !------------------------------------------------------------------- if (trim(bgc_data_type)=='clim' .AND. tr_bgc_Nit) then - ! call read_clim_data (readm, 0, ixm, month, ixp, & + ! call read_clim_data (readm, 0, ixm, mmonth, ixp, & ! nit_file, nit_data, & ! field_loc_center, field_type_scalar) fieldname = 'nitrate' - call read_clim_data_nc (readm, 0, ixm, month, ixp, & + call read_clim_data_nc (readm, 0, ixm, mmonth, ixp, & nit_file, fieldname, nit_data, & field_loc_center, field_type_scalar) call interpolate_data (nit_data, nitdat) @@ -584,7 +584,7 @@ end subroutine faero_default subroutine faero_data - use ice_calendar, only: month, mday, istep, sec + use ice_calendar, only: mmonth, mday, istep, msec use ice_domain_size, only: max_blocks use ice_blocks, only: nx_block, ny_block use ice_flux_bgc, only: faero_atm @@ -625,12 +625,12 @@ subroutine faero_data !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = 99 ! other two points will be used if (mday < midmonth) ixp = 99 @@ -647,23 +647,23 @@ subroutine faero_data ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. ! aero_file = trim(atm_data_dir)//'faero.nc' aero_file = '/usr/projects/climate/eclare/DATA/gx1v3/faero.nc' fieldname='faero_atm001' - call read_clim_data_nc (readm, 0, ixm, month, ixp, & + call read_clim_data_nc (readm, 0, ixm, mmonth, ixp, & aero_file, fieldname, aero1_data, & field_loc_center, field_type_scalar) fieldname='faero_atm002' - call read_clim_data_nc (readm, 0, ixm, month, ixp, & + call read_clim_data_nc (readm, 0, ixm, mmonth, ixp, & aero_file, fieldname, aero2_data, & field_loc_center, field_type_scalar) fieldname='faero_atm003' - call read_clim_data_nc (readm, 0, ixm, month, ixp, & + call read_clim_data_nc (readm, 0, ixm, mmonth, ixp, & aero_file, fieldname, aero3_data, & field_loc_center, field_type_scalar) @@ -727,12 +727,12 @@ subroutine fzaero_data !------------------------------------------------------------------- midmonth = 15 ! data is given on 15th of every month -! midmonth = fix(p5 * real(daymo(month))) ! exact middle +! midmonth = fix(p5 * real(daymo(mmonth))) ! exact middle ! Compute record numbers for surrounding months maxrec = 12 - ixm = mod(month+maxrec-2,maxrec) + 1 - ixp = mod(month, maxrec) + 1 + ixm = mod(mmonth+maxrec-2,maxrec) + 1 + ixp = mod(mmonth, maxrec) + 1 if (mday >= midmonth) ixm = -99 ! other two points will be used if (mday < midmonth) ixp = -99 @@ -749,14 +749,14 @@ subroutine fzaero_data ! Read 2 monthly values readm = .false. - if (istep==1 .or. (mday==midmonth .and. sec==0)) readm = .true. + if (istep==1 .or. (mday==midmonth .and. msec==0)) readm = .true. ! aero_file = trim(atm_data_dir)//'faero.nc' ! Cam5 monthly total black carbon deposition on the gx1 grid" aero_file = '/usr/projects/climate/njeffery/DATA/CAM/Hailong_Wang/Cam5_bc_monthly_popgrid.nc' fieldname='bcd' - call read_clim_data_nc (readm, 0, ixm, month, ixp, & + call read_clim_data_nc (readm, 0, ixm, mmonth, ixp, & aero_file, fieldname, aero_data, & field_loc_center, field_type_scalar) diff --git a/cicecore/cicedynB/general/ice_init.F90 b/cicecore/cicedynB/general/ice_init.F90 index 91c0d10b1..46b614c47 100644 --- a/cicecore/cicedynB/general/ice_init.F90 +++ b/cicecore/cicedynB/general/ice_init.F90 @@ -64,10 +64,11 @@ subroutine input_data n_iso, n_aero, n_zaero, n_algae, & n_doc, n_dic, n_don, n_fed, n_fep, & max_nstrm - use ice_calendar, only: year_init, istep0, histfreq, histfreq_n, & + use ice_calendar, only: year_init, month_init, day_init, sec_init, & + istep0, histfreq, histfreq_n, & dumpfreq, dumpfreq_n, diagfreq, & npt, dt, ndtd, days_per_year, use_leap_years, & - write_ic, dump_last + write_ic, dump_last, npt_unit use ice_arrays_column, only: oceanmixed_ice use ice_restart_column, only: restart_age, restart_FY, restart_lvl, & restart_pond_cesm, restart_pond_lvl, restart_pond_topo, restart_aero, & @@ -154,7 +155,7 @@ subroutine input_data !----------------------------------------------------------------- namelist /setup_nml/ & - days_per_year, use_leap_years, year_init, istep0, & + days_per_year, use_leap_years, istep0, npt_unit, & dt, npt, ndtd, numin, & runtype, runid, bfbflag, numax, & ice_ic, restart, restart_dir, restart_file, & @@ -165,6 +166,7 @@ subroutine input_data dbug, histfreq, histfreq_n, hist_avg, & history_dir, history_file, history_precision, cpl_bgc, & conserv_check, & + year_init, month_init, day_init, sec_init, & write_ic, incond_dir, incond_file, version_name namelist /grid_nml/ & @@ -250,6 +252,9 @@ subroutine input_data days_per_year = 365 ! number of days in a year use_leap_years= .false.! if true, use leap years (Feb 29) year_init = 0 ! initial year + month_init = 1 ! initial month + day_init = 1 ! initial day + sec_init = 0 ! initial second istep0 = 0 ! no. of steps taken in previous integrations, ! real (dumped) or imagined (to set calendar) #ifndef CESMCOUPLED @@ -258,6 +263,7 @@ subroutine input_data numin = 11 ! min allowed unit number numax = 99 ! max allowed unit number npt = 99999 ! total number of time steps (dt) + npt_unit = '1' ! units of npt 'y', 'm', 'd', 's', '1' diagfreq = 24 ! how often diag output is written print_points = .false. ! if true, print point data print_global = .true. ! if true, print global diagnostic data @@ -584,9 +590,13 @@ subroutine input_data call broadcast_scalar(days_per_year, master_task) call broadcast_scalar(use_leap_years, master_task) call broadcast_scalar(year_init, master_task) + call broadcast_scalar(month_init, master_task) + call broadcast_scalar(day_init, master_task) + call broadcast_scalar(sec_init, master_task) call broadcast_scalar(istep0, master_task) call broadcast_scalar(dt, master_task) call broadcast_scalar(npt, master_task) + call broadcast_scalar(npt_unit, master_task) call broadcast_scalar(diagfreq, master_task) call broadcast_scalar(print_points, master_task) call broadcast_scalar(print_global, master_task) @@ -1621,7 +1631,11 @@ subroutine input_data write(nu_diag,1031) ' runid = ', trim(runid) write(nu_diag,1031) ' runtype = ', trim(runtype) write(nu_diag,1021) ' year_init = ', year_init + write(nu_diag,1021) ' month_init = ', month_init + write(nu_diag,1021) ' day_init = ', day_init + write(nu_diag,1021) ' sec_init = ', sec_init write(nu_diag,1021) ' istep0 = ', istep0 + write(nu_diag,1031) ' npt_unit = ', trim(npt_unit) write(nu_diag,1021) ' npt = ', npt write(nu_diag,1021) ' diagfreq = ', diagfreq write(nu_diag,1011) ' print_global = ', print_global diff --git a/cicecore/cicedynB/general/ice_step_mod.F90 b/cicecore/cicedynB/general/ice_step_mod.F90 index b21908e77..38dc328f2 100644 --- a/cicecore/cicedynB/general/ice_step_mod.F90 +++ b/cicecore/cicedynB/general/ice_step_mod.F90 @@ -1026,7 +1026,7 @@ subroutine step_radiation (dt, iblk) kaer_tab, waer_tab, gaer_tab, kaer_bc_tab, waer_bc_tab, & gaer_bc_tab, bcenh, swgrid, igrid use ice_blocks, only: block, get_block - use ice_calendar, only: calendar_type, days_per_year, nextsw_cday, yday, sec + use ice_calendar, only: calendar_type, days_per_year, nextsw_cday, yday, msec use ice_domain, only: blocks_ice use ice_domain_size, only: ncat, n_aero, nilyr, nslyr, n_zaero, n_algae, nblyr use ice_flux, only: swvdr, swvdf, swidr, swidf, coszen, fsnow @@ -1145,7 +1145,7 @@ subroutine step_radiation (dt, iblk) calendar_type=calendar_type, & days_per_year=days_per_year, & nextsw_cday=nextsw_cday, yday=yday, & - sec=sec, & + sec=msec, & kaer_tab=kaer_tab, kaer_bc_tab=kaer_bc_tab(:,:), & waer_tab=waer_tab, waer_bc_tab=waer_bc_tab(:,:), & gaer_tab=gaer_tab, gaer_bc_tab=gaer_bc_tab(:,:), & diff --git a/cicecore/cicedynB/infrastructure/comm/mpi/ice_timers.F90 b/cicecore/cicedynB/infrastructure/comm/mpi/ice_timers.F90 index 6f9c8b0c6..046cf9336 100644 --- a/cicecore/cicedynB/infrastructure/comm/mpi/ice_timers.F90 +++ b/cicecore/cicedynB/infrastructure/comm/mpi/ice_timers.F90 @@ -60,6 +60,7 @@ module ice_timers #endif timer_bound, &! boundary updates timer_bgc, &! biogeochemistry + timer_forcing, &! forcing timer_evp_1d, &! timer only loop timer_evp_2d ! timer including conversion 1d/2d ! timer_tmp ! for temporary timings @@ -179,6 +180,7 @@ subroutine init_ice_timers call get_ice_timer(timer_hist, 'History ',nblocks,distrb_info%nprocs) call get_ice_timer(timer_bound, 'Bound', nblocks,distrb_info%nprocs) call get_ice_timer(timer_bgc, 'BGC', nblocks,distrb_info%nprocs) + call get_ice_timer(timer_forcing, 'Forcing', nblocks,distrb_info%nprocs) #if (defined CESMCOUPLED) call get_ice_timer(timer_cplrecv, 'Cpl-recv', nblocks,distrb_info%nprocs) call get_ice_timer(timer_rcvsnd, 'Rcv->Snd', nblocks,distrb_info%nprocs) diff --git a/cicecore/cicedynB/infrastructure/comm/serial/ice_timers.F90 b/cicecore/cicedynB/infrastructure/comm/serial/ice_timers.F90 index 3074c1dc9..4599de42e 100644 --- a/cicecore/cicedynB/infrastructure/comm/serial/ice_timers.F90 +++ b/cicecore/cicedynB/infrastructure/comm/serial/ice_timers.F90 @@ -52,6 +52,7 @@ module ice_timers timer_hist, &! diagnostics/history timer_bound, &! boundary updates timer_bgc, &! biogeochemistry + timer_forcing, &! forcing timer_evp_1d, &! timer only loop timer_evp_2d ! timer including conversion 1d/2d ! timer_tmp ! for temporary timings @@ -193,8 +194,9 @@ subroutine init_ice_timers call get_ice_timer(timer_hist, 'History ',nblocks,distrb_info%nprocs) call get_ice_timer(timer_bound, 'Bound', nblocks,distrb_info%nprocs) call get_ice_timer(timer_bgc, 'BGC', nblocks,distrb_info%nprocs) - call get_ice_timer(timer_evp_1d, '1d-evp', nblocks,distrb_info%nprocs) - call get_ice_timer(timer_evp_2d, '2d-evp', nblocks,distrb_info%nprocs) + call get_ice_timer(timer_forcing, 'Forcing', nblocks,distrb_info%nprocs) + call get_ice_timer(timer_evp_1d, '1d-evp', nblocks,distrb_info%nprocs) + call get_ice_timer(timer_evp_2d, '2d-evp', nblocks,distrb_info%nprocs) ! call get_ice_timer(timer_tmp, ' ',nblocks,distrb_info%nprocs) !----------------------------------------------------------------------- diff --git a/cicecore/cicedynB/infrastructure/ice_domain.F90 b/cicecore/cicedynB/infrastructure/ice_domain.F90 index cc57ea585..f34d1967e 100644 --- a/cicecore/cicedynB/infrastructure/ice_domain.F90 +++ b/cicecore/cicedynB/infrastructure/ice_domain.F90 @@ -449,7 +449,6 @@ subroutine init_domain_distribution(KMTG,ULATG) if (my_task == master_task) then ! cannot use ice_read_write due to circular dependency #ifdef USE_NETCDF - write(nu_diag,*) 'read ',trim(distribution_wght_file),minval(wght),maxval(wght) status = nf90_open(distribution_wght_file, NF90_NOWRITE, fid) if (status /= nf90_noerr) then call abort_ice (subname//'ERROR: Cannot open '//trim(distribution_wght_file)) @@ -457,6 +456,7 @@ subroutine init_domain_distribution(KMTG,ULATG) status = nf90_inq_varid(fid, 'wght', varid) status = nf90_get_var(fid, varid, wght) status = nf90_close(fid) + write(nu_diag,*) 'read ',trim(distribution_wght_file),minval(wght),maxval(wght) #else call abort_ice(subname//'ERROR: USE_NETCDF cpp not defined', & file=__FILE__, line=__LINE__) diff --git a/cicecore/cicedynB/infrastructure/ice_grid.F90 b/cicecore/cicedynB/infrastructure/ice_grid.F90 index a354efb6b..d50cf5fa1 100644 --- a/cicecore/cicedynB/infrastructure/ice_grid.F90 +++ b/cicecore/cicedynB/infrastructure/ice_grid.F90 @@ -2400,7 +2400,7 @@ subroutine get_bathymetry do iblk = 1, nblocks do j = 1, ny_block do i = 1, nx_block - k = kmt(i,j,iblk) + k = min(nint(kmt(i,j,iblk)),nlevel) if (k > puny) bathymetry(i,j,iblk) = depth(k) enddo enddo diff --git a/cicecore/cicedynB/infrastructure/ice_restart_driver.F90 b/cicecore/cicedynB/infrastructure/ice_restart_driver.F90 index 5a6c79503..1a5681b38 100644 --- a/cicecore/cicedynB/infrastructure/ice_restart_driver.F90 +++ b/cicecore/cicedynB/infrastructure/ice_restart_driver.F90 @@ -197,7 +197,7 @@ subroutine restartfile (ice_ic) use ice_boundary, only: ice_HaloUpdate_stress use ice_blocks, only: nghost, nx_block, ny_block - use ice_calendar, only: istep0, npt + use ice_calendar, only: istep0, npt, calendar use ice_communicate, only: my_task, master_task use ice_domain, only: nblocks, halo_info use ice_domain_size, only: nilyr, nslyr, ncat, & @@ -244,6 +244,7 @@ subroutine restartfile (ice_ic) file=__FILE__, line=__LINE__) call init_restart_read(ice_ic) + call calendar() diag = .true. @@ -529,7 +530,8 @@ subroutine restartfile_v4 (ice_ic) use ice_broadcast, only: broadcast_scalar use ice_blocks, only: nghost, nx_block, ny_block - use ice_calendar, only: istep0, istep1, time, time_forc, calendar, npt + use ice_calendar, only: istep0, istep1, timesecs, calendar, npt, & + set_date_from_timesecs use ice_communicate, only: my_task, master_task use ice_domain, only: nblocks, distrb_info use ice_domain_size, only: nilyr, nslyr, ncat, nx_global, ny_global, & @@ -571,6 +573,9 @@ subroutine restartfile_v4 (ice_ic) real (kind=dbl_kind), dimension(:,:), allocatable :: & work_g1, work_g2 + real (kind=dbl_kind) :: & + time_forc ! historic, now local + character(len=*), parameter :: subname = '(restartfile_v4)' call icepack_query_tracer_sizes(ntrcr_out=ntrcr) @@ -602,14 +607,15 @@ subroutine restartfile_v4 (ice_ic) if (use_restart_time) then if (my_task == master_task) then - read (nu_restart) istep0,time,time_forc - write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc + read (nu_restart) istep0,timesecs,time_forc + write(nu_diag,*) 'Restart read at istep=',istep0,timesecs endif call broadcast_scalar(istep0,master_task) istep1 = istep0 - call broadcast_scalar(time,master_task) - call broadcast_scalar(time_forc,master_task) - call calendar(time) + call broadcast_scalar(timesecs,master_task) +! call broadcast_scalar(time_forc,master_task) + call set_date_from_timesecs(timesecs) + call calendar() else diff --git a/cicecore/cicedynB/infrastructure/io/io_binary/ice_restart.F90 b/cicecore/cicedynB/infrastructure/io/io_binary/ice_restart.F90 index b1a2d026b..91d57ea48 100644 --- a/cicecore/cicedynB/infrastructure/io/io_binary/ice_restart.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_binary/ice_restart.F90 @@ -31,6 +31,8 @@ module ice_restart public :: init_restart_write, init_restart_read, & read_restart_field, write_restart_field, final_restart + real(kind=dbl_kind) :: time_forc = -99. ! historic now local + !======================================================================= contains @@ -42,7 +44,8 @@ module ice_restart subroutine init_restart_read(ice_ic) - use ice_calendar, only: istep0, istep1, time, time_forc, npt, nyr + use ice_calendar, only: istep0, istep1, timesecs, npt, myear, & + set_date_from_timesecs use ice_communicate, only: my_task, master_task use ice_dyn_shared, only: kdyn use ice_read_write, only: ice_open, ice_open_ext @@ -105,17 +108,18 @@ subroutine init_restart_read(ice_ic) call ice_open(nu_restart,trim(filename),0) endif if (use_restart_time) then - read (nu_restart) istep0,time,time_forc,nyr + read (nu_restart) istep0,timesecs,time_forc,myear else read (nu_restart) iignore,rignore,rignore ! use namelist values endif - write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc + write(nu_diag,*) 'Restart read at istep=',istep0,timesecs endif call broadcast_scalar(istep0,master_task) - call broadcast_scalar(time,master_task) + call broadcast_scalar(timesecs,master_task) call broadcast_scalar(time_forc,master_task) - call broadcast_scalar(nyr,master_task) + call broadcast_scalar(myear,master_task) + call set_date_from_timesecs(timesecs) istep1 = istep0 @@ -375,8 +379,8 @@ end subroutine init_restart_read subroutine init_restart_write(filename_spec) - use ice_calendar, only: sec, month, mday, nyr, istep1, & - time, time_forc, year_init + use ice_calendar, only: msec, mmonth, mday, myear, istep1, & + timesecs use ice_communicate, only: my_task, master_task use ice_dyn_shared, only: kdyn use ice_read_write, only: ice_open, ice_open_ext @@ -391,8 +395,7 @@ subroutine init_restart_write(filename_spec) tr_pond_topo, tr_pond_lvl, tr_brine integer (kind=int_kind) :: & - nbtrcr, & ! number of bgc tracers - iyear, imonth, iday ! year, month, day + nbtrcr ! number of bgc tracers character(len=char_len_long) :: filename @@ -414,14 +417,10 @@ subroutine init_restart_write(filename_spec) 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 + myear,'-',mmonth,'-',mday,'-',msec end if ! write pointer (path/file) @@ -434,7 +433,7 @@ subroutine init_restart_write(filename_spec) else call ice_open(nu_dump,filename,0) endif - write(nu_dump) istep1,time,time_forc,nyr + write(nu_dump) istep1,timesecs,time_forc,myear write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -445,7 +444,7 @@ subroutine init_restart_write(filename_spec) 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)),'.eap.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_eap,filename,0) @@ -454,7 +453,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_eap) istep1,time,time_forc + write(nu_dump_eap) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -465,7 +464,7 @@ subroutine init_restart_write(filename_spec) 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)),'.fsd.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_fsd,filename,0) @@ -474,7 +473,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_fsd) istep1,time,time_forc + write(nu_dump_fsd) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -485,7 +484,7 @@ subroutine init_restart_write(filename_spec) 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)),'.FY.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_FY,filename,0) @@ -494,7 +493,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_FY) istep1,time,time_forc + write(nu_dump_FY) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -505,7 +504,7 @@ subroutine init_restart_write(filename_spec) 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)),'.iage.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_age,filename,0) @@ -514,7 +513,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_age) istep1,time,time_forc + write(nu_dump_age) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -525,7 +524,7 @@ subroutine init_restart_write(filename_spec) 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)),'.lvl.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_lvl,filename,0) @@ -534,7 +533,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_lvl) istep1,time,time_forc + write(nu_dump_lvl) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -545,7 +544,7 @@ subroutine init_restart_write(filename_spec) 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)),'.pond_cesm.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_pond,filename,0) @@ -554,7 +553,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_pond) istep1,time,time_forc + write(nu_dump_pond) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -565,7 +564,7 @@ subroutine init_restart_write(filename_spec) 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)),'.pond_lvl.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_pond,filename,0) @@ -574,7 +573,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_pond) istep1,time,time_forc + write(nu_dump_pond) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -585,7 +584,7 @@ subroutine init_restart_write(filename_spec) 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)),'.pond_topo.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_pond,filename,0) @@ -594,7 +593,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_pond) istep1,time,time_forc + write(nu_dump_pond) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -605,7 +604,7 @@ subroutine init_restart_write(filename_spec) 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)),'.brine.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_hbrine,filename,0) @@ -614,7 +613,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_hbrine) istep1,time,time_forc + write(nu_dump_hbrine) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -625,7 +624,7 @@ subroutine init_restart_write(filename_spec) 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)),'.bgc.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_bgc,filename,0) @@ -634,7 +633,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_bgc) istep1,time,time_forc + write(nu_dump_bgc) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif endif @@ -644,7 +643,7 @@ subroutine init_restart_write(filename_spec) 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)),'.iso.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_iso,filename,0) @@ -653,7 +652,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_iso) istep1,time,time_forc + write(nu_dump_iso) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -664,7 +663,7 @@ subroutine init_restart_write(filename_spec) 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)),'.aero.', & - iyear,'-',month,'-',mday,'-',sec + myear,'-',mmonth,'-',mday,'-',msec if (restart_ext) then call ice_open_ext(nu_dump_aero,filename,0) @@ -673,7 +672,7 @@ subroutine init_restart_write(filename_spec) endif if (my_task == master_task) then - write(nu_dump_aero) istep1,time,time_forc + write(nu_dump_aero) istep1,timesecs,time_forc write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) endif @@ -803,7 +802,7 @@ end subroutine write_restart_field subroutine final_restart() - use ice_calendar, only: istep1, time, time_forc + use ice_calendar, only: istep1, timesecs use ice_communicate, only: my_task, master_task logical (kind=log_kind) :: & @@ -843,7 +842,7 @@ subroutine final_restart() if (solve_zsal .or. nbtrcr > 0) & close(nu_dump_bgc) - write(nu_diag,*) 'Restart read/written ',istep1,time,time_forc + write(nu_diag,*) 'Restart read/written ',istep1,timesecs endif end subroutine final_restart diff --git a/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 b/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 index b3024302e..9c6b30ee1 100644 --- a/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_history_write.F90 @@ -47,8 +47,9 @@ subroutine ice_write_hist (ns) use ice_arrays_column, only: hin_max, floe_rad_c use ice_blocks, only: nx_block, ny_block use ice_broadcast, only: broadcast_scalar - use ice_calendar, only: time, sec, idate, idate0, write_ic, & - histfreq, dayyr, days_per_year, use_leap_years + use ice_calendar, only: msec, timesecs, idate, idate0, write_ic, & + histfreq, days_per_year, use_leap_years, dayyr, & + year_init, month_init, day_init use ice_communicate, only: my_task, master_task use ice_domain, only: distrb_info use ice_domain_size, only: nx_global, ny_global, max_nstrm, max_blocks @@ -80,7 +81,6 @@ subroutine ice_write_hist (ns) integer (kind=int_kind), dimension(5) :: dimidcz integer (kind=int_kind), dimension(3) :: dimid_nverts integer (kind=int_kind), dimension(6) :: dimidex -! real (kind=real_kind) :: ltime real (kind=dbl_kind) :: ltime2 character (char_len) :: title character (char_len_long) :: ncfile(max_nstrm) @@ -133,8 +133,7 @@ subroutine ice_write_hist (ns) if (my_task == master_task) then -! ltime=time/int(secday) - ltime2=time/int(secday) + ltime2 = timesecs/secday call construct_filename(ncfile(ns),'nc',ns) @@ -1038,9 +1037,9 @@ subroutine ice_write_hist (ns) 'ERROR: global attribute source') if (use_leap_years) then - write(title,'(a,i3,a)') 'This year has ',int(dayyr),' days' + write(title,'(a,i3,a)') 'This year has ',dayyr,' days' else - write(title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days' + write(title,'(a,i3,a)') 'All years have exactly ',dayyr,' days' endif status = nf90_put_att(ncid,nf90_global,'comment',title) if (status /= nf90_noerr) call abort_ice(subname// & @@ -1051,7 +1050,7 @@ subroutine ice_write_hist (ns) if (status /= nf90_noerr) call abort_ice(subname// & 'ERROR: global attribute date1') - write(title,'(a,i6)') 'seconds elapsed into model date: ',sec + write(title,'(a,i6)') 'seconds elapsed into model date: ',msec status = nf90_put_att(ncid,nf90_global,'comment3',title) if (status /= nf90_noerr) call abort_ice(subname// & 'ERROR: global attribute date2') @@ -1091,7 +1090,6 @@ subroutine ice_write_hist (ns) status = nf90_inq_varid(ncid,'time',varid) if (status /= nf90_noerr) call abort_ice(subname// & 'ERROR: getting time varid') -!sgl status = nf90_put_var(ncid,varid,ltime) status = nf90_put_var(ncid,varid,ltime2) if (status /= nf90_noerr) call abort_ice(subname// & 'ERROR: writing time variable') diff --git a/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_restart.F90 b/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_restart.F90 index 53c7dac60..e744caf09 100644 --- a/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_restart.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_netcdf/ice_restart.F90 @@ -42,8 +42,8 @@ module ice_restart subroutine init_restart_read(ice_ic) - use ice_calendar, only: sec, month, mday, nyr, istep0, istep1, & - time, time_forc, npt + use ice_calendar, only: msec, mmonth, mday, myear, & + istep0, istep1, npt use ice_communicate, only: my_task, master_task character(len=char_len_long), intent(in), optional :: ice_ic @@ -53,7 +53,7 @@ subroutine init_restart_read(ice_ic) character(len=char_len_long) :: & filename, filename0 - integer (kind=int_kind) :: status + integer (kind=int_kind) :: status, status1 character(len=*), parameter :: subname = '(init_restart_read)' @@ -79,24 +79,36 @@ subroutine init_restart_read(ice_ic) '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) + status1 = nf90_noerr + status = nf90_get_att(ncid, nf90_global, 'istep1', istep0) + if (status /= nf90_noerr) status1 = status +! 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, 'myear', myear) + if (status /= nf90_noerr) status = nf90_get_att(ncid, nf90_global, 'nyr', myear) + if (status /= nf90_noerr) status1 = status + status = nf90_get_att(ncid, nf90_global, 'mmonth', mmonth) + if (status /= nf90_noerr) status = nf90_get_att(ncid, nf90_global, 'month', mmonth) + if (status /= nf90_noerr) status1 = status status = nf90_get_att(ncid, nf90_global, 'mday', mday) - status = nf90_get_att(ncid, nf90_global, 'sec', sec) - endif + if (status /= nf90_noerr) status1 = status + status = nf90_get_att(ncid, nf90_global, 'msec', msec) + if (status /= nf90_noerr) status = nf90_get_att(ncid, nf90_global, 'sec', msec) + if (status /= nf90_noerr) status1 = status + if (status1 /= nf90_noerr) call abort_ice(subname// & + 'ERROR: reading restart time '//trim(filename)) endif ! use namelist values if use_restart_time = F endif call broadcast_scalar(istep0,master_task) - call broadcast_scalar(time,master_task) - call broadcast_scalar(time_forc,master_task) - call broadcast_scalar(nyr,master_task) - +! call broadcast_scalar(time,master_task) + call broadcast_scalar(myear,master_task) + call broadcast_scalar(mmonth,master_task) + call broadcast_scalar(mday,master_task) + call broadcast_scalar(msec,master_task) +! call broadcast_scalar(time_forc,master_task) + istep1 = istep0 ! if runid is bering then need to correct npt for istep0 @@ -118,8 +130,7 @@ end subroutine init_restart_read 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_calendar, only: msec, mmonth, mday, myear, istep1 use ice_communicate, only: my_task, master_task use ice_domain_size, only: nx_global, ny_global, ncat, nilyr, nslyr, & n_iso, n_aero, nblyr, n_zaero, n_algae, n_doc, & @@ -145,7 +156,6 @@ subroutine init_restart_write(filename_spec) integer (kind=int_kind) :: & k, n, & ! index nx, ny, & ! global array size - iyear, & ! year nbtrcr ! number of bgc tracers character(len=char_len_long) :: filename @@ -186,12 +196,10 @@ subroutine init_restart_write(filename_spec) if (present(filename_spec)) then filename = trim(filename_spec) else - iyear = nyr + year_init - 1 - 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 + myear,'-',mmonth,'-',mday,'-',msec end if ! write pointer (path/file) @@ -208,12 +216,12 @@ subroutine init_restart_write(filename_spec) '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,'time',time) +! status = nf90_put_att(ncid,nf90_global,'time_forc',time_forc) + status = nf90_put_att(ncid,nf90_global,'myear',myear) + status = nf90_put_att(ncid,nf90_global,'mmonth',mmonth) status = nf90_put_att(ncid,nf90_global,'mday',mday) - status = nf90_put_att(ncid,nf90_global,'sec',sec) + status = nf90_put_att(ncid,nf90_global,'msec',msec) nx = nx_global ny = ny_global @@ -795,7 +803,7 @@ end subroutine write_restart_field subroutine final_restart() - use ice_calendar, only: istep1, time, time_forc + use ice_calendar, only: istep1, idate use ice_communicate, only: my_task, master_task integer (kind=int_kind) :: status @@ -806,7 +814,7 @@ subroutine final_restart() status = nf90_close(ncid) if (my_task == master_task) & - write(nu_diag,*) 'Restart read/written ',istep1,time,time_forc + write(nu_diag,*) 'Restart read/written ',istep1,idate #else call abort_ice(subname//'ERROR: USE_NETCDF cpp not defined', & diff --git a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 index 7e16f2591..72a1ed97f 100644 --- a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 @@ -41,8 +41,8 @@ subroutine ice_write_hist (ns) use ice_blocks, only: nx_block, ny_block use ice_broadcast, only: broadcast_scalar - use ice_calendar, only: time, sec, idate, idate0, write_ic, & - histfreq, dayyr, days_per_year, use_leap_years + use ice_calendar, only: msec, timesecs, idate, idate0, write_ic, & + histfreq, days_per_year, use_leap_years, dayyr use ice_communicate, only: my_task, master_task use ice_constants, only: c0, c360, spval, spval_dbl use ice_domain, only: distrb_info, nblocks @@ -76,7 +76,6 @@ subroutine ice_write_hist (ns) character (char_len_long) :: ncfile(max_nstrm) integer (kind=int_kind) :: iotype - integer (kind=int_kind) :: iyear, imonth, iday integer (kind=int_kind) :: icategory,ind,i_aice,boundid character (char_len) :: start_time,current_date,current_time @@ -176,8 +175,8 @@ subroutine ice_write_hist (ns) call ice_pio_initdecomp(ndim3=nzslyr, ndim4=ncat_hist, iodesc=iodesc4ds) call ice_pio_initdecomp(ndim3=nfsd_hist, ndim4=ncat_hist, iodesc=iodesc4df) - ltime2 = time/int(secday) - ltime = real(time/int(secday),kind=real_kind) + ltime2 = timesecs/secday + ltime = real(timesecs/secday,kind=real_kind) ! option of turning on double precision history files lprecision = pio_real @@ -861,16 +860,16 @@ subroutine ice_write_hist (ns) status = pio_put_att(File,pio_global,'source',trim(title)) if (use_leap_years) then - write(title,'(a,i3,a)') 'This year has ',int(dayyr),' days' + write(title,'(a,i3,a)') 'This year has ',dayyr,' days' else - write(title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days' + write(title,'(a,i3,a)') 'All years have exactly ',dayyr,' days' endif status = pio_put_att(File,pio_global,'comment',trim(title)) write(title,'(a,i8.8)') 'File written on model date ',idate status = pio_put_att(File,pio_global,'comment2',trim(title)) - write(title,'(a,i6)') 'seconds elapsed into model date: ',sec + write(title,'(a,i6)') 'seconds elapsed into model date: ',msec status = pio_put_att(File,pio_global,'comment3',trim(title)) title = 'CF-1.0' diff --git a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 index eb703abcd..12d5d8e71 100644 --- a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 +++ b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 @@ -41,8 +41,8 @@ module ice_restart subroutine init_restart_read(ice_ic) - use ice_calendar, only: istep0, istep1, time, time_forc, nyr, month, & - mday, sec, npt + use ice_calendar, only: istep0, istep1, myear, mmonth, & + mday, msec, npt use ice_communicate, only: my_task, master_task use ice_domain_size, only: ncat use ice_read_write, only: ice_open @@ -54,7 +54,7 @@ subroutine init_restart_read(ice_ic) character(len=char_len_long) :: & filename, filename0 - integer (kind=int_kind) :: status + integer (kind=int_kind) :: status, status1 integer (kind=int_kind) :: iotype @@ -87,28 +87,40 @@ subroutine init_restart_read(ice_ic) call ice_pio_initdecomp(ndim3=ncat , iodesc=iodesc3d_ncat,remap=.true.) if (use_restart_time) then - status = pio_get_att(File, pio_global, 'istep1', istep0) - status = pio_get_att(File, pio_global, 'time', time) - status = pio_get_att(File, pio_global, 'time_forc', time_forc) - call pio_seterrorhandling(File, PIO_BCAST_ERROR) - status = pio_get_att(File, pio_global, 'nyr', nyr) - call pio_seterrorhandling(File, PIO_INTERNAL_ERROR) - if (status == PIO_noerr) then - status = pio_get_att(File, pio_global, 'month', month) + status1 = PIO_noerr + status = pio_get_att(File, pio_global, 'istep1', istep0) +! status = pio_get_att(File, pio_global, 'time', time) +! status = pio_get_att(File, pio_global, 'time_forc', time_forc) + call pio_seterrorhandling(File, PIO_BCAST_ERROR) + status = pio_get_att(File, pio_global, 'myear', myear) + if (status /= PIO_noerr) status = pio_get_att(File, pio_global, 'nyr', myear) + if (status /= PIO_noerr) status1 = status + status = pio_get_att(File, pio_global, 'mmonth', mmonth) + if (status /= PIO_noerr) status = pio_get_att(File, pio_global, 'month', mmonth) + if (status /= PIO_noerr) status1 = status status = pio_get_att(File, pio_global, 'mday', mday) - status = pio_get_att(File, pio_global, 'sec', sec) - endif + if (status /= PIO_noerr) status1 = status + status = pio_get_att(File, pio_global, 'msec', msec) + if (status /= PIO_noerr) status = pio_get_att(File, pio_global, 'sec', msec) + if (status /= PIO_noerr) status1 = status + if (status1 /= PIO_noerr) & + call abort_ice(subname//"ERROR: reading restart time ") + call pio_seterrorhandling(File, PIO_INTERNAL_ERROR) endif ! use namelist values if use_restart_time = F ! endif if (my_task == master_task) then - write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc + write(nu_diag,*) 'Restart read at istep=',istep0,myear,mmonth,mday,msec endif call broadcast_scalar(istep0,master_task) - call broadcast_scalar(time,master_task) - call broadcast_scalar(time_forc,master_task) - call broadcast_scalar(nyr,master_task) + call broadcast_scalar(myear,master_task) + call broadcast_scalar(mmonth,master_task) + call broadcast_scalar(mday,master_task) + call broadcast_scalar(msec,master_task) +! call broadcast_scalar(time,master_task) +! call broadcast_scalar(time_forc,master_task) + call broadcast_scalar(myear,master_task) istep1 = istep0 @@ -126,8 +138,7 @@ end subroutine init_restart_read subroutine init_restart_write(filename_spec) - use ice_calendar, only: sec, month, mday, nyr, istep1, & - time, time_forc, year_init + use ice_calendar, only: msec, mmonth, mday, myear, istep1 use ice_communicate, only: my_task, master_task use ice_domain_size, only: nx_global, ny_global, ncat, nilyr, nslyr, & n_iso, n_aero, nblyr, n_zaero, n_algae, n_doc, & @@ -155,9 +166,6 @@ subroutine init_restart_write(filename_spec) ! local variables - integer (kind=int_kind) :: & - iyear, imonth, iday ! year, month, day - character(len=char_len_long) :: filename integer (kind=int_kind) :: dimid_ni, dimid_nj, dimid_ncat, & @@ -196,14 +204,10 @@ subroutine init_restart_write(filename_spec) 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 + myear,'-',mmonth,'-',mday,'-',msec end if if (restart_format(1:3) /= 'bin') filename = trim(filename) // '.nc' @@ -224,12 +228,12 @@ subroutine init_restart_write(filename_spec) clobber=.true., cdf64=lcdf64, iotype=iotype) status = pio_put_att(File,pio_global,'istep1',istep1) - status = pio_put_att(File,pio_global,'time',time) - status = pio_put_att(File,pio_global,'time_forc',time_forc) - status = pio_put_att(File,pio_global,'nyr',nyr) - status = pio_put_att(File,pio_global,'month',month) +! status = pio_put_att(File,pio_global,'time',time) +! status = pio_put_att(File,pio_global,'time_forc',time_forc) + status = pio_put_att(File,pio_global,'myear',myear) + status = pio_put_att(File,pio_global,'mmonth',mmonth) status = pio_put_att(File,pio_global,'mday',mday) - status = pio_put_att(File,pio_global,'sec',sec) + status = pio_put_att(File,pio_global,'msec',msec) status = pio_def_dim(File,'ni',nx_global,dimid_ni) status = pio_def_dim(File,'nj',ny_global,dimid_nj) @@ -702,7 +706,7 @@ subroutine read_restart_field(nu,nrec,work,atype,vname,ndim3,diag, & status = pio_inq_varid(File,trim(vname),vardesc) - if (status /= 0) then + if (status /= PIO_noerr) then call abort_ice(subname//"ERROR: CICE restart? Missing variable: "//trim(vname)) endif @@ -854,7 +858,7 @@ end subroutine write_restart_field subroutine final_restart() - use ice_calendar, only: istep1, time, time_forc + use ice_calendar, only: istep1, idate, msec use ice_communicate, only: my_task, master_task character(len=*), parameter :: subname = '(final_restart)' @@ -864,7 +868,7 @@ subroutine final_restart() call pio_closefile(File) if (my_task == master_task) & - write(nu_diag,*) 'Restart read/written ',istep1,time,time_forc + write(nu_diag,*) 'Restart read/written ',istep1,idate,msec end subroutine final_restart diff --git a/cicecore/drivers/mct/cesm1/CICE_InitMod.F90 b/cicecore/drivers/mct/cesm1/CICE_InitMod.F90 index da745d965..c3de87f68 100644 --- a/cicecore/drivers/mct/cesm1/CICE_InitMod.F90 +++ b/cicecore/drivers/mct/cesm1/CICE_InitMod.F90 @@ -64,7 +64,7 @@ subroutine cice_init(mpicom_ice) floe_binwidth, c_fsd_range use ice_state, only: alloc_state use ice_flux_bgc, only: alloc_flux_bgc - use ice_calendar, only: dt, dt_dyn, time, istep, istep1, write_ic, & + use ice_calendar, only: dt, dt_dyn, istep, istep1, write_ic, & init_calendar, calendar use ice_communicate, only: init_communicate, my_task, master_task use ice_diagnostics, only: init_diags @@ -156,7 +156,7 @@ subroutine cice_init(mpicom_ice) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) - call calendar(time) ! determine the initial date + call calendar ! determine the initial date call init_forcing_ocn(dt) ! initialize sss and sst from data call init_state ! initialize the ice state @@ -233,7 +233,7 @@ subroutine init_restart use ice_arrays_column, only: dhsn use ice_blocks, only: nx_block, ny_block - use ice_calendar, only: time, calendar + use ice_calendar, only: calendar use ice_constants, only: c0 use ice_domain, only: nblocks use ice_domain_size, only: ncat, n_iso, n_aero, nfsd @@ -295,7 +295,7 @@ subroutine init_restart if (trim(runtype) == 'continue') then ! start from core restart file call restartfile() ! given by pointer in ice_in - call calendar(time) ! update time parameters + call calendar ! update time parameters if (kdyn == 2) call read_restart_eap ! EAP else if (restart) then ! ice_ic = core restart file call restartfile (ice_ic) ! or 'default' or 'none' diff --git a/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 b/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 index d53014b7b..cdcfb6175 100644 --- a/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 +++ b/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 @@ -44,7 +44,7 @@ module CICE_RunMod subroutine CICE_Run - use ice_calendar, only: istep, istep1, time, dt, stop_now, calendar + use ice_calendar, only: istep, istep1, dt, stop_now, calendar, advance_timestep use ice_forcing, only: get_forcing_atmo, get_forcing_ocn, & get_wave_spec use ice_forcing_bgc, only: get_forcing_bgc, get_atm_bgc, & @@ -81,12 +81,14 @@ subroutine CICE_Run ! call ice_step - istep = istep + 1 ! update time step counters - istep1 = istep1 + 1 - time = time + dt ! determine the time and date +! istep = istep + 1 ! update time step counters +! istep1 = istep1 + 1 +! time = time + dt ! determine the time and date ! call calendar(time) ! at the end of the timestep + call advance_timestep() ! advance timestep and update calendar data + call ice_timer_start(timer_couple) ! atm/ocn coupling ! for standalone @@ -108,7 +110,7 @@ subroutine CICE_Run call init_flux_atm ! Initialize atmosphere fluxes sent to coupler call init_flux_ocn ! initialize ocean fluxes sent to coupler - call calendar(time) ! at the end of the timestep +! call calendar(time) ! at the end of the timestep call ice_timer_stop(timer_couple) ! atm/ocn coupling diff --git a/cicecore/drivers/mct/cesm1/ice_comp_esmf.F90 b/cicecore/drivers/mct/cesm1/ice_comp_esmf.F90 index 8ae80abdc..414c6bded 100644 --- a/cicecore/drivers/mct/cesm1/ice_comp_esmf.F90 +++ b/cicecore/drivers/mct/cesm1/ice_comp_esmf.F90 @@ -48,7 +48,7 @@ module ice_comp_esmf use ice_constants, only : c0, c1, spval_dbl, rad_to_deg, radius, secday use ice_communicate, only : my_task, master_task, MPI_COMM_ICE use ice_calendar, only : istep, istep1, force_restart_now, write_ic,& - idate, idate0, mday, time, month, daycal, & + idate, idate0, mday, time, month, & sec, dt, dt_dyn, calendar, & calendar_type, nextsw_cday, days_per_year, & nyr, new_year, time2sec, year_init @@ -183,7 +183,6 @@ subroutine ice_init_esmf(comp, import_state, export_state, EClock, rc) integer :: shrlogunit,shrloglev ! old values integer :: iam,ierr integer :: lbnum - integer :: daycal(13) !number of cumulative days per month integer :: nleaps ! number of leap days before current year integer :: mpicom_loc, mpicom_vm, gsize integer :: nfields diff --git a/cicecore/drivers/mct/cesm1/ice_comp_mct.F90 b/cicecore/drivers/mct/cesm1/ice_comp_mct.F90 index 7162d6397..b29d62daa 100644 --- a/cicecore/drivers/mct/cesm1/ice_comp_mct.F90 +++ b/cicecore/drivers/mct/cesm1/ice_comp_mct.F90 @@ -47,10 +47,9 @@ module ice_comp_mct use ice_constants, only : ice_init_constants use ice_communicate, only : my_task, master_task, MPI_COMM_ICE use ice_calendar, only : istep, istep1, force_restart_now, write_ic,& - idate, idate0, mday, time, month, daycal, & + idate, idate0, mday, month, nyr, & sec, dt, dt_dyn, calendar, & - calendar_type, nextsw_cday, days_per_year, & - nyr, new_year, time2sec, year_init + calendar_type, nextsw_cday, days_per_year use ice_timers use ice_kinds_mod, only : int_kind, dbl_kind, char_len_long, log_kind @@ -151,13 +150,11 @@ subroutine ice_init_mct( EClock, cdata_i, x2i_i, i2x_i, NLFilename ) integer :: curr_tod ! Current time of day (s) integer :: ref_ymd ! Reference date (YYYYMMDD) integer :: ref_tod ! reference time of day (s) - integer :: iyear ! yyyy integer :: nyrp ! yyyy integer :: dtime ! time step integer :: shrlogunit,shrloglev ! old values integer :: iam,ierr integer :: lbnum - integer :: daycal(13) !number of cumulative days per month integer :: nleaps ! number of leap days before current year integer :: mpicom_loc ! temporary mpicom logical (kind=log_kind) :: atm_aero, tr_aero, tr_zaero @@ -302,10 +299,9 @@ subroutine ice_init_mct( EClock, cdata_i, x2i_i, i2x_i, NLFilename ) ! - on restart run ! - istep0, time and time_forc are read from restart file ! - istep1 is set to istep0 - ! - idate is determined from time via the call to calendar (see below) + ! - date information is determined from restart ! - on initial run - ! - iyear, month and mday obtained from sync clock - ! - time determined from iyear, month and mday + ! - nyr, month, mday, sec obtained from sync clock ! - istep0 and istep1 are set to 0 call seq_timemgr_EClockGetData(EClock, & @@ -335,37 +331,26 @@ subroutine ice_init_mct( EClock, cdata_i, x2i_i, i2x_i, NLFilename ) idate0 = curr_ymd idate = curr_ymd -! idate0 = curr_ymd - (year_init*10000) -! idate = curr_ymd - (year_init*10000) - if (idate < 0) then - write(nu_diag,*) trim(subname),' ERROR curr_ymd,year_init =',curr_ymd,year_init + write(nu_diag,*) trim(subname),' ERROR curr_ymd =',curr_ymd write(nu_diag,*) trim(subname),' ERROR idate lt zero',idate call shr_sys_abort(subname//' :: ERROR idate lt zero') endif - iyear = (idate/10000) ! integer year of basedate - month = (idate-iyear*10000)/100 ! integer month of basedate - mday = idate-iyear*10000-month*100 ! day of month of basedate + nyr = (idate/10000) ! integer year of basedate + month = (idate-nyr*10000)/100 ! integer month of basedate + mday = idate-nyr*10000-month*100 ! day of month of basedate + sec = start_tod ! seconds if (my_task == master_task) then write(nu_diag,*) trim(subname),' curr_ymd = ',curr_ymd - write(nu_diag,*) trim(subname),' cice year_init = ',year_init write(nu_diag,*) trim(subname),' cice start date = ',idate - write(nu_diag,*) trim(subname),' cice start ymds = ',iyear,month,mday,start_tod + write(nu_diag,*) trim(subname),' cice start ymds = ',nyr,month,mday,start_tod endif - if (calendar_type /= "GREGORIAN") then - call time2sec(iyear-year_init,month,mday,time) - else - call time2sec(iyear-(year_init-1),month,mday,time) - endif - - time = time+start_tod - call shr_sys_flush(nu_diag) end if - call calendar(time) ! update calendar info + call calendar ! update calendar info if (write_ic) call accum_hist(dt) ! write initial conditions !--------------------------------------------------------------------------- diff --git a/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 b/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 index 644ef72fa..f47a68392 100644 --- a/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 +++ b/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 @@ -44,7 +44,7 @@ module CICE_RunMod subroutine CICE_Run - use ice_calendar, only: istep, istep1, time, dt, stop_now, calendar + use ice_calendar, only: istep, istep1, time, dt, calendar use ice_forcing, only: get_forcing_atmo, get_forcing_ocn, & get_wave_spec use ice_forcing_bgc, only: get_forcing_bgc, get_atm_bgc, & diff --git a/cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90 b/cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90 index da3d95369..8bba18470 100644 --- a/cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90 +++ b/cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90 @@ -29,7 +29,7 @@ module ice_comp_nuopc use ice_grid , only : tlon, tlat, hm, tarea, ULON, ULAT use ice_communicate , only : init_communicate, my_task, master_task, mpi_comm_ice use ice_calendar , only : force_restart_now, write_ic - use ice_calendar , only : idate, mday, time, month, daycal, time2sec, year_init + use ice_calendar , only : idate, mday, time, month, time2sec, year_init use ice_calendar , only : sec, dt, calendar, calendar_type, nextsw_cday, istep use ice_kinds_mod , only : dbl_kind, int_kind, char_len, char_len_long use ice_scam , only : scmlat, scmlon, single_column diff --git a/cicecore/drivers/standalone/cice/CICE_InitMod.F90 b/cicecore/drivers/standalone/cice/CICE_InitMod.F90 index 8b507740d..60f71fa8a 100644 --- a/cicecore/drivers/standalone/cice/CICE_InitMod.F90 +++ b/cicecore/drivers/standalone/cice/CICE_InitMod.F90 @@ -64,8 +64,8 @@ subroutine cice_init floe_binwidth, c_fsd_range use ice_state, only: alloc_state use ice_flux_bgc, only: alloc_flux_bgc - use ice_calendar, only: dt, dt_dyn, time, istep, istep1, write_ic, & - init_calendar, calendar + use ice_calendar, only: dt, dt_dyn, istep, istep1, write_ic, & + init_calendar, advance_timestep, calc_timesteps use ice_communicate, only: init_communicate, my_task, master_task use ice_diagnostics, only: init_diags use ice_domain, only: init_domain_blocks @@ -156,8 +156,6 @@ subroutine cice_init if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) -! call calendar(time) ! determine the initial date - call init_forcing_ocn(dt) ! initialize sss and sst from data call init_state ! initialize the ice state call init_transport ! initialize horizontal transport @@ -175,6 +173,7 @@ subroutine cice_init call init_diags ! initialize diagnostic output points call init_history_therm ! initialize thermo history variables call init_history_dyn ! initialize dynamic history variables + call calc_timesteps ! update timestep counter if not using npt_unit="1" call icepack_query_tracer_flags(tr_aero_out=tr_aero, tr_zaero_out=tr_zaero) call icepack_query_tracer_flags(tr_iso_out=tr_iso) @@ -191,10 +190,12 @@ subroutine cice_init if (trim(runtype) == 'continue' .or. restart) & call init_shortwave ! initialize radiative transfer - istep = istep + 1 ! update time step counters - istep1 = istep1 + 1 - time = time + dt ! determine the time and date - call calendar(time) ! at the end of the first timestep +! tcraig, use advance_timestep here +! istep = istep + 1 ! update time step counters +! istep1 = istep1 + 1 +! time = time + dt ! determine the time and date +! call calendar(time) ! at the end of the first timestep + call advance_timestep() !-------------------------------------------------------------------- ! coupler communication or forcing data initialization @@ -231,7 +232,7 @@ subroutine init_restart use ice_arrays_column, only: dhsn use ice_blocks, only: nx_block, ny_block - use ice_calendar, only: time, calendar + use ice_calendar, only: calendar use ice_constants, only: c0 use ice_domain, only: nblocks use ice_domain_size, only: ncat, n_iso, n_aero, nfsd @@ -293,7 +294,7 @@ subroutine init_restart if (trim(runtype) == 'continue') then ! start from core restart file call restartfile() ! given by pointer in ice_in - call calendar(time) ! update time parameters + call calendar() ! update time parameters if (kdyn == 2) call read_restart_eap ! EAP else if (restart) then ! ice_ic = core restart file call restartfile (ice_ic) ! or 'default' or 'none' diff --git a/cicecore/drivers/standalone/cice/CICE_RunMod.F90 b/cicecore/drivers/standalone/cice/CICE_RunMod.F90 index 9f6f42f28..72acff992 100644 --- a/cicecore/drivers/standalone/cice/CICE_RunMod.F90 +++ b/cicecore/drivers/standalone/cice/CICE_RunMod.F90 @@ -43,7 +43,7 @@ module CICE_RunMod subroutine CICE_Run - use ice_calendar, only: istep, istep1, time, dt, stop_now, calendar + use ice_calendar, only: istep, istep1, dt, stop_now, advance_timestep use ice_forcing, only: get_forcing_atmo, get_forcing_ocn, & get_wave_spec use ice_forcing_bgc, only: get_forcing_bgc, get_atm_bgc, & @@ -82,11 +82,12 @@ subroutine CICE_Run call ice_step - istep = istep + 1 ! update time step counters - istep1 = istep1 + 1 - time = time + dt ! determine the time and date - - call calendar(time) ! at the end of the timestep +! tcraig, use advance_timestep now +! istep = istep + 1 ! update time step counters +! istep1 = istep1 + 1 +! time = time + dt ! determine the time and date +! call calendar(time) ! at the end of the timestep + call advance_timestep() ! advance time #ifndef CICE_IN_NEMO if (stop_now >= 1) exit timeLoop diff --git a/cicecore/drivers/unittest/calchk/calchk.F90 b/cicecore/drivers/unittest/calchk/calchk.F90 new file mode 100644 index 000000000..ba740501f --- /dev/null +++ b/cicecore/drivers/unittest/calchk/calchk.F90 @@ -0,0 +1,588 @@ + + program calchk + + use ice_kinds_mod, only: int_kind, dbl_kind + use ice_calendar, only: myear, mmonth, mday, msec + use ice_calendar, only: year_init, month_init, day_init, sec_init + use ice_calendar, only: dt, ndtd, istep0, diagfreq, npt, npt_unit + use ice_calendar, only: months_per_year, daymo, timesecs, seconds_per_day + use ice_calendar, only: use_leap_years, days_per_year + use ice_calendar, only: compute_elapsed_days + use ice_calendar, only: update_date, calc_timesteps + use ice_calendar, only: init_calendar, calendar + use ice_calendar, only: set_date_from_timesecs + use ice_calendar, only: calendar_date2time, calendar_time2date + use ice_calendar, only: compute_calendar_data + implicit none + + integer(kind=int_kind) :: yearmax + integer(kind=int_kind) :: nday,nptc + integer(kind=int_kind) :: n,m,ny,nm,nd,nf1,nf2,xadd,nfa,nfb,nfc,ns1,ns2 + integer(kind=int_kind) :: yi,mi,di,si + integer(kind=int_kind) :: dyear,dmon,dday,dsec + integer(kind=int_kind) :: fyear,fmon,fday,fsec + character(len=32) :: calstr,unitstr,signstr + integer (kind=int_kind) :: tdaymo (months_per_year) ! days per month + integer (kind=int_kind) :: tdaycal(months_per_year+1) ! day count per month + integer (kind=int_kind) :: tdayyr ! days in year + + integer(kind=int_kind), parameter :: ntests = 8 + character(len=8) :: errorflag0,errorflag(1:ntests),errorflagtmp + character(len=32) :: testname(ntests) + integer(kind=int_kind) :: yearv(ntests),monv(ntests),dayv(ntests),secv(ntests),ndayv(ntests) ! computed values + integer(kind=int_kind) :: yearc(ntests),monc(ntests),dayc(ntests),secc(ntests),ndayc(ntests) ! correct results + real(kind=dbl_kind) :: timesecsv(ntests),timesecsc(ntests) + character(len=*), parameter :: & + passflag = 'PASS', & + failflag = 'FAIL' + + write(6,*) ' ' + write(6,*) 'Running CALCHK' + write(6,*) ' ' + + errorflag0 = passflag + errorflag(:) = passflag + testname(:) = '' + testname(1) = 'compute_elapsed_days' + testname(2) = 'set_date_from_timesecs' + testname(3) = 'calendar advance' + testname(4) = 'date2time time2date' + testname(5) = 'big add/sub update_date' + testname(6) = 'small add/sub update_date' + testname(7) = 'special checks' + testname(8) = 'calc_timesteps' + + ndtd = 1 + + ! test yearmax years from year 0 + yearmax = 1000 +! yearmax = 100000 + + ! test 3 calendars + do n = 1,3 + + errorflag(:) = passflag + + if (n == 1) then + use_leap_years = .false. + days_per_year = 365 + calstr = 'noleap' + elseif (n == 2) then + use_leap_years = .false. + days_per_year = 360 + calstr = '360day' + elseif (n == 3) then + use_leap_years = .true. + days_per_year = 365 + calstr = 'gregorian' + endif + + istep0 = 1000 + year_init = 0 + month_init = 1 + day_init = 1 + sec_init = 0 + myear = -1 + mmonth = -1 + mday = -1 + dt = 86400._dbl_kind + diagfreq = 99999999 + call init_calendar() + + !----------------- + ! This test makes sure compute_elapsed_days works for different calendars + ! and multiple years. This also checks that the timesecs value computed + ! in calendar and passed into set_date_from_timesecs returns the correct date. + ! In test1, nday should increment 1 day each loop and the final number + ! of days is known for 1000 and 100000 years (precomputed) + ! In test2, set_date_from_timesecs will reset myear, mmonth, mday, msec + !----------------- + + ndayc(1) = -1 ! prior day + do ny = 0,yearmax + do nm = 1,months_per_year + do nd = 1,daymo(nm) + + errorflagtmp = passflag + yearv(1) = ny + monv(1) = nm + dayv(1) = nd + secv(1) = 0 + + ! check days increment by 1 + ndayv(1) = compute_elapsed_days(yearv(1),monv(1),dayv(1)) + if (ndayv(1) - ndayc(1) /= 1) then + errorflagtmp = failflag + errorflag(1) = failflag + write(6,*) 'ERROR1: did not increment one day',yearv(1),monv(1),dayv(1),ndayv(1) + endif + + ! establish internal date and update internal calendar including timesecs + myear = yearv(1) + mmonth = monv(1) + mday = dayv(1) + msec = secv(1) + call calendar() + timesecsv(1) = timesecs + + ! check set_date_from_timesecs + yearc(2) = myear + monc(2) = mmonth + dayc(2) = mday + secc(2) = msec + timesecsc(2) = timesecs + ndayc(2) = ndayv(1) + myear = -1 + mmonth = -1 + mday = -1 + msec = -1 + timesecs = -1 + call set_date_from_timesecs(timesecsc(2)) + if (myear /= yearc(2) .or. mmonth /= monc(2) .or. mday /= dayc(2) .or. msec /= secc(2) .or. timesecs /= timesecsc(2)) then + errorflagtmp = failflag + errorflag(2) = failflag + write(6,*) 'ERROR2: timesecs error' + write(6,1001) 'e2',ndayc(2),yearc(2),'-',monc(2),'-',dayc(2),':',secc(2),' timesecs = ',timesecsc(2) + endif + if (errorflagtmp /= passflag .or. & + ndayv(1) <= 10 .or. mod(ndayv(1),yearmax*10) == 0 .or. & + (yearv(1) == yearmax .and. monv(1) == months_per_year)) then + write(6,1001) ' CHECK1: ',ndayv(1),yearv(1) ,'-',monv(1),'-',dayv(1),':',secv(1) ,' timesecs = ',timesecsv(1) + endif + ndayc(1) = ndayv(1) + enddo + enddo + enddo + + ! check total number of days run in yearmax years + if (yearmax == 1000) then + if (n == 1) then + ndayc(1) = 365364 + elseif (n == 2) then + ndayc(1) = 360359 + elseif (n == 3) then + ndayc(1) = 365607 + endif + if (ndayv(1) /= ndayc(1)) then + errorflag(1) = failflag + write(6,*) 'ERROR1a: final nday incorrect', ndayv(1), ndayc(1) + endif + endif + + ! check total number of days run in yearmax years + if (yearmax == 100000) then + if (n == 1) then + ndayc(1) = 36500364 + elseif (n == 2) then + ndayc(1) = 36000359 + elseif (n == 3) then + ndayc(1) = 36524615 + endif + if (ndayv(1) /= ndayc(1)) then + errorflag(1) = failflag + write(6,*) 'ERROR1a: final nday incorrect', ndayv(1), ndayc(1) + endif + endif + + !----------------- + ! check adding arbitrary amounts to each date unit and see if calendar reconciles properly + ! then subtract same arbitrary amounts in reverse order and make sure it ends at original value + !----------------- + + yearv(1) = 1000 + monv(1) = 1 + dayv(1) = 1 + secv(1) = 0 + myear = yearv(1) + mmonth = monv(1) + mday = dayv(1) + msec = secv(1) + call calendar() + nday = compute_elapsed_days(myear,mmonth,mday) + dyear = 0 + dmon = 0 + dday = 0 + dsec = 0 + do nfa = 1,-1,-2 + write(6,*) ' ' + write(6,1001) ' CHECK3: ',nday,myear ,'-',mmonth ,'-',mday ,':',msec ,' timesecs = ',timesecs + do nfb = 1,10 + do nfc = 1,4 + if (nfa == 1) then + nf1 = nfb + nf2 = nfc + signstr = 'Add' + elseif (nfa == -1) then + nf1 = 11-nfb + nf2 = 5-nfc + signstr = 'Sub' + endif + fyear = 0 + fmon = 0 + fday = 0 + fsec = 0 + if (nf2 == 1) then + xadd = nf1*nf1 + unitstr = 'years' + myear = myear + nfa*xadd + if (nfa == 1) dyear = dyear + nfa*xadd + fyear = nfa*xadd + elseif (nf2 == 2) then + xadd = nf1*nf1 + unitstr = 'months' + mmonth = mmonth + nfa*xadd + if (nfa == 1) dmon = dmon + nfa*xadd + fmon = nfa*xadd + elseif (nf2 == 3) then + xadd = nf1*nf1*nf1*nf1 + unitstr = 'days' + mday = mday + nfa*xadd + if (nfa == 1) dday = dday + nfa*xadd + fday = nfa*xadd + elseif (nf2 == 4) then + xadd = nf1*nf1*nf1*nf1*nf1*nf1*nf1 + unitstr = 'seconds' + msec = msec + nfa*xadd + if (nfa == 1) dsec = dsec + nfa*xadd + fsec = nfa*xadd + endif + call calendar() + nday = compute_elapsed_days(myear,mmonth,mday) + write(6,1002) ' CHECK3: '//trim(signstr)//' ',xadd,trim(unitstr) + write(6,1001) ' CHECK3: ',nday,myear ,'-',mmonth ,'-',mday ,':',msec ,' timesecs = ',timesecs + + !----------------- + ! This checks update_date add and subtract to make sure the original value is returned + !----------------- + + yearc(6) = myear + monc(6) = mmonth + dayc(6) = mday + secc(6) = msec + timesecsc(6) = timesecs + yearv(6) = yearc(6) + monv(6) = monc(6) + dayv(6) = dayc(6) + secv(6) = secc(6) + call update_date(yearv(6),monv(6),dayv(6),secv(6),fyear,fmon,fday,fsec) + write(6,1001) ' CHECK6: ',-1,yearv(6),'-',monv(6),'-',dayv(6),':',secv(6) + if (yearc(6) == yearv(6) .and. monc(6) == monv(6) .and. dayc(6) == dayv(6) .and. secc(6) == secv(6) .and. timesecsc(6) == timesecsv(6)) then + errorflag(6) = failflag + write(6,*) ' ' + write(6,*) 'ERROR6a: update date error' + write(6,1001) 'e6',nday,yearv(6),'-',monv(6),'-',dayv(6),':',secv(6),' timesecs = ',timesecsv(6) + write(6,1001) ' ',nday,yearc(6),'-',monc(6),'-',dayc(6),':',secc(6),' timesecs = ',timesecsc(6) + write(6,*) ' ',fyear,fmon,fday,fsec + write(6,*) ' ' + endif + call update_date(yearv(6),monv(6),dayv(6),secv(6),-fyear,-fmon,-fday,-fsec) + call calendar_date2time(yearc(6),monc(6),dayc(6),secc(6),timesecsv(6)) + if (yearc(6) /= yearv(6) .or. monc(6) /= monv(6) .or. dayc(6) /= dayv(6) .or. secc(6) /= secv(6) .or. timesecsc(6) /= timesecsv(6)) then + errorflag(6) = failflag + write(6,*) ' ' + write(6,*) 'ERROR6b: update date error' + write(6,1001) 'e6',nday,yearv(6),'-',monv(6),'-',dayv(6),':',secv(6),' timesecs = ',timesecsv(6) + write(6,1001) ' ',nday,yearc(6),'-',monc(6),'-',dayc(6),':',secc(6),' timesecs = ',timesecsc(6) + write(6,*) ' ',fyear,fmon,fday,fsec + write(6,*) ' ' + endif + + !----------------- + ! This checks date2time and time2date leveraging the pseudo random dates + ! plus various reference settings. Different reference dates means + ! timesecs won't match, so don't check them. + !----------------- + + yi = myear/2 + mi = max(mmonth/2,1) + di = max(mday*7/8,1) + si = max(msec*7/8,1) + yearc(4) = myear + monc(4) = mmonth + dayc(4) = mday + secc(4) = msec + timesecsc(4) = timesecs + yearv(4) = -1 + monv(4) = -1 + dayv(4) = -1 + secv(4) = -1 + timesecsv(4) = -1 + call calendar_date2time(yearc(4),monc(4),dayc(4),secc(4),timesecsv(4),yi,mi,di,si) + call calendar_time2date(timesecsv(4),yearv(4),monv(4),dayv(4),secv(4),yi,mi,di,si) + write(6,*) 'CHECK4: ',timesecsv(4) + if (yearc(4) /= yearv(4) .or. monc(4) /= monv(4) .or. dayc(4) /= dayv(4) .or. secc(4) /= secv(4)) then + errorflag(4) = failflag + write(6,*) ' ' + write(6,*) 'ERROR4: date2time time2date error' + write(6,1001) 'e4',nday,yearv(4),'-',monv(4),'-',dayv(4),':',secv(4),' timesecs = ',timesecsv(4) + write(6,1001) ' ',nday,yearc(4),'-',monc(4),'-',dayc(4),':',secc(4),' timesecs = ',timesecsc(4) + write(6,*) ' ' + endif + + enddo + enddo + + yearv(3) = myear + monv(3) = mmonth + dayv(3) = mday + secv(3) = msec + timesecsv(3) = timesecs + if (nfa == 1) then + if (n == 1) then + yearc(3) = 1487 + monc(3) = 1 + dayc(3) = 21 + secc(3) = 22825 + ndayc(3) = 542775 + elseif (n == 2) then + yearc(3) = 1488 + monc(3) = 1 + dayc(3) = 13 + secc(3) = 22825 + ndayc(3) = 535692 + elseif (n == 3) then + yearc(3) = 1487 + monc(3) = 1 + dayc(3) = 5 + secc(3) = 22825 + ndayc(3) = 543120 + endif + elseif (nfa == -1) then + yearc(3) = yearv(1) + monc(3) = monv(1) + dayc(3) = dayv(1) + secc(3) = secv(1) + if (n == 1) then + ndayc(3) = 365000 + elseif (n == 2) then + ndayc(3) = 360000 + elseif (n == 3) then + ndayc(3) = 365243 + endif + endif + + ! check answers + if (yearv(3) /= yearc(3) .or. monv(3) /= monc(3) .or. dayv(3) /= dayc(3) .or. secv(3) /= secc(3)) then + errorflag(3) = failflag + write(6,*) ' ' + write(6,*) 'ERROR3: calendar advance error' + write(6,1001) 'e3',nday,yearc(3),'-',monc(3),'-',dayc(3),':',secc(3),' timesecs = ',timesecsc(3) + write(6,1001) ' ',nday,yearv(3),'-',monv(3),'-',dayv(3),':',secv(3),' timesecs = ',timesecsv(3) + write(6,*) ' ' + endif + enddo + + write(6,*) ' ' + yearv(1) = 1000 + monv(1) = 1 + dayv(1) = 1 + secv(1) = 0 + yearv(5) = yearv(1) + monv(5) = monv(1) + dayv(5) = dayv(1) + secv(5) = secv(1) + write(6,1001) ' CHECK5a: ',-1,yearv(5) ,'-',monv(5) ,'-',dayv(5) ,':',secv(5) + write(6,1002) ' Add ',dyear,'years' + write(6,1002) ' Add ',dmon,'months' + write(6,1002) ' Add ',dday,'days' + write(6,1002) ' Add ',dsec,'seconds' + call update_date(yearv(5),monv(5),dayv(5),secv(5),dyear,dmon,dday,dsec) + write(6,1001) ' CHECK5a: ',-1,yearv(5) ,'-',monv(5) ,'-',dayv(5) ,':',secv(5) + write(6,*) ' ' + + ! correct answers + if (n == 1) then + yearc(5) = 1487 + monc(5) = 1 + dayc(5) = 24 + secc(5) = 22825 + ndayc(5) = 542775 + elseif (n == 2) then + yearc(5) = 1488 + monc(5) = 1 + dayc(5) = 13 + secc(5) = 22825 + ndayc(5) = 535692 + elseif (n == 3) then + yearc(5) = 1487 + monc(5) = 1 + dayc(5) = 7 + secc(5) = 22825 + ndayc(5) = 543120 + endif + + ! check answers + if (yearv(5) /= yearc(5) .or. monv(5) /= monc(5) .or. dayv(5) /= dayc(5) .or. secv(5) /= secc(5)) then + errorflag(5) = failflag + write(6,*) ' ' + write(6,*) 'ERROR5a: calendar advance error' + write(6,1001) 'e5',nday,yearc(5),'-',monc(5),'-',dayc(5),':',secc(5),' timesecs = ',timesecs + write(6,1001) ' ',nday,yearv(5),'-',monv(5),'-',dayv(5),':',secv(5),' timesecs = ',timesecs + write(6,*) ' ' + endif + + write(6,1001) ' CHECK5b: ',-1,yearv(5) ,'-',monv(5) ,'-',dayv(5) ,':',secv(5) + write(6,1002) ' Sub ',dyear,'years' + write(6,1002) ' Sub ',dmon,'months' + write(6,1002) ' Sub ',dday,'days' + write(6,1002) ' Sub ',dsec,'seconds' + call update_date(yearv(5),monv(5),dayv(5),secv(5),-dyear,-dmon,-dday,-dsec) + write(6,1001) ' CHECK5b: ',-1,yearv(5) ,'-',monv(5) ,'-',dayv(5) ,':',secv(5) + + ! correct answers + yearc(5) = yearv(1) + monc(5) = monv(1) + dayc(5) = dayv(1) + secc(5) = secv(1) + if (yearv(5) /= yearc(5) .or. monv(5) /= monc(5) .or. dayv(5) /= dayc(5) .or. secv(5) /= secc(5)) then + errorflag(5) = failflag + write(6,*) ' ' + write(6,*) 'ERROR5b: calendar advance error' + write(6,1001) 'e5',nday,yearc(5),'-',monc(5),'-',dayc(5),':',secc(5),' timesecs = ',timesecs + write(6,1001) ' ',nday,yearv(5),'-',monv(5),'-',dayv(5),':',secv(5),' timesecs = ',timesecs + write(6,*) ' ' + endif + + !------------------------- + ! Special checks: + ! Add a month to the last day of each month + ! Check date2time for seconds + !------------------------- + + write(6,*) ' ' + do ny = 1,5 + do nm = 1, months_per_year + if (ny == 1) yearv(7) = 1900 + if (ny == 2) yearv(7) = 1999 + if (ny == 3) yearv(7) = 2000 + if (ny == 4) yearv(7) = 2004 + if (ny == 5) yearv(7) = 2005 + call compute_calendar_data(yearv(7),tdaymo,tdaycal,tdayyr) + monv(7) = nm + dayv(7) = tdaymo(nm) + secv(7) = 0 + if (tdaymo(mod(nm,months_per_year)+1) >= tdaymo(nm)) then + monc(7) = mod(nm,months_per_year)+1 + dayc(7) = dayv(7) + else + monc(7) = mod(nm+1,months_per_year)+1 + dayc(7) = tdaymo(nm) - tdaymo(mod(nm,months_per_year)+1) + endif + yearc(7) = yearv(7) + if (monc(7) < monv(7)) yearc(7) = yearv(7) + 1 + secc(7) = secv(7) + call update_date(yearv(7),monv(7),dayv(7),secv(7),dmon=1) + write(6,1001) ' CHECK7a:',1,yearv(7),'-',monv(7),'-',dayv(7),':',secv(7) + if (yearv(7) /= yearc(7) .or. monv(7) /= monc(7) .or. dayv(7) /= dayc(7) .or. secv(7) /= secc(7)) then + errorflag(7) = failflag + write(6,*) ' ' + write(6,*) 'ERROR7a: add 1 month to end of month error' + write(6,1001) 'e7',-1,yearc(7),'-',monc(7),'-',dayc(7),':',secc(7) + write(6,1001) ' ',-1,yearv(7),'-',monv(7),'-',dayv(7),':',secv(7) + write(6,*) ' ' + endif + enddo + enddo + + do ns1 = 0,seconds_per_day,seconds_per_day/4 + do ns2 = 0,seconds_per_day,seconds_per_day/4 + yearv(7) = 2002 + monv(7) = 3 + call compute_calendar_data(yearv(7),tdaymo,tdaycal,tdayyr) + dayv(7) = tdaymo(monv(7)) + call calendar_date2time(yearv(7),monv(7),dayv(7),ns2,timesecsv(7),yearv(7),monv(7),dayv(7),ns1) + write(6,*) 'CHECK7b:',ns1,ns2,timesecsv(7) + if (timesecsv(7) /= ns2-ns1) then + errorflag(7) = failflag + write(6,*) ' ' + write(6,*) 'ERROR7b: sec diff same date error' + write(6,*) ' ',ns1,ns2,timesecsv(7),ns2-ns1 + write(6,*) ' ' + endif + call calendar_date2time(yearv(7),monv(7)+1,1,ns2,timesecsv(7),yearv(7),monv(7),dayv(7),ns1) + write(6,*) 'CHECK7c:',ns1,ns2,timesecsv(7) + if (timesecsv(7) /= ns2-ns1+seconds_per_day) then + errorflag(7) = failflag + write(6,*) ' ' + write(6,*) 'ERROR7c: sec diff next day error' + write(6,*) ' ',ns1,ns2,timesecsv(7),ns2-ns1+seconds_per_day + write(6,*) ' ' + endif + enddo + enddo + + !------------------------- + ! calc_timesteps + !------------------------- + + myear = 2000 + mmonth = 2 + mday = 1 + msec = 0 + do nf1 = 1,6 + npt = 10 + dt = 3600._dbl_kind + + if (nf1 == 1) then + npt_unit = '1' + nptc = 10 + endif + if (nf1 == 2) then + npt_unit = 's' + npt = 36000. + nptc = 10 + endif + if (nf1 == 3) then + npt_unit = 'h' + nptc = 10 + endif + if (nf1 == 4) then + npt_unit = 'd' + nptc = 240 + endif + if (nf1 == 5) then + npt_unit = 'm' + if (n == 1) nptc = 7272 + if (n == 2) nptc = 7200 + if (n == 3) nptc = 7296 + endif + if (nf1 == 6) then + npt_unit = 'y' + if (n == 1) nptc = 87600 + if (n == 2) nptc = 86400 + if (n == 3) nptc = 87672 + endif + call calc_timesteps() + write(6,*) 'CHECK8:',npt + if (npt /= nptc) then + errorflag(8) = failflag + write(6,*) 'ERROR8: npt error',npt,nptc + endif + enddo + + !------------------------- + ! write test results + !------------------------- + + write(6,*) ' ' + write(6,*) 'Test Results: ',yearmax,' years' + do m = 1,ntests + write(6,*) trim(errorflag(m))," ... ",trim(calstr)," ",trim(testname(m)) + if (errorflag(m) == failflag) errorflag0=failflag + enddo + write(6,*) ' ' + + enddo ! do n + + 1001 format(a,i10,1x,i7.4,a,i2.2,a,i2.2,a,i5.5,a,e23.16) + 1002 format(a,i10,1x,a) + + write(6,*) ' ' + if (errorflag0 == passflag) then + write(6,*) 'CALCHK COMPLETED SUCCESSFULLY' + else + write(6,*) 'CALCHK FAILED' + endif + + end program + diff --git a/cicecore/drivers/unittest/helloworld/helloworld.F90 b/cicecore/drivers/unittest/helloworld/helloworld.F90 new file mode 100644 index 000000000..651436bea --- /dev/null +++ b/cicecore/drivers/unittest/helloworld/helloworld.F90 @@ -0,0 +1,8 @@ + + program hello_world + + write(6,*) 'hello_world' + write(6,*) 'COMPLETED SUCCESSFULLY' + + end program + diff --git a/cicecore/shared/ice_calendar.F90 b/cicecore/shared/ice_calendar.F90 index e7107f42a..622bd8000 100644 --- a/cicecore/shared/ice_calendar.F90 +++ b/cicecore/shared/ice_calendar.F90 @@ -2,7 +2,7 @@ ! Calendar routines for managing time ! -! authors: Elizabeth C. Hunke, LANL +! Authors: Elizabeth C. Hunke, LANL ! Tony Craig, NCAR ! Craig MacLachlan, UK Met Office ! @@ -10,10 +10,12 @@ ! Converted to free form source (F90). ! 2010 CM : Fixed support for Gregorian calendar: subroutines ! sec2time, time2sec and set_calendar added. +! 2020 TC : Significant refactor to move away from time as prognostic module ice_calendar use ice_kinds_mod + use ice_communicate, only: my_task, master_task use ice_constants, only: c0, c1, c100, c30, c360, c365, c3600, & c4, c400 use ice_domain_size, only: max_nstrm @@ -25,78 +27,88 @@ module ice_calendar implicit none private - public :: init_calendar, calendar, time2sec, sec2time, hc_jday + ! INTERFACES - integer (kind=int_kind), public :: & - days_per_year , & ! number of days in one year - daymo(12) , & ! number of days in each month - daycal(13) ! day number at end of month + public :: init_calendar ! initialize calendar + public :: calc_timesteps ! initialize number of timesteps (after namelist and restart are read) + public :: advance_timestep ! advance model 1 timestep and update calendar + public :: calendar ! update model internal calendar/time information + public :: set_date_from_timesecs ! set model date from time in seconds + ! (relative to init date) + ! needed for binary restarts + public :: hc_jday ! converts "calendar" date to HYCOM julian day - ! 360-day year data - integer (kind=int_kind) :: & - daymo360(12) , & ! number of days in each month - daycal360(13) ! day number at end of month - data daymo360 / 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30/ - data daycal360/ 0,30, 60, 90,120,150,180,210,240,270,300,330,360/ + ! semi-private, only used directly by unit tester + public :: compute_elapsed_days ! compute elapsed days since 0000-01-01 + public :: compute_days_between ! compute elapsed days between two dates + public :: update_date ! input date and delta date, compute new date + public :: calendar_date2time ! convert date to time relative to init date + public :: calendar_time2date ! convert time to date relative to init date + public :: compute_calendar_data ! compute info about calendar for a given year - ! 365-day year data - integer (kind=int_kind) :: & - daymo365(12) , & ! number of days in each month - daycal365(13) ! day number at end of month - data daymo365 / 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31/ - data daycal365/ 0,31, 59, 90,120,151,181,212,243,273,304,334,365/ + ! private functions + private :: set_calendar ! sets model calendar type (noleap, etc) - ! 366-day year data (leap year) - integer (kind=int_kind) :: & - daymo366(12) , & ! number of days in each month - daycal366(13) ! day number at end of month - data daymo366 / 31, 29, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31/ - data daycal366/ 0,31, 60, 91,121,152,182,213,244,274,305,335,366/ + ! PUBLIC + + character(len=*), public, parameter :: & + ice_calendar_gregorian = 'Gregorian', & ! calendar name, actually proleptic gregorian here + ice_calendar_noleap = 'NO_LEAP', & ! 365 day per year calendar + ice_calendar_360day = '360day' ! 360 day calendar with 30 days per month - real (kind=dbl_kind), parameter :: & - days_per_4c = 146097.0_dbl_kind, & - days_per_c = 36524.0_dbl_kind, & - days_per_4y = 1461.0_dbl_kind, & - days_per_y = 365.0_dbl_kind + integer (kind=int_kind), public, parameter :: & + months_per_year = 12, & ! months per year + hours_per_day = 24 ! hours per day integer (kind=int_kind), public :: & - istep , & ! local step counter for time loop - istep0 , & ! counter, number of steps taken in previous run + seconds_per_day , & ! seconds per day + days_per_year , & ! number of days in one year + daymo(months_per_year), & ! number of days in each month + daycal(months_per_year+1) ! accumulated days in year to end of prior month + + integer (kind=int_kind), public :: & + ! step counters + istep , & ! local step counter for current run in time loop + istep0 , & ! counter, number of steps at start of run istep1 , & ! counter, number of steps at current timestep + ! basic model time variables + myear , & ! year number + mmonth , & ! month number, 1 to months_per_year mday , & ! day of the month - hour , & ! hour of the day - month , & ! month number, 1 to 12 - monthp , & ! last month + msec , & ! elapsed seconds into date + ! initial time year_init, & ! initial year - nyr , & ! year number + month_init,& ! initial month + day_init, & ! initial day of month + sec_init , & ! initial seconds + ! other stuff idate , & ! date (yyyymmdd) - idate0 , & ! initial date (yyyymmdd) - sec , & ! elapsed seconds into date + idate0 , & ! initial date (yyyymmdd), associated with year_init, month_init, day_init + dayyr , & ! number of days in the current year npt , & ! total number of time steps (dt) + npt0 , & ! original npt value in npt0_unit ndtd , & ! number of dynamics subcycles: dt_dyn=dt/ndtd stop_now , & ! if 1, end program execution write_restart, & ! if 1, write restart now diagfreq , & ! diagnostic output frequency (10 = once per 10 dt) dumpfreq_n , & ! restart output frequency (10 = once per 10 d,m,y) nstreams , & ! number of history output streams - histfreq_n(max_nstrm) ! history output frequency + histfreq_n(max_nstrm) ! history output frequency + + logical (kind=log_kind), public :: & + new_year , & ! new year = .true. + new_month , & ! new month = .true. + new_day , & ! new day = .true. + new_hour ! new hour = .true. real (kind=dbl_kind), public :: & dt , & ! thermodynamics timestep (s) dt_dyn , & ! dynamics/transport/ridging timestep (s) - time , & ! total elapsed time (s) - time_forc , & ! time of last forcing update (s) + timesecs , & ! total elapsed time (s) yday , & ! day of the year - tday , & ! absolute day number - dayyr , & ! number of days per year - nextsw_cday , & ! julian day of next shortwave calculation - basis_seconds ! Seconds since calendar zero + nextsw_cday ! julian day of next shortwave calculation logical (kind=log_kind), public :: & - new_year , & ! new year = .true. - new_month , & ! new month = .true. - new_day , & ! new day = .true. - new_hour , & ! new hour = .true. use_leap_years , & ! use leap year functionality if true write_ic , & ! write initial condition now dump_last , & ! write restart file on last time step @@ -104,6 +116,8 @@ module ice_calendar write_history(max_nstrm) ! write history now character (len=1), public :: & + npt_unit, & ! run length unit, 'y', 'm', 'd', 'h', 's', '1' + npt0_unit, & ! original run length unit, 'y', 'm', 'd', 'h', 's', '1' histfreq(max_nstrm), & ! history output frequency, 'y','m','d','h','1' dumpfreq ! restart frequency, 'y','m','d' @@ -111,17 +125,33 @@ module ice_calendar calendar_type ! differentiates Gregorian from other calendars ! default = ' ' + ! PRIVATE + + integer (kind=int_kind) :: & + hour ! hour of the day + + ! 360-day year data + integer (kind=int_kind) :: & + daymo360(months_per_year) ! number of days in each month + data daymo360 / 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30/ + + ! 365-day year data + integer (kind=int_kind) :: & + daymo365(months_per_year) ! number of days in each month + data daymo365 / 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31/ + + ! 366-day year data (leap year) + integer (kind=int_kind) :: & + daymo366(months_per_year) ! number of days in each month + data daymo366 / 31, 29, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31/ + + !======================================================================= contains !======================================================================= - ! Initialize calendar variables -! -! authors: Elizabeth C. Hunke, LANL -! Tony Craig, NCAR -! Craig MacLachlan, UK Met Office subroutine init_calendar @@ -134,99 +164,172 @@ subroutine init_calendar if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) + seconds_per_day = nint(secday) + if ((abs(real(seconds_per_day,kind=dbl_kind)/secday)-1.0_dbl_kind) > 1.0e-7) then + write(nu_diag,*) trim(subname),' ERROR secday should basically be an integer',secday + call abort_ice(subname//'ERROR: improper secday') + endif + istep = 0 ! local timestep number - time=istep0*dt ! s - yday=c0 ! absolute day number - mday=0 ! day of the month - month=0 ! month - nyr=0 ! year - idate=00000101 ! date - sec=0 ! seconds into date + myear=year_init ! year + mmonth=month_init ! month + mday=day_init ! day of the month + msec=sec_init ! seconds into date istep1 = istep0 ! number of steps at current timestep ! real (dumped) or imagined (use to set calendar) + idate0 = (myear)*10000 + mmonth*100 + mday ! date (yyyymmdd) stop_now = 0 ! end program execution if stop_now=1 dt_dyn = dt/real(ndtd,kind=dbl_kind) ! dynamics et al timestep force_restart_now = .false. - ! Check that the number of days per year is set correctly when using - ! leap years. If not, set days_per_year correctly and warn the user. - if (use_leap_years .and. days_per_year /= 365) then - days_per_year = 365 - write(nu_diag,*) 'Warning: days_per_year has been set to 365', & - ' because use_leap_years = .true.' - end if - #ifdef CESMCOUPLED ! calendar_type set by coupling #else - calendar_type = ' ' - if (use_leap_years .and. days_per_year == 365) calendar_type = 'Gregorian' -#endif - - dayyr = real(days_per_year, kind=dbl_kind) - if (days_per_year == 360) then - daymo = daymo360 - daycal = daycal360 - elseif (days_per_year == 365) then - daymo = daymo365 - daycal = daycal365 - else - call abort_ice(subname//'ERROR: days_per_year must be 360 or 365') + calendar_type = '' + if (use_leap_years) then + if (days_per_year == 365) then + calendar_type = trim(ice_calendar_gregorian) + else + call abort_ice(subname//'ERROR: use_leap_years is true, must set days_per_year to 365') + endif + else + if (days_per_year == 365) then + calendar_type = trim(ice_calendar_noleap) + elseif (days_per_year == 360) then + calendar_type = trim(ice_calendar_360day) + else + call abort_ice(subname//'ERROR: days_per_year only 365 or 360 supported') + endif endif +#endif - ! Get the time in seconds from calendar zero to start of initial year - call time2sec(year_init,1,1,basis_seconds) + call set_calendar(myear) + call calendar() - ! determine initial date (assumes namelist year_init, istep0 unchanged) - sec = mod(time,secday) ! elapsed seconds into date at - ! end of dt - tday = (time-sec)/secday + c1 ! absolute day number + end subroutine init_calendar - ! Convert the current timestep into a calendar date - call sec2time(nyr,month,mday,basis_seconds+sec) +!======================================================================= +! Initialize timestep counter +! This converts npt_unit and npt to a number of timesteps stored in npt +! npt0 and npt0_unit remember the original values +! It is safe to call this more than once, but it should be called only after +! the initial model run date is known (from namelist or restart) and before +! the first timestep - yday = mday + daycal(month) ! day of the year - nyr = nyr - year_init + 1 ! year number + subroutine calc_timesteps - idate0 = (nyr+year_init-1)*10000 + month*100 + mday ! date (yyyymmdd) + real (kind=dbl_kind) :: secday ! seconds per day + real (kind=dbl_kind) :: dtimesecs ! time in seconds of run + integer (kind=int_kind) :: yeare,monthe,daye,sece ! time at end of run + character(len=*),parameter :: subname='(calc_timesteps)' - end subroutine init_calendar + call icepack_query_parameters(secday_out=secday) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) -!======================================================================= + yeare = myear + monthe = mmonth + daye = mday + sece = msec + npt0 = npt + npt0_unit = npt_unit + + if (npt_unit == 'y') then + call update_date(yeare,monthe,daye,sece,dyear=npt) + call calendar_date2time(yeare,monthe,daye,sece,dtimesecs,myear,mmonth,mday,msec) + elseif (npt_unit == 'm') then + call update_date(yeare,monthe,daye,sece,dmon=npt) + call calendar_date2time(yeare,monthe,daye,sece,dtimesecs,myear,mmonth,mday,msec) + elseif (npt_unit == 'd') then + dtimesecs = real(npt,kind=dbl_kind)*secday + call update_date(yeare,monthe,daye,sece,dday=npt) + elseif (npt_unit == 'h') then + dtimesecs = real(npt,kind=dbl_kind)*secday/real(hours_per_day,kind=dbl_kind) + call update_date(yeare,monthe,daye,sece,dsec=nint(dtimesecs)) + elseif (npt_unit == 's') then + call update_date(yeare,monthe,daye,sece,dsec=npt) + dtimesecs = real(npt,kind=dbl_kind) + elseif (npt_unit == '1') then + dtimesecs = dt*real(npt,kind=dbl_kind) + call update_date(yeare,monthe,daye,sece,dsec=nint(dtimesecs)) + else + write(nu_diag,*) trim(subname),' ERROR invalid npt_unit = ',trim(npt_unit) + call abort_ice(subname//'ERROR: invalid npt_unit') + endif + + npt = nint(dtimesecs/dt) + npt_unit = '1' + + if (my_task == master_task) then + write(nu_diag,*) ' ' + write(nu_diag,'(1x,2a,i9,a,f13.2)') subname,' modified npt from ',npt0,' '//trim(npt0_unit)//' with dt= ',dt + write(nu_diag,'(1x,2a,i9,a,f13.2)') subname,' to ',npt ,' '//trim(npt_unit )//' with dt= ',dt + write(nu_diag,'(1x,2a,i6.4,a,i2.2,a,i2.2,a,i5.5)') subname,' start time is',myear,'-',mmonth,'-',mday,':',msec + write(nu_diag,'(1x,2a,i6.4,a,i2.2,a,i2.2,a,i5.5)') subname,' end time is',yeare,'-',monthe,'-',daye,':',sece + write(nu_diag,*) ' ' + endif + + ! check that npt is very close to an integer + if ((abs(real(npt,kind=dbl_kind)*dt/dtimesecs)-1.0_dbl_kind) > 1.0e-7) then + write(nu_diag,*) trim(subname),' ERROR dt and npt not consistent',npt,dt + call abort_ice(subname//'ERROR: improper npt') + endif + end subroutine calc_timesteps + +!======================================================================= ! Determine the date at the end of the time step -! -! authors: Elizabeth C. Hunke, LANL -! Tony Craig, NCAR -! Craig MacLachlan, UK Met Office - subroutine calendar(ttime) + subroutine advance_timestep() - use ice_communicate, only: my_task, master_task + ! local variables - real (kind=dbl_kind), intent(in) :: & - ttime ! time variable + integer(kind=int_kind) :: & + idt ! integer dt + character(len=*),parameter :: subname='(advance_timestep)' + + if (trim(npt_unit) /= '1') then + write(nu_diag,*) trim(subname),' ERROR npt_unit should be converted to timesteps by now ',trim(npt_unit) + write(nu_diag,*) trim(subname),' ERROR you may need to call calc_timesteps to convert from other units' + call abort_ice(subname//'ERROR: npt_unit incorrect') + endif + + istep = istep + 1 + istep1 = istep1 + 1 + idt = nint(dt) + ! dt is historically a real but it should be an integer + ! make sure dt is very close to an integer + if ((abs(real(idt,kind=dbl_kind)/dt)-1.0_dbl_kind) > 1.0e-7) then + write(nu_diag,*) trim(subname),' ERROR dt error, needs to be integer number of seconds, dt=',dt + call abort_ice(subname//'ERROR: improper dt') + endif + msec = msec + idt + call calendar() + + end subroutine advance_timestep + +!======================================================================= +! Update the calendar and time manager info + + subroutine calendar() + +! real (kind=dbl_kind), intent(in), optional :: & +! ttime ! time variable ! local variables integer (kind=int_kind) :: & ns , & ! loop index - nyrp,mdayp,hourp , & ! previous year, day, hour + yearp,monthp,dayp,hourp , & ! previous year, month, day, hour elapsed_days , & ! since beginning this run elapsed_months , & ! since beginning this run - elapsed_hours , & ! since beginning this run - month0 - real (kind=dbl_kind) :: secday ! seconds per day + elapsed_hours ! since beginning this run character(len=*),parameter :: subname='(calendar)' - call icepack_query_parameters(secday_out=secday) - call icepack_warnings_flush(nu_diag) - if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & - file=__FILE__, line=__LINE__) - - nyrp=nyr - monthp=month - mdayp=mday + yearp=myear + monthp=mmonth + dayp=mday hourp=hour new_year=.false. new_month=.false. @@ -235,308 +338,575 @@ subroutine calendar(ttime) write_history(:)=.false. write_restart=0 - sec = mod(ttime,secday) ! elapsed seconds into date at - ! end of dt - tday = (ttime-sec)/secday + c1 ! absolute day number + call update_date(myear,mmonth,mday,msec) + call set_calendar(myear) - ! Deterime the current date from the timestep - call sec2time(nyr,month,mday,basis_seconds+ttime) + idate = (myear)*10000 + mmonth*100 + mday ! date (yyyymmdd) + yday = daycal(mmonth) + mday ! day of the year + elapsed_months = (myear - year_init)*months_per_year + mmonth - month_init + elapsed_days = compute_days_between(year_init,month_init,day_init,myear,mmonth,mday) + elapsed_hours = elapsed_days * hours_per_day + call calendar_date2time(myear,mmonth,mday,msec,timesecs) - yday = mday + daycal(month) ! day of the year - nyr = nyr - year_init + 1 ! year number - - hour = int((ttime)/c3600) + c1 ! hour - - month0 = int((idate0 - int(idate0 / 10000) * 10000) / 100) - - elapsed_months = (nyr - 1)*12 + (month - month0) - elapsed_days = int((istep * dt) / secday) - elapsed_hours = int(ttime/3600) - - idate = (nyr+year_init-1)*10000 + month*100 + mday ! date (yyyymmdd) + !--- compute other stuff #ifndef CESMCOUPLED if (istep >= npt+1) stop_now = 1 if (istep == npt .and. dump_last) write_restart = 1 ! last timestep #endif - if (nyr /= nyrp) new_year = .true. - if (month /= monthp) new_month = .true. - if (mday /= mdayp) new_day = .true. - if (hour /= hourp) new_hour = .true. + if (myear /= yearp) new_year = .true. + if (mmonth /= monthp) new_month = .true. + if (mday /= dayp) new_day = .true. + if (hour /= hourp) new_hour = .true. + ! History writing flags do ns = 1, nstreams - if (histfreq(ns)=='1' .and. histfreq_n(ns)/=0) then - if (mod(istep1, histfreq_n(ns))==0) & - write_history(ns)=.true. - endif + + select case (histfreq(ns)) + case ("y", "Y") + if (new_year .and. histfreq_n(ns)/=0) then + if (mod(myear, histfreq_n(ns))==0) & + write_history(ns) = .true. + endif + case ("m", "M") + if (new_month .and. histfreq_n(ns)/=0) then + if (mod(elapsed_months,histfreq_n(ns))==0) & + write_history(ns) = .true. + endif + case ("d", "D") + if (new_day .and. histfreq_n(ns)/=0) then + if (mod(elapsed_days,histfreq_n(ns))==0) & + write_history(ns) = .true. + endif + case ("h", "H") + if (new_hour .and. histfreq_n(ns)/=0) then + if (mod(elapsed_hours,histfreq_n(ns))==0) & + write_history(ns) = .true. + endif + case ("1") + if (histfreq_n(ns)/=0) then + if (mod(istep1, histfreq_n(ns))==0) & + write_history(ns)=.true. + endif + end select + enddo - if (dumpfreq == '1') then + ! Restart writing flag + + select case (dumpfreq) + case ("y", "Y") + if (new_year .and. mod(myear, dumpfreq_n)==0) & + write_restart = 1 + case ("m", "M") + if (new_month .and. mod(elapsed_months,dumpfreq_n)==0) & + write_restart = 1 + case ("d", "D") + if (new_day .and. mod(elapsed_days, dumpfreq_n)==0) & + write_restart = 1 + case ("h", "H") + if (new_hour .and. mod(elapsed_hours, dumpfreq_n)==0) & + write_restart = 1 + case ("1") if (mod(istep1, dumpfreq_n)==0) & write_restart = 1 - endif + end select - if (istep > 1) then + if (force_restart_now) write_restart = 1 - do ns = 1, nstreams - - select case (histfreq(ns)) - case ("y", "Y") - if (new_year .and. histfreq_n(ns)/=0) then - if (mod(nyr, histfreq_n(ns))==0) & - write_history(ns) = .true. - endif - case ("m", "M") - if (new_month .and. histfreq_n(ns)/=0) then - if (mod(elapsed_months,histfreq_n(ns))==0) & - write_history(ns) = .true. - endif - case ("d", "D") - if (new_day .and. histfreq_n(ns)/=0) then - if (mod(elapsed_days,histfreq_n(ns))==0) & - write_history(ns) = .true. - endif - case ("h", "H") - if (new_hour .and. histfreq_n(ns)/=0) then - if (mod(elapsed_hours,histfreq_n(ns))==0) & - write_history(ns) = .true. - endif - end select - - enddo ! nstreams - - select case (dumpfreq) - case ("y", "Y") - if (new_year .and. mod(nyr, dumpfreq_n)==0) & - write_restart = 1 - case ("m", "M") - if (new_month .and. mod(elapsed_months,dumpfreq_n)==0) & - write_restart = 1 - case ("d", "D") - if (new_day .and. mod(elapsed_days, dumpfreq_n)==0) & - write_restart = 1 - case ("h", "H") - if (new_hour .and. mod(elapsed_hours, dumpfreq_n)==0) & - write_restart = 1 - end select - - if (force_restart_now) write_restart = 1 - - endif ! istep > 1 - - if (my_task == master_task .and. mod(istep,diagfreq) == 0 & + if (my_task == master_task .and. mod(istep1,diagfreq) == 0 & .and. stop_now /= 1) then write(nu_diag,*) ' ' write(nu_diag,'(a7,i10,4x,a6,i10,4x,a4,i10)') & - 'istep1:', istep1, 'idate:', idate, 'sec:', sec + 'istep1:', istep1, 'idate:', idate, 'sec:', msec endif end subroutine calendar !======================================================================= +! Set the model calendar data for year -! Convert the date to seconds since calendar zero. -! ** This is based on the UM routine TIME2SEC ** -! -! authors: Craig MacLachlan, UK Met Office + subroutine set_calendar(year) - subroutine time2sec(year,month,day,tsec) + integer (kind=int_kind), intent(in) :: year ! current year - integer (kind=int_kind), intent(in) :: year ! year - integer (kind=int_kind), intent(in) :: month ! month - integer (kind=int_kind), intent(in) :: day ! year - real (kind=dbl_kind), intent(out) :: tsec ! seconds since calendar zero + ! Internal variable + character(len=*),parameter :: subname='(set_calendar)' - ! local variables + call compute_calendar_data(year,daymo,daycal,dayyr) + + end subroutine set_calendar + +!======================================================================= +! Add and reconcile date +! delta time arguments are optional + + subroutine update_date(ayear,amon,aday,asec,dyear,dmon,dday,dsec) - real (kind=dbl_kind) :: days_since_calz ! days since calendar zero - real (kind=dbl_kind) :: secday ! seconds per day - integer (kind=int_kind) :: years_since_calz ! days since calendar zero - character(len=*),parameter :: subname='(time2sec)' + integer (kind=int_kind), intent(inout) :: ayear, amon, aday, asec ! year, month, day, sec + integer (kind=int_kind), intent(in), optional :: dyear, dmon, dday, dsec ! delta year, month, day, sec + + ! local variables + integer (kind=int_kind) :: tdaymo (months_per_year) ! days per month + integer (kind=int_kind) :: tdaycal(months_per_year+1) ! day count per month + integer (kind=int_kind) :: tdayyr ! days in year + real (kind=dbl_kind) :: secday ! seconds per day + integer (kind=int_kind) :: isecday ! seconds per day + integer (kind=int_kind) :: delta + character(len=*),parameter :: subname='(update_date)' call icepack_query_parameters(secday_out=secday) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) + isecday = nint(secday) + + ! order matters. think about adding 1 month and 10 days to the 25th of a month + ! what is the right order? + ! will add all deltas then reconcile years then months then days then seconds + + if (present(dyear)) ayear = ayear + dyear + if (present(dmon)) amon = amon + dmon + if (present(dday)) aday = aday + dday + if (present(dsec)) asec = asec + dsec + + ! adjust negative data first + ! reconcile months - years + do while (amon <= 0) + delta = int((abs(amon))/months_per_year) + 1 + ayear = ayear - delta + amon = amon + delta*months_per_year + enddo + call compute_calendar_data(ayear,tdaymo,tdaycal,tdayyr) + + ! reconcile days - months - years + do while (aday <= 0) + amon = amon - 1 + do while (amon <= 0) + delta = int((abs(amon))/months_per_year) + 1 + ayear = ayear - delta + amon = amon + delta*months_per_year + call compute_calendar_data(ayear,tdaymo,tdaycal,tdayyr) + enddo + aday = aday + tdaymo(amon) + enddo - if (dayyr == 360) then - days_since_calz = c360*year + c30*(month-1) + day - c1 - tsec = secday * days_since_calz + ! reconcile seconds - days - months - years + if (asec < 0) then + delta = int(abs(asec)/isecday) + 1 + aday = aday - delta + asec = asec + delta*isecday + endif + do while (aday <= 0) + amon = amon - 1 + do while (amon <= 0) + delta = int((abs(amon))/months_per_year) + 1 + ayear = ayear - delta + amon = amon + delta*months_per_year + call compute_calendar_data(ayear,tdaymo,tdaycal,tdayyr) + enddo + aday = aday + tdaymo(amon) + enddo - else - - if (use_leap_years) then + ! check for negative data + if (ayear < 0 .or. amon <= 0 .or. aday <= 0 .or. asec < 0) then + write(nu_diag,*) trim(subname),' ERROR in dateA, ',ayear,amon,aday,asec + call abort_ice(subname//'ERROR: in date') + endif - call set_calendar(year) + ! reconcile months - years + do while (amon > months_per_year) + delta = int((amon-1)/months_per_year) + ayear = ayear + delta + amon = amon - delta*months_per_year + enddo + call compute_calendar_data(ayear,tdaymo,tdaycal,tdayyr) + + ! reconcile days - months - years + do while (aday > tdaymo(amon)) + aday = aday - tdaymo(amon) + amon = amon + 1 + do while (amon > months_per_year) + delta = int((amon-1)/months_per_year) + ayear = ayear + delta + amon = amon - delta*months_per_year + call compute_calendar_data(ayear,tdaymo,tdaycal,tdayyr) + enddo + enddo + + ! reconcile seconds - days - months - years + if (asec >= isecday) then + delta = int(asec/isecday) + aday = aday + delta + asec = asec - delta*isecday + endif + do while (aday > tdaymo(amon)) + aday = aday - tdaymo(amon) + amon = amon + 1 + do while (amon > months_per_year) + delta = int((amon-1)/months_per_year) + ayear = ayear + delta + amon = amon - delta*months_per_year + call compute_calendar_data(ayear,tdaymo,tdaycal,tdayyr) + enddo + enddo - ! Add on the days from this year - days_since_calz = day + daycal(month) - c1 + ! check for negative data, just in case + if (ayear < 0 .or. amon <= 0 .or. aday <= 0 .or. asec < 0) then + write(nu_diag,*) trim(subname),' ERROR in dateB, ',ayear,amon,aday,asec + call abort_ice(subname//'ERROR: in date') + endif - ! Subtract a year because we only want to count whole years - years_since_calz = year - 1 - - ! Add days from preceeding years - days_since_calz = days_since_calz & - + int(years_since_calz/c400)*days_per_4c - years_since_calz = years_since_calz & - - int(years_since_calz/c400)*400 + end subroutine update_date - days_since_calz = days_since_calz & - + int(years_since_calz/c100)*days_per_c - years_since_calz = years_since_calz & - - int(years_since_calz/c100)*100 +!======================================================================= - days_since_calz = days_since_calz & - + int(years_since_calz/c4)*days_per_4y - years_since_calz = years_since_calz & - - int(years_since_calz/c4)*4 +! Set internal calendar date from timesecs input +! Needed for binary restarts where only timesecs is on the restart file - days_since_calz = days_since_calz & - + years_since_calz*days_per_y + subroutine set_date_from_timesecs(ttimesecs) - tsec = secday * days_since_calz - - else ! Using fixed 365-day calendar - - days_since_calz = c365*year + daycal365(month) + day - c1 - tsec = secday * days_since_calz + real (kind=dbl_kind), intent(in) :: ttimesecs ! seconds since init date - end if + ! Internal variable + character(len=*),parameter :: subname='(set_date_from_timesecs)' - end if + timesecs = ttimesecs + call calendar_time2date(ttimesecs,myear,mmonth,mday,msec,year_init,month_init,day_init,sec_init) - end subroutine time2sec + end subroutine set_date_from_timesecs !======================================================================= +! Compute elapsed days from year0,month0,day0 to year1,month1,day1 +! Same day results in 0 elapsed days -! Convert the time in seconds since calendar zero to a date. -! -! authors: Craig MacLachlan, UK Met Office + integer function compute_days_between(year0,month0,day0,year1,month1,day1) - subroutine sec2time(year,month,day,tsec) + integer (kind=int_kind), intent(in) :: year0 ! start year + integer (kind=int_kind), intent(in) :: month0 ! start month + integer (kind=int_kind), intent(in) :: day0 ! start day + integer (kind=int_kind), intent(in) :: year1 ! end year + integer (kind=int_kind), intent(in) :: month1 ! end month + integer (kind=int_kind), intent(in) :: day1 ! end day - integer (kind=int_kind), intent(out) :: year ! year - integer (kind=int_kind), intent(out) :: month ! month - integer (kind=int_kind), intent(out) :: day ! year - real (kind=dbl_kind), intent(in) :: tsec ! seconds since calendar zero + ! Internal variable + logical (kind=log_kind) :: isleap ! Leap year logical + integer (kind=int_kind) :: nday0, nday1 + character(len=*),parameter :: subname='(compute_days_between)' - ! local variables + nday0 = compute_elapsed_days(year0,month0,day0) + nday1 = compute_elapsed_days(year1,month1,day1) - real (kind=dbl_kind) :: days_since_calz ! days since calendar zero - real (kind=dbl_kind) :: secday ! seconds per day - integer (kind=int_kind) :: k ! counter - character(len=*),parameter :: subname='(sec2time)' + compute_days_between = nday1 - nday0 - call icepack_query_parameters(secday_out=secday) - call icepack_warnings_flush(nu_diag) - if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & - file=__FILE__, line=__LINE__) + end function compute_days_between + +!======================================================================= +! compute calendar data based on year + + subroutine compute_calendar_data(ayear,adaymo,adaycal,adayyr) + + integer (kind=int_kind), intent(in) :: ayear ! year + integer (kind=int_kind), intent(out) :: adaymo(:) ! days per month + integer (kind=int_kind), intent(out) :: adaycal(:) ! day count per month + integer (kind=int_kind), intent(out) :: adayyr ! days per year + + ! Internal variable + logical (kind=log_kind) :: isleap ! Leap year logical + integer (kind=int_kind) :: n + character(len=*),parameter :: subname='(compute_calendar_data)' + + if (ayear < 0) then + write(nu_diag,*) trim(subname),' ERROR in ayear = ',ayear + call abort_ice(subname//'ERROR: in ayear') + endif - days_since_calz = int(tsec/secday) + if (size(adaymo) /= months_per_year .or. & + size(adaycal) /= months_per_year+1 ) then + call abort_ice(subname//'ERROR: in argument sizes') + endif + + if (trim(calendar_type) == trim(ice_calendar_gregorian)) then - if (dayyr == 360) then + isleap = .false. ! not a leap year + if (mod(ayear, 4) == 0) isleap = .true. + if (mod(ayear,100) == 0) isleap = .false. + if (mod(ayear,400) == 0) isleap = .true. - year = int(days_since_calz/c360) - month = mod(int(days_since_calz/c30),12) + 1 - day = mod(int(days_since_calz),30) + 1 + if (isleap) then + adaymo = daymo366 + else + adaymo = daymo365 + endif + elseif (trim(calendar_type) == trim(ice_calendar_360day)) then + adaymo = daymo360 else + adaymo = daymo365 + endif - if (use_leap_years) then - - year = int(days_since_calz/days_per_4c)*400 - days_since_calz = days_since_calz & - - int(days_since_calz/days_per_4c)*days_per_4c - - if (days_since_calz == 4*days_per_c) then - year = year + 400 - days_since_calz = days_per_y + 1 - else - year = year + int(days_since_calz/days_per_c)*100 - days_since_calz = days_since_calz & - - int(days_since_calz/days_per_c)*days_per_c - - year = year + int(days_since_calz/days_per_4y)*4 - days_since_calz = days_since_calz & - - int(days_since_calz/days_per_4y)*days_per_4y - - if (days_since_calz == 4*days_per_y) then - year = year + 4 - days_since_calz = days_per_y + 1 - else - year = year + int(days_since_calz/days_per_y) + 1 - days_since_calz = days_since_calz & - - int(days_since_calz/days_per_y)*days_per_y + c1 - endif - endif + adaycal(1) = 0 + do n = 1, months_per_year + adaycal(n+1) = adaycal(n) + adaymo(n) + enddo + adayyr=adaycal(months_per_year+1) - ! Ensure the calendar variables are correct for this year. - call set_calendar(year) + end subroutine compute_calendar_data - ! Calculate the month - month = 1 - do k = 1, 12 - if (days_since_calz > daycal(k)) month = k - enddo +!======================================================================= +! Compute elapsed days from 0000-01-01 to year1,month1,day1 +! 0000-01-01 is 0 elapsed days - ! Calculate the day of the month - day = days_since_calz - daycal(month) - - else ! Using fixed 365-day calendar - - year = int(days_since_calz/c365) - days_since_calz = days_since_calz - year*365 + 1 - - ! Calculate the month - month = 1 - do k = 1, 12 - if (days_since_calz > daycal365(k)) month = k - enddo + integer function compute_elapsed_days(ayear,amonth,aday) - ! Calculate the day of the month - day = days_since_calz - daycal365(month) + integer (kind=int_kind), intent(in) :: ayear ! year + integer (kind=int_kind), intent(in) :: amonth ! month + integer (kind=int_kind), intent(in) :: aday ! day - end if + ! Internal variable + integer (kind=int_kind) :: ced_nday, n + integer (kind=int_kind) :: lyear,lmonth,lday,lsec + integer (kind=int_kind) :: tdaymo (months_per_year) ! days per month + integer (kind=int_kind) :: tdaycal(months_per_year+1) ! day count per month + integer (kind=int_kind) :: tdayyr ! days in year + character(len=*),parameter :: subname='(compute_elapsed_days)' + + ! use 0000-01-01 as base, year 0 is a leap year + ! this must be implemented consistent with set_calendar + + lyear = ayear + lmonth = amonth + lday = aday + lsec = 0 + + if (lyear < 0 .or. lmonth <= 0 .or. lday <= 0) then + write(nu_diag,*) trim(subname),' ERROR for year,month,day = ',lyear,lmonth,lday + call abort_ice(subname//'ERROR: illegal date') + elseif (lmonth > months_per_year) then + call update_date(lyear,lmonth,lday,lsec) + endif - end if + ! compute days from year 0000-01-01 to year-01-01 + ! don't loop thru years for performance reasons + if (trim(calendar_type) == trim(ice_calendar_gregorian)) then + if (lyear == 0) then + ced_nday = 0 + else + ced_nday = lyear * 365 + 1 + (lyear-1)/4 - (lyear-1)/100 + (lyear-1)/400 + endif + else + ced_nday = lyear * daycal(months_per_year+1) + endif - end subroutine sec2time + ! now compute days in this year + call compute_calendar_data(lyear,tdaymo,tdaycal,tdayyr) + + do n = 1, lmonth-1 + ced_nday = ced_nday + tdaymo(n) + enddo + + if (lday <= tdaymo(lmonth)) then + ced_nday = ced_nday + lday - 1 + else + write(nu_diag,*) trim(subname),' ERROR for year,month,day = ',ayear,amonth,aday + call abort_ice(subname//'ERROR: illegal day in month') + endif + + compute_elapsed_days = ced_nday + + end function compute_elapsed_days !======================================================================= +! Compute time in seconds from input calendar date +! relative to year_init, month_init, day_init, sec_init unless _ref values passed in +! For santity, must pass all four ref values or none -! Set the "days per month", "days per year", etc variables for the -! current year. -! -! authors: Craig MacLachlan, UK Met Office + subroutine calendar_date2time(ayear,amon,aday,asec,atimesecs,year_ref,mon_ref,day_ref,sec_ref) - subroutine set_calendar(year) + integer(kind=int_kind), intent(in) :: & + ayear,amon,aday,asec ! year, month, day, sec of ttimesecs + real (kind=dbl_kind), intent(out) :: atimesecs ! seconds since init date + integer(kind=int_kind), intent(in), optional :: & + year_ref,mon_ref,day_ref,sec_ref ! year, month, day, sec reference time - integer (kind=int_kind), intent(in) :: year ! current year + ! Internal variable + real (kind=dbl_kind) :: secday + integer (kind=int_kind) :: elapsed_days ! since beginning this run + integer (kind=int_kind) :: lyear_ref,lmon_ref,lday_ref,lsec_ref ! local reference year, month, day, sec + integer (kind=int_kind) :: cnt + character(len=*),parameter :: subname='(calendar_date2time)' + + ! set reference date and check that 0 or 4 optional arguments are passed + cnt = 0 + if (present(year_ref)) then + lyear_ref = year_ref + cnt = cnt + 1 + else + lyear_ref = year_init + endif + if (present(mon_ref)) then + lmon_ref = mon_ref + cnt = cnt + 1 + else + lmon_ref = month_init + endif + if (present(day_ref)) then + lday_ref = day_ref + cnt = cnt + 1 + else + lday_ref = day_init + endif + if (present(sec_ref)) then + lsec_ref = sec_ref + cnt = cnt + 1 + else + lsec_ref = sec_init + endif + if (cnt /= 0 .and. cnt /= 4) then + write(nu_diag,*) trim(subname),' ERROR in ref args, must pass 0 or 4 ' + call abort_ice(subname//'ERROR: in ref args, must pass 0 or 4') + endif + + call icepack_query_parameters(secday_out=secday) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + elapsed_days = compute_days_between(lyear_ref,lmon_ref,lday_ref,ayear,amon,aday) + atimesecs = real(elapsed_days,kind=dbl_kind)*secday + & + real(asec,kind=dbl_kind) - real(lsec_ref,kind=dbl_kind) + + end subroutine calendar_date2time + +!======================================================================= +! Compute calendar date from input time in seconds +! relative to year_init, month_init, day_init, sec_init or ref data if passed. +! For sanity, require all four or no ref values. +! Implemented to minimize accumulating errors and avoid overflows +! and perform well. + + subroutine calendar_time2date(atimesecs,ayear,amon,aday,asec,year_ref,mon_ref,day_ref,sec_ref) + + real (kind=dbl_kind), intent(in) :: atimesecs ! seconds since init date + integer(kind=int_kind), intent(out) :: & + ayear,amon,aday,asec ! year, month, day, sec of timesecs + integer(kind=int_kind), intent(in), optional :: & + year_ref,mon_ref,day_ref,sec_ref ! year, month, day, sec reference time ! Internal variable - logical (kind=log_kind) :: isleap ! Leap year logical - character(len=*),parameter :: subname='(set_calendar)' + integer (kind=int_kind) :: ndays + integer (kind=int_kind) :: tyear, tmon, tday, tsec ! temporaries + integer (kind=int_kind) :: tdaymo (months_per_year) ! days per month + integer (kind=int_kind) :: tdaycal(months_per_year+1) ! day count per month + integer (kind=int_kind) :: tdayyr ! days in year + real (kind=dbl_kind) :: secday, rdays, ltimesecs + integer (kind=int_kind) :: lyear_ref,lmon_ref,lday_ref,lsec_ref ! local reference year, month, day, sec + integer (kind=int_kind) :: cnt + character(len=*),parameter :: subname='(calendar_time2date)' + + call icepack_query_parameters(secday_out=secday) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) - isleap = .false. ! not a leap year - if (mod(year, 4) == 0) isleap = .true. - if (mod(year,100) == 0) isleap = .false. - if (mod(year,400) == 0) isleap = .true. - - ! Ensure the calendar is set correctly - if (isleap) then - daycal = daycal366 - daymo = daymo366 - dayyr=real(daycal(13), kind=dbl_kind) - days_per_year=int(dayyr) + ! we could allow negative atimesecs, but this shouldn't be needed + if (atimesecs < 0._dbl_kind) then + write(nu_diag,*) trim(subname),' ERROR in atimesecs ',atimesecs + call abort_ice(subname//'ERROR: in atimesecs') + endif + + ! set reference date and check that 0 or 4 optional arguments are passed + cnt = 0 + if (present(year_ref)) then + lyear_ref = year_ref + cnt = cnt + 1 + else + lyear_ref = year_init + endif + if (present(mon_ref)) then + lmon_ref = mon_ref + cnt = cnt + 1 else - daycal = daycal365 - daymo = daymo365 - dayyr=real(daycal(13), kind=dbl_kind) - days_per_year=int(dayyr) + lmon_ref = month_init + endif + if (present(day_ref)) then + lday_ref = day_ref + cnt = cnt + 1 + else + lday_ref = day_init + endif + if (present(sec_ref)) then + lsec_ref = sec_ref + cnt = cnt + 1 + else + lsec_ref = sec_init + endif + if (cnt /= 0 .and. cnt /= 4) then + write(nu_diag,*) trim(subname),' ERROR in ref args, must pass 0 or 4 ' + call abort_ice(subname//'ERROR: in ref args, must pass 0 or 4') + endif + +! ------------------------------------------------------------------- +! tcraig, this is risky because atimesecs is real and could be very large +! ayear = lyear_ref +! amon = lmon_ref +! aday = lday_ref +! asec = lsec_ref +! +! call update_date(ayear,amon,aday,asec,dsec=nint(atimesecs)) +! return +! ------------------------------------------------------------------- + + ! initial guess + tyear = lyear_ref + tmon = 1 + tday = 1 + tsec = 0 + + ! add initial seconds to timesecs and treat lsec_ref as zero + ltimesecs = atimesecs + real(lsec_ref,kind=dbl_kind) + + ! first estimate of tyear + call compute_calendar_data(tyear,tdaymo,tdaycal,tdayyr) + rdays = ltimesecs/secday + tyear = tyear + int(rdays)/tdayyr + + ! reduce estimate of tyear if ndays > rdays + ndays = compute_days_between(lyear_ref,lmon_ref,lday_ref,tyear,tmon,tday) + if (ndays > int(rdays)) then + tyear = tyear - (ndays - int(rdays))/tdayyr - 1 + ndays = compute_days_between(lyear_ref,lmon_ref,lday_ref,tyear,tmon,tday) + endif + call compute_calendar_data(tyear,tdaymo,tdaycal,tdayyr) + + ! compute residual days, switch to integers, compute date + rdays = ltimesecs/secday + tday = int(rdays) - ndays + 1 + + do while (tday > tdaymo(tmon)) + tday = tday - tdaymo(tmon) + tmon = tmon + 1 + do while (tmon > months_per_year) + tmon = tmon - months_per_year + tyear = tyear + 1 + call compute_calendar_data(tyear,tdaymo,tdaycal,tdayyr) + enddo + enddo + + ndays = compute_days_between(lyear_ref,lmon_ref,lday_ref,tyear,tmon,tday) + tsec = int(ltimesecs - real(ndays,kind=dbl_kind)*secday) + if (tsec > secday) then + write(nu_diag,*) trim(subname),' ERROR in seconds, ',tyear,tmon,tday,tsec + call abort_ice(subname//'ERROR: in seconds') endif - end subroutine set_calendar + ayear = tyear + amon = tmon + aday = tday + asec = tsec + + end subroutine calendar_time2date !======================================================================= @@ -551,16 +921,22 @@ real(kind=dbl_kind) function hc_jday(iyear,imm,idd,ihour) real(kind=dbl_kind) :: dtime integer(kind=int_kind) :: iyear,iyr,imm,idd,idoy,ihr integer(kind=int_kind), optional :: ihour + integer (kind=int_kind) :: n if (present(ihour)) then !----------------- ! yyyy mm dd HH !----------------- iyr=iyear-1901 + dtime = floor(365.25_dbl_kind*iyr)*c1 + idd*c1 + ihour/24._dbl_kind if (mod(iyr,4)==3) then - dtime = floor(365.25_dbl_kind*iyr)*c1 + daycal366(imm)*c1 + idd*c1 + ihour/24._dbl_kind + do n = 1,imm-1 + dtime = dtime + daymo366(n) + enddo else - dtime = floor(365.25_dbl_kind*iyr)*c1 + daycal365(imm)*c1 + idd*c1 + ihour/24._dbl_kind + do n = 1,imm-1 + dtime = dtime + daymo365(n) + enddo endif else diff --git a/cicecore/shared/ice_init_column.F90 b/cicecore/shared/ice_init_column.F90 index b3937c0cd..e9e668a9d 100644 --- a/cicecore/shared/ice_init_column.F90 +++ b/cicecore/shared/ice_init_column.F90 @@ -188,7 +188,7 @@ subroutine init_shortwave swgrid, igrid use ice_blocks, only: block, get_block, nx_block, ny_block use ice_calendar, only: dt, calendar_type, & - days_per_year, nextsw_cday, yday, sec + days_per_year, nextsw_cday, yday, msec use ice_diagnostics, only: npnt, print_points, pmloc, piloc, pjloc use ice_domain, only: nblocks, blocks_ice use ice_flux, only: alvdf, alidf, alvdr, alidr, & @@ -356,7 +356,7 @@ subroutine init_shortwave calendar_type=calendar_type, & days_per_year=days_per_year, & nextsw_cday=nextsw_cday, yday=yday, & - sec=sec, & + sec=msec, & kaer_tab=kaer_tab, kaer_bc_tab=kaer_bc_tab(:,:), & waer_tab=waer_tab, waer_bc_tab=waer_bc_tab(:,:), & gaer_tab=gaer_tab, gaer_bc_tab=gaer_bc_tab(:,:), & diff --git a/configuration/scripts/Makefile b/configuration/scripts/Makefile index 7b39d5c8d..e0b7799d6 100644 --- a/configuration/scripts/Makefile +++ b/configuration/scripts/Makefile @@ -75,7 +75,7 @@ AR := ar .SUFFIXES: .SUFFIXES: .F90 .F .c .o -.PHONY: all cice libcice targets target db_files db_flags clean realclean +.PHONY: all cice libcice targets target db_files db_flags clean realclean helloworld calchk all: $(EXEC) cice: $(EXEC) @@ -92,7 +92,9 @@ cice: $(EXEC) targets: @echo " " - @echo "Supported Makefile Targets are: cice, libcice, makdep, depends, clean, realclean, targets, db_files, db_flags" + @echo "Supported Makefile Targets are: cice, libcice, makdep, depends, clean, realclean" + @echo " Diagnostics: targets, db_files, db_flags" + @echo " Unit Tests : helloworld, calchk" target: targets db_files: @@ -134,6 +136,20 @@ $(DEPGEN): $(OBJS_DEPGEN) @ echo "Building makdep" $(SCC) -o $@ $(CFLAGS_HOST) $< +#------------------------------------------------------------------------------- +# unit tests +#------------------------------------------------------------------------------- + +# this builds all dependent source code automatically even though only a subset might actually be used +# this is no different than the cice target and in fact the binary is called cice +# it exists just to create separation as needed for unit tests +calchk: $(EXEC) + +# this builds just a subset of source code specified explicitly and requires a separate target +HWOBJS := helloworld.o +helloworld: $(HWOBJS) + $(LD) -o $(EXEC) $(LDFLAGS) $(HWOBJS) $(ULIBS) $(SLIBS) + #------------------------------------------------------------------------------- # build rules: MACFILE, cmd-line, or env vars must provide the needed macros #------------------------------------------------------------------------------- diff --git a/configuration/scripts/cice.build b/configuration/scripts/cice.build index b9aed44fe..d75d74253 100755 --- a/configuration/scripts/cice.build +++ b/configuration/scripts/cice.build @@ -142,6 +142,10 @@ if !($?ICE_MACHINE_BLDTHRDS) then set ICE_MACHINE_BLDTHRDS = 1 endif +if (${directmake} == 0) then + set target = ${ICE_TARGET} +endif + if (${directmake} == 1) then echo "make ${target}" ${ICE_MACHINE_MAKE} -j ${ICE_MACHINE_BLDTHRDS} VPFILE=Filepath EXEC=${ICE_RUNDIR}/cice \ @@ -185,12 +189,12 @@ if (${quiet} == "true") then echo " quiet mode on... patience" ${ICE_MACHINE_MAKE} -j ${ICE_MACHINE_BLDTHRDS} VPFILE=Filepath EXEC=${ICE_RUNDIR}/cice \ -f ${ICE_CASEDIR}/Makefile MACFILE=${ICE_CASEDIR}/Macros.${ICE_MACHCOMP} \ - DEPFILE=${ICE_CASEDIR}/makdep.c cice >& ${ICE_BLDLOG_FILE} + DEPFILE=${ICE_CASEDIR}/makdep.c ${target} >& ${ICE_BLDLOG_FILE} set bldstat = ${status} else ${ICE_MACHINE_MAKE} -j ${ICE_MACHINE_BLDTHRDS} VPFILE=Filepath EXEC=${ICE_RUNDIR}/cice \ -f ${ICE_CASEDIR}/Makefile MACFILE=${ICE_CASEDIR}/Macros.${ICE_MACHCOMP} \ - DEPFILE=${ICE_CASEDIR}/makdep.c cice |& tee ${ICE_BLDLOG_FILE} + DEPFILE=${ICE_CASEDIR}/makdep.c ${target} |& tee ${ICE_BLDLOG_FILE} set bldstat = ${status} endif diff --git a/configuration/scripts/cice.run.setup.csh b/configuration/scripts/cice.run.setup.csh index 901671a36..ea8efeb03 100755 --- a/configuration/scripts/cice.run.setup.csh +++ b/configuration/scripts/cice.run.setup.csh @@ -95,9 +95,15 @@ if ( \$status == 0 ) then echo "CICE run completed successfully" echo "\`date\` \${0}: CICE run completed successfully" >> \${ICE_CASEDIR}/README.case else - echo "CICE run did NOT complete" - echo "\`date\` \${0}: CICE run did NOT complete" >> \${ICE_CASEDIR}/README.case - exit -1 + grep 'COMPLETED SUCCESSFULLY' \${checkfile} + if ( \$status == 0 ) then + echo "Run completed successfully" + echo "\`date\` \${0}: Run completed successfully" >> \${ICE_CASEDIR}/README.case + else + echo "CICE run did NOT complete" + echo "\`date\` \${0}: CICE run did NOT complete" >> \${ICE_CASEDIR}/README.case + exit -1 + endif endif if ( \${diagtype} == 0) then diff --git a/configuration/scripts/cice.settings b/configuration/scripts/cice.settings index d2653a29d..3bd85f5f9 100755 --- a/configuration/scripts/cice.settings +++ b/configuration/scripts/cice.settings @@ -13,6 +13,7 @@ setenv ICE_RSTDIR ${ICE_RUNDIR}/restart setenv ICE_HSTDIR ${ICE_RUNDIR}/history setenv ICE_LOGDIR ${ICE_CASEDIR}/logs setenv ICE_DRVOPT standalone/cice +setenv ICE_TARGET cice setenv ICE_IOTYPE netcdf # binary, netcdf, pio1, pio2 setenv ICE_CLEANBUILD true setenv ICE_CPPDEFS "" diff --git a/configuration/scripts/cice_decomp.csh b/configuration/scripts/cice_decomp.csh index b20f8d129..d481da854 100755 --- a/configuration/scripts/cice_decomp.csh +++ b/configuration/scripts/cice_decomp.csh @@ -98,6 +98,12 @@ else if (${grid} == 'tx1') then set blckx = 10; set blcky = 10 endif +# this is for unit testing +else if (${grid} == 'none') then + set nxglob = 1 + set nyglob = 1 + set blckx = 1; set blcky = 1 + else echo "${0:t}: ERROR unknown grid ${grid}" exit -9 diff --git a/configuration/scripts/ice_in b/configuration/scripts/ice_in index b1fa561a8..f66cb4e69 100644 --- a/configuration/scripts/ice_in +++ b/configuration/scripts/ice_in @@ -4,6 +4,7 @@ year_init = 1997 istep0 = 0 dt = 3600.0 + npt_unit = '1' npt = 24 ndtd = 1 runtype = 'initial' diff --git a/configuration/scripts/options/set_env.calchk b/configuration/scripts/options/set_env.calchk new file mode 100644 index 000000000..7dfe9612e --- /dev/null +++ b/configuration/scripts/options/set_env.calchk @@ -0,0 +1,2 @@ +setenv ICE_DRVOPT unittest/calchk +setenv ICE_TARGET calchk diff --git a/configuration/scripts/options/set_env.helloworld b/configuration/scripts/options/set_env.helloworld new file mode 100644 index 000000000..60587fb91 --- /dev/null +++ b/configuration/scripts/options/set_env.helloworld @@ -0,0 +1,2 @@ +setenv ICE_DRVOPT unittest/helloworld +setenv ICE_TARGET helloworld diff --git a/configuration/scripts/options/set_nml.run10day b/configuration/scripts/options/set_nml.run10day index deae3e993..05160c42d 100644 --- a/configuration/scripts/options/set_nml.run10day +++ b/configuration/scripts/options/set_nml.run10day @@ -1,4 +1,5 @@ -npt = 240 +npt_unit = 'd' +npt = 10 dumpfreq = 'd' dumpfreq_n = 10 histfreq = 'd','x','x','x','x' diff --git a/configuration/scripts/options/set_nml.run1day b/configuration/scripts/options/set_nml.run1day index d7b70f973..a4ed751d5 100644 --- a/configuration/scripts/options/set_nml.run1day +++ b/configuration/scripts/options/set_nml.run1day @@ -1,4 +1,5 @@ -npt = 24 +npt_unit = 'd' +npt = 1 dumpfreq = 'd' dumpfreq_n = 1 diag_type = 'stdout' diff --git a/configuration/scripts/options/set_nml.run1year b/configuration/scripts/options/set_nml.run1year index 9a5baadfd..4e481870c 100644 --- a/configuration/scripts/options/set_nml.run1year +++ b/configuration/scripts/options/set_nml.run1year @@ -1,4 +1,5 @@ -npt = 8760 +npt_unit = 'y' +npt = 1 dumpfreq = 'm' dumpfreq_n = 12 diagfreq = 24 diff --git a/configuration/scripts/options/set_nml.run2day b/configuration/scripts/options/set_nml.run2day index 8129d59f6..60ece02f0 100644 --- a/configuration/scripts/options/set_nml.run2day +++ b/configuration/scripts/options/set_nml.run2day @@ -1,4 +1,5 @@ -npt = 48 +npt_unit = 'd' +npt = 2 dumpfreq = 'd' dumpfreq_n = 2 histfreq = 'd','x','x','x','x' diff --git a/configuration/scripts/options/set_nml.run3day b/configuration/scripts/options/set_nml.run3day index 1fbf7a115..1a839468e 100644 --- a/configuration/scripts/options/set_nml.run3day +++ b/configuration/scripts/options/set_nml.run3day @@ -1,4 +1,5 @@ -npt = 72 +npt_unit = 'd' +npt = 3 dumpfreq = 'd' dumpfreq_n = 2 diag_type = 'stdout' diff --git a/configuration/scripts/options/set_nml.run3dt b/configuration/scripts/options/set_nml.run3dt index 102a19d80..4ff27ce22 100644 --- a/configuration/scripts/options/set_nml.run3dt +++ b/configuration/scripts/options/set_nml.run3dt @@ -1,3 +1,4 @@ +npt_unit = '1' npt = 3 dump_last = .true. histfreq = '1','x','x','x','x' diff --git a/configuration/scripts/options/set_nml.run5day b/configuration/scripts/options/set_nml.run5day index 4113c48e6..88d498a89 100644 --- a/configuration/scripts/options/set_nml.run5day +++ b/configuration/scripts/options/set_nml.run5day @@ -1,4 +1,5 @@ -npt = 120 +npt_unit = 'd' +npt = 5 dumpfreq = 'd' dumpfreq_n = 5 histfreq = 'd','x','x','x','x' diff --git a/configuration/scripts/options/set_nml.run60day b/configuration/scripts/options/set_nml.run60day index 01fd59504..96f6dea1c 100644 --- a/configuration/scripts/options/set_nml.run60day +++ b/configuration/scripts/options/set_nml.run60day @@ -1,4 +1,5 @@ -npt = 1440 +npt_unit = 'd' +npt = 60 dumpfreq = 'd' dumpfreq_n = 30 histfreq = 'd','x','x','x','x' diff --git a/configuration/scripts/options/set_nml.run90day b/configuration/scripts/options/set_nml.run90day index 06db1a3d8..34d31166f 100644 --- a/configuration/scripts/options/set_nml.run90day +++ b/configuration/scripts/options/set_nml.run90day @@ -1,4 +1,5 @@ -npt = 2160 +npt_unit = 'd' +npt = 90 dumpfreq = 'd' dumpfreq_n = 30 histfreq = 'd','x','x','x','x' diff --git a/configuration/scripts/options/test_nml.restart1 b/configuration/scripts/options/test_nml.restart1 index 82f934720..6ab0bd88b 100644 --- a/configuration/scripts/options/test_nml.restart1 +++ b/configuration/scripts/options/test_nml.restart1 @@ -1,4 +1,5 @@ -npt = 240 +npt = 10 +npt_unit = 'd' dumpfreq = 'd' dumpfreq_n = 5 runtype = 'initial' diff --git a/configuration/scripts/options/test_nml.restart2 b/configuration/scripts/options/test_nml.restart2 index 4ae10c5a6..c12887eb0 100644 --- a/configuration/scripts/options/test_nml.restart2 +++ b/configuration/scripts/options/test_nml.restart2 @@ -1,4 +1,5 @@ -npt = 120 +npt = 5 +npt_unit = 'd' dumpfreq = 'd' dumpfreq_n = 5 runtype = 'continue' diff --git a/configuration/scripts/tests/base_suite.ts b/configuration/scripts/tests/base_suite.ts old mode 100755 new mode 100644 diff --git a/configuration/scripts/tests/baseline.script b/configuration/scripts/tests/baseline.script index a1ab4e055..9f2722785 100644 --- a/configuration/scripts/tests/baseline.script +++ b/configuration/scripts/tests/baseline.script @@ -20,38 +20,88 @@ endif # Baseline comparing run if (${ICE_BASECOM} != ${ICE_SPVAL}) then - set test_dir = ${ICE_RUNDIR}/restart - set base_dir = ${ICE_BASELINE}/${ICE_BASECOM}/${ICE_TESTNAME}/restart - - set baseline_log = `ls -1t ${ICE_BASELINE}/${ICE_BASECOM}/${ICE_TESTNAME}/cice.runlog* | head -1` set btimeloop = -1 set bdynamics = -1 set bcolumn = -1 - if (${baseline_log} != "" ) then - set btimeloop = `grep TimeLoop ${baseline_log} | grep Timer | cut -c 22-32` - set bdynamics = `grep Dynamics ${baseline_log} | grep Timer | cut -c 22-32` - set bcolumn = `grep Column ${baseline_log} | grep Timer | cut -c 22-32` - if (${btimeloop} == "") set btimeloop = -1 - if (${bdynamics} == "") set bdynamics = -1 - if (${bcolumn} == "") set bcolumn = -1 - endif - echo "" - echo "Regression Compare Mode:" - echo "base_dir: ${base_dir}" - echo "test_dir: ${test_dir}" + if (${ICE_TEST} == "unittest") then + set test_file = `ls -1t ${ICE_RUNDIR}/cice.runlog* | head -1` + set base_file = `ls -1t ${ICE_BASELINE}/${ICE_BASECOM}/${ICE_TESTNAME}/cice.runlog* | head -1` + + echo "" + echo "bfb Log Compare Mode:" + echo "base_file: ${base_file}" + echo "test_file: ${test_file}" + + ${ICE_CASEDIR}/casescripts/comparelog.csh ${base_file} ${test_file} notcicefile + set bfbstatus = $status + + else + + set test_dir = ${ICE_RUNDIR}/restart + set base_dir = ${ICE_BASELINE}/${ICE_BASECOM}/${ICE_TESTNAME}/restart + + set baseline_log = `ls -1t ${ICE_BASELINE}/${ICE_BASECOM}/${ICE_TESTNAME}/cice.runlog* | head -1` + if (${baseline_log} != "" ) then + set btimeloop = `grep TimeLoop ${baseline_log} | grep Timer | cut -c 22-32` + set bdynamics = `grep Dynamics ${baseline_log} | grep Timer | cut -c 22-32` + set bcolumn = `grep Column ${baseline_log} | grep Timer | cut -c 22-32` + if (${btimeloop} == "") set btimeloop = -1 + if (${bdynamics} == "") set bdynamics = -1 + if (${bcolumn} == "") set bcolumn = -1 + endif + + echo "" + echo "Regression Compare Mode:" + echo "base_dir: ${base_dir}" + echo "test_dir: ${test_dir}" + + ${ICE_CASEDIR}/casescripts/comparebfb.csh ${base_dir} ${test_dir} + set bfbstatus = $status + + if ( ${bfbstatus} != 0 ) then + + set test_file = `ls -1t ${ICE_RUNDIR}/cice.runlog* | head -1` + set base_file = `ls -1t ${ICE_BASELINE}/${ICE_BASECOM}/${ICE_TESTNAME}/cice.runlog* | head -1` + + echo "" + echo "bfb Log Compare Mode:" + echo "base_file: ${base_file}" + echo "test_file: ${test_file}" + + ${ICE_CASEDIR}/casescripts/comparelog.csh ${base_file} ${test_file} + set logstatus = $status + + if ( ${logstatus} == 0 ) then + echo "PASS ${ICE_TESTNAME} complog ${ICE_BASECOM}" >> ${ICE_CASEDIR}/test_output + echo "Regression baseline and test dataset may be the same" + else if ( ${logstatus} == 1 ) then + echo "FAIL ${ICE_TESTNAME} complog ${ICE_BASECOM} different-data" >> ${ICE_CASEDIR}/test_output + echo "Regression baseline and test dataset are not the same" + else if ( ${logstatus} == 2 ) then + echo "MISS ${ICE_TESTNAME} complog ${ICE_BASECOM} missing-data" >> ${ICE_CASEDIR}/test_output + echo "Missing data" + else + echo "FAIL ${ICE_TESTNAME} complog ${ICE_BASECOM} usage-error" >> ${ICE_CASEDIR}/test_output + echo "Regression baseline and test dataset error in usage" + endif + + endif + + endif - ${ICE_CASEDIR}/casescripts/comparebfb.csh ${base_dir} ${test_dir} - set bfbstatus = $status if ( ${bfbstatus} == 0 ) then echo "PASS ${ICE_TESTNAME} compare ${ICE_BASECOM} ${btimeloop} ${bdynamics} ${bcolumn}" >> ${ICE_CASEDIR}/test_output echo "Regression baseline and test dataset are identical" + else if ( ${bfbstatus} == 1 ) then + echo "FAIL ${ICE_TESTNAME} compare ${ICE_BASECOM} ${btimeloop} ${bdynamics} ${bcolumn} different-data" >> ${ICE_CASEDIR}/test_output + echo "Regression baseline and test dataset are different" else if ( ${bfbstatus} == 2 ) then echo "MISS ${ICE_TESTNAME} compare ${ICE_BASECOM} ${btimeloop} ${bdynamics} ${bcolumn} missing-data" >> ${ICE_CASEDIR}/test_output echo "Missing data" else - echo "FAIL ${ICE_TESTNAME} compare ${ICE_BASECOM} ${btimeloop} ${bdynamics} ${bcolumn} different-data" >> ${ICE_CASEDIR}/test_output - echo "Regression baseline and test dataset are different" + echo "FAIL ${ICE_TESTNAME} compare ${ICE_BASECOM} ${btimeloop} ${bdynamics} ${bcolumn} usage-error" >> ${ICE_CASEDIR}/test_output + echo "Regression baseline and test dataset error in usage" endif endif @@ -88,12 +138,15 @@ if (${ICE_BFBCOMP} != ${ICE_SPVAL}) then if (${bfbstatus} == 0) then echo "PASS ${ICE_TESTNAME} bfbcomp ${ICE_BFBCOMP}" >> ${ICE_CASEDIR}/test_output echo "bfb baseline and test dataset are identical" + else if (${bfbstatus} == 1) then + echo "FAIL ${ICE_TESTNAME} bfbcomp ${ICE_BFBCOMP} different-data" >> ${ICE_CASEDIR}/test_output + echo "bfbcomp and test dataset are different" else if (${bfbstatus} == 2) then echo "MISS ${ICE_TESTNAME} bfbcomp ${ICE_BFBCOMP} missing-data" >> ${ICE_CASEDIR}/test_output echo "Missing data" else - echo "FAIL ${ICE_TESTNAME} bfbcomp ${ICE_BFBCOMP} different-data" >> ${ICE_CASEDIR}/test_output - echo "bfbcomp and test dataset are different" + echo "FAIL ${ICE_TESTNAME} bfbcomp ${ICE_BFBCOMP} usage-error" >> ${ICE_CASEDIR}/test_output + echo "bfbcomp and test dataset usage error" endif endif diff --git a/configuration/scripts/tests/comparelog.csh b/configuration/scripts/tests/comparelog.csh index 8c1ff3a3c..d9e4a7a89 100755 --- a/configuration/scripts/tests/comparelog.csh +++ b/configuration/scripts/tests/comparelog.csh @@ -3,8 +3,9 @@ # Compare prognostic output in two log files #----------------------------------------------------------- -# usage: comparelog.csh base_file test_file +# usage: comparelog.csh base_file test_file [notcicefile] # does diff of two files +# optional 3rd argument indicates the file is not a cice file so diff entire thing # # Return Codes (depends on quality of error checking) # 0 = pass @@ -13,13 +14,26 @@ # 9 = error set filearg = 0 +set cicefile = 0 +set notcicefile = "notcicefile" if ( $#argv == 2 ) then + set cicefile = 1 set filearg = 1 set base_data = $argv[1] set test_data = $argv[2] -else + if ("$argv[1]" == "${notcicefile}") set filearg = 0 + if ("$argv[2]" == "${notcicefile}") set filearg = 0 +else if ( $#argv == 3 ) then + set cicefile = 0 + set filearg = 1 + set base_data = $argv[1] + set test_data = $argv[2] + if ("$argv[3]" != "${notcicefile}") set filearg = 0 +endif + +if (${filearg} == 0) then echo "Error in ${0}" - echo "Usage: ${0} " + echo "Usage: ${0} [notcicefile]" echo " does diff of two files" exit 9 endif @@ -28,7 +42,7 @@ set failure = 0 set base_out = "comparelog_base_out_file.log" set test_out = "comparelog_test_out_file.log" -if ($filearg == 1) then +if (${filearg} == 1) then echo "base_data: $base_data" echo "test_data: $test_data" if ( -f ${base_data} && -f ${test_data}) then @@ -38,12 +52,18 @@ if ($filearg == 1) then else touch ${base_out} - cat ${base_data} | grep -A 99999999 istep1: | grep -e istep1: -e = | grep -iv "min, max, sum" >&! ${base_out} touch ${test_out} - cat ${test_data} | grep -A 99999999 istep1: | grep -e istep1: -e = | grep -iv "min, max, sum" >&! ${test_out} + + if (${cicefile} == 1) then + cat ${base_data} | grep -A 99999999 "total ice area (km^2)" | grep -e istep1: -e = | grep -iv "min, max, sum" >&! ${base_out} + cat ${test_data} | grep -A 99999999 "total ice area (km^2)" | grep -e istep1: -e = | grep -iv "min, max, sum" >&! ${test_out} + else + cp -f ${base_data} ${base_out} + cp -f ${test_data} ${test_out} + endif set basenum = `cat ${base_out} | wc -l` - set testnum = `cat ${base_out} | wc -l` + set testnum = `cat ${test_out} | wc -l` set filediff = `diff -w ${base_out} ${test_out} | wc -l` if (${basenum} > 0 && ${testnum} > 0) then diff --git a/configuration/scripts/tests/io_suite.ts b/configuration/scripts/tests/io_suite.ts old mode 100755 new mode 100644 diff --git a/configuration/scripts/tests/quick_suite.ts b/configuration/scripts/tests/quick_suite.ts index 9384f0333..48646673d 100644 --- a/configuration/scripts/tests/quick_suite.ts +++ b/configuration/scripts/tests/quick_suite.ts @@ -2,5 +2,5 @@ smoke gx3 8x2 diag1,run5day smoke gx3 1x1 diag1,run1day restart gbox128 8x1 diag1 -restart gx3 4x2 debug,diag1,run5day +restart gx3 4x2 debug,diag1 smoke gx3 4x1 diag1,run5day,thread smoke_gx3_8x2_diag1_run5day diff --git a/configuration/scripts/tests/test_unittest.script b/configuration/scripts/tests/test_unittest.script new file mode 100644 index 000000000..0fcd148a6 --- /dev/null +++ b/configuration/scripts/tests/test_unittest.script @@ -0,0 +1,24 @@ + +#---------------------------------------------------- +# Run the CICE model +# cice.run returns -1 if run did not complete successfully + +./cice.run +set res="$status" + +set log_file = `ls -t1 ${ICE_RUNDIR}/cice.runlog* | head -1` + +mv -f ${ICE_CASEDIR}/test_output ${ICE_CASEDIR}/test_output.prev +cat ${ICE_CASEDIR}/test_output.prev | grep -iv "${ICE_TESTNAME} run" >! ${ICE_CASEDIR}/test_output +mv -f ${ICE_CASEDIR}/test_output ${ICE_CASEDIR}/test_output.prev +cat ${ICE_CASEDIR}/test_output.prev | grep -iv "${ICE_TESTNAME} test" >! ${ICE_CASEDIR}/test_output +rm -f ${ICE_CASEDIR}/test_output.prev + +set grade = FAIL +if ( $res == 0 ) then + set grade = PASS +endif + +echo "$grade ${ICE_TESTNAME} run " >> ${ICE_CASEDIR}/test_output +echo "$grade ${ICE_TESTNAME} test " >> ${ICE_CASEDIR}/test_output + diff --git a/configuration/scripts/tests/unittest_suite.ts b/configuration/scripts/tests/unittest_suite.ts new file mode 100644 index 000000000..2e9dcc7cf --- /dev/null +++ b/configuration/scripts/tests/unittest_suite.ts @@ -0,0 +1,4 @@ +# Test Grid PEs Sets BFB-compare +unittest gx3 1x1 helloworld +unittest gx3 1x1 calchk + diff --git a/doc/source/cice_index.rst b/doc/source/cice_index.rst index 59ddc4122..43b9cf58a 100644 --- a/doc/source/cice_index.rst +++ b/doc/source/cice_index.rst @@ -139,6 +139,7 @@ either Celsius or Kelvin units). "daymo", "number of days in one month", "" "daycal", "day number at end of month", "" "days_per_year", ":math:`\bullet` number of days in one year", "365" + "day_init", ":math:`\bullet` the initial day of the month", "" "dbl_kind", "definition of double precision", "selected_real_kind(13)" "dbug", ":math:`\bullet` write extra diagnostics", ".false." "Delta", "function of strain rates (see Section :ref:`dynam`)", "1/s" @@ -258,9 +259,9 @@ either Celsius or Kelvin units). "fswthru_idr", "near IR direct shortwave penetrating to ocean", "W/m\ :math:`^2`" "fswthru_idf", "near IR diffuse shortwave penetrating to ocean", "W/m\ :math:`^2`" "fswthru_ai", "grid-box-mean shortwave penetrating to ocean (fswthru)", "W/m\ :math:`^2`" - "fyear", "current data year", "" - "fyear_final", "last data year", "" - "fyear_init", ":math:`\bullet` initial data year", "" + "fyear", "current forcing data year", "" + "fyear_final", "last forcing data year", "" + "fyear_init", ":math:`\bullet` initial forcing data year", "" "**G**", "", "" "gravit", "gravitational acceleration", "9.80616 m/s\ :math:`^2`" "grid_file", ":math:`\bullet` input file for grid info", "" @@ -378,18 +379,21 @@ either Celsius or Kelvin units). "max_blocks", "maximum number of blocks per processor", "" "max_ntrcr", "maximum number of tracers available", "5" "maxraft", "maximum thickness of ice that rafts", "1. m" - "mday", "day of the month", "" + "mday", "model day of the month", "" "meltb", "basal ice melt", "m" "meltl", "lateral ice melt", "m" "melts", "snow melt", "m" "meltt", "top ice melt", "m" "min_salin", "threshold for brine pockets", "0.1 ppt" "mlt_onset", "day of year that surface melt begins", "" - "month", "the month number", "" + "mmonth", "model month number", "" "monthp", "previous month number", "" + "month_init", ":math:`\bullet` the initial month", "" "mps_to_cmpdy", "m per s to cm per day conversion", "8.64\ :math:`\times`\ 10\ :math:`^6`" + "msec", "model seconds elasped into day", "" "mtask", "local processor number that writes debugging data", "" "mu_rdg", ":math:`\bullet` e-folding scale of ridged ice", "" + "myear", "model year", "" "my_task", "task ID for the current processor", "" "**N**", "", "" "n_aero", "number of aerosol species", "" @@ -416,7 +420,8 @@ either Celsius or Kelvin units). "nlt_bgc_[chem]", "ocean sources and sinks for biogeochemistry", "" "nml_filename", "namelist file name", "" "nprocs", ":math:`\bullet` total number of processors", "" - "npt", ":math:`\bullet` total number of time steps (dt)", "" + "npt", ":math:`\bullet` total run length values associate with npt_unit", "" + "npt_unit", "units of the run length, number set by npt", "" "ns_boundary_type", ":math:`\bullet` type of north-south boundary condition", "" "nslyr", "number of snow layers in each category", "" "nspint", "number of solar spectral intervals", "" @@ -443,7 +448,6 @@ either Celsius or Kelvin units). "nvarz", "number of category, vertical grid fields written to history", "" "nx(y)_block", "total number of gridpoints on block in x(y) direction", "" "nx(y)_global", "number of physical gridpoints in x(y) direction, global domain", "" - "nyr", "year number", "" "**O**", "", "" "ocean_bio", "concentrations of bgc constituents in the ocean", "" "oceanmixed_file", ":math:`\bullet` data file containing ocean forcing data", "" @@ -555,8 +559,8 @@ either Celsius or Kelvin units). "scale_factor", "scaling factor for shortwave radiation components", "" "seabed_stress", "if true, calculate seabed stress", "F" "seabed_stress_method", "method for calculating seabed stress (‘LKD’ or ‘probabilistic’)", "LKD" - "sec", "seconds elasped into idate", "" "secday", "number of seconds in a day", "86400." + "sec_init", ":math:`\bullet` the initial second", "" "shcoef", "transfer coefficient for sensible heat", "" "shear", "strain rate II component", "1/s" "shlat", "southern latitude of artificial mask edge", "30\ :math:`^\circ`\ N" @@ -596,12 +600,11 @@ either Celsius or Kelvin units). "tarear", "1/tarea", "1/m\ :math:`^2`" "tareas", "area of southern hemisphere T-cells", "m\ :math:`^2`" "tcstr", "string identifying T grid for history variables", "" - "tday", "absolute day number", "" "Tf", "freezing temperature", "C" "Tffresh", "freezing temp of fresh ice", "273.15 K" "tfrz_option", ":math:`\bullet` form of ocean freezing temperature", "" "thinS", "minimum ice thickness for brine tracer", "" - "time", "total elapsed time", "s" + "timesecs", "total elapsed time in seconds", "s" "time_beg", "beginning time for history averages", "" "time_bounds", "beginning and ending time for history averages", "" "time_end", "ending time for history averages", "" @@ -681,7 +684,7 @@ either Celsius or Kelvin units). "**X**", "", "" "**Y**", "", "" "ycycle", ":math:`\bullet` number of years in forcing data cycle", "" - "yday", "day of the year", "" + "yday", "day of the year, computed in the model calendar", "" "yield_curve", "type of yield curve", "ellipse" "yieldstress11(12, 22)", "yield stress tensor components", "" "year_init", ":math:`\bullet` the initial year", "" diff --git a/doc/source/developer_guide/dg_dynamics.rst b/doc/source/developer_guide/dg_dynamics.rst index c94d47b35..47b54bde2 100644 --- a/doc/source/developer_guide/dg_dynamics.rst +++ b/doc/source/developer_guide/dg_dynamics.rst @@ -90,12 +90,9 @@ Time Manager Time manager data is module data in **cicecore/shared/ice_calendar.F90**. Much of the time manager data is public and operated on during the model timestepping. The model timestepping actually takes -place in the **CICE_RunMod.F90** file which is part of the driver code and tends to look like this:: +place in the **CICE_RunMod.F90** file which is part of the driver code. - call ice_step - istep = istep + 1 ! update time step counters - istep1 = istep1 + 1 - time = time + dt ! determine the time and date +The time manager was updated in early 2021. Additional information about the time manager can be found here, :ref:`timemanagerplus` diff --git a/doc/source/user_guide/ug_case_settings.rst b/doc/source/user_guide/ug_case_settings.rst index ccf7f0356..86240e904 100644 --- a/doc/source/user_guide/ug_case_settings.rst +++ b/doc/source/user_guide/ug_case_settings.rst @@ -78,6 +78,7 @@ can be modified as needed. "ICE_HSTDIR", "string", "unused", "${ICE_RUNDIR}/history" "ICE_LOGDIR", "string", "log directory", "${ICE_CASEDIR}/logs" "ICE_DRVOPT", "string", "unused", "standalone/cice" + "ICE_TARGET", "string", "build target", "set by cice.setup" "ICE_IOTYPE", "string", "I/O format", "set by cice.setup" " ", "netcdf", "serial netCDF" " ", "pio", "parallel netCDF" @@ -143,6 +144,7 @@ setup_nml "``conserv_check``", "logical", "check conservation", "``.false.``" "``cpl_bgc``", "logical", "couple bgc thru driver", "``.false.``" "``days_per_year``", "integer", "number of days in a model year", "365" + "``day_init``", "integer", "the initial day of the month if not using restart", "1" "``dbug``", "logical", "write extra diagnostics", "``.false.``" "``diagfreq``", "integer", "frequency of diagnostic output in timesteps", "24" "``diag_type``", "``stdout``", "write diagnostic output to stdout", "``stdout``" @@ -178,8 +180,15 @@ setup_nml "``latpnt``", "real", "latitude of (2) diagnostic points", "90.0,-65.0" "``lcdf64``", "logical", "use 64-bit netcdf format", "``.false.``" "``lonpnt``", "real", "longitude of (2) diagnostic points", "0.0,-45.0" + "``month_init``", "integer", "the initial month if not using restart", "1" "``ndtd``", "integer", "number of dynamics/advection/ridging/steps per thermo timestep", "1" - "``npt``", "integer", "total number of time steps to take", "99999" + "``npt``", "integer", "total number of npt_units to run the model", "99999" + "``npt_unit``", "``d``", "run ``npt`` days", "1" + "", "``h``", "run ``npt`` hours", "" + "", "``m``", "run ``npt`` months", "" + "", "``s``", "run ``npt`` seconds", "" + "", "``y``", "run ``npt`` years", "" + "", "``1``", "run ``npt`` timesteps", "" "``numin``", "integer", "minimum internal IO unit number", "11" "``numax``", "integer", "maximum internal IO unit number", "99" "``pointer_file``", "string", "restart pointer filename", "'ice.restart_file'" @@ -194,6 +203,7 @@ setup_nml "``runid``", "string", "label for run (currently CESM only)", "'unknown'" "``runtype``", "``continue``", "restart using ``pointer_file``", "``initial``" "", "``initial``", "start from ``ice_ic``", "" + "``sec_init``", "integer", "the initial second if not using restart", "0" "``use_leap_years``", "logical", "include leap days", "``.false.``" "``use_restart_time``", "logical", "set initial date using restart file", "``.true.``" "``version_name``", "string", "model version", "'unknown_version_name'" @@ -355,6 +365,8 @@ dynamics_nml "", "", "", "" "``advection``", "``remap``", "linear remapping advection scheme", "``remap``" "", "``upwind``", "donor cell advection", "" + "``algo_nonlin``", "``anderson``", "use nonlinear anderson algorithm for implicit solver", "picard" + "", "``picard``", "use picard algorithm", "" "``alphab``", "real", ":math:`\alpha_{b}` factor in :cite:`Lemieux16`", "20.0" "``arlx``", "real", "revised_evp value", "300.0" "``brlx``", "real", "revised_evp value", "300.0" @@ -394,11 +406,11 @@ dynamics_nml "``monitor_pgmres``", "logical", "write velocity norm at each PGMRES iteration", "``.false.``" "``mu_rdg``", "real", "e-folding scale of ridged ice for ``krdg_partic`` = 1 in m^0.5", "3.0" "``ndte``", "integer", "number of EVP subcycles", "120" - "``ortho_type``", "``mgs``", "Use modified Gram-Shchmidt in FGMRES solver", "``mgs``" - "", "``cgs``", "Use classical Gram-Shchmidt in FGMRES solver", "" - "``precond``", "``pgmres``", "Use GMRES as preconditioner for FGMRES solver", "``pgmres``" - "", "``diag``", "Use Jacobi preconditioner for the FGMRES solver", "" + "``ortho_type``", "``cgs``", "Use classical Gram-Shchmidt in FGMRES solver", "``mgs``" + "", "``mgs``", "Use modified Gram-Shchmidt in FGMRES solver", "" + "``precond``", "``diag``", "Use Jacobi preconditioner for the FGMRES solver", "``pgmres``" "", "``ident``", "Don't use a preconditioner for the FGMRES solver", "" + "", "``pgmres``", "Use GMRES as preconditioner for FGMRES solver", "" "``Pstar``", "real", "constant in Hibler strength formula (N/m\ :math:`^2`)", "2.75e4" "``reltol_nonlin``", "real", "relative tolerance for nonlinear solver", "1e-8" "``reltol_fgmres``", "real", "relative tolerance for FGMRES solver", "1e-2" @@ -411,6 +423,7 @@ dynamics_nml "", "``geostropic``", "computed from ocean velocity", "" "``threshold_hw``", "real", "Max water depth for grounding (see :cite:`Amundrud04`)", "30." "``yield_curve``", "``ellipse``", "elliptical yield curve", "``ellipse``" + "``use_mean_vrel``", "logical", "Use mean of two previous iterations for vrel in VP", "``.true.``" "", "", "", "" shortwave_nml @@ -525,6 +538,7 @@ forcing_nml "``restart_coszen``", "logical", "read/write coszen in restart files", "``.false.``" "``restore_ocn``", "logical", "restore sst to data", "``.false.``" "``restore_ice``", "logical", "restore ice state along lateral boundaries", "``.false.``" + "``rotate_wind``", "logical", "rotate wind from east/north to computation grid", "``.true.``" "``tfrz_option``", "``linear_salt``", "linear functino of salinity (ktherm=1)", "``mushy``" "", "``minus1p8``", "constant ocean freezing temperature (:math:`-1.8^{\circ} C`)", "" "", "``mushy``", "matches mushy-layer thermo (ktherm=2)", "" diff --git a/doc/source/user_guide/ug_implementation.rst b/doc/source/user_guide/ug_implementation.rst index cbfe37b0c..88afdcf52 100644 --- a/doc/source/user_guide/ug_implementation.rst +++ b/doc/source/user_guide/ug_implementation.rst @@ -529,12 +529,72 @@ schemes and the aerosol tracers, and the level-ice pond parameterization additionally requires the level-ice tracers. +.. _timemanagerplus: + +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +Time Manager and Initialization +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +The time manager is an important piece of the CICE model. + +.. _timemanager: + +**************************** +Time Manager +**************************** + +The primary prognostic variables in the time manager are ``myear``, +``mmonth``, ``mday``, and ``msec``. These are integers and identify +the current model year, month, day, and second respectively. +The model timestep is ``dt`` with units of seconds. See :ref:`parameters` +for additional information about choosing an appropriate timestep. +The internal variables ``istep``, ``istep0``, and ``istep1`` keep +track of the number of timesteps. ``istep`` is the counter for +the current run and is set to 0 at the start of each run. ``istep0`` +is the step count at the start of a long multi-restart run, and +``istep1`` is the step count of a long multi-restart run. + +In general, the time manager should be advanced by calling +*advance\_timestep*. This subroutine in **ice\_calendar.F90** +automatically advances the model time by ``dt``. It also advances +the istep numbers and calls subroutine *calendar* to update +additional calendar data. + +The namelist variable ``use_restart_time`` specifies whether to +use the time and step numbers saved on a restart file or whether +to set the initial model time to the namelist values defined by +``year_init``, ``month_init``, ``day_init``, and ``sec_init``. +Normally, ``use_restart_time`` is set to false on the initial run +and then set to true on subsequent restart runs of the same +case to allow time to advance thereafter. More information about +the restart capability can be found here, :ref:`restartfiles`. + +The time manager was updated in early 2021. The standalone model +was modified, and some tests were done in a coupled framework after +modifications to the high level coupling interface. For some coupled models, the +coupling interface may need to be updated when updating CICE with the new time manager. +In particular, the old prognostic variable ``time`` no longer exists in CICE, +``year_init`` only defines the model initial year, and +the calendar subroutine is called without any arguments. One can +set the namelist variables ``year_init``, ``month_init``, ``day_init``, +``sec_init``, and ``dt`` in conjuction with ``days_per_year`` and +``use_leap_years`` to initialize the model date, timestep, and calendar. +To overwrite the default/namelist settings in the coupling layer, +set the **ice\_calendar.F90** variables ``myear``, ``mmonth``, ``mday``, +``msec`` and ``dt`` after the namelists have been read. Subroutine +*calendar* should then be called to update all the calendar data. +Finally, subroutine *advance\_timestep* should be used to advance +the model time manager. It advances the step numbers, advances +time by ``dt``, and updates the calendar data. The older method +of manually advancing the steps and adding ``dt`` to ``time`` should +be deprecated. + .. _init: -~~~~~~~~~~~~~~~~~~~~~~~~~~~ -Initialization and coupling -~~~~~~~~~~~~~~~~~~~~~~~~~~~ +**************************** +Initialization and Restarts +**************************** The ice model’s parameters and variables are initialized in several steps. Many constants and physical parameters are set in @@ -612,9 +672,9 @@ reset to ‘none.’ .. _parameters: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +********************************** Choosing an appropriate time step -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +********************************** The time step is chosen based on stability of the transport component (both horizontal and in thickness space) and on resolution of the @@ -705,6 +765,8 @@ the problem, and ``brlx`` represents the effective subcycling Model output ~~~~~~~~~~~~ +There are a number of model output streams and formats. + .. _history: ************* @@ -759,7 +821,8 @@ is now a character string corresponding to ``histfreq`` or ‘x’ for none. files, no matter what the frequency is.) If there are no namelist flags with a given ``histfreq`` value, or if an element of ``histfreq_n`` is 0, then no file will be written at that frequency. The output period can be -discerned from the filenames. +discerned from the filenames. Because all history is average output, it's +not possible to write instananeous output at any frequency except every timestep. For example, in the namelist: @@ -784,6 +847,14 @@ as long as for a single frequency. If you only want monthly output, the most efficient setting is ``histfreq`` = ’m’,’x’,’x’,’x’,’x’. The code counts the number of desired streams (``nstreams``) based on ``histfreq``. +There is no restart capability built into the history implementation. If the +model stops in the middle of a history accumulation period, that data is lost +on restart, and the accumulation is zeroed out at startup. That means the +dump frequency (see :ref:`restartfiles`) and history frequency need to be +somewhat coordinated. For +example, if monthly history files are requested, the dump frequency should be +set to an integer number of months. + The history variable names must be unique for netCDF, so in cases where a variable is written at more than one frequency, the variable name is appended with the frequency in files after the first one. In the example @@ -908,6 +979,8 @@ The timers use *MPI\_WTIME* for parallel runs and the F90 intrinsic | 16 | BGC | biogeochemistry | +--------------+-------------+----------------------------------------------------+ +.. _restartfiles: + ************* Restart files ************* @@ -937,7 +1010,8 @@ Additional namelist flags provide further control of restart behavior. of a run when it is otherwise not scheduled to occur. The flag ``use_restart_time`` enables the user to choose to use the model date provided in the restart files. If ``use_restart_time`` = false then the -initial model date stamp is determined from the namelist parameters. +initial model date stamp is determined from the namelist parameters, +``year_init``, ``month_init``, ``day_init``, and ``sec_init``.. lcdf64 = true sets 64-bit netCDF output, allowing larger file sizes. Routines for gathering, scattering and (unformatted) reading and writing @@ -957,5 +1031,6 @@ initialized with no ice. The gx3 case was run for 1 year using the 1997 forcing data provided with the code. The gx1 case was run for 20 years, so that the date of restart in the file is 1978-01-01. Note that the restart dates provided in the restart files can be overridden using the -namelist variables ``use_restart_time``, ``year_init`` and ``istep0``. The +namelist variables ``use_restart_time``, ``year_init``, ``month_init``, +``day_init``, and ``sec_init``. The forcing time can also be overridden using ``fyear_init``. diff --git a/doc/source/user_guide/ug_testing.rst b/doc/source/user_guide/ug_testing.rst index 61aa1c05f..4eb7b6795 100644 --- a/doc/source/user_guide/ug_testing.rst +++ b/doc/source/user_guide/ug_testing.rst @@ -625,6 +625,49 @@ Test Suite Examples The setenv syntax is for csh/tcsh. In bash, the syntax would be SUITE_BUILD=true. +.. _unittesting: + +Unit Testing +--------------- + +Unit testing is supported in the CICE scripts. Unit tests are implemented +via a distinct top level driver that tests CICE model features explicitly. +These drivers can be found in **cicecore/drivers/unittest/**. In addition, +there are some script files that also support the unit testing. + +The unit tests build and run very much like the standard CICE model. +A case is created and model output is saved to the case logs directory. +Unit tests can be run as part of a test suite and the output is +compared against an earlier set of output using a simple diff of the +log files. + +For example, to run the existing calendar unit test as a case, + +.. code-block:: bash + + ./cice.setup -m onyx -e intel --case calchk01 -p 1x1 -s calchk + cd calchk01 + ./cice.build + ./cice.submit + +Or to run the existing calendar unit test as a test, + +.. code-block:: bash + + ./cice.setup -m onyx -e intel --test unittest -p 1x1 --testid cc01 -s calchk --bgen cice.cc01 + cd onyx_intel_unittest_gx3_1x1_calchk.cc01/ + ./cice.build + ./cice.submit + +To create a new unit test, add a new driver in **cicecore/driver/unittest**. +The directory name should be the name of the test. +Then create the appropriate set_nml or set_env files for the new unittest name +in **configuration/scripts/options**. In particular, **ICE_DRVOPT** and +**ICE_TARGET** need to be defined in a set_env file. Finally, edit +**configuration/scripts/Makefile** and create a target for the unit test. +The unit tests calchk or helloworld can be used as examples. + + .. _testreporting: Test Reporting