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.

Sign up for GitHub

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

Jump to bottom

Hamocc hybrid coord2 #179

Merged
merged 10 commits into from
Oct 3, 2022
Merged

Hamocc hybrid coord2 #179

Show file tree
Hide file tree
Changes from 9 commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
2961edb
Change mixed layer depth index "kmle" from parameter to 2D field.
TomasTorsvik Mar 31, 2022
7675858
Remove redundant comments in carchm.F90
TomasTorsvik Sep 28, 2022
cbf7682
Restructure cyano.F90 with outer k-loop
TomasTorsvik Sep 28, 2022
3665d94
Replace default "kmle" with "kmle_static" parameter
TomasTorsvik Sep 28, 2022
3742fc0
Set kmle from CVMix
TomasTorsvik Sep 28, 2022
f3444ba
Make hOBL varable avilable for iHAMOCC
TomasTorsvik Sep 28, 2022
eb9b960
Fix problem with "spval" exceeding integer range
TomasTorsvik Sep 28, 2022
0a1ca6e
Revert back to original format for mo_apply_rivin.F90
TomasTorsvik Sep 29, 2022
60e5eb7
Remove explicit dependence on vcoord_type_tag in iHAMOCC
TomasTorsvik Sep 29, 2022
cf29c31
Re-introduce indices on kmle in no_apply_rivin.F90
TomasTorsvik Sep 29, 2022
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
10 changes: 4 additions & 6 deletions hamocc/beleg_vars.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -220,18 +220,16 @@ SUBROUTINE BELEG_VARS(kpaufr,kpie,kpje,kpke,kbnd,pddpo,prho,omask, &

! Initialise preformed tracers in the mixed layer; note that the
! whole field has been initialised to zero above
DO k=1,kmle
DO j=1,kpje
DO i=1,kpie
IF(omask(i,j) .GT. 0.5) THEN
ocetra(i,j,k,iprefo2) =ocetra(i,j,k,ioxygen)
ocetra(i,j,k,iprefpo4)=ocetra(i,j,k,iphosph)
ocetra(i,j,k,iprefalk)=ocetra(i,j,k,ialkali)
ocetra(i,j,k,iprefdic)=ocetra(i,j,k,isco212)
ocetra(i,j,1:kmle(i,j),iprefo2) = ocetra(i,j,1:kmle(i,j),ioxygen)
ocetra(i,j,1:kmle(i,j),iprefpo4) = ocetra(i,j,1:kmle(i,j),iphosph)
ocetra(i,j,1:kmle(i,j),iprefalk) = ocetra(i,j,1:kmle(i,j),ialkali)
ocetra(i,j,1:kmle(i,j),iprefdic) = ocetra(i,j,1:kmle(i,j),isco212)
ENDIF
ENDDO
ENDDO
ENDDO


! Initial values for sediment
Expand Down
11 changes: 5 additions & 6 deletions hamocc/carchm.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -96,11 +96,11 @@ SUBROUTINE CARCHM(kpie,kpje,kpke,kbnd, &
!**********************************************************************
use mo_carbch, only: atm,atmflx,co2fxd,co2fxu,co2star,co3,hi,keqb,kwco2sol,ocetra,omegaa,omegac,pco2d,satn2o,satoxy
use mo_chemcon, only: al1,al2,al3,al4,an0,an1,an2,an3,an4,an5,an6,atn2o,bl1,bl2,bl3,calcon,ox0,ox1,ox2,ox3,ox4,ox5,ox6, &
& oxyco,tzero
use mo_control_bgc, only: dtbgc
& oxyco,tzero
use mo_control_bgc, only: dtbgc
use mo_param1_bgc, only: ialkali,iatmo2,iatmco2,iatmdms,iatmn2,iatmn2o,ian2o,icalc,idicsat,idms,igasnit,ioxygen,iphosph, &
& isco212,isilica
use mo_vgrid, only: dp_min,kbo,ptiestu
& isco212,isilica
use mo_vgrid, only: dp_min,kmle,kbo,ptiestu

#ifdef BROMO
use mo_param1_bgc, only: iatmbromo,ibromo
Expand Down Expand Up @@ -390,8 +390,7 @@ SUBROUTINE CARCHM(kpie,kpje,kpke,kbnd, &
ta = ocetra(i,j,k,ialkali) / rrho
CALL carchm_solve_DICsat(s,atco2*rpp0,ta,sit,pt,Kh,K1,K2,Kb,Kw,Ks1,Kf, &
Ksi,K1p,K2p,K3p,tc_sat,niter)
ocetra(i,j,k, idicsat)=tc_sat * rrho ! convert mol/kg to kmol/m^3
ocetra(i,j,k+1,idicsat)=tc_sat * rrho ! k+1 = the rest of the mixed layer
ocetra(i,j,1:kmle(i,j),idicsat) = tc_sat * rrho ! convert mol/kg to kmlo/m^3

#ifdef cisonew
! Ocean-Atmosphere fluxes for carbon isotopes
Expand Down
87 changes: 41 additions & 46 deletions hamocc/cyano.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@
! along with BLOM. If not, see https://www.gnu.org/licenses/.


SUBROUTINE CYANO(kpie,kpje,kpke,kbnd,pddpo,omask,ptho)
SUBROUTINE CYANO(kpie,kpje,kpke,kbnd,pddpo,omask,ptho)
!**********************************************************************
!
!**** *CYANO* - .
Expand Down Expand Up @@ -61,74 +61,69 @@ SUBROUTINE CYANO(kpie,kpje,kpke,kbnd,pddpo,omask,ptho)
! .
!**********************************************************************

use mo_carbch, only: ocetra
use mo_biomod, only: bluefix,intnfix,rnit,tf0,tf1,tf2,tff
use mo_param1_bgc, only: ialkali,iano3,igasnit,iphosph,ioxygen
use mo_vgrid, only: kmle
use mo_carbch, only: ocetra
use mo_biomod, only: bluefix,intnfix,rnit,tf0,tf1,tf2,tff
use mo_param1_bgc, only: ialkali,iano3,igasnit,iphosph,ioxygen
use mo_vgrid, only: kmle
#ifdef natDIC
use mo_param1_bgc, only: inatalkali
use mo_param1_bgc, only: inatalkali
#endif

implicit none
implicit none

INTEGER, intent(in) :: kpie,kpje,kpke,kbnd
REAL, intent(in) :: pddpo(kpie,kpje,kpke)
REAL, intent(in) :: omask(kpie,kpje)
REAL, intent(in) :: ptho(1-kbnd:kpie+kbnd,1-kbnd:kpje+kbnd,kpke)
INTEGER, intent(in) :: kpie,kpje,kpke,kbnd
REAL, intent(in) :: pddpo(kpie,kpje,kpke)
REAL, intent(in) :: omask(kpie,kpje)
REAL, intent(in) :: ptho(1-kbnd:kpie+kbnd,1-kbnd:kpje+kbnd,kpke)

! Local variables
INTEGER :: i,j,k
REAL :: oldocetra,dano3
REAL :: ttemp,nfixtfac
! Local variables
INTEGER :: i,j,k
REAL :: oldocetra,dano3
REAL :: ttemp,nfixtfac

intnfix(:,:)=0.0

intnfix(:,:)=0.0

!
! N-fixation by cyano bacteria (followed by remineralisation and nitrification),
! N-fixation by cyano bacteria (followed by remineralisation and nitrification),
! it is assumed here that this process is limited to the mixed layer
!
DO k=1,kmle
!$OMP PARALLEL DO PRIVATE(i,oldocetra,dano3,ttemp,nfixtfac)
DO j=1,kpje
DO i=1,kpie
IF(omask(i,j).gt.0.5) THEN
IF(ocetra(i,j,k,iano3).LT.(rnit*ocetra(i,j,k,iphosph))) THEN
DO j=1,kpje
DO i=1,kpie
IF(omask(i,j).gt.0.5) THEN
DO k=1,kmle(i,j)
IF(ocetra(i,j,k,iano3).LT.(rnit*ocetra(i,j,k,iphosph))) THEN

oldocetra = ocetra(i,j,k,iano3)
ttemp = min(40.,max(-3.,ptho(i,j,k)))
oldocetra = ocetra(i,j,k,iano3)
ttemp = min(40.,max(-3.,ptho(i,j,k)))

! Temperature dependence of nitrogen fixation, Kriest and Oschlies 2015.
nfixtfac = MAX(0.0,tf2*ttemp*ttemp + tf1*ttemp + tf0)/tff
nfixtfac = MAX(0.0,tf2*ttemp*ttemp + tf1*ttemp + tf0)/tff

ocetra(i,j,k,iano3)=ocetra(i,j,k,iano3)*(1-bluefix*nfixtfac) &
& +bluefix*nfixtfac*rnit*ocetra(i,j,k,iphosph)
ocetra(i,j,k,iano3)=ocetra(i,j,k,iano3)*(1-bluefix*nfixtfac) &
& + bluefix*nfixtfac*rnit*ocetra(i,j,k,iphosph)

dano3=ocetra(i,j,k,iano3)-oldocetra
dano3=ocetra(i,j,k,iano3)-oldocetra

ocetra(i,j,k,igasnit)=ocetra(i,j,k,igasnit)-dano3*(1./2.)
ocetra(i,j,k,igasnit)=ocetra(i,j,k,igasnit)-dano3*(1./2.)

! Note: to fix one mole N2 requires: N2+H2O+y*O2 = 2* HNO3 <-> y=2.5 mole O2.
! I.e., to release one mole HNO3 = H+ + NO3- requires 1.25 mole O2
ocetra(i,j,k,ioxygen)=ocetra(i,j,k,ioxygen)-dano3*1.25
ocetra(i,j,k,ioxygen)=ocetra(i,j,k,ioxygen)-dano3*1.25

! Nitrogen fixation followed by remineralisation and nitrification decreases
! alkalinity by 1 mole per mole nitrogen fixed (Wolf-Gladrow et al. 2007)
ocetra(i,j,k,ialkali)=ocetra(i,j,k,ialkali)-dano3
ocetra(i,j,k,ialkali)=ocetra(i,j,k,ialkali)-dano3
#ifdef natDIC
ocetra(i,j,k,inatalkali)=ocetra(i,j,k,inatalkali)-dano3
ocetra(i,j,k,inatalkali)=ocetra(i,j,k,inatalkali)-dano3
#endif

intnfix(i,j) = intnfix(i,j) + &
& (ocetra(i,j,k,iano3)-oldocetra)*pddpo(i,j,k)

ENDIF
ENDIF
ENDDO
ENDDO
!$OMP END PARALLEL DO
ENDDO

intnfix(i,j) = intnfix(i,j) + &
& (ocetra(i,j,k,iano3)-oldocetra)*pddpo(i,j,k)

ENDIF
ENDDO
ENDIF
ENDDO
ENDDO
TomasTorsvik marked this conversation as resolved.
Show resolved Hide resolved

RETURN
END
END SUBROUTINE CYANO
2 changes: 1 addition & 1 deletion hamocc/mo_apply_rivin.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -120,7 +120,7 @@ subroutine apply_rivin(kpie,kpje,kpke,pddpo,omask,rivin)

! Distribute riverine inputs over the model mixed layer
volij = 0.
DO k=1,kmle
DO k=1,kmle(i,j)
volij=volij+pddpo(i,j,k)
ENDDO

Expand Down
9 changes: 8 additions & 1 deletion hamocc/mo_intfcblom.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -244,14 +244,16 @@ subroutine blom2hamocc(m,n,mm,nn)
!******************************************************************************
!
use mod_constants, only: onem
use mod_xc, only: ii,jdm,jj,kdm,kk,ifp,isp,ilp,idm
use mod_xc, only: ii,jdm,jj,kdm,kk,ifp,isp,ilp,idm
use mod_grid, only: scpx,scpy
use mod_state, only: dp,temp,saln
use mod_eos, only: rho,p_alpha
use mod_difest, only: hOBL
use mod_tracers, only: ntrbgc,itrbgc,trc
use mo_param1_bgc, only: ks,nsedtra,npowtra,natm
use mo_carbch, only: ocetra,atm
use mo_sedmnt, only: sedlay,powtra,sedhpl,burial
use mo_vgrid, only: kmle, kmle_static

implicit none

Expand Down Expand Up @@ -292,6 +294,11 @@ subroutine blom2hamocc(m,n,mm,nn)
! --- - dimension of grid box in meters
bgc_dx(i,j) = scpx(i,j)/1.e2
bgc_dy(i,j) = scpy(i,j)/1.e2
!
! --- - index of level above OBL depth
! --- isopycninc coords: hOBL(i,j) = hOBL_static = 3. => kmle(i,j) = 2
! --- hybrid coords: hOBL defined according to cvmix_kpp_compute_kOBL_depth
kmle(i,j) = nint(hOBL(i,j))-1
enddo
enddo
!$OMP END PARALLEL DO
Expand Down
19 changes: 16 additions & 3 deletions hamocc/mo_vgrid.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -53,16 +53,18 @@ module mo_vgrid
!******************************************************************************
implicit none

INTEGER, PARAMETER :: kmle = 2 ! k-end index for layers that
! represent the mixed layer in BLOM
INTEGER, PARAMETER :: kmle_static = 2 ! k-end index for layers that
! represent the mixed layer in BLOM.
! Default value used for isopycnic coordinates.
REAL, PARAMETER :: dp_ez = 100.0 ! depth of euphotic zone
REAL, PARAMETER :: dp_min = 1.0E-12 ! min layer thickness layers thinner
REAL, PARAMETER :: dp_min = 1.0E-12 ! min layer thickness layers thinner
! than this are ignored by HAMOCC
REAL, PARAMETER :: dp_min_sink = 1.0 ! min layer thickness for sinking (layers thinner than
! this are ignored and set to the concentration of the
! layer above). Note that the bottom layer index kbo(i,j)
! is defined as the lowermost layer thicker than dp_min_sink.

INTEGER, DIMENSION(:,:), ALLOCATABLE :: kmle
INTEGER, DIMENSION(:,:), ALLOCATABLE :: kbo
INTEGER, DIMENSION(:,:), ALLOCATABLE :: kwrbioz
INTEGER, DIMENSION(:,:), ALLOCATABLE :: k0100,k0500,k1000,k2000,k4000
Expand Down Expand Up @@ -263,6 +265,17 @@ subroutine alloc_mem_vgrid(kpie,kpje,kpke)
ptiestw(:,:,:) = 0.0


IF(mnproc.eq.1) THEN
WRITE(io_stdo_bgc,*)'Memory allocation for variable kmle ...'
WRITE(io_stdo_bgc,*)'First dimension : ',kpie
WRITE(io_stdo_bgc,*)'Second dimension : ',kpje
ENDIF

ALLOCATE(kmle(kpie,kpje),stat=errstat)
if(errstat.ne.0) stop 'not enough memory kmle'
kmle(:,:) = kmle_static


IF(mnproc.eq.1) THEN
WRITE(io_stdo_bgc,*)'Memory allocation for variable kbo ...'
WRITE(io_stdo_bgc,*)'First dimension : ',kpie
Expand Down
4 changes: 2 additions & 2 deletions hamocc/ocprod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -938,7 +938,7 @@ subroutine ocprod(kpie,kpje,kpke,kbnd,pdlxp,pdlyp,pddpo,omask,ptho,pi_ph)
! Set minimum particle number to nmldmin in the mixed layer. This is to prevent
! very small values of nos (and asscociated high sinking speed if there is mass)
! in high latitudes during winter
if ( k <= kmle ) then
if ( k <= kmle(i,j) ) then
ocetra(i,j,k,inos) = MAX(nmldmin,ocetra(i,j,k,inos))
endif

Expand Down Expand Up @@ -974,7 +974,7 @@ subroutine ocprod(kpie,kpje,kpke,kbnd,pdlxp,pdlyp,pddpo,omask,ptho,pi_ph)

! As a first step, assume that shear in the mixed layer is high and
! zero below.
if ( k <= kmle ) then
if ( k <= kmle(i,j) ) then
fshear = fsh
else
fshear = 0.
Expand Down
38 changes: 17 additions & 21 deletions hamocc/preftrc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@
! along with BLOM. If not, see https://www.gnu.org/licenses/.


SUBROUTINE PREFTRC(kpie,kpje,omask)
SUBROUTINE PREFTRC(kpie,kpje,omask)
!****************************************************************
!
!**** *PREFTRC* - update preformed tracers in the mixed layer.
Expand All @@ -43,31 +43,27 @@ SUBROUTINE PREFTRC(kpie,kpje,omask)
!
!**************************************************************************

use mo_carbch, only: ocetra
use mo_param1_bgc, only: ialkali,ioxygen,iphosph,iprefalk,iprefdic,iprefo2,iprefpo4,isco212
use mo_vgrid, only: kmle
use mo_carbch, only: ocetra
use mo_param1_bgc, only: ialkali,ioxygen,iphosph,iprefalk,iprefdic,iprefo2,iprefpo4,isco212
use mo_vgrid, only: kmle

implicit none
implicit none

INTEGER :: kpie,kpje
REAL :: omask(kpie,kpje)
INTEGER :: kpie,kpje
REAL :: omask(kpie,kpje)

INTEGER :: i,j,k
INTEGER :: i,j

do k=1,kmle
!$OMP PARALLEL DO PRIVATE(i)
do j=1,kpje
do i=1,kpie
do j=1,kpje
do i=1,kpie
if (omask(i,j) .gt. 0.5 ) then
ocetra(i,j,k,iprefo2) =ocetra(i,j,k,ioxygen)
ocetra(i,j,k,iprefpo4)=ocetra(i,j,k,iphosph)
ocetra(i,j,k,iprefalk)=ocetra(i,j,k,ialkali)
ocetra(i,j,k,iprefdic)=ocetra(i,j,k,isco212)
ocetra(i,j,1:kmle(i,j),iprefo2) = ocetra(i,j,1:kmle(i,j),ioxygen)
ocetra(i,j,1:kmle(i,j),iprefpo4) = ocetra(i,j,1:kmle(i,j),iphosph)
ocetra(i,j,1:kmle(i,j),iprefalk) = ocetra(i,j,1:kmle(i,j),ialkali)
ocetra(i,j,1:kmle(i,j),iprefdic) = ocetra(i,j,1:kmle(i,j),isco212)
endif
enddo
enddo
!$OMP END PARALLEL DO
enddo
enddo
enddo


END SUBROUTINE PREFTRC
END SUBROUTINE PREFTRC
10 changes: 8 additions & 2 deletions phy/mod_difest.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -78,13 +78,19 @@ module mod_difest
implicit none
c
private
c
c Initialize hOBL with hOBL_static = 3. for consistency with bulk
c mixed layer formulation in iHAMOCC: kmle = nint(hOBL) - 1 = 2
real, PARAMETER :: hOBL_static = 3.
c
real, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy,kdm+1) ::
. rig
real, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy,kdm) ::
. du2l,drhol,up,vp
real, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy) ::
. OBLdepth
real, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy) ::
. hOBL
integer, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy,kdm) ::
. mskv,msku
integer, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy) ::
Expand Down Expand Up @@ -174,7 +180,7 @@ module mod_difest
. cs=98.96,minOBLdepth=1.0)
c
public :: OBLdepth, inivar_difest, init_difest, difest_isobml,
. difest_lateral_hybrid, difest_vertical_hybrid
. difest_lateral_hybrid, difest_vertical_hybrid, hOBL
c
contains
c
Expand Down Expand Up @@ -203,6 +209,7 @@ subroutine inivar_difest
enddo
do i=1-nbdy,ii+nbdy
OBLdepth(i,j)=spval
hOBL(i,j) = hOBL_static
enddo
enddo
c$OMP END PARALLEL DO
Expand Down Expand Up @@ -904,7 +911,6 @@ subroutine difest_vertical_hyb(m,n,mm,nn,k1m,k1n)
real :: Simmons_coeff, zBottomMinusOffset
real :: bl1, bl2, bl3, bl4
integer ki, ksfc, ktmp, kOBL, kn1
real, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy) :: hOBL
c
surf_layer_ext = 0.1
bl1 = 8e-5
Expand Down