Skip to content

Commit

Permalink
Alternative congelation formulation following Plante et al. 2024 (CIC…
Browse files Browse the repository at this point in the history 
…E-Consortium#494)

Add congel_freeze namelist and add 'one-step' option to the default 'two-step' option.

Namelist flag congel_freeze chooses which formulation to use. The original is ‘two-step’, since only the mushy boundary moves in the first step and the phase change happens in the next step. Plante et al. (‘one-step’) moves the boundary and performs the phase change in a single time step.  Maintain 'two-step' as default for now.
  • Loading branch information
@eclare108213 @dabail10
eclare108213 authored and dabail10 committed Oct 25, 2024
1 parent 65987dd commit 0e7e470
Showing 1 changed file with 26 additions and 14 deletions.
40 changes: 26 additions & 14 deletions columnphysics/icepack_therm_vertical.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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
Expand Down Expand Up @@ -86,7 +87,7 @@ subroutine thermo_vertical (nilyr, nslyr, &
fsnow, fpond, &
fbot, Tbot, &
Tsnice, sss, &
rsnw, &
sst, rsnw, &
lhcoef, shcoef, &
fswsfc, fswint, &
Sswabs, Iswabs, &
Expand Down Expand Up @@ -168,6 +169,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
Expand Down Expand Up @@ -413,6 +415,7 @@ subroutine thermo_vertical (nilyr, nslyr, &
congel, snoice, &
mlt_onset, frz_onset, &
zSin, sss, &
sst, &
dsnow, rsnw)
if (icepack_warnings_aborted(subname)) return

Expand Down Expand Up @@ -1094,6 +1097,7 @@ subroutine thickness_changes (nilyr, nslyr, &
congel, snoice, &
mlt_onset, frz_onset,&
zSin, sss, &
sst, &
dsnow, rsnw)

integer (kind=int_kind), intent(in) :: &
Expand Down Expand Up @@ -1163,6 +1167,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
Expand Down Expand Up @@ -1212,9 +1217,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
Expand Down Expand Up @@ -1335,15 +1340,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

Expand Down Expand Up @@ -2760,7 +2772,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), &
Expand Down

0 comments on commit 0e7e470

Please sign in to comment.