Skip to content

Commit

Permalink
Update icepack_step_therm2, make fsd arguments optional (#440)
Browse files Browse the repository at this point in the history
* Update icepack_step_therm2, make fsd arguments optional

* Update interface documentation

* Update user guide, change references from master to main

* update optional arguments in lower level routines for fsd step_therm2

* update documentation about optional arguments

* update documentation about optional arguments
  • Loading branch information
apcraig authored Jul 6, 2023
1 parent 31e85f2 commit d024340
Show file tree
Hide file tree
Showing 5 changed files with 184 additions and 146 deletions.
190 changes: 112 additions & 78 deletions columnphysics/icepack_therm_itd.F90
Original file line number Diff line number Diff line change
Expand Up @@ -927,15 +927,7 @@ subroutine lateral_melt (dt, ncat, &
meltl , & ! lateral ice melt (m/step-->cm/day)
fzsal ! salt flux from zsalinity (kg/m2/s)

real (kind=dbl_kind), dimension (:), intent(in) :: &
floe_rad_c , & ! fsd size bin centre in m (radius)
floe_binwidth ! fsd size bin width in m (radius)

real (kind=dbl_kind), dimension (:), intent(out) :: &
d_afsd_latm ! change in fsd due to lateral melt (m)

real (kind=dbl_kind), dimension(nbtrcr), &
intent(inout) :: &
real (kind=dbl_kind), dimension(nbtrcr), intent(inout) :: &
flux_bio ! biology tracer flux from layer bgc (mmol/m^2/s)

real (kind=dbl_kind), dimension(:), intent(inout) :: &
Expand All @@ -944,6 +936,13 @@ subroutine lateral_melt (dt, ncat, &
real (kind=dbl_kind), dimension(:), intent(inout) :: &
fiso_ocn ! isotope flux to ocean (kg/m^2/s)

real (kind=dbl_kind), dimension (:), intent(in) :: &
floe_rad_c , & ! fsd size bin centre in m (radius)
floe_binwidth ! fsd size bin width in m (radius)

real (kind=dbl_kind), dimension (:), intent(out) :: &
d_afsd_latm ! change in fsd due to lateral melt (m)

! local variables

integer (kind=int_kind) :: &
Expand All @@ -970,15 +969,14 @@ subroutine lateral_melt (dt, ncat, &
delta_an , & ! change in the ITD
rsiden ! delta_an/aicen

real (kind=dbl_kind), dimension (nfsd,ncat) :: &
real (kind=dbl_kind), dimension (:,:), allocatable :: &
afsdn , & ! floe size distribution tracer
afsdn_init ! initial value

real (kind=dbl_kind), dimension (nfsd) :: &
df_flx, & ! finite difference for FSD
afsd_tmp, d_afsd_tmp

real (kind=dbl_kind), dimension(nfsd+1) :: &
real (kind=dbl_kind), dimension (:), allocatable :: &
df_flx , & ! finite difference for FSD
afsd_tmp , & !
d_afsd_tmp, & !
f_flx !

real (kind=dbl_kind) :: &
Expand All @@ -1002,19 +1000,24 @@ subroutine lateral_melt (dt, ncat, &
G_radialn = c0
delta_an = c0
rsiden = c0
afsdn = c1
afsdn_init = c0
df_flx = c0
f_flx = c0

if (tr_fsd) then
call icepack_cleanup_fsd (ncat, nfsd, trcrn(nt_fsd:nt_fsd+nfsd-1,:))
if (icepack_warnings_aborted(subname)) return

afsdn = trcrn(nt_fsd:nt_fsd+nfsd-1,:)
aicen_init = aicen
afsdn_init = afsdn ! for diagnostics
d_afsd_latm(:) = c0
allocate(afsdn(nfsd,ncat))
allocate(afsdn_init(nfsd,ncat))
allocate(df_flx(nfsd))
allocate(afsd_tmp(nfsd))
allocate(d_afsd_tmp(nfsd))
allocate(f_flx(nfsd+1))

aicen_init = aicen
afsdn = trcrn(nt_fsd:nt_fsd+nfsd-1,:)
afsdn_init = afsdn ! for diagnostics
df_flx = c0
d_afsd_latm = c0
f_flx = c0
end if

if (tr_fsd .and. wlat > puny) then
Expand Down Expand Up @@ -1257,6 +1260,14 @@ subroutine lateral_melt (dt, ncat, &
+ afsdn(k,n)*aicen(n) - afsdn_init(k,n)*aicen_init(n)
end do
end do

deallocate(afsdn)
deallocate(afsdn_init)
deallocate(df_flx)
deallocate(afsd_tmp)
deallocate(d_afsd_tmp)
deallocate(f_flx)

end if

end subroutine lateral_melt
Expand Down Expand Up @@ -1403,7 +1414,7 @@ subroutine add_new_ice (ncat, nilyr, &
real (kind=dbl_kind), intent(in) :: &
wave_sig_ht ! significant height of waves globally (m)

real (kind=dbl_kind), dimension(:), intent(in) :: &
real (kind=dbl_kind), dimension(:), intent(in) :: &
wave_spectrum ! ocean surface wave spectrum, E(f) (m^2 s)

real(kind=dbl_kind), dimension(:), intent(in) :: &
Expand All @@ -1414,11 +1425,6 @@ subroutine add_new_ice (ncat, nilyr, &
floe_rad_c , & ! fsd size bin centre in m (radius)
floe_binwidth ! fsd size bin width in m (radius)

real (kind=dbl_kind), dimension(ncat) :: & ! for now
! change in thickness distribution (area)
d_an_latg , & ! due to fsd lateral growth
d_an_newi ! new ice formation

real (kind=dbl_kind), dimension(:), intent(out) :: &
! change in thickness distribution (area)
d_afsd_latg , & ! due to fsd lateral growth
Expand Down Expand Up @@ -1468,12 +1474,17 @@ subroutine add_new_ice (ncat, nilyr, &
vicen_init ! volume per unit area of ice (m)

! floe size distribution
real (kind=dbl_kind), dimension (nfsd,ncat) :: &
real (kind=dbl_kind), dimension (:,:), allocatable :: &
afsdn ! floe size distribution tracer (originally areal_mfstd_init)

! real (kind=dbl_kind), dimension (nfsd) :: &
! afsd , & ! fsd tracer for each thickness category

real (kind=dbl_kind), dimension(ncat) :: & ! for now
! change in thickness distribution (area)
d_an_latg , & ! due to fsd lateral growth
d_an_newi ! new ice formation

real (kind=dbl_kind), dimension (ncat) :: &
d_an_tot, & ! change in the ITD due to lateral growth and new ice
area2 ! area after lateral growth and before new ice formation
Expand All @@ -1497,7 +1508,6 @@ subroutine add_new_ice (ncat, nilyr, &
hsurp = c0
hi0new = c0
ai0new = c0
afsdn(:,:) = c0
d_an_latg(:) = c0
d_an_tot(:) = c0
d_an_newi(:) = c0
Expand All @@ -1520,6 +1530,8 @@ subroutine add_new_ice (ncat, nilyr, &
endif

if (tr_fsd) then
allocate(afsdn(nfsd,ncat))
afsdn(:,:) = c0
call icepack_cleanup_fsd (ncat, nfsd, trcrn(nt_fsd:nt_fsd+nfsd-1,:))
if (icepack_warnings_aborted(subname)) return
endif
Expand Down Expand Up @@ -1636,7 +1648,7 @@ subroutine add_new_ice (ncat, nilyr, &

if (vi0new > c0) then

if (tr_fsd) & ! lateral growth of existing ice
if (tr_fsd) then ! lateral growth of existing ice
! calculate change in conc due to lateral growth
! update vi0new, without change to afsdn or aicen
call fsd_lateral_growth (ncat, nfsd, &
Expand All @@ -1647,8 +1659,8 @@ subroutine add_new_ice (ncat, nilyr, &
lead_area, latsurf_area, &
G_radial, d_an_latg, &
tot_latg)

if (icepack_warnings_aborted(subname)) return
if (icepack_warnings_aborted(subname)) return
endif

ai0mod = aice0
! separate frazil ice growth from lateral ice growth
Expand Down Expand Up @@ -1822,7 +1834,7 @@ subroutine add_new_ice (ncat, nilyr, &
trcrn(nt_FY,n) = min(trcrn(nt_FY,n), c1)
endif

if (tr_fsd) & ! evolve the floe size distribution
if (tr_fsd) then ! evolve the floe size distribution
! both new frazil ice and lateral growth
call fsd_add_new_ice (ncat, n, nfsd, &
dt, ai0new, &
Expand All @@ -1837,8 +1849,8 @@ subroutine add_new_ice (ncat, nilyr, &
d_afsd_newi, &
afsdn, aicen_init, &
aicen, trcrn)

if (icepack_warnings_aborted(subname)) return
if (icepack_warnings_aborted(subname)) return
endif

if (vicen(n) > puny) then
if (tr_iage) &
Expand Down Expand Up @@ -1938,7 +1950,7 @@ subroutine add_new_ice (ncat, nilyr, &
!-----------------------------------------------------------------
! Biogeochemistry
!-----------------------------------------------------------------
if (tr_brine .or. nbtrcr > 0) &
if (tr_brine .or. nbtrcr > 0) then
call add_new_ice_bgc(dt, nblyr, &
ncat, nilyr, nltrcr, &
bgrid, cgrid, igrid, &
Expand All @@ -1948,6 +1960,11 @@ subroutine add_new_ice (ncat, nilyr, &
nbtrcr, sss, ocean_bio,&
flux_bio, hsurp)
if (icepack_warnings_aborted(subname)) return
endif

if (tr_fsd) then
deallocate(afsdn)
endif

end subroutine add_new_ice

Expand Down Expand Up @@ -1996,14 +2013,16 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &

integer (kind=int_kind), intent(in) :: &
ncat , & ! number of thickness categories
nfsd , & ! number of floe size categories
nltrcr , & ! number of zbgc tracers
nblyr , & ! number of bio layers
nilyr , & ! number of ice layers
nslyr ! number of snow layers

integer (kind=int_kind), intent(in), optional :: &
nfsd ! number of floe size categories

logical (kind=log_kind), intent(in) :: &
update_ocn_f ! if true, update fresh water and salt fluxes
update_ocn_f ! if true, update fresh water and salt fluxes

real (kind=dbl_kind), dimension(0:ncat), intent(inout) :: &
hin_max ! category boundaries (m)
Expand All @@ -2013,42 +2032,27 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &
Tf , & ! freezing temperature (C)
sss , & ! sea surface salinity (ppt)
rside , & ! fraction of ice that melts laterally
frzmlt , & ! freezing/melting potential (W/m^2)
wave_sig_ht ! significant height of waves in ice (m)

real (kind=dbl_kind), intent(in), optional :: &
wlat ! lateral melt rate (m/s)

real (kind=dbl_kind), dimension(:), intent(in) :: &
wave_spectrum ! ocean surface wave spectrum E(f) (m^2 s)

real(kind=dbl_kind), dimension(:), intent(in) :: &
wavefreq, & ! wave frequencies (s^-1)
dwavefreq ! wave frequency bin widths (s^-1)

real (kind=dbl_kind), dimension (:), intent(in) :: &
floe_rad_c , & ! fsd size bin centre in m (radius)
floe_binwidth ! fsd size bin width in m (radius)
frzmlt ! freezing/melting potential (W/m^2)

integer (kind=int_kind), dimension (:), intent(in) :: &
trcr_depend, & ! = 0 for aicen tracers, 1 for vicen, 2 for vsnon
n_trcr_strata ! number of underlying tracer layers

real (kind=dbl_kind), dimension (:,:), intent(in) :: &
trcr_base ! = 0 or 1 depending on tracer dependency
! argument 2: (1) aice, (2) vice, (3) vsno
trcr_base ! = 0 or 1 depending on tracer dependency
! argument 2: (1) aice, (2) vice, (3) vsno

integer (kind=int_kind), dimension (:,:), intent(in) :: &
nt_strata ! indices of underlying tracer layers
nt_strata ! indices of underlying tracer layers

real (kind=dbl_kind), dimension (nblyr+2), intent(in) :: &
bgrid ! biology nondimensional vertical grid points
bgrid ! biology nondimensional vertical grid points

real (kind=dbl_kind), dimension (nblyr+1), intent(in) :: &
igrid ! biology vertical interface points
igrid ! biology vertical interface points

real (kind=dbl_kind), dimension (nilyr+1), intent(in) :: &
cgrid ! CICE vertical coordinate
cgrid ! CICE vertical coordinate

real (kind=dbl_kind), dimension(:), intent(in) :: &
salinz , & ! initial salinity profile
Expand All @@ -2068,6 +2072,9 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &
frazil , & ! frazil ice growth (m/step-->cm/day)
frazil_diag ! frazil ice growth diagnostic (m/step-->cm/day)

real (kind=dbl_kind), intent(in), optional :: &
wlat ! lateral melt rate (m/s)

real (kind=dbl_kind), dimension(:), intent(inout) :: &
aicen_init,& ! initial concentration of ice
vicen_init,& ! initial volume per unit area of ice (m)
Expand All @@ -2081,14 +2088,7 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &
trcrn ! tracers

logical (kind=log_kind), dimension(:), intent(inout) :: &
first_ice ! true until ice forms

real (kind=dbl_kind), dimension(:), intent(out) :: &
! change in floe size distribution (area)
d_afsd_latg , & ! due to fsd lateral growth
d_afsd_newi , & ! new ice formation
d_afsd_latm , & ! lateral melt
d_afsd_weld ! welding
first_ice ! true until ice forms

real (kind=dbl_kind), intent(inout), optional :: &
frz_onset ! day of year that freezing begins (congel or frazil)
Expand All @@ -2098,12 +2098,34 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &

! water isotopes
real (kind=dbl_kind), dimension(:), intent(inout), optional :: &
fiso_ocn ! isotope flux to ocean (kg/m^2/s)
fiso_ocn ! isotope flux to ocean (kg/m^2/s)

real (kind=dbl_kind), intent(in), optional :: &
HDO_ocn , & ! ocean concentration of HDO (kg/kg)
H2_16O_ocn , & ! ocean concentration of H2_16O (kg/kg)
H2_18O_ocn ! ocean concentration of H2_18O (kg/kg)

real (kind=dbl_kind), intent(in), optional :: &
wave_sig_ht ! significant height of waves in ice (m)

real (kind=dbl_kind), dimension(:), intent(in), optional :: &
wave_spectrum ! ocean surface wave spectrum E(f) (m^2 s)

real(kind=dbl_kind), dimension(:), intent(in), optional :: &
wavefreq, & ! wave frequencies (s^-1)
dwavefreq ! wave frequency bin widths (s^-1)

real (kind=dbl_kind), dimension(:), intent(out), optional :: &
! change in floe size distribution (area)
d_afsd_latg, & ! due to fsd lateral growth
d_afsd_newi, & ! new ice formation
d_afsd_latm, & ! lateral melt
d_afsd_weld ! welding

real (kind=dbl_kind), dimension (:), intent(in), optional :: &
floe_rad_c, & ! fsd size bin centre in m (radius)
floe_binwidth ! fsd size bin width in m (radius)

!autodocument_end

! local variables
Expand All @@ -2129,7 +2151,18 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &
endif
endif
if (tr_fsd) then
if (.not.(present(wlat))) then
if (.not.(present(nfsd) .and. &
present(wlat) .and. &
present(wave_sig_ht) .and. &
present(wave_spectrum) .and. &
present(wavefreq) .and. &
present(dwavefreq) .and. &
present(d_afsd_latg) .and. &
present(d_afsd_newi) .and. &
present(d_afsd_latm) .and. &
present(d_afsd_weld) .and. &
present(floe_rad_c) .and. &
present(floe_binwidth))) then
call icepack_warnings_add(subname//' error in FSD arguments, tr_fsd=T')
call icepack_warnings_setabort(.true.,__FILE__,__LINE__)
return
Expand Down Expand Up @@ -2244,12 +2277,13 @@ subroutine icepack_step_therm2 (dt, ncat, nltrcr, &
if (icepack_warnings_aborted(subname)) return

! Floe welding during freezing conditions
if (tr_fsd) &
call fsd_weld_thermo (ncat, nfsd, &
dt, frzmlt, &
aicen, trcrn, &
d_afsd_weld)
if (icepack_warnings_aborted(subname)) return
if (tr_fsd) then
call fsd_weld_thermo (ncat, nfsd, &
dt, frzmlt, &
aicen, trcrn, &
d_afsd_weld)
if (icepack_warnings_aborted(subname)) return
endif

!-----------------------------------------------------------------
! For the special case of a single category, adjust the area and
Expand Down
Loading

0 comments on commit d024340

Please sign in to comment.