Skip to content

Commit

Permalink
Merge pull request #1 from apcraig/iso1_interface
Browse files Browse the repository at this point in the history 
Clean up Icepack interface for isotopes
  • Loading branch information
@eclare108213
eclare108213 authored Mar 10, 2020
2 parents d3f5108 + d704ff8 commit a87cc48
Show file tree
Hide file tree
Showing 15 changed files with 575 additions and 289 deletions.
102 changes: 57 additions & 45 deletions columnphysics/icepack_atmo.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,8 @@ module icepack_atmo
use icepack_parameters, only: pih, dragio, rhoi, rhos, rhow
use icepack_parameters, only: atmbndy, calc_strair, formdrag
use icepack_parameters, only: highfreq, natmiter
use icepack_tracers, only: 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 @@ -61,13 +63,11 @@ subroutine atmo_boundary_layer (sfctype, &
lhcoef, shcoef, &
Cdn_atm, &
Cdn_atm_ratio_n, &
n_iso, &
Qa_iso, Qref_iso, &
iso_flag, &
uvel, vvel, &
Uref )

use icepack_tracers, only: tr_iso

character (len=3), intent(in) :: &
sfctype ! ice or ocean

Expand All @@ -79,9 +79,6 @@ subroutine atmo_boundary_layer (sfctype, &
integer (kind=int_kind), intent(in) :: &
natmiter ! number of iterations for boundary layer calculations

integer (kind=int_kind), optional, intent(in) :: &
n_iso ! number of isotopes

real (kind=dbl_kind), intent(in) :: &
Tsf , & ! surface temperature of ice or ocean
potT , & ! air potential temperature (K)
Expand Down Expand Up @@ -111,11 +108,14 @@ subroutine atmo_boundary_layer (sfctype, &
shcoef , & ! transfer coefficient for sensible heat
lhcoef ! transfer coefficient for latent heat

logical (kind=log_kind), intent(in), optional :: &
iso_flag ! flag to trigger iso calculations

real (kind=dbl_kind), intent(in), optional, dimension(:) :: &
Qa_iso ! specific isotopic humidity (kg/kg)
Qa_iso ! specific isotopic humidity (kg/kg)

real (kind=dbl_kind), intent(inout), optional, dimension(:) :: &
Qref_iso ! reference specific isotopic humidity (kg/kg)
Qref_iso ! reference specific isotopic humidity (kg/kg)

real (kind=dbl_kind), intent(in) :: &
uvel , & ! x-direction ice speed (m/s)
Expand Down Expand Up @@ -166,10 +166,18 @@ subroutine atmo_boundary_layer (sfctype, &
psixh ! stability function at zlvl (heat and water)

real (kind=dbl_kind), parameter :: &
zTrf = c2 ! reference height for air temp (m)
zTrf = c2 ! reference height for air temp (m)

logical (kind=log_kind) :: &
l_iso_flag ! local flag to trigger iso calculations

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

l_iso_flag = .false.
if (present(iso_flag)) then
l_iso_flag = iso_flag
endif

al2 = log(zref/zTrf)

!------------------------------------------------------------
Expand Down Expand Up @@ -370,8 +378,8 @@ subroutine atmo_boundary_layer (sfctype, &
Uref = vmag * rd / rdn
endif

if (present(Qref_iso) .and. present(Qa_iso) &
.and. present(n_iso)) then
if (l_iso_flag) then
if (present(Qref_iso) .and. present(Qa_iso)) then
Qref_iso(:) = c0
if (tr_iso) then
do n = 1, n_iso
Expand All @@ -380,6 +388,11 @@ subroutine atmo_boundary_layer (sfctype, &
Qref_iso(n) = Qa_iso(n) - ratio*delq*fac
enddo
endif
else
call icepack_warnings_add(subname//' l_iso_flag true but optional arrays missing')
call icepack_warnings_setabort(.true.,__FILE__,__LINE__)
return
endif
endif

end subroutine atmo_boundary_layer
Expand Down Expand Up @@ -494,7 +507,7 @@ end subroutine atmo_boundary_const
!=======================================================================

! Neutral drag coefficients for ocean and atmosphere also compute the
! intermediate necessary variables ridge height, distance, floe size
! intermediate necessary variables ridge height, distance, floe size
! based upon Tsamados et al. (2014), JPO, DOI: 10.1175/JPO-D-13-0215.1.
! Places where the code varies from the paper are commented.
!
Expand Down Expand Up @@ -600,7 +613,7 @@ subroutine neutral_drag_coeffs (apnd, hpnd, &
ctecaf, & ! constante
ctecwf, & ! constante
sca, & ! wind attenuation function
scw, & ! ocean attenuation function
scw, & ! ocean attenuation function
lp, & ! pond length (m)
ctecar, &
ctecwk, &
Expand Down Expand Up @@ -724,7 +737,7 @@ subroutine neutral_drag_coeffs (apnd, hpnd, &

!------------------------------------------------------------
! Skin drag (atmo)
!------------------------------------------------------------
!------------------------------------------------------------

Cdn_atm_skin = csa*(c1 - mrdg*tmp1/distrdg)
Cdn_atm_skin = max(min(Cdn_atm_skin,camax),c0)
Expand All @@ -751,7 +764,7 @@ subroutine neutral_drag_coeffs (apnd, hpnd, &

!------------------------------------------------------------
! Skin drag bottom ice (ocean)
!------------------------------------------------------------
!------------------------------------------------------------

Cdn_ocn_skin = csw * (c1 - mrdgo*tmp1/dkeel)
Cdn_ocn_skin = max(min(Cdn_ocn_skin,cwmax), c0)
Expand Down Expand Up @@ -837,13 +850,10 @@ subroutine icepack_atm_boundary(sfctype, &
lhcoef, shcoef, &
Cdn_atm, &
Cdn_atm_ratio_n, &
n_iso, &
Qa_iso, Qref_iso, &
uvel, vvel, &
Uref)

use icepack_tracers, only: tr_iso

character (len=3), intent(in) :: &
sfctype ! ice or ocean

Expand All @@ -857,9 +867,6 @@ subroutine icepack_atm_boundary(sfctype, &
Qa , & ! specific humidity (kg/kg)
rhoa ! air density (kg/m^3)

integer (kind=int_kind), optional, intent(in) :: &
n_iso ! number of isotopes

real (kind=dbl_kind), intent(inout) :: &
Cdn_atm , & ! neutral drag coefficient
Cdn_atm_ratio_n ! ratio drag coeff / neutral drag coeff
Expand Down Expand Up @@ -895,33 +902,37 @@ subroutine icepack_atm_boundary(sfctype, &
! local variables

real (kind=dbl_kind) :: &
worku, workv, workr

integer (kind=int_kind) :: &
niso
l_uvel, l_vvel, l_Uref

real (kind=dbl_kind), dimension(:), allocatable :: &
Qaiso, Qriso
l_Qa_iso, l_Qref_iso ! local copies of Qa_iso, Qref_iso

logical (kind=log_kind) :: &
iso_flag ! flag to turn on iso calcs in other subroutines

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

worku = c0
workv = c0
workr = c0
l_uvel = c0
l_vvel = c0
l_Uref = c0
if (present(uvel)) then
worku = uvel
l_uvel = uvel
endif
if (present(vvel)) then
workv = vvel
l_vvel = vvel
endif
if (present(n_iso) .and. present(Qa_iso) .and. &
present(Qref_iso)) then
allocate(Qaiso(n_iso),Qriso(n_iso))
Qaiso = Qa_iso
Qriso = Qref_iso
niso = n_iso
if (present(Qa_iso) .and. present(Qref_iso)) then
iso_flag = .true.
allocate(l_Qa_iso(size(Qa_iso,dim=1)))
allocate(l_Qref_iso(size(Qref_iso,dim=1)))
l_Qa_iso = Qa_iso
l_Qref_iso = Qref_iso
else
allocate(Qaiso(1),Qriso(1))
iso_flag = .false.
allocate(l_Qa_iso(1))
allocate(l_Qref_iso(1))
l_Qa_iso = c0
l_Qref_iso = c0
endif

Cdn_atm_ratio_n = c1
Expand Down Expand Up @@ -950,22 +961,23 @@ subroutine icepack_atm_boundary(sfctype, &
lhcoef, shcoef, &
Cdn_atm, &
Cdn_atm_ratio_n, &
n_iso=niso, &
Qa_iso=Qaiso, Qref_iso=Qriso,&
uvel=worku, vvel=workv, &
Uref=workr)
iso_flag = iso_flag, &
Qa_iso=l_Qa_iso, &
Qref_iso=l_Qref_iso, &
uvel=l_uvel, vvel=l_vvel, &
Uref=l_Uref)
if (icepack_warnings_aborted(subname)) return
endif ! atmbndy

if (present(Uref)) then
Uref = workr
Uref = l_Uref
endif

if (present(Qref_iso)) then
Qref_iso = Qriso
Qref_iso = l_Qref_iso
endif

deallocate(Qaiso,Qriso)
deallocate(l_Qa_iso,l_Qref_iso)

end subroutine icepack_atm_boundary

Expand Down
19 changes: 11 additions & 8 deletions columnphysics/icepack_flux.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@ module icepack_flux
use icepack_parameters, only: c1, emissivity
use icepack_warnings, only: warnstr, icepack_warnings_add
use icepack_warnings, only: icepack_warnings_setabort, icepack_warnings_aborted
use icepack_tracers, only: tr_iso

implicit none
private
Expand Down Expand Up @@ -162,14 +163,16 @@ subroutine merge_fluxes (aicen, &
Qref = Qref + Qrefn * aicen

! Isotopes
if (present(Qrefn_iso) .and. present(Qref_iso)) then
Qref_iso (:) = Qref_iso (:) + Qrefn_iso (:) * aicen
endif
if (present(fiso_ocnn) .and. present(fiso_ocn)) then
fiso_ocn (:) = fiso_ocn (:) + fiso_ocnn (:) * aicen
endif
if (present(fiso_evapn) .and. present(fiso_evap)) then
fiso_evap(:) = fiso_evap(:) + fiso_evapn(:) * aicen
if (tr_iso) then
if (present(Qrefn_iso) .and. present(Qref_iso)) then
Qref_iso (:) = Qref_iso (:) + Qrefn_iso (:) * aicen
endif
if (present(fiso_ocnn) .and. present(fiso_ocn)) then
fiso_ocn (:) = fiso_ocn (:) + fiso_ocnn (:) * aicen
endif
if (present(fiso_evapn) .and. present(fiso_evap)) then
fiso_evap(:) = fiso_evap(:) + fiso_evapn(:) * aicen
endif
endif

! ocean fluxes
Expand Down
Loading

0 comments on commit a87cc48

Please sign in to comment.