Update to CICE-Consortium/main (2024-05-01)

.github/workflows/test-cice.yml

@@ -46,7 +46,6 @@ jobs:
  run: |
  sudo xcode-select -r
  sudo xcode-select -s /Library/Developer/CommandLineTools
- sudo ln -s /Library/Developer/CommandLineTools/SDKs/MacOSX.sdk/usr/include/* /usr/local/include/
  echo "xcrun --show-sdk-path: $(xcrun --show-sdk-path)"
  echo "xcode-select -p: $(xcode-select -p)"
  - name: system info
@@ -55,7 +54,7 @@ jobs:
  type wget
  type curl
  type csh
- echo "readlink \$(which csh): $(python -c 'import os, sys; print os.path.realpath(sys.argv[1])' $(which csh))"
+ echo "readlink \$(which csh): $(python -c 'import os, sys; print(os.path.realpath(sys.argv[1]))' $(which csh))"
  echo "csh --version: $(csh --version)"
  echo "uname -a: $(uname -a)"
  echo "sw_vers: $(sw_vers)"

cicecore/cicedyn/analysis/ice_history.F90

@@ -362,6 +362,7 @@ subroutine init_hist (dt)
  f_sidmasslat = 'mxxxx'
  f_sndmasssnf = 'mxxxx'
  f_sndmassmelt = 'mxxxx'
+ f_sndmassdyn = 'mxxxx'
  f_siflswdtop = 'mxxxx'
  f_siflswutop = 'mxxxx'
  f_siflswdbot = 'mxxxx'
@@ -402,6 +403,11 @@ subroutine init_hist (dt)
  f_siu = f_CMIP
  f_siv = f_CMIP
  f_sispeed = f_CMIP
+ f_sndmasssubl = f_CMIP
+ f_sndmasssnf = f_CMIP
+ f_sndmassmelt = f_CMIP
+ f_sndmassdyn = f_CMIP
+ f_sidmasssi = f_CMIP
  endif
 
  if (grid_ice == 'CD' .or. grid_ice == 'C') then
@@ -447,6 +453,14 @@ subroutine init_hist (dt)
  if (f_Tsnz (1:1) /= 'x') f_VGRDs = .true.
  if (tr_fsd) f_NFSD = .true.
 
+ call broadcast_scalar (f_tlon, master_task)
+ call broadcast_scalar (f_tlat, master_task)
+ call broadcast_scalar (f_ulon, master_task)
+ call broadcast_scalar (f_ulat, master_task)
+ call broadcast_scalar (f_nlon, master_task)
+ call broadcast_scalar (f_nlat, master_task)
+ call broadcast_scalar (f_elon, master_task)
+ call broadcast_scalar (f_elat, master_task)
  call broadcast_scalar (f_tmask, master_task)
  call broadcast_scalar (f_umask, master_task)
  call broadcast_scalar (f_nmask, master_task)
@@ -646,6 +660,7 @@ subroutine init_hist (dt)
  call broadcast_scalar (f_sidmasslat, master_task)
  call broadcast_scalar (f_sndmasssnf, master_task)
  call broadcast_scalar (f_sndmassmelt, master_task)
+ call broadcast_scalar (f_sndmassdyn, master_task)
  call broadcast_scalar (f_siflswdtop, master_task)
  call broadcast_scalar (f_siflswutop, master_task)
  call broadcast_scalar (f_siflswdbot, master_task)
@@ -1640,7 +1655,7 @@ subroutine init_hist (dt)
  "none", c1, c0, &
  ns1, f_sidmassevapsubl)
 
- call define_hist_field(n_sndmasssubl,"sndmassubl","kg m-2 s-1",tstr2D, tcstr, &
+ call define_hist_field(n_sndmasssubl,"sndmasssubl","kg m-2 s-1",tstr2D, tcstr, &
  "snow mass change from evaporation and sublimation", &
  "none", c1, c0, &
  ns1, f_sndmasssubl)
@@ -1670,6 +1685,11 @@ subroutine init_hist (dt)
  "none", c1, c0, &
  ns1, f_sndmassmelt)
 
+ call define_hist_field(n_sndmassdyn,"sndmassdyn","kg m-2 s-1",tstr2D, tcstr, &
+ "snow mass change from dynamics ridging", &
+ "none", c1, c0, &
+ ns1, f_sndmassdyn)
+
  call define_hist_field(n_siflswdtop,"siflswdtop","W/m2",tstr2D, tcstr, &
  "down shortwave flux over sea ice", &
  "positive downward", c1, c0, &
@@ -1973,6 +1993,21 @@ subroutine init_hist (dt)
  ! floe size distribution
  call init_hist_fsd_4Df
 
+ !-----------------------------------------------------------------
+ ! fill icoord array with namelist values
+ !-----------------------------------------------------------------
+
+ icoord=.true.
+
+ icoord(n_tlon ) = f_tlon
+ icoord(n_tlat ) = f_tlat
+ icoord(n_ulon ) = f_ulon
+ icoord(n_ulat ) = f_ulat
+ icoord(n_nlon ) = f_nlon
+ icoord(n_nlat ) = f_nlat
+ icoord(n_elon ) = f_elon
+ icoord(n_elat ) = f_elat
+
  !-----------------------------------------------------------------
  ! fill igrd array with namelist values
  !-----------------------------------------------------------------
@@ -2137,7 +2172,7 @@ subroutine accum_hist (dt)
  taubxN, taubyN, strocnxN, strocnyN, &
  strairxE, strairyE, strtltxE, strtltyE, strintxE, strintyE, &
  taubxE, taubyE, strocnxE, strocnyE, &
- fmU, fmN, fmE, daidtt, dvidtt, daidtd, dvidtd, fsurf, &
+ fmU, fmN, fmE, daidtt, dvidtt, daidtd, dvidtd, dvsdtd, fsurf, &
  fcondtop, fcondbot, fsurfn, fcondtopn, flatn, fsensn, albcnt, snwcnt, &
  stressp_1, stressm_1, stress12_1, &
  stresspT, stressmT, stress12T, &
@@ -3045,7 +3080,7 @@ subroutine accum_hist (dt)
  do j = jlo, jhi
  do i = ilo, ihi
  if (aice(i,j,iblk) > puny) then
- worka(i,j) = evaps(i,j,iblk)*rhos
+ worka(i,j) = evaps(i,j,iblk)
  endif
  enddo
  enddo
@@ -3057,7 +3092,7 @@ subroutine accum_hist (dt)
  do j = jlo, jhi
  do i = ilo, ihi
  if (aice(i,j,iblk) > puny) then
- worka(i,j) = aice(i,j,iblk)*fsnow(i,j,iblk)*rhos
+ worka(i,j) = aice(i,j,iblk)*fsnow(i,j,iblk)
  endif
  enddo
  enddo
@@ -3076,6 +3111,18 @@ subroutine accum_hist (dt)
  call accum_hist_field(n_sndmassmelt, iblk, worka(:,:), a2D)
  endif
 
+ if (f_sndmassdyn(1:1) /= 'x') then
+ worka(:,:) = c0
+ do j = jlo, jhi
+ do i = ilo, ihi
+ if (aice(i,j,iblk) > puny) then
+ worka(i,j) = dvsdtd(i,j,iblk)*rhos
+ endif
+ enddo
+ enddo
+ call accum_hist_field(n_sndmassdyn, iblk, worka(:,:), a2D)
+ endif
+
  if (f_siflswdtop(1:1) /= 'x') then
  worka(:,:) = c0
  do j = jlo, jhi

cicecore/cicedyn/analysis/ice_history_shared.F90

@@ -57,7 +57,7 @@ module ice_history_shared
  history_rearranger ! history file rearranger, box or subset for pio
 
  character (len=char_len), public :: &
- hist_suffix(max_nstrm) ! appended to 'h' in filename when not 'x'
+ hist_suffix(max_nstrm) ! appended to history_file in filename
 
  integer (kind=int_kind), public :: &
  history_iotasks , & ! iotasks, root, stride defines io pes for pio
@@ -131,6 +131,7 @@ module ice_history_shared
  avail_hist_fields(max_avail_hist_fields)
 
  integer (kind=int_kind), parameter, public :: &
+ ncoord = 8 , & ! number of coordinate variables: TLON, TLAT, ULON, ULAT, NLON, NLAT, ELON, ELAT
  nvar_grd = 21 , & ! number of grid fields that can be written
  ! excluding grid vertices
  nvar_grdz = 6 ! number of category/vertical grid fields written
@@ -165,6 +166,7 @@ module ice_history_shared
  avgct(max_nstrm) ! average sample counter
 
  logical (kind=log_kind), public :: &
+ icoord(ncoord) , & ! true if coord field is written to output file
  igrd (nvar_grd), & ! true if grid field is written to output file
  igrdz(nvar_grdz) ! true if category/vertical grid field is written
 
@@ -194,6 +196,10 @@ module ice_history_shared
  !---------------------------------------------------------------
 
  logical (kind=log_kind), public :: &
+ f_tlon = .true., f_tlat = .true., &
+ f_ulon = .true., f_ulat = .true., &
+ f_nlon = .true., f_nlat = .true., &
+ f_elon = .true., f_elat = .true., &
  f_tmask = .true., f_umask = .true., &
  f_nmask = .true., f_emask = .true., &
  f_blkmask = .true., &
@@ -308,6 +314,7 @@ module ice_history_shared
  f_sidmasslat = 'x', &
  f_sndmasssnf = 'x', &
  f_sndmassmelt = 'x', &
+ f_sndmassdyn = 'x', &
  f_siflswdtop = 'x', &
  f_siflswutop = 'x', &
  f_siflswdbot = 'x', &
@@ -362,6 +369,10 @@ module ice_history_shared
  !---------------------------------------------------------------
 
  namelist / icefields_nml / &
+ f_tlon , f_tlat , &
+ f_ulon , f_ulat , &
+ f_nlon , f_nlat , &
+ f_elon , f_elat , &
  f_tmask , f_umask , &
  f_nmask , f_emask , &
  f_blkmask , &
@@ -475,6 +486,7 @@ module ice_history_shared
  f_sidmasslat, &
  f_sndmasssnf, &
  f_sndmassmelt, &
+ f_sndmassdyn, &
  f_siflswdtop, &
  f_siflswutop, &
  f_siflswdbot, &
@@ -529,6 +541,15 @@ module ice_history_shared
  !---------------------------------------------------------------
 
  integer (kind=int_kind), parameter, public :: &
+ n_tlon = 1, &
+ n_tlat = 2, &
+ n_ulon = 3, &
+ n_ulat = 4, &
+ n_nlon = 5, &
+ n_nlat = 6, &
+ n_elon = 7, &
+ n_elat = 8, &
+
  n_tmask = 1, &
  n_umask = 2, &
  n_nmask = 3, &
@@ -665,6 +686,7 @@ module ice_history_shared
  n_sidmasslat, &
  n_sndmasssnf, &
  n_sndmassmelt, &
+ n_sndmassdyn, &
  n_siflswdtop, &
  n_siflswutop, &
  n_siflswdbot, &
@@ -735,18 +757,22 @@ subroutine construct_filename(ncfile,suffix,ns)
  dt
  use ice_restart_shared, only: lenstr
 
- character (char_len_long), intent(inout) :: ncfile
- character (len=2), intent(in) :: suffix
+ character (len=*), intent(inout) :: ncfile
+ character (len=*), intent(in) :: suffix
  integer (kind=int_kind), intent(in) :: ns
 
  integer (kind=int_kind) :: iyear, imonth, iday, isec
- character (len=1) :: cstream
+ integer (kind=int_kind) :: n
+ character (len=char_len) :: cstream
+ character (len=char_len_long), save :: ncfile_last(max_nstrm) = 'UnDefineD'
  character(len=*), parameter :: subname = '(construct_filename)'
 
  iyear = myear
  imonth = mmonth
  iday = mday
  isec = int(msec - dt,int_kind)
+ cstream = ''
+ if (hist_suffix(ns) /= 'x') cstream = hist_suffix(ns)
 
  ! construct filename
  if (write_ic) then
@@ -771,9 +797,6 @@ subroutine construct_filename(ncfile,suffix,ns)
  endif
  endif
 
- cstream = ''
- if (hist_suffix(ns) /= 'x') cstream = hist_suffix(ns)
-
  if (hist_avg(ns)) then ! write averaged data
  if (histfreq(ns) == '1' .and. histfreq_n(ns) == 1) then ! timestep
  write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') &
@@ -809,6 +832,25 @@ subroutine construct_filename(ncfile,suffix,ns)
 
  endif
 
+ ! Check whether the filename is already in use.
+ ! Same filename in multiple history streams leads to files being overwritten (not good).
+ ! The current filename convention means we just have to check latest filename, 
+ ! not all filenames ever generated because of use of current model date/time in filename.
+
+ ! write(nu_diag,'(2a,i2,1x,a)') subname, 'debug ncfile= ',ns,trim(ncfile)
+ do n = 1,max_nstrm
+ ! write(nu_diag,'(2a,i2,1x,a)') subname, 'debug nfile_last= ',n,trim(ncfile_last(n))
+ if (ncfile == ncfile_last(n)) then
+ write(nu_diag,*) subname,' history stream = ',ns
+ write(nu_diag,*) subname,' history filename = ',trim(ncfile)
+ write(nu_diag,*) subname,' filename in use for stream ',n
+ write(nu_diag,*) subname,' filename for stream ',trim(ncfile_last(n))
+ write(nu_diag,*) subname,' Use namelist hist_suffix so history filenames are unique'
+ call abort_ice(subname//' ERROR: history filename already used for another history stream '//trim(ncfile))
+ endif
+ enddo
+ ncfile_last(ns) = ncfile
+
  end subroutine construct_filename
 
 !=======================================================================
@@ -869,7 +911,7 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, &
  if(present(mask_ice_free_points)) l_mask_ice_free_points = mask_ice_free_points
 
  if (histfreq(ns) == 'x') then
- call abort_ice(subname//'ERROR: define_hist_fields has histfreq x')
+ call abort_ice(subname//' ERROR: define_hist_fields has histfreq x')
  endif
 
  if (ns == 1) id(:) = 0
@@ -879,7 +921,7 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, &
  if (vhistfreq(ns1:ns1) == histfreq(ns)) then
 
  if (ns1 > 1 .and. index(vhistfreq(1:ns1-1),'x') /= 0) then
- call abort_ice(subname//'ERROR: history frequency variable f_' // vname // ' can''t contain ''x'' along with active frequencies')
+ call abort_ice(subname//' ERROR: history frequency variable f_' // vname // ' can''t contain ''x'' along with active frequencies')
  endif
 
  num_avail_hist_fields_tot = num_avail_hist_fields_tot + 1
@@ -909,7 +951,7 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, &
  write(nu_diag,*) subname,' num_avail_hist_fields_tot = ',num_avail_hist_fields_tot
  write(nu_diag,*) subname,' max_avail_hist_fields = ',max_avail_hist_fields
  endif
- call abort_ice(subname//'ERROR: Need in computation of max_avail_hist_fields')
+ call abort_ice(subname//' ERROR: Need in computation of max_avail_hist_fields')
  endif
 
  if (num_avail_hist_fields_tot /= &
@@ -925,7 +967,7 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, &
  if (my_task == master_task) then
  write(nu_diag,*) subname,' num_avail_hist_fields_tot = ',num_avail_hist_fields_tot
  endif
- call abort_ice(subname//'ERROR: in num_avail_hist_fields')
+ call abort_ice(subname//' ERROR: in num_avail_hist_fields')
  endif
 
  id(ns) = num_avail_hist_fields_tot

cicecore/cicedyn/general/ice_flux.F90

@@ -107,6 +107,7 @@ module ice_flux
  strintyE, & ! divergence of internal ice stress, y at E points (N/m^2)
  daidtd , & ! ice area tendency due to transport (1/s)
  dvidtd , & ! ice volume tendency due to transport (m/s)
+ dvsdtd , & ! snow volume tendency due to transport (m/s)
  dagedtd , & ! ice age tendency due to transport (s/s)
  dardg1dt, & ! rate of area loss by ridging ice (1/s)
  dardg2dt, & ! rate of area gain by new ridges (1/s)
@@ -319,6 +320,7 @@ module ice_flux
  dsnow, & ! change in snow thickness (m/step-->cm/day)
  daidtt, & ! ice area tendency thermo. (s^-1)
  dvidtt, & ! ice volume tendency thermo. (m/s)
+ dvsdtt, & ! snow volume tendency thermo. (m/s)
  dagedtt,& ! ice age tendency thermo. (s/s)
  mlt_onset, &! day of year that sfc melting begins
  frz_onset, &! day of year that freezing begins (congel or frazil)
@@ -419,6 +421,7 @@ subroutine alloc_flux
  strintyU (nx_block,ny_block,max_blocks), & ! divergence of internal ice stress, y (N/m^2)
  daidtd (nx_block,ny_block,max_blocks), & ! ice area tendency due to transport (1/s)
  dvidtd (nx_block,ny_block,max_blocks), & ! ice volume tendency due to transport (m/s)
+ dvsdtd (nx_block,ny_block,max_blocks), & ! snow volume tendency due to transport (m/s)
  dagedtd (nx_block,ny_block,max_blocks), & ! ice age tendency due to transport (s/s)
  dardg1dt (nx_block,ny_block,max_blocks), & ! rate of area loss by ridging ice (1/s)
  dardg2dt (nx_block,ny_block,max_blocks), & ! rate of area gain by new ridges (1/s)
@@ -530,6 +533,7 @@ subroutine alloc_flux
  dsnow (nx_block,ny_block,max_blocks), & ! change in snow thickness (m/step-->cm/day)
  daidtt (nx_block,ny_block,max_blocks), & ! ice area tendency thermo. (s^-1)
  dvidtt (nx_block,ny_block,max_blocks), & ! ice volume tendency thermo. (m/s)
+ dvsdtt (nx_block,ny_block,max_blocks), & ! snow volume tendency thermo. (m/s)
  dagedtt (nx_block,ny_block,max_blocks), & ! ice age tendency thermo. (s/s)
  mlt_onset (nx_block,ny_block,max_blocks), & ! day of year that sfc melting begins
  frz_onset (nx_block,ny_block,max_blocks), & ! day of year that freezing begins (congel or frazil)
@@ -918,7 +922,7 @@ end subroutine init_flux_ocn
 
  subroutine init_history_therm
 
- use ice_state, only: aice, vice, trcr
+ use ice_state, only: aice, vice, vsno, trcr
  use ice_arrays_column, only: &
  hfreebd, hdraft, hridge, distrdg, hkeel, dkeel, lfloe, dfloe, &
  Cdn_atm_skin, Cdn_atm_floe, Cdn_atm_pond, Cdn_atm_rdg, &
@@ -965,6 +969,7 @@ subroutine init_history_therm
  meltl (:,:,:) = c0
  daidtt (:,:,:) = aice(:,:,:) ! temporary initial area
  dvidtt (:,:,:) = vice(:,:,:) ! temporary initial volume
+ dvsdtt (:,:,:) = vsno(:,:,:) ! temporary initial volume
  if (tr_iage) then
  dagedtt(:,:,:) = trcr(:,:,nt_iage,:) ! temporary initial age
  else
@@ -1022,7 +1027,7 @@ end subroutine init_history_therm
 
  subroutine init_history_dyn
 
- use ice_state, only: aice, vice, trcr, strength, divu, shear, vort
+ use ice_state, only: aice, vice, vsno, trcr, strength, divu, shear, vort
  use ice_grid, only: grid_ice
 
  logical (kind=log_kind) :: &
@@ -1061,6 +1066,7 @@ subroutine init_history_dyn
  opening (:,:,:) = c0
  daidtd (:,:,:) = aice(:,:,:) ! temporary initial area
  dvidtd (:,:,:) = vice(:,:,:) ! temporary initial volume
+ dvsdtd (:,:,:) = vsno(:,:,:) ! temporary initial volume
  if (tr_iage) &
  dagedtd (:,:,:) = trcr(:,:,nt_iage,:) ! temporary initial age
  fmU (:,:,:) = c0