Update gx1 initial condition, Icepack, tests, and version. Fix hmix d…

…efault value (#586)

* Update gx1 ic

- update set_nml.gx1prod to match current production system
- add apr 1 test case for gx1
- update landice tests, use gx1
- fix hmix default #585
- delete old code in ice_forcing.F90 (should have been done in earlier PR)

* update gx1coreii initial condition

* update icepack

* add gx1 debug test, expected to fail

* Update version number to 6.2.0

* Update gx3 and gx1 input filenames
Update CI wget implementation #588
Create variable in ice_forcing.F90 called mixed_layer_depth_default and use it instead of c20
Fix bug in hourly output, created in time manager update #589
Set start year for all runs to 2005 and turn leap years on by default
Update documentation in history output section to add information about hist_avg namelist, noted in #566

* update ic filenames, tests, and documentation

* update color links in test results wiki

* update hmix initialization

.github/workflows/test-cice.yml

@@ -17,7 +17,7 @@ on:
 
 defaults:
   run:
-    shell: /bin/csh {0}
+    shell: /bin/csh -e {0}
 
 jobs:
   build:
@@ -104,9 +104,9 @@ jobs:
       - name: download input data
         run: |
           cd $HOME/cice-dirs/input
-          wget https://zenodo.org/record/3728358/files/CICE_data_gx3_grid_ic-20200320.tar.gz && tar xvfz CICE_data_gx3_grid_ic-20200320.tar.gz
-          wget https://zenodo.org/record/3728362/files/CICE_data_gx3_forcing_NCAR_bulk-20200320.tar.gz && tar xvfz CICE_data_gx3_forcing_NCAR_bulk-20200320.tar.gz
-          wget https://zenodo.org/record/3728364/files/CICE_data_gx3_forcing_JRA55-20200320.tar.gz && tar xvfz CICE_data_gx3_forcing_JRA55-20200320.tar.gz
+          wget --progress=dot:giga https://zenodo.org/record/3728358/files/CICE_data_gx3_grid_ic-20200320.tar.gz && tar xvfz CICE_data_gx3_grid_ic-20200320.tar.gz
+          wget --progress=dot:giga https://zenodo.org/record/3728362/files/CICE_data_gx3_forcing_NCAR_bulk-20200320.tar.gz && tar xvfz CICE_data_gx3_forcing_NCAR_bulk-20200320.tar.gz
+          wget --progress=dot:giga https://zenodo.org/record/3728364/files/CICE_data_gx3_forcing_JRA55-20200320.tar.gz && tar xvfz CICE_data_gx3_forcing_JRA55-20200320.tar.gz
           pwd
           ls -alR
 #      - name: run case

cicecore/cicedynB/general/ice_forcing.F90

@@ -169,6 +169,9 @@ module ice_forcing
 
       ! PRIVATE:
 
+      real (dbl_kind), parameter :: &
+         mixed_layer_depth_default = c20  ! default mixed layer depth in m
+
       logical (kind=log_kind), parameter :: &
          forcing_debug = .false.   ! local debug flag
 
@@ -2668,247 +2671,6 @@ subroutine JRA55_data
 
       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
-      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, use_leap_years
-
-      integer (kind=int_kind), intent(in) :: &
-          yr              ! current forcing year
-
-      integer (kind=int_kind) :: & 
-          ncid        , & ! netcdf file id
-          i, j, n1, iblk, &
-          yrp         , & ! year after yr in forcing cycle
-          recnum      , & ! record number
-          maxrec      , & ! maximum record number
-          recslot     , & ! spline slot for current record
-          dataloc         ! = 1 for data located in middle of time interval
-                          ! = 2 for date located at end of time interval
-
-      real (kind=dbl_kind) :: &
-          sec3hr          , & ! number of seconds in 3 hours
-          secday          , & ! number of seconds in day
-          eps, tt         , & ! interpolation coeff calc
-          Tffresh         , &
-          vmin, vmax
-
-      logical (kind=log_kind) :: debug_n_d = .false.
-
-      character (char_len_long) :: uwind_file_old
-      character(len=64) :: fieldname !netcdf field name
-      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)
-      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 (debug_n_d .and. my_task == master_task) then
-         write (nu_diag,*) subname,'recnum',recnum
-         write (nu_diag,*) subname,'maxrec',maxrec
-         write (nu_diag,*) subname,'days_per_year', days_per_year
-      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
-      ! Both states and fluxes have 1st record defined as 00z Jan 1
-      ! interpolate states, do not interpolate fluxes
-      ! fluxes are held constant from [init period, end period)
-      !-------------------------------------------------------------------
-      ! File is NETCDF with winds in NORTH and EAST direction
-      ! file variable names are:
-      ! glbrad   (shortwave W/m^2)
-      ! dlwsfc   (longwave W/m^2)
-      ! wndewd   (eastward wind m/s)
-      ! wndnwd   (northward wind m/s)
-      ! airtmp   (air temperature K)
-      ! spchmd   (specific humidity kg/kg)
-      ! ttlpcp   (precipitation kg/m s-1)
-      !-------------------------------------------------------------------
-
-      uwind_file_old = uwind_file
-      call file_year(uwind_file,yr)
-      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
-
-         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
-               yrp = fyear_init + mod(myear,ycycle)  ! next year
-               recnum = 1
-               call file_year(uwind_file,yrp)
-               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 (debug_n_d .and. my_task == master_task) then
-            write(nu_diag,*) subname,' read recnum = ',recnum,n1
-         endif
-
-         fieldname = 'airtmp'
-         call ice_read_nc(ncid,recnum,fieldname,Tair_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-
-         fieldname = 'wndewd'
-         call ice_read_nc(ncid,recnum,fieldname,uatm_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-
-         fieldname = 'wndnwd'
-         call ice_read_nc(ncid,recnum,fieldname,vatm_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-
-         fieldname = 'spchmd'
-         call ice_read_nc(ncid,recnum,fieldname,Qa_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-
-       ! only read one timestep for fluxes, 3 hr average, no interpolation
-       if (n1 == 1) then
-         fieldname = 'glbrad'
-         call ice_read_nc(ncid,recnum,fieldname,fsw_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-
-         fieldname = 'dlwsfc'
-         call ice_read_nc(ncid,recnum,fieldname,flw_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-
-         fieldname = 'ttlpcp'
-         call ice_read_nc(ncid,recnum,fieldname,fsnow_data(:,:,n1,:),debug_n_d, &
-              field_loc=field_loc_center, &
-              field_type=field_type_scalar)
-       endif
-
-      enddo
-
-      call ice_close_nc(ncid)
-
-      ! reset uwind_file to original year
-      call file_year(uwind_file,yr)
-
-      ! 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 (debug_n_d .and. my_task == master_task) then
-         write(nu_diag,*) subname,' 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
-      ! 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 (debug_n_d .or. dbug) then
-           if (my_task.eq.master_task) &
-               write (nu_diag,*) subname,'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,'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,'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,'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,'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,'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,'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,'Qa',vmin,vmax
-
-        endif                   ! dbug
-
-      end subroutine Jra55_data_old
-
 !=======================================================================
 !
 ! AOMIP shortwave forcing
@@ -4403,7 +4165,7 @@ subroutine ocn_data_ncar(dt)
             do iblk = 1, nblocks
               if (hm(i,j,iblk) == c1) then
                 if (n == 2) sss    (i,j,iblk) = work1(i,j,iblk)
-                if (n == 3) hmix   (i,j,iblk) = work1(i,j,iblk)
+                if (n == 3) hmix   (i,j,iblk) = max(mixed_layer_depth_default,work1(i,j,iblk))
                 if (n == 4) uocn   (i,j,iblk) = work1(i,j,iblk)
                 if (n == 5) vocn   (i,j,iblk) = work1(i,j,iblk)
                 if (n == 6) ss_tltx(i,j,iblk) = work1(i,j,iblk)
@@ -4518,7 +4280,7 @@ subroutine ocn_data_oned
       ss_tlty(:,:,:) = c0
       frzmlt (:,:,:) = c0              ! freezing/melting potential (W/m^2)
       qdp    (:,:,:) = c0              ! deep ocean heat flux (W/m^2)
-      hmix   (:,:,:) = c20             ! ocean mixed layer depth
+      hmix   (:,:,:) = mixed_layer_depth_default   ! ocean mixed layer depth
 
       end subroutine ocn_data_oned
 

cicecore/shared/ice_calendar.F90

@@ -62,6 +62,7 @@ module ice_calendar
 
       integer (kind=int_kind), public :: &
          seconds_per_day       , & ! seconds per day
+         seconds_per_hour      , & ! seconds per hour
          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
@@ -169,12 +170,18 @@ subroutine init_calendar
          write(nu_diag,*) trim(subname),' ERROR secday should basically be an integer',secday
          call abort_ice(subname//'ERROR: improper secday')
       endif
+      seconds_per_hour = nint(secday/real(hours_per_day,kind=dbl_kind))
+      if (abs(seconds_per_hour*hours_per_day - seconds_per_day) > 0) then
+         write(nu_diag,*) trim(subname),' ERROR seconds per day and hours per day inconsistent'
+         call abort_ice(subname//'ERROR: improper seconds_per_hour')
+      endif
 
       istep = 0         ! local timestep number
       myear=year_init   ! year
       mmonth=month_init ! month
       mday=day_init     ! day of the month
       msec=sec_init     ! seconds into date
+      hour=0            ! computed in calendar, but needs some reasonable initial value
       istep1 = istep0   ! number of steps at current timestep
                         ! real (dumped) or imagined (use to set calendar)
       idate0 = (myear)*10000 + mmonth*100 + mday ! date (yyyymmdd) 
@@ -343,6 +350,7 @@ subroutine calendar()
 
       idate = (myear)*10000 + mmonth*100 + mday ! date (yyyymmdd) 
       yday = daycal(mmonth) + mday            ! day of the year
+      hour = (msec+1)/(seconds_per_hour)
       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

cicecore/version.txt

@@ -1 +1 @@
-CICE 6.1.4
+CICE 6.2.0

configuration/scripts/ice_in

@@ -1,7 +1,10 @@
 &setup_nml
     days_per_year  = 365
-    use_leap_years = .false.
-    year_init      = 1997
+    use_leap_years = .true.
+    year_init      = 2005
+    month_init     = 1
+    day_init       = 1
+    sec_init       = 0
     istep0         = 0
     dt             = 3600.0
     npt_unit       = '1'
@@ -217,7 +220,7 @@
     bgc_data_type   = 'default'
     fe_data_type    = 'default'
     ice_data_type   = 'default'
-    fyear_init      = 1997
+    fyear_init      = 2005
     ycycle          = 1
     atm_data_format = 'bin'
     atm_data_dir    = '/glade/u/home/tcraig/cice_data/'

configuration/scripts/options/set_nml.box2001

@@ -1,4 +1,5 @@
 days_per_year = 360
+use_leap_years = .false.
 npt = 240
 ice_ic = 'default'
 restart = .false.

configuration/scripts/options/set_nml.boxdyn

@@ -2,6 +2,7 @@ nilyr =	1
 ice_ic = 'default'
 restart = .false.
 days_per_year = 360
+use_leap_years = .false.
 npt            = 72
 dumpfreq       = 'd'
 dumpfreq_n     = 2

configuration/scripts/options/set_nml.gx1

@@ -3,7 +3,7 @@ runtype = 'initial'
 year_init      = 2005
 use_leap_years = .true.
 use_restart_time = .false.
-ice_ic  = 'ICE_MACHINE_INPUTDATA/CICE_data/ic/gx1/iced_gx1_v5.nc'
+ice_ic  = 'ICE_MACHINE_INPUTDATA/CICE_data/ic/gx1/iced_gx1_v6.2005-01-01.nc'
 grid_format = 'bin'
 grid_type   = 'displaced_pole'
 grid_file   = 'ICE_MACHINE_INPUTDATA/CICE_data/grid/gx1/grid_gx1.bin'
@@ -17,5 +17,5 @@ atm_data_format = 'nc'
 atm_data_type   = 'JRA55_gx1'
 atm_data_dir    = 'ICE_MACHINE_INPUTDATA/CICE_data/forcing/gx1/JRA55'
 precip_units    = 'mks'
-ocn_data_dir    = 'ICE_MACHINE_INPUTDATA/CICE_data/forcing/gx1/COREII'
+ocn_data_dir    = 'ICE_MACHINE_INPUTDATA/CICE_data/forcing/gx1/CESM/MONTHLY'
 bgc_data_dir    = 'ICE_MACHINE_INPUTDATA/CICE_data/forcing/gx1/WOA/MONTHLY'

configuration/scripts/options/set_nml.gx1apr

@@ -0,0 +1,5 @@
+year_init   = 2005
+month_init  = 4
+day_init    = 1
+sec_init    = 0
+ice_ic  = 'ICE_MACHINE_INPUTDATA/CICE_data/ic/gx1/iced_gx1_v6.2005-04-01.nc'