diff --git a/columnphysics/icepack_parameters.F90 b/columnphysics/icepack_parameters.F90 index f875409b5..e76b23b09 100644 --- a/columnphysics/icepack_parameters.F90 +++ b/columnphysics/icepack_parameters.F90 @@ -161,6 +161,11 @@ module icepack_parameters ! only for use with tr_aero or tr_zaero conserv_check = .false. ! if true, do conservations checks and abort + character(len=char_len), public :: & + congel_freeze = 'two-step' ! congelation computation + ! 'two-step' = original formulation + ! 'one-step' = Plante et al, The Cryosphere, 2024 + character(len=char_len), public :: & tfrz_option = 'mushy' ! form of ocean freezing temperature ! 'minus1p8' = -1.8 C @@ -465,7 +470,7 @@ subroutine icepack_init_parameters( & atmbndy_in, calc_strair_in, formdrag_in, highfreq_in, natmiter_in, & atmiter_conv_in, calc_dragio_in, & tfrz_option_in, kitd_in, kcatbound_in, hs0_in, frzpnd_in, & - saltflux_option_in, & + saltflux_option_in, congel_freeze_in, & floeshape_in, wave_spec_in, wave_spec_type_in, nfreq_in, & dpscale_in, rfracmin_in, rfracmax_in, pndaspect_in, hs1_in, hp1_in, & bgc_flux_type_in, z_tracers_in, scale_bgc_in, solve_zbgc_in, & @@ -575,15 +580,20 @@ subroutine icepack_init_parameters( & phi_i_mushy_in ! liquid fraction of congelation ice character(len=*), intent(in), optional :: & - tfrz_option_in ! form of ocean freezing temperature - ! 'minus1p8' = -1.8 C - ! 'linear_salt' = -depressT * sss - ! 'mushy' conforms with ktherm=2 + congel_freeze_in ! congelation computation + ! 'two-step' = original formulation + ! 'one-step' = Plante et al, The Cryosphere, 2024 character(len=*), intent(in), optional :: & - saltflux_option_in ! Salt flux computation - ! 'constant' reference value of ice_ref_salinity - ! 'prognostic' prognostic salt flux + tfrz_option_in ! form of ocean freezing temperature + ! 'minus1p8' = -1.8 C + ! 'linear_salt' = -depressT * sss + ! 'mushy' conforms with ktherm=2 + + character(len=*), intent(in), optional :: & + saltflux_option_in ! Salt flux computation + ! 'constant' reference value of ice_ref_salinity + ! 'prognostic' prognostic salt flux !----------------------------------------------------------------------- ! Parameters for radiation @@ -949,6 +959,7 @@ subroutine icepack_init_parameters( & if (present(highfreq_in) ) highfreq = highfreq_in if (present(natmiter_in) ) natmiter = natmiter_in if (present(atmiter_conv_in) ) atmiter_conv = atmiter_conv_in + if (present(congel_freeze_in) ) congel_freeze = congel_freeze_in if (present(tfrz_option_in) ) tfrz_option = tfrz_option_in if (present(saltflux_option_in) ) saltflux_option = saltflux_option_in if (present(kitd_in) ) kitd = kitd_in @@ -1192,7 +1203,7 @@ subroutine icepack_query_parameters( & atmbndy_out, calc_strair_out, formdrag_out, highfreq_out, natmiter_out, & atmiter_conv_out, calc_dragio_out, & tfrz_option_out, kitd_out, kcatbound_out, hs0_out, frzpnd_out, & - saltflux_option_out, & + saltflux_option_out, congel_freeze_out, & floeshape_out, wave_spec_out, wave_spec_type_out, nfreq_out, & dpscale_out, rfracmin_out, rfracmax_out, pndaspect_out, hs1_out, hp1_out, & bgc_flux_type_out, z_tracers_out, scale_bgc_out, solve_zbgc_out, & @@ -1311,16 +1322,21 @@ subroutine icepack_query_parameters( & phi_i_mushy_out ! liquid fraction of congelation ice character(len=*), intent(out), optional :: & - tfrz_option_out ! form of ocean freezing temperature - ! 'minus1p8' = -1.8 C - ! 'constant' = Tocnfrz - ! 'linear_salt' = -depressT * sss - ! 'mushy' conforms with ktherm=2 + congel_freeze_out ! congelation computation + ! 'two-step' = original formulation + ! 'one-step' = Plante et al, The Cryosphere, 2024 + + character(len=*), intent(out), optional :: & + tfrz_option_out ! form of ocean freezing temperature + ! 'minus1p8' = -1.8 C + ! 'constant' = Tocnfrz + ! 'linear_salt' = -depressT * sss + ! 'mushy' conforms with ktherm=2 character(len=*), intent(out), optional :: & - saltflux_option_out ! Salt flux computation - ! 'constant' reference value of ice_ref_salinity - ! 'prognostic' prognostic salt flux + saltflux_option_out ! Salt flux computation + ! 'constant' reference value of ice_ref_salinity + ! 'prognostic' prognostic salt flux !----------------------------------------------------------------------- @@ -1719,6 +1735,7 @@ subroutine icepack_query_parameters( & if (present(highfreq_out) ) highfreq_out = highfreq if (present(natmiter_out) ) natmiter_out = natmiter if (present(atmiter_conv_out) ) atmiter_conv_out = atmiter_conv + if (present(congel_freeze_out) ) congel_freeze_out = congel_freeze if (present(tfrz_option_out) ) tfrz_option_out = tfrz_option if (present(saltflux_option_out) ) saltflux_option_out = saltflux_option if (present(kitd_out) ) kitd_out = kitd @@ -1926,6 +1943,7 @@ subroutine icepack_write_parameters(iounit) write(iounit,*) " highfreq = ", highfreq write(iounit,*) " natmiter = ", natmiter write(iounit,*) " atmiter_conv = ", atmiter_conv + write(iounit,*) " congel_freeze = ", trim(congel_freeze) write(iounit,*) " tfrz_option= ", trim(tfrz_option) write(iounit,*) " saltflux_option = ", trim(saltflux_option) write(iounit,*) " kitd = ", kitd diff --git a/columnphysics/icepack_therm_vertical.F90 b/columnphysics/icepack_therm_vertical.F90 index 7cd9f2533..233e38546 100644 --- a/columnphysics/icepack_therm_vertical.F90 +++ b/columnphysics/icepack_therm_vertical.F90 @@ -27,7 +27,7 @@ module icepack_therm_vertical use icepack_parameters, only: ustar_min, fbot_xfer_type, formdrag, calc_strair use icepack_parameters, only: rfracmin, rfracmax, dpscale, frzpnd, snwgrain, snwlvlfac use icepack_parameters, only: phi_i_mushy, floeshape, floediam, use_smliq_pnd, snwredist - use icepack_parameters, only: saltflux_option + use icepack_parameters, only: saltflux_option, congel_freeze use icepack_parameters, only: icepack_chkoptargflag use icepack_tracers, only: tr_iage, tr_FY, tr_aero, tr_pond, tr_fsd, tr_iso @@ -46,6 +46,7 @@ module icepack_therm_vertical use icepack_mushy_physics, only: icepack_mushy_temperature_mush use icepack_mushy_physics, only: liquidus_temperature_mush use icepack_mushy_physics, only: icepack_enthalpy_mush, enthalpy_of_melting + use icepack_mushy_physics, only: enthalpy_mush_liquid_fraction, enthalpy_brine use icepack_aerosol, only: update_aerosol use icepack_isotope, only: update_isotope @@ -85,7 +86,7 @@ subroutine thermo_vertical (nilyr, nslyr, & fsnow, fpond, & fbot, Tbot, & Tsnice, sss, & - rsnw, & + sst, rsnw, & lhcoef, shcoef, & fswsfc, fswint, & Sswabs, Iswabs, & @@ -164,6 +165,7 @@ subroutine thermo_vertical (nilyr, nslyr, & real (kind=dbl_kind), intent(in) :: & fbot , & ! ice-ocean heat flux at bottom surface (W/m^2) Tbot , & ! ice bottom surface temperature (deg C) + sst , & ! sea surface temperature (C) sss ! ocean salinity ! coupler fluxes to atmosphere @@ -403,6 +405,7 @@ subroutine thermo_vertical (nilyr, nslyr, & congel, snoice, & mlt_onset, frz_onset, & zSin, sss, & + sst, & dsnow, rsnw) if (icepack_warnings_aborted(subname)) return @@ -1013,6 +1016,7 @@ subroutine thickness_changes (nilyr, nslyr, & congel, snoice, & mlt_onset, frz_onset,& zSin, sss, & + sst, & dsnow, rsnw) integer (kind=int_kind), intent(in) :: & @@ -1082,6 +1086,7 @@ subroutine thickness_changes (nilyr, nslyr, & zSin ! ice layer salinity (ppt) real (kind=dbl_kind), intent(in) :: & + sst , & ! sea surface temperature (C) sss ! ocean salinity (PSU) ! local variables @@ -1131,9 +1136,9 @@ subroutine thickness_changes (nilyr, nslyr, & qmlt ! enthalpy of melted ice (J m-3) = zero in BL99 formulation real (kind=dbl_kind) :: & - qbotm , & - qbotp , & - qbot0 , & + qbotm , & ! enthalpy of newly formed congelation ice + qbotp , & ! enthalpy needed to grow new congelation ice (includes ocean enthalpy) + qbotw , & ! enthalpy transferred to ocean during congelation freezing mass , & ! total mass from snow density tracers (kg/m^2) massi , & ! ice mass change factor tmp1 ! temporary scalar @@ -1254,15 +1259,22 @@ subroutine thickness_changes (nilyr, nslyr, & if (ktherm == 2) then - qbotm = icepack_enthalpy_mush(Tbot, sss) - qbotp = -Lfresh * rhoi * (c1 - phi_i_mushy) - qbot0 = qbotm - qbotp - - dhi = ebot_gro / qbotp ! dhi > 0 - + if (congel_freeze == 'one-step') then + ! Plante et al., The Cryosphere, 18, 1685-1708, 2024 + qbotm = enthalpy_mush_liquid_fraction(Tbot, phi_i_mushy) + qbotw = enthalpy_brine(sst) + qbotp = qbotm - qbotw + dhi = ebot_gro / qbotp ! dhi > 0 + hstot = dzi(nilyr)*zSin(nilyr) + dhi*sss*phi_i_mushy + else ! two-step + qbotm = icepack_enthalpy_mush(Tbot, sss) + qbotp = -Lfresh * rhoi * (c1 - phi_i_mushy) + qbotw = qbotm - qbotp + dhi = ebot_gro / qbotp ! dhi > 0 + hstot = dzi(nilyr)*zSin(nilyr) + dhi*sss + endif hqtot = dzi(nilyr)*zqin(nilyr) + dhi*qbotm - hstot = dzi(nilyr)*zSin(nilyr) + dhi*sss - emlt_ocn = emlt_ocn - qbot0 * dhi + emlt_ocn = emlt_ocn - qbotw * dhi else @@ -2653,7 +2665,7 @@ subroutine icepack_step_therm1(dt, ncat, nilyr, nslyr, & fsnow=fsnow, fpond=fpond, & fbot=fbot, Tbot=Tbot, & Tsnice=Tsnice, sss=sss, & - rsnw=rsnw, & + sst=sst, rsnw=rsnw, & lhcoef=lhcoef, shcoef=shcoef, & fswsfc=fswsfcn (n), fswint=fswintn (n), & Sswabs=Sswabsn(:,n), Iswabs=Iswabsn (:,n), & diff --git a/configuration/driver/icedrv_init.F90 b/configuration/driver/icedrv_init.F90 index c1e738ce6..b64284104 100644 --- a/configuration/driver/icedrv_init.F90 +++ b/configuration/driver/icedrv_init.F90 @@ -99,7 +99,7 @@ subroutine input_data natmiter, kitd, kcatbound character (len=char_len) :: shortwave, albedo_type, conduct, fbot_xfer_type, & - cpl_frazil, tfrz_option, saltflux_option, & + cpl_frazil, congel_freeze, tfrz_option, saltflux_option, & frzpnd, atmbndy, wave_spec_type, snwredist, snw_aging_table logical (kind=log_kind) :: sw_redist, use_smliq_pnd, snwgrain, update_ocn_f @@ -174,12 +174,12 @@ subroutine input_data update_ocn_f, l_mpond_fresh, ustar_min, & fbot_xfer_type, oceanmixed_ice, emissivity, & formdrag, highfreq, natmiter, & - atmiter_conv, calc_dragio, & + atmiter_conv, calc_dragio, congel_freeze, & tfrz_option, saltflux_option, ice_ref_salinity, & default_season, wave_spec_type, cpl_frazil, & precip_units, fyear_init, ycycle, & atm_data_type, ocn_data_type, bgc_data_type, & - lateral_flux_type, & + lateral_flux_type, & atm_data_file, ocn_data_file, bgc_data_file, & ice_data_file, & atm_data_format, ocn_data_format, bgc_data_format, & @@ -225,6 +225,7 @@ subroutine input_data dSdt_slow_mode_out=dSdt_slow_mode, & phi_c_slow_mode_out=phi_c_slow_mode, Tliquidus_max_out=Tliquidus_max, & phi_i_mushy_out=phi_i_mushy, conserv_check_out=conserv_check, & + congel_freeze_out=congel_freeze, & tfrz_option_out=tfrz_option, saltflux_option_out=saltflux_option, & ice_ref_salinity_out=ice_ref_salinity, kalg_out=kalg, & fbot_xfer_type_out=fbot_xfer_type, puny_out=puny, & @@ -791,6 +792,7 @@ subroutine input_data write(nu_diag,1005) ' hi_min = ', hi_min write(nu_diag,1030) ' fbot_xfer_type = ', trim(fbot_xfer_type) write(nu_diag,1010) ' oceanmixed_ice = ', oceanmixed_ice + write(nu_diag,1030) ' congel_freeze = ', trim(congel_freeze) write(nu_diag,1030) ' tfrz_option = ', trim(tfrz_option) write(nu_diag,*) ' saltflux_option = ', & trim(saltflux_option) @@ -979,6 +981,7 @@ subroutine input_data dSdt_slow_mode_in=dSdt_slow_mode, & phi_c_slow_mode_in=phi_c_slow_mode, Tliquidus_max_in=Tliquidus_max, & phi_i_mushy_in=phi_i_mushy, conserv_check_in=conserv_check, & + congel_freeze_in=congel_freeze, & tfrz_option_in=tfrz_option, saltflux_option_in=saltflux_option, & ice_ref_salinity_in=ice_ref_salinity, kalg_in=kalg, & fbot_xfer_type_in=fbot_xfer_type, & diff --git a/configuration/scripts/icepack_in b/configuration/scripts/icepack_in index bbf7a5973..4ac97ad81 100644 --- a/configuration/scripts/icepack_in +++ b/configuration/scripts/icepack_in @@ -113,6 +113,7 @@ cpl_frazil = 'fresh_ice_correction' update_ocn_f = .false. l_mpond_fresh = .false. + congel_freeze = 'two-step' tfrz_option = 'mushy' ice_ref_salinity = 4.0 saltflux_option = 'constant' diff --git a/configuration/scripts/machines/env.chicoma_intel b/configuration/scripts/machines/env.chicoma_intel index 182280f95..cb7f12679 100755 --- a/configuration/scripts/machines/env.chicoma_intel +++ b/configuration/scripts/machines/env.chicoma_intel @@ -63,7 +63,7 @@ setenv ICE_MACHINE_WKDIR /lustre/scratch5/$user/ICEPACK_RUNS setenv ICE_MACHINE_INPUTDATA /usr/projects/climate/eclare/DATA/Consortium setenv ICE_MACHINE_BASELINE /lustre/scratch5/$user/ICEPACK_BASELINE setenv ICE_MACHINE_SUBMIT "sbatch " -setenv ICE_MACHINE_ACCT t23_cice +setenv ICE_MACHINE_ACCT t24_cice setenv ICE_MACHINE_QUEUE "debug" setenv ICE_MACHINE_TPNODE 128 # tasks per node setenv ICE_MACHINE_BLDTHRDS 12 diff --git a/configuration/scripts/options/set_nml.congel b/configuration/scripts/options/set_nml.congel new file mode 100644 index 000000000..a71cda48b --- /dev/null +++ b/configuration/scripts/options/set_nml.congel @@ -0,0 +1 @@ + congel_freeze = 'one-step' diff --git a/configuration/scripts/tests/base_suite.ts b/configuration/scripts/tests/base_suite.ts index 43e405dab..3a68d7dd7 100644 --- a/configuration/scripts/tests/base_suite.ts +++ b/configuration/scripts/tests/base_suite.ts @@ -23,6 +23,7 @@ smoke col 1x1 debug,run1year,fluxopenw smoke col 1x1 debug,run1year,dyn,fluxopenw smoke col 1x1 debug,run1year,debug,snicartest smoke col 1x1 debug,run1year,debug,snicar +smoke col 1x1 debug,run1year,debug,congel restart col 1x1 debug restart col 1x1 diag1 restart col 1x1 pondlvl @@ -43,4 +44,5 @@ restart col 1x1 modal restart col 1x1 fluxopenw restart col 1x1 snicartest restart col 1x1 snicar +restart col 1x1 congel diff --git a/doc/source/master_list.bib b/doc/source/master_list.bib index 7ffb0a841..445bf758c 100644 --- a/doc/source/master_list.bib +++ b/doc/source/master_list.bib @@ -32,7 +32,7 @@ @string{JGR @string{JGRO = {J. Geophys. Res. Oceans}} @string{JGRB = {J. Geophys. Res. Biogeo.}} @string{JGRA = {J. Geophys. Res. Atmos.}} -@string{JCT = {J. Comput. Phys.}} +@string{TC = {The Cryosphere}} @string{QJRMS = {Quart. J. Royal Met. Soc.}} @string{GRL = {Geophys. Res. Lett.}} @string{JAS = {J. Atmos. Sci.}} @@ -244,7 +244,7 @@ @Article{Maykut95 @Article{Murray96 author = "R.J. Murray", title = "{Explicit generation of orthogonal grids for ocean models}", - journal = JCT, + journal = JCP, year = {1996}, volume = {126}, pages = {251-273}, @@ -805,13 +805,23 @@ @article{Roberts19 @article{Dang19, author = {Dang, C. and Zender, C. S. and Flanner, M. G.}, title = {Intercomparison and improvement of two-stream shortwave radiative transfer schemes in Earth system models for a unified treatment of cryospheric surfaces}, -journal = {The Cryosphere}, +journal = {TC}, volume = 13, pages = {2325--2343}, doi = {10.5194/tc-13-2325-2019}, year = {2019} } +@article{Plante24, +author = {Plante, M. and J.-F. Lemieux and L. B. Tremblay and A. Tivy and J. Angnatok and F. Roy and G. Smith and F. Dupont and A. K. Turner}, +title = {Using Icepack to reproduce ice mass balance buoy observations in landfast ice: improvements from the mushy-layer thermodynamcis}, +journal = {TC}, +volume = 18, +pages = {1685--1708}, +doi = {10.5194/tc-18-1685-2024}, +year = {20124} +} + % ********************************************** % For new entries, see example entry in BIB_TEMPLATE.txt % ********************************************** diff --git a/doc/source/science_guide/sg_thermo.rst b/doc/source/science_guide/sg_thermo.rst index d12b90cfe..1577a7e5e 100755 --- a/doc/source/science_guide/sg_thermo.rst +++ b/doc/source/science_guide/sg_thermo.rst @@ -1737,9 +1737,21 @@ conductive heat flux at the bottom surface: If ice is melting at the bottom surface, :math:`q` in Equation :eq:`bottom-melting` is the enthalpy of the bottom ice layer. If -ice is growing, :math:`q` is the enthalpy of new ice with temperature -:math:`T_f` and salinity :math:`S_{max}` (``ktherm`` = 1) or ocean surface -salinity (``ktherm`` = 2). This ice is added to the bottom layer. +ice is growing, :math:`q` is the enthalpy of new ice added to the bottom layer with +temperature :math:`T_f` and salinity :math:`S_{max}` (``ktherm`` = 1) for the Bitz99 +:cite:`Bitz99` formulation. The original mushy thermodynamics +(``ktherm`` = 2, ``congel_freeze`` = 'two-step') formed new ice in two steps, first +moving the lower ice boundary into the ocean to form a mushy layer with an initial +liquid fraction :math:`phi_{init}`, then freezing the new ice in the next step. In +this case, the freezing temperature is calculated using the ocean surface salinity +:math:`SSS`. A second mushy thermo method +(``ktherm`` = 2, ``congel_freeze`` = 'one-step') freezes the ice immediately using +the brine salinity :math:`phi_{init} * SSS` :cite:`Plante24`. +In the two-step method, enthalpy not +yet used to freeze the ice is returned to the ocean, causing frazil ice to form +instead. Together, the congelation and frazil ice add up to a similar total amount +of new ice, but the differing freezing mechanism complicates comparisons with observational +data and there is an unnecessary lag in ice-ocean coupling. In general, frazil ice formed in the ocean is added to the thinnest ice category. The new ice is grown in the open water area of the grid cell diff --git a/doc/source/user_guide/ug_case_settings.rst b/doc/source/user_guide/ug_case_settings.rst index 0b7e4fa31..cb5e5cb67 100644 --- a/doc/source/user_guide/ug_case_settings.rst +++ b/doc/source/user_guide/ug_case_settings.rst @@ -344,6 +344,8 @@ forcing_nml "``calc_strair``", "``.false.``", "read wind stress and speed from files", "``.true.``" "", "``.true.``", "calculate wind stress and speed", "" "``calc_Tsfc``", "logical", "calculate surface temperature", "``.true.``" + "``congel_freeze``", "``one-step``", "immediately freeze congelation ice", "``two-step``" + "", "``two-step``", "delayed freezing of congelation ice", "" "``cpl_frazil``", "``external``", "frazil water/salt fluxes are handled outside of Icepack", "``fresh_ice_correction``" "", "``fresh_ice_correction``", "correct fresh-ice frazil water/salt fluxes for mushy physics", "" "", "``internal``", "send full frazil water/salt fluxes for mushy physics", ""