Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Clean up Icepack interfaces #38

Merged
merged 3 commits into from
Jun 24, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
19 changes: 3 additions & 16 deletions columnphysics/icepack_aerosol.F90
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ module icepack_aerosol
use icepack_kinds
use icepack_parameters, only: c0, c1, c2, p5, puny, rhoi, rhos, hs_min
use icepack_parameters, only: hi_ssl, hs_ssl, hs_ssl_min
use icepack_tracers, only: max_aero
use icepack_tracers, only: max_aero, nilyr, nslyr, nblyr, ntrcr, nbtrcr, n_aero
use icepack_warnings, only: warnstr, icepack_warnings_add
use icepack_warnings, only: icepack_warnings_setabort, icepack_warnings_aborted

Expand All @@ -31,8 +31,6 @@ module icepack_aerosol
! Called from icepack_step_therm1 when tr_aero=T (not used for zbgc tracers)

subroutine update_aerosol(dt, &
nilyr, nslyr, &
n_aero, &
meltt, melts, &
meltb, congel, &
snoice, &
Expand All @@ -43,9 +41,6 @@ subroutine update_aerosol(dt, &
vicen, vsnon, aicen, &
faero_atm, faero_ocn)

integer (kind=int_kind), intent(in) :: &
nilyr, nslyr, n_aero

real (kind=dbl_kind), intent(in) :: &
dt, & ! time step
meltt, & ! thermodynamic melt/growth rates
Expand Down Expand Up @@ -428,12 +423,10 @@ end subroutine update_aerosol
! Aerosol in snow for vertical biogeochemistry with mushy thermodynamics
! Called from icepack_algae.F90 when z_tracers=T (replaces update_aerosol)

subroutine update_snow_bgc (dt, nblyr, &
nslyr, &
subroutine update_snow_bgc(dt, &
meltt, melts, &
meltb, congel, &
snoice, nbtrcr, &
fsnow, ntrcr, &
snoice, fsnow, &
trcrn, bio_index, &
aice_old, zbgc_snow, &
vice_old, vsno_old, &
Expand All @@ -442,12 +435,6 @@ subroutine update_snow_bgc (dt, nblyr, &
zbgc_atm, flux_bio, &
bio_index_o)

integer (kind=int_kind), intent(in) :: &
nbtrcr, & ! number of distinct snow tracers
nblyr, & ! number of bio layers
nslyr, & ! number of snow layers
ntrcr ! number of tracers

integer (kind=int_kind), dimension (nbtrcr), intent(in) :: &
bio_index, &
bio_index_o ! provides index of scavenging (kscavz) data array
Expand Down
216 changes: 55 additions & 161 deletions columnphysics/icepack_algae.F90

Large diffs are not rendered by default.

7 changes: 2 additions & 5 deletions columnphysics/icepack_atmo.F90
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ module icepack_atmo
use icepack_parameters, only: pih, dragio, rhoi, rhos, rhow
use icepack_parameters, only: atmbndy, calc_strair, formdrag
use icepack_parameters, only: icepack_chkoptargflag
use icepack_tracers, only: n_iso
use icepack_tracers, only: ncat, n_iso
use icepack_tracers, only: tr_iso
use icepack_warnings, only: warnstr, icepack_warnings_add
use icepack_warnings, only: icepack_warnings_setabort, icepack_warnings_aborted
Expand Down Expand Up @@ -521,13 +521,10 @@ subroutine neutral_drag_coeffs (apnd, hpnd, &
hdraft, hridge, &
distrdg, hkeel, &
dkeel, lfloe, &
dfloe, ncat)
dfloe)

use icepack_tracers, only: tr_pond

integer (kind=int_kind), intent(in) :: &
ncat

real (kind=dbl_kind), dimension (:), intent(in) :: &
apnd ,& ! melt pond fraction of sea ice
hpnd ,& ! mean melt pond depth over sea ice
Expand Down
151 changes: 62 additions & 89 deletions columnphysics/icepack_brine.F90
Original file line number Diff line number Diff line change
Expand Up @@ -11,15 +11,15 @@ module icepack_brine
use icepack_parameters, only: gravit, rhoi, rhow, rhos, depressT
use icepack_parameters, only: salt_loss, min_salin, rhosi
use icepack_parameters, only: dts_b, l_sk
use icepack_tracers, only: ntrcr, nt_qice, nt_sice
use icepack_tracers, only: nilyr, nblyr, ntrcr, nt_qice, nt_sice
use icepack_tracers, only: nt_Tsfc
use icepack_zbgc_shared, only: k_o, exp_h, Dm, Ra_c, viscos_dynamic, thinS
use icepack_zbgc_shared, only: bgrid, cgrid, igrid, swgrid, icgrid
use icepack_zbgc_shared, only: remap_zbgc
use icepack_warnings, only: warnstr, icepack_warnings_add
use icepack_warnings, only: icepack_warnings_setabort, icepack_warnings_aborted

use icepack_mushy_physics, only: icepack_mushy_temperature_mush, icepack_mushy_liquid_fraction
use icepack_therm_shared, only: calculate_Tin_from_qin

implicit none

Expand Down Expand Up @@ -124,29 +124,14 @@ end subroutine preflushing_changes
! NOTE: This subroutine uses thermosaline_vertical output to compute
! average ice permeability and the surface ice porosity

subroutine compute_microS_mushy (nilyr, nblyr, &
bgrid, cgrid, igrid, &
trcrn, hice_old, hbr_old, &
subroutine compute_microS_mushy (trcrn, hice_old, hbr_old, &
sss, sst, bTin, &
iTin, bphin, &
kperm, bphi_min, &
kperm, bphi_min, &
bSin, brine_sal, brine_rho, &
iphin, ibrine_rho, ibrine_sal, &
iDin, iSin)

integer (kind=int_kind), intent(in) :: &
nilyr , & ! number of ice layers
nblyr ! number of bio layers

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

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

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

real (kind=dbl_kind), intent(in) :: &
hice_old , & ! previous timestep ice height (m)
sss , & ! ocean salinity (ppt)
Expand Down Expand Up @@ -269,17 +254,15 @@ subroutine compute_microS_mushy (nilyr, nblyr, &
! Define ice multiphase structure
!-----------------------------------------------------------------

call prepare_hbrine (nblyr, &
bSin, bTin, iTin, &
call prepare_hbrine (bSin, bTin, iTin, &
brine_sal, brine_rho, &
ibrine_sal, ibrine_rho, &
bphin, iphin, &
kperm, bphi_min, &
igrid, sss, iSin)
sss, iSin)
if (icepack_warnings_aborted(subname)) return

call calculate_drho(nblyr, igrid, bgrid, &
brine_rho, ibrine_rho, drho)
call calculate_drho(brine_rho, ibrine_rho, drho)
if (icepack_warnings_aborted(subname)) return

do k= 2, nblyr+1
Expand All @@ -294,21 +277,16 @@ end subroutine compute_microS_mushy

!=======================================================================

subroutine prepare_hbrine (nblyr, &
bSin, bTin, iTin, &
subroutine prepare_hbrine (bSin, bTin, iTin, &
brine_sal, brine_rho, &
ibrine_sal, ibrine_rho, &
bphin, iphin, &
kperm, bphi_min, &
i_grid, sss, iSin)

integer (kind=int_kind), intent(in) :: &
nblyr ! number of bio layers
sss, iSin)

real (kind=dbl_kind), dimension (:), intent(in) :: &
bSin , & ! salinity of ice layers on bio grid (ppt)
bTin , & ! temperature of ice layers on bio grid for history (C)
i_grid ! biology grid interface points
bTin ! temperature of ice layers on bio grid for history (C)

real (kind=dbl_kind), dimension (:), intent(inout) :: &
brine_sal , & ! equilibrium brine salinity (ppt)
Expand Down Expand Up @@ -350,7 +328,7 @@ subroutine prepare_hbrine (nblyr, &
if (k == 1) then
igrm = 0
else
igrm = i_grid(k) - i_grid(k-1)
igrm = igrid(k) - igrid(k-1)
endif

brine_sal(k) = a1*bTin(k) &
Expand Down Expand Up @@ -389,9 +367,9 @@ subroutine prepare_hbrine (nblyr, &
kperm = k_min
endif

igrp = i_grid(k+1) - i_grid(k )
igrm = i_grid(k ) - i_grid(k-1)
rigr = c1 / (i_grid(k+1)-i_grid(k-1))
igrp = igrid(k+1) - igrid(k )
igrm = igrid(k ) - igrid(k-1)
rigr = c1 / (igrid(k+1)-igrid(k-1))

ibrine_sal(k) = (brine_sal(k+1)*igrp + brine_sal(k)*igrm) * rigr
ibrine_rho(k) = (brine_rho(k+1)*igrp + brine_rho(k)*igrm) * rigr
Expand Down Expand Up @@ -422,12 +400,12 @@ end subroutine prepare_hbrine
! ice. This volume fraction may be > 1 in which case there is brine
! above the ice surface (ponds).

subroutine update_hbrine (meltt, &
subroutine update_hbrine (meltt, &
melts, dt, &
hin, hsn, &
hin_old, hbr, &
hbr_old, &
fbri, &
hbr_old, &
fbri, &
dhS_top, dhS_bottom, &
dh_top_chl, dh_bot_chl, &
kperm, bphi_min, &
Expand Down Expand Up @@ -557,17 +535,7 @@ end subroutine update_hbrine
! Find density difference about interface grid points
! for gravity drainage parameterization

subroutine calculate_drho (nblyr, i_grid, b_grid, &
brine_rho, ibrine_rho, drho)

integer (kind=int_kind), intent(in) :: &
nblyr ! number of bio layers

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

real (kind=dbl_kind), dimension (nblyr+1), intent(in) :: &
i_grid ! biology grid interface points
subroutine calculate_drho (brine_rho, ibrine_rho, drho)

real (kind=dbl_kind), dimension (nblyr+2), intent(in) :: &
brine_rho ! Internal brine density (kg/m^3)
Expand Down Expand Up @@ -602,29 +570,29 @@ subroutine calculate_drho (nblyr, i_grid, b_grid, &
rho_2b(:) = c0
drho (:) = c0 ! surface is snow or atmosphere

do k = 1, nblyr+1 ! i_grid values
do k = 1, nblyr+1 ! igrid values

!----------------------------------------------
! h_avg(k) = i_grid(k)
! Calculate rho_a(k), ie average rho above i_grid(k)
! h_avg(k) = igrid(k)
! Calculate rho_a(k), ie average rho above igrid(k)
! first part is good
!----------------------------------------------

if (k == 2) then
rho_a(2) = (brine_rho(2)*b_grid(2) &
rho_a(2) = (brine_rho(2)*bgrid(2) &
+ (ibrine_rho(2) + brine_rho(2)) &
* p5*(i_grid(2)-b_grid(2)) )/i_grid(2)
* p5*(igrid(2)-bgrid(2)) )/igrid(2)
rho_b(2) = brine_rho(2)

elseif (k > 2 .AND. k < nblyr+1) then
rho_a(k) = (rho_a(k-1)*i_grid(k-1) + (ibrine_rho(k-1) + brine_rho(k)) &
* p5*(b_grid(k)-i_grid(k-1)) + (ibrine_rho(k ) + brine_rho(k)) &
* p5*(i_grid(k)-b_grid(k)))/i_grid(k)
rho_a(k) = (rho_a(k-1)*igrid(k-1) + (ibrine_rho(k-1) + brine_rho(k)) &
* p5*(bgrid(k)-igrid(k-1)) + (ibrine_rho(k ) + brine_rho(k)) &
* p5*(igrid(k)-bgrid(k)))/igrid(k)
rho_b(k) = brine_rho(k)
else
rho_a(nblyr+1) = (rho_a(nblyr)*i_grid(nblyr) + (ibrine_rho(nblyr) + &
brine_rho(nblyr+1))*p5*(b_grid(nblyr+1)-i_grid(nblyr)) + &
brine_rho(nblyr+1)*(i_grid(nblyr+1)-b_grid(nblyr+1)))/i_grid(nblyr+1)
rho_a(nblyr+1) = (rho_a(nblyr)*igrid(nblyr) + (ibrine_rho(nblyr) + &
brine_rho(nblyr+1))*p5*(bgrid(nblyr+1)-igrid(nblyr)) + &
brine_rho(nblyr+1)*(igrid(nblyr+1)-bgrid(nblyr+1)))/igrid(nblyr+1)
rho_a(1) = brine_rho(2) !for k == 1 use grid point value
rho_b(nblyr+1) = brine_rho(nblyr+1)
rho_b(1) = brine_rho(2)
Expand All @@ -636,10 +604,10 @@ subroutine calculate_drho (nblyr, i_grid, b_grid, &
! Calculate average above and below k rho_2a
!----------------------------------------------

do k = 1, nblyr+1 !i_grid values
do k = 1, nblyr+1 !igrid values
if (k == 1) then
rho_2a(1) = (rho_a(1)*b_grid(2) + p5*(brine_rho(2) + ibrine_rho(2)) &
* (i_grid(2)-b_grid(2)))/i_grid(2)
rho_2a(1) = (rho_a(1)*bgrid(2) + p5*(brine_rho(2) + ibrine_rho(2)) &
* (igrid(2)-bgrid(2)))/igrid(2)
rho_2b(1) = brine_rho(2)
else
mstop = 2*(k-1) + 1
Expand All @@ -653,7 +621,7 @@ subroutine calculate_drho (nblyr, i_grid, b_grid, &
endif

do mm = mstart,mstop
rho_2a(k) =(rho_a(nblyr+1) + rhow*(c2*i_grid(k)-c1))*p5/i_grid(k)
rho_2a(k) =(rho_a(nblyr+1) + rhow*(c2*igrid(k)-c1))*p5/igrid(k)
enddo
rho_2b(k) = brine_rho(k+1)
endif
Expand All @@ -667,26 +635,22 @@ end subroutine calculate_drho
!autodocument_start icepack_init_hbrine
! Initialize brine height tracer

subroutine icepack_init_hbrine(bgrid, igrid, cgrid, &
icgrid, swgrid, nblyr, nilyr, phi_snow)

integer (kind=int_kind), intent(in) :: &
nilyr, & ! number of ice layers
nblyr ! number of bio layers
subroutine icepack_init_hbrine(bgrid_out, igrid_out, cgrid_out, &
icgrid_out, swgrid_out, phi_snow)

real (kind=dbl_kind), intent(inout) :: &
real (kind=dbl_kind), optional, intent(inout) :: &
phi_snow ! porosity at the ice-snow interface

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

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

real (kind=dbl_kind), dimension (nilyr+1), intent(out) :: &
cgrid , & ! CICE vertical coordinate
icgrid , & ! interface grid for CICE (shortwave variable)
swgrid ! grid for ice tracers used in dEdd scheme
real (kind=dbl_kind), optional, dimension (:), intent(out) :: &
cgrid_out , & ! CICE vertical coordinate
icgrid_out , & ! interface grid for CICE (shortwave variable)
swgrid_out ! grid for ice tracers used in dEdd scheme

!autodocument_end

Expand All @@ -700,8 +664,17 @@ subroutine icepack_init_hbrine(bgrid, igrid, cgrid, &

character(len=*),parameter :: subname='(icepack_init_hbrine)'

!-----------------------------------------------------------------

if (phi_snow .le. c0) phi_snow = c1-rhos/rhoi
if (present(phi_snow)) then
if (phi_snow .le. c0) phi_snow = c1-rhos/rhoi
endif

allocate(bgrid (nblyr+2))
allocate(igrid (nblyr+1))
allocate(cgrid (nilyr+1))
allocate(icgrid(nilyr+1))
allocate(swgrid(nilyr+1))

!-----------------------------------------------------------------
! Calculate bio gridn: 0 to 1 corresponds to ice top to bottom
Expand Down Expand Up @@ -757,21 +730,21 @@ subroutine icepack_init_hbrine(bgrid, igrid, cgrid, &
swgrid(k) = zspace * (real(k,kind=dbl_kind)-c1p5)
enddo

if (present( bgrid_out)) bgrid_out=bgrid
if (present( cgrid_out)) cgrid_out=cgrid
if (present( igrid_out)) igrid_out=igrid
if (present(icgrid_out)) icgrid_out=icgrid
if (present(swgrid_out)) swgrid_out=swgrid

end subroutine icepack_init_hbrine

!=======================================================================
!autodocument_start icepack_init_zsalinity
! **DEPRECATED**, all code removed
! Interface provided for backwards compatibility

subroutine icepack_init_zsalinity(nblyr,ntrcr_o, Rayleigh_criteria, &
Rayleigh_real, trcrn_bgc, nt_bgc_S, ncat, sss)

integer (kind=int_kind), intent(in) :: &
nblyr , & ! number of biolayers
ntrcr_o, & ! number of non bio tracers
ncat , & ! number of categories
nt_bgc_S ! zsalinity index
subroutine icepack_init_zsalinity(Rayleigh_criteria, &
Rayleigh_real, trcrn_bgc, sss)

logical (kind=log_kind), intent(inout) :: &
Rayleigh_criteria
Expand Down
Loading
Loading