Update insert_sic feature

cicecore/cicedyn/general/ice_init.F90

@@ -84,7 +84,7 @@ subroutine input_data
           restart, restart_ext, restart_coszen, use_restart_time, &
           runtype, restart_file, restart_dir, runid, pointer_file, &
           restart_format, restart_rearranger, restart_iotasks, restart_root, &
-          restart_stride, restart_deflate, restart_chunksize, insert_sic
+          restart_stride, restart_deflate, restart_chunksize, restart_mod
       use ice_history_shared, only: &
           history_precision, hist_avg, history_format, history_file, incond_file, &
           history_dir, incond_dir, version_name, history_rearranger, &
@@ -196,7 +196,7 @@ subroutine input_data
         ice_ic,         restart,        restart_dir,     restart_file,  &
         restart_ext,    use_restart_time, restart_format, lcdf64,       &
         restart_root,   restart_stride, restart_iotasks, restart_rearranger, &
-        restart_deflate, restart_chunksize,                             &
+        restart_deflate, restart_chunksize, restart_mod,                &
         pointer_file,   dumpfreq,       dumpfreq_n,      dump_last,     &
         diagfreq,       diag_type,      diag_file,       history_format,&
         history_root,   history_stride, history_iotasks, history_rearranger, &
@@ -294,7 +294,7 @@ subroutine input_data
         atm_data_dir,   ocn_data_dir,    bgc_data_dir,                  &
         atm_data_format, ocn_data_format, rotate_wind,                  &
         oceanmixed_file, atm_data_version,semi_implicit_Tsfc,           &
-        vapor_flux_correction, insert_sic
+        vapor_flux_correction
 
       !-----------------------------------------------------------------
       ! default values
@@ -573,7 +573,7 @@ subroutine input_data
       restore_ocn     = .false.   ! restore sst if true
       trestore        = 90        ! restoring timescale, days (0 instantaneous)
       restore_ice     = .false.   ! restore ice state on grid edges if true
-      insert_sic      = .false.   ! if true, on restart update concentration from a file
+      restart_mod     = 'none'    ! restart modification option
       debug_forcing   = .false.   ! true writes diagnostics for input forcing
 
       latpnt(1) =  90._dbl_kind   ! latitude of diagnostic point 1 (deg)
@@ -992,6 +992,7 @@ subroutine input_data
       call broadcast_scalar(restart_rearranger,   master_task)
       call broadcast_scalar(restart_deflate,      master_task)
       call broadcast_array(restart_chunksize,     master_task)
+      call broadcast_scalar(restart_mod,          master_task)
       call broadcast_scalar(lcdf64,               master_task)
       call broadcast_scalar(pointer_file,         master_task)
       call broadcast_scalar(ice_ic,               master_task)
@@ -1171,7 +1172,6 @@ subroutine input_data
       call broadcast_scalar(restore_ocn,          master_task)
       call broadcast_scalar(trestore,             master_task)
       call broadcast_scalar(restore_ice,          master_task)
-      call broadcast_scalar(insert_sic,           master_task)
       call broadcast_scalar(debug_forcing,        master_task)
       call broadcast_array (latpnt(1:2),          master_task)
       call broadcast_array (lonpnt(1:2),          master_task)
@@ -2622,6 +2622,7 @@ subroutine input_data
          write(nu_diag,1011) ' restart          = ', restart
          write(nu_diag,1031) ' restart_dir      = ', trim(restart_dir)
          write(nu_diag,1011) ' restart_ext      = ', restart_ext
+         write(nu_diag,1031) ' restart_mod      = ', trim(restart_mod)
          write(nu_diag,1011) ' restart_coszen   = ', restart_coszen
          write(nu_diag,1031) ' restart_format   = ', trim(restart_format)
          write(nu_diag,1021) ' restart_deflate  = ', restart_deflate
@@ -2685,7 +2686,6 @@ subroutine input_data
          write(nu_diag,1011) ' restore_ice      = ', restore_ice
          if (restore_ice .or. restore_ocn) &
          write(nu_diag,1021) ' trestore         = ', trestore
-         write(nu_diag,1011) ' insert_sic       = ', insert_sic
 
          write(nu_diag,*) ' '
          write(nu_diag,'(a31,2f8.2)') 'Diagnostic point 1: lat, lon =', &

cicecore/cicedyn/infrastructure/ice_restart_driver.F90

@@ -25,7 +25,7 @@ module ice_restart_driver
           field_loc_Eface, field_loc_Nface, &
           field_type_scalar, field_type_vector
       use ice_restart_shared, only: restart_dir, pointer_file, &
-          runid, use_restart_time, lenstr, restart_coszen, insert_sic
+          runid, use_restart_time, lenstr, restart_coszen, restart_mod
       use ice_restart
       use ice_exit, only: abort_ice
       use ice_fileunits, only: nu_diag, nu_rst_pointer, nu_restart, nu_dump
@@ -299,7 +299,7 @@ subroutine restartfile (ice_ic)
           aice0, aicen, vicen, vsnon, trcrn, aice_init, uvel, vvel, &
           uvelE, vvelE, uvelN, vvelN, &
           trcr_base, nt_strata, n_trcr_strata
-      use icepack_itd, only: cleanup_itd  !for insert_sic
+      use icepack_itd, only: cleanup_itd  !for restart_mod
       use ice_arrays_column, only: first_ice, hin_max
       use ice_flux, only: fpond, fresh, fsalt, fhocn
       use ice_flux_bgc, only: faero_ocn, fiso_ocn, flux_bio
@@ -729,12 +729,21 @@ subroutine restartfile (ice_ic)
       !-----------------------------------------------------------------
       ! update concentration from a file
       !-----------------------------------------------------------------
-      if (insert_sic) then
+      if (restart_mod /= "none") then
 
-         if (my_task == master_task) &
-              write(nu_diag,*) "Inserting sic"
+         select case (trim(restart_mod))
+
+         case('adjust_aice')
+            call direct_adjust_aice
+
+         case('adjust_aice_test')
+            call direct_adjust_aice(test=.true.)
 
-         call direct_insert_sic
+         case default
+            call abort_ice(subname//'ERROR: unsupported restart_mod='//trim(restart_mod), &
+               file=__FILE__, line=__LINE__)
+
+         end select
 
          !-----------------------------------------------------------------
          ! Ensure ice is binned in correct categories
@@ -785,7 +794,7 @@ subroutine restartfile (ice_ic)
             enddo     ! j
          enddo        ! iblk
 
-      endif !insert_sic
+      endif !restart_mod
 
       end subroutine restartfile
 
@@ -1173,14 +1182,14 @@ end subroutine restartfile_v4
 !
 !  Alan J. Wallcraft, COAPS/FSU, Nov 2024
 
-      subroutine direct_insert_sic
+      subroutine direct_adjust_aice(test)
 
         use ice_blocks, only: nghost, nx_block, ny_block
         use ice_domain, only: nblocks
         use ice_domain_size, only: nilyr, nslyr, ncat, max_blocks
         use ice_grid, only: tmask
         use ice_communicate, only: my_task, master_task
-        use ice_constants, only: c0, c1, c4, p5, p2, p1, p01, p001, &
+        use ice_constants, only: c0, c1, c4, c100, p5, p2, p1, p01, p001, &
              field_loc_center, field_loc_NEcorner, &
              field_type_scalar, field_type_vector
         use ice_fileunits, only: nu_diag
@@ -1197,6 +1206,9 @@ subroutine direct_insert_sic
         ! use icepack_mushy_physics, only: enthalpy_mush
         use icepack_intfc, only: icepack_init_trcr
 
+        logical(kind=log_kind), optional, intent(in) :: &
+             test         ! use internally generated aice
+
         ! --- local variables
         real(kind=dbl_kind) :: &
              q     ,    & ! scale factor
@@ -1213,17 +1225,21 @@ subroutine direct_insert_sic
              vsnon_old, & ! old value of snow volume to check when adding ice
              Tsfc         ! surface temp.
         integer (kind=int_kind) :: &
-             fid             ! file id for netCDF file
+             fid          ! file id for netCDF file
         integer (kind=int_kind) :: &
-             i, j, k, n, iblk  ! counting  indices
+             i, j, k, n, iblk ! counting  indices
         logical (kind=log_kind) :: &
-             diag  ! diagnostic output
+             diag,      & ! diagnostic output
+             ltest        ! local value of test argument
         real (kind=dbl_kind), dimension (nx_block,ny_block,max_blocks) :: &
              work1
         real (kind=dbl_kind), dimension(nilyr) :: &
-             qin             ! ice enthalpy (J/m3)
+             qin          ! ice enthalpy (J/m3)
         real (kind=dbl_kind), dimension(nslyr) :: &
-             qsn             ! snow enthalpy (J/m3)
+             qsn          ! snow enthalpy (J/m3)
+        character(len=char_len_long) :: &
+             aice_filename, &! filename to read in
+             aice_fldname    ! fieldname to read in
 
         ! parameters from icepack
         real (kind=dbl_kind) :: &
@@ -1232,9 +1248,15 @@ subroutine direct_insert_sic
         integer (kind=int_kind) :: &
              nt_Tsfc, nt_sice, nt_qice, nt_qsno, &
              ktherm
-        character(len=*), parameter :: subname = '(direct_insert_sic)'
+        character(len=*), parameter :: subname = '(direct_adjust_aice)'
 
         diag = .true.
+        ltest = .false.
+        if (present(test)) then
+           ltest = test
+        endif
+        aice_filename = trim(restart_dir)//'sic.nc'
+        aice_fldname = 'sic'
 
         ! get parameters from icepack
         call icepack_query_parameters( &
@@ -1258,15 +1280,22 @@ subroutine direct_insert_sic
         if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
              file=__FILE__, line=__LINE__)
 
-        if (my_task == master_task) then
-           write(nu_diag,*) "Direct_insert_sic from:"
-           write(nu_diag,*) trim(restart_dir)//'/sic.nc'
-        endif !master_task
-        call ice_open_nc(trim(restart_dir)//'/sic.nc', fid)
-        call ice_read_nc(fid,1,'sic',work1,diag, &
-                         field_loc=field_loc_center, &
-                         field_type=field_type_scalar)
-        call ice_close_nc(fid)
+        if (ltest) then
+           if (my_task == master_task) then
+              write(nu_diag,*) subname//" direct_adjust_aice rounding to nearest 1/100th"
+           endif
+           work1 = nint(aice*c100)/c100   ! round to nearest 1/100th
+        else
+           if (my_task == master_task) then
+              write(nu_diag,*) subname//" direct_adjust_aice from "//trim(aice_filename)
+           endif
+
+           call ice_open_nc(trim(aice_filename), fid)
+           call ice_read_nc(fid,1,trim(aice_fldname),work1,diag, &
+                            field_loc=field_loc_center, &
+                            field_type=field_type_scalar)
+           call ice_close_nc(fid)
+        endif
 
         edge_om = p2  ! nominal ice edge zone
         diff_om = p1  ! allowed model vs obs difference
@@ -1442,7 +1471,7 @@ subroutine direct_insert_sic
            enddo                ! i
         enddo                   ! iblk
 
-      end subroutine direct_insert_sic
+      end subroutine direct_adjust_aice
 
 !=======================================================================
 

cicecore/shared/ice_restart_shared.F90

@@ -14,9 +14,6 @@ module ice_restart_shared
          restart_coszen, &   ! if true, read/write coszen
          use_restart_time ! if true, use time written in core restart file
 
-     logical(kind=log_kind), public :: &
-         insert_sic        ! if true, update restart concentation from a file
-
       character (len=char_len), public :: &
          runtype           ! initial, continue, hybrid, branch
 
@@ -30,6 +27,7 @@ module ice_restart_shared
 
       character (len=char_len), public :: &
          restart_format      , & ! format of restart files 'nc'
+         restart_mod         , & ! restart modification option, "none", "adjust_aice"
          restart_rearranger      ! restart file rearranger, box or subset for pio
 
       integer (kind=int_kind), public :: &

configuration/scripts/ice_in

@@ -21,6 +21,7 @@
     restart_stride = -99
     restart_deflate = 0
     restart_chunksize = 0, 0
+    restart_mod    = 'none'
     lcdf64         = .false.
     numin          = 21
     numax          = 89

configuration/scripts/options/set_nml.restaicetest

@@ -0,0 +1,2 @@
+restart_mod = "adjust_aice_test"
+

configuration/scripts/tests/base_suite.ts

@@ -86,3 +86,4 @@ restart        gx3     4x4        diag1,gx3ncarbulk,short
 smoke          gx3     4x1        calcdragio
 restart        gx3     4x2        atmbndyconstant
 restart        gx3     4x2        atmbndymixed
+smoke          gx3    12x2        diag1,run5day,restaicetest,debug

doc/source/cice_index.rst

@@ -595,6 +595,7 @@ section :ref:`tabnamelist`.
    "restart_file", "restart file prefix", ""
    "restart_format", "restart file format", ""
    "restart_iotasks", "restart output total number of tasks used", ""
+   "restart_mod", "restart modification mode", ""
    "restart_rearranger", "restart output io rearranger method", ""
    "restart_root", "restart output io root task id", ""
    "restart_stride", "restart output io task stride", ""

doc/source/developer_guide/dg_assim.rst

@@ -0,0 +1,37 @@
+:tocdepth: 3
+
+.. _dataassimilation:
+
+Data Assimilation
+======================
+
+Data assimilation (DA) is the scientific process of combining external
+data with numerical model forecasts.  There are several ways this can
+be done including by adjusting the model initial conditions (internally
+or externally) or adjusting the model solution as it evolves in time.
+Various data assimilation options are being introduced in CICE and are
+described below.
+
+.. _restartmod:
+
+Data Assimilation on restart
+------------------------------------
+
+The namelist variable, ``restart_mod``, specifies the restart DA mode.  
+By default, this namelist value is set to ``none`` which disables the feature.  
+The current active options are ``adjust_aice`` and ``adjust_aice_test``.
+
+With ``adjust_aice`` and ``adjust_aice_test``, the category averaged aice 
+value is modified at restart to specified values using the method implemented in 
+**cicecore/cicedyn/infrastructure/ice_restart_driver.F90** subroutine
+**direct_adjust_aice**.  This method adjusts aice, vice, vsno, qice, and
+sice in all categories to be consistent with the category average aice
+specified.  It also adjusts several thermodynamic variables such as 
+temperature and salinity (see :cite:`Posey15`).  
+``adjust_aice`` reads in a sea ice concentration
+field from an external file.  The field is currently hardwired to 'sic' and the
+file is currently hardwired to 'sic.nc'.  The field must be on the model grid.
+``adjust_aice_test`` modifies the
+aice field read on restart internally.  The current implementation rounds
+the aice values read at restart to the nearest 1/100th.  This mode exists
+primarily to test the feature.

doc/source/developer_guide/index.rst

@@ -16,6 +16,7 @@ Developer Guide
    dg_infra.rst
    dg_driver.rst
    dg_forcing.rst
+   dg_assim.rst
    dg_icepack.rst
    dg_scripts.rst
    dg_tools.rst

doc/source/master_list.bib

@@ -1041,7 +1041,7 @@ @incollection{Arakawa77
 
 @article{Horvat15,
    author = "C. Horvat and E. Tziperman",
-   journal = {The Cryosphere},
+   journal = TC,
    number = {6},
    pages = {2119-2134},
    title = "{A prognostic model of the sea-ice floe size and thickness distribution}",
@@ -1113,6 +1113,16 @@ @Article{Tsujino18
   pages = {79-139},
   url = {http://dx.doi.org/10.1016/j.ocemod.2018.07.002}
 }
+
+@Article{Posey15,
+  author = "P.G. Posey and E.J. Metzger and A.J. Wallcraft and D.A. Hebert and R.A. Allard and O.M. Smedstad and M.W. Phelps and F. Fetterer and J.S. Stewart and W.N. Meier and S.R. Helfrich",
+  title = "{Improving Arctic sea ice edge forecasts by assimilating high horizontal resolution sea ice concentration data into the US Navy's ice forecast system}",
+  journal = TC,
+  year = {2015},
+  volume = {9},
+  pages = {1735-1745},
+  url = {https://doi.org/10.5194/tc-9-1735-2015}
+}
 %  **********************************************
 %  For new entries, see example entry in BIB_TEMPLATE.txt
 %  **********************************************

doc/source/user_guide/ug_case_settings.rst

@@ -266,6 +266,9 @@ setup_nml
    "", "``pnetcdf2``", "write restart files with pnetcdf cdf2 (netcdf3-64bit-offset) format", ""
    "", "``pnetcdf5``", "write restart files with pnetcdf cdf5 (netcdf3-64bit-data) format", ""
    "``restart_iotasks``", "integer", "pe io tasks for restart output with restart_root and restart_stride (PIO only), -99=internal default", "-99"
+   "``restart_mod``", "``adjust_aice``", "adjust aice on restart read from file", "none"
+   "", "``adjust_aice_test``", "adjust aice on restart read rounding", ""
+   "", "``none``", "no modification of restart at read", ""
    "``restart_rearranger``", "``box``", "box io rearranger option for restart output (PIO only)", "default"
    "", "``default``", "internal default io rearranger option for restart output", ""
    "", "``subset``", "subset io rearranger option for restart output", ""