Skip to content

Commit

Permalink
Merge latest master into feature-hamocc_beyond-CMIP6 branch (NorESMhu…
Browse files Browse the repository at this point in the history 
…b#232)

* Dynamic mapping of pore water tracers to ocean tracers (NorESMhub#192)

* Initial restructuring of sediment-related tracer declaration and initialization

* Introducing mapping function

* Remove unncessary comments

* Fixed diagnostics bug and updated index naming

* Added initial support for NUOPC driver.

* Lon-lat variable sediment porosity (NorESMhub#189)

Introducing a static 3D sediment porosity field that can be optionally read in with effects on molecular pore water diffusion and shifting.

* Added wave forcing fields.

* Renamed folder for MCT driver.

* Moved MCT specific file from drivers/cpl_share/ to drivers/mct/.

* Rename drivers/mct/mod_swtfrz.F to drivers/mct/mod_swtfrz.F90.

* Rewrite to drivers/mct/mod_swtfrz.F90 to free format Fortran.

* Remove redundant definition of kOBL.

* Redefine kOBL, cast as integer

* Fixing variable sediment porosity - field initialization in case of `sedbypass=true` (NorESMhub#198)

* Removing bodensed -  Initialization of sediment parameters and fields now in mo_sedmnt

* This is the first commit of MKS units. All variables in the subroutines have been converted to MKS [meter, kg, seconds] instead of CGS [cm, gram, seconds] with necessary coefficients. The default option which is CGS reproduce old results. The new option MKS cannot reproduce because of machine precision.

* Hamocc hybrid coord2 (NorESMhub#179)

Make the surface mixed layer depth fractional index `hOBL` available for use in iHAMOCC, and adjust the internal iHAMOCC index `kmle` according to `hOBL`. Default value `kmle = 2` is retained for consistency with isopycnic coordinates.

* BLOM CIME cpp updates to run in NorESM

* bug fixes for the CGS MKS conversion

* cesm thermal forcing bug fixes for reproducibility

* BLOM MKS update to export winds into the CESM using proper units.

* input values in ocn_in case is updated for mks setup

* default cgsmks value changed

* Initialize some variables in the k-epsilon model.

* Fix porosity read (NorESMhub#201)

* Fixing the reading of variable porosity input field in preparation for the NorESM 2.0.6 release

Cherry-picked from private Ncycleprivate branch 0d56930e2fdd62caba964d375b57304942568926

* Provide number of layers (3rd dim) via ks and not hard-coded

* minor clean-up

* Correct unit of diagnostic variable dp_trc.

* Made conservation and checksum diagnostics selectable by namelist options (default off).

* pCO2, Piston velocity and solubility output (NorESMhub#202)

* add pCO2m (moist), CO2 piston velocity and solubility output - caution: kwco2 piston velocity now really holds only piston velocity (and not times solubility)

* Bugfix pnetcdf (NorESMhub#208)

* Add variables used by PNETCDF to explicit use staements.

* Move implicit none statments

* update explicit use statement for pnetcdf

* fixed units and renamed calcium burial to CaCO3 burial (NorESMhub#212)

Fixed sediment clay units.

* - Made the "fuk95" configuration work with MKS units.
- Removed "CGS" CPP flag.
- Changed some unit conversion factors from variables to parameters.
- Introduced rho0 (= 1/alpha0) parameter.
- Updated copyright statements.

* Correct unit conversion of mixed layer depth to pressure.

* Updated NorESM coupling scripts for the use of MKS units.

* Fixed check of unit system when building as NorESM component.

* Add option for surface pH output (NorESMhub#221)

* Remove unused parameters in wrt* subroutine calls in ncout_hamocc.F90

* Import get_bgc_namelist only in subroutine where it is needed. (NorESMhub#225)

* fix missing ' (NorESMhub#228)

Fixing a missing ' that only showed up when using `cisonew`

---------

Co-authored-by: Mats Bentsen <mben@norceresearch.no>
Co-authored-by: Tomas Torsvik <tomas.torsvik.work@gmail.com>
Co-authored-by: Mehmet Ilicak <milicak@itu.edu.tr>
Co-authored-by: Tomas Torsvik <43031053+TomasTorsvik@users.noreply.github.com>
Co-authored-by: Jörg Schwinger <jorg.schwinger@norceresearch.no>
  • Loading branch information
@jmaerz @matsbn
@TomasTorsvik @milicak @JorgSchwinger
6 people authored Feb 9, 2023
1 parent 26e65d6 commit c0ae6c2
Show file tree
Hide file tree
Showing 43 changed files with 828 additions and 557 deletions.
21 changes: 15 additions & 6 deletions ben02/mod_ben02.F
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
! ------------------------------------------------------------------------------
! Copyright (C) 2002-2021 Mats Bentsen, Mehmet Ilicak
! Copyright (C) 2002-2022 Mats Bentsen, Mehmet Ilicak
!
! This file is part of BLOM.
!
Expand All @@ -26,7 +26,7 @@ module mod_ben02
c
use mod_types, only: i2, r4
use mod_config, only: expcnf
use mod_constants, only: t0deg, spval
use mod_constants, only: t0deg, spval, L_mks2cgs
use mod_calendar, only: date_offset, calendar_noerr,
. calendar_errstr
use mod_time, only: date, calendar, nday_in_year, nday_of_year,
Expand Down Expand Up @@ -183,10 +183,17 @@ module mod_ben02
. atm_cswa_era ! short-wave radiation adjustment factor
! (NCEP)
c
#ifdef MKS
data atm_ice_csmt_ncep,atm_rnf_csmt_ncep /2.e10,1.e9/,
. atm_crnf_ncep,atm_cswa_ncep /0.82073,0.88340/,
. atm_ice_csmt_era,atm_rnf_csmt_era /0.0,1.e9/,
. atm_crnf_era,atm_cswa_era /0.7234,0.9721/
#else
data atm_ice_csmt_ncep,atm_rnf_csmt_ncep /2.e14,1.e13/,
. atm_crnf_ncep,atm_cswa_ncep /0.82073,0.88340/,
. atm_ice_csmt_era,atm_rnf_csmt_era /0.0,1.e13/,
. atm_crnf_era,atm_cswa_era /0.7234,0.9721/
#endif
c
real ::
. zu, ! measurement height of wind [m]
Expand Down Expand Up @@ -2090,11 +2097,13 @@ subroutine inifrc_ben02clim
integer, dimension(1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy,12) :: smtmsk
real dx2,dy2,prc_sum,eva_sum,rnf_sum,swa_sum,lwa_sum,lht_sum,
. sht_sum,fwf_fac,dangle,garea,le,albedo,fac,swa_ave,lwa_ave,
. lht_ave,sht_ave,crnf,cswa
. lht_ave,sht_ave,crnf,cswa,A_cgs2mks
real*4 rw4
integer i,j,k,l,il,jl
integer*2 rn2,ri2,rj2
c
A_cgs2mks=1./(L_mks2cgs**2)
c
c --- Allocate memory for additional monthly forcing fields.
allocate(taud (1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy,12),
. tauxd (1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy,12),
Expand Down Expand Up @@ -2775,7 +2784,7 @@ subroutine inifrc_ben02clim
do k=1,12
do l=1,isp(j)
do i=max(1,ifp(j,l)),min(ii,ilp(j,l))
garea=scp2(i,j)*1.e-4 ! [m^2]
garea=scp2(i,j)*A_cgs2mks ! [m^2]
c
c --- ----- freshwater fluxes [m/s]
util1(i,j)=util1(i,j)+precip(i,j,k)*fwf_fac*garea
Expand Down Expand Up @@ -2819,7 +2828,7 @@ subroutine inifrc_ben02clim
do j=1,jj
do l=1,isp(j)
do i=max(1,ifp(j,l)),min(ii,ilp(j,l))
garea=scp2(i,j)*1.e-4 ! [m^2]
garea=scp2(i,j)*A_cgs2mks ! [m^2]
c
c --- ----- heat fluxes
albedo=albs_f*ricclm(i,j,k)+albw(i,j)*(1.-ricclm(i,j,k))
Expand All @@ -2838,7 +2847,7 @@ subroutine inifrc_ben02clim
call xcsum(lht_sum,util3,ip)
call xcsum(sht_sum,util4,ip)
c
fac=1.e4/(12.*area)
fac=(L_mks2cgs*L_mks2cgs)/(12.*area)
swa_ave=swa_sum*fac
lwa_ave=lwa_sum*fac
lht_ave=lht_sum*fac
Expand Down
7 changes: 4 additions & 3 deletions ben02/sfcstr_ben02.F
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
! ------------------------------------------------------------------------------
! Copyright (C) 2004-2020 Mats Bentsen
! Copyright (C) 2004-2022 Mats Bentsen, Mehmet Ilicak
!
! This file is part of BLOM.
!
Expand All @@ -24,6 +24,7 @@ subroutine sfcstr_ben02(m,n,mm,nn,k1m,k1n)
c --- ------------------------------------------------------------------
c
use mod_xc
use mod_constants, only: P_mks2cgs
use mod_forcing, only: ztx, mty, taux, tauy
use mod_seaice, only: ficem, hicem, tauxice, tauyice
use mod_checksum, only: csdiag, chksummsk
Expand All @@ -44,14 +45,14 @@ subroutine sfcstr_ben02(m,n,mm,nn,k1m,k1n)
do i=max(1,ifu(j,l)),min(ii,ilu(j,l))
facice=(ficem(i,j)+ficem(i-1,j))
. *min(2.,hicem(i,j)+hicem(i-1,j))*.25
taux(i,j)=10.*(ztx(i,j)*(1.-facice)+tauxice(i,j)*facice)
taux(i,j)=P_mks2cgs*(ztx(i,j)*(1.-facice)+tauxice(i,j)*facice)
enddo
enddo
do l=1,isv(j)
do i=max(1,ifv(j,l)),min(ii,ilv(j,l))
facice=(ficem(i,j)+ficem(i,j-1))
. *min(2.,hicem(i,j)+hicem(i,j-1))*.25
tauy(i,j)=10.*(mty(i,j)*(1.-facice)+tauyice(i,j)*facice)
tauy(i,j)=P_mks2cgs*(mty(i,j)*(1.-facice)+tauyice(i,j)*facice)
enddo
enddo
enddo
Expand Down
57 changes: 31 additions & 26 deletions ben02/thermf_ben02.F
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
! ------------------------------------------------------------------------------
! Copyright (C) 2002-2021 Mats Bentsen
! Copyright (C) 2002-2022 Mats Bentsen, Mehmet Ilicak
!
! This file is part of BLOM.
!
Expand All @@ -21,7 +21,8 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c
c --- NERSC version of thermf.
c
use mod_constants, only: spcifh, t0deg, epsil, onem
use mod_constants, only: spcifh, t0deg, alpha0, epsilt, onem,
. g2kg, kg2g, L_mks2cgs, M_mks2cgs
use mod_time, only: nday_in_year, nday_of_year, nstep,
. nstep_in_day, baclin,
. xmi, l1mi, l2mi, l3mi, l4mi, l5mi
Expand Down Expand Up @@ -71,7 +72,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
. fice,hice,hsnw,tsrf,fice0,hice0,hsnw0,qsww,qnsw,tice,albi,
. tsmlt,albi_h,qswi,dh,qsnwf,fcond,qdamp,qsmlt,qo2i,qbot,swfac,
. dtml,q,volice,df,dvi,dvs,fwflx,sstc,rice,trxflx,sssc,srxflx,
. totsfl,totwfl,sflxc,totsrp,totsrn
. totsfl,totwfl,sflxc,totsrp,totsrn,A_cgs2mks
#ifdef TRC
integer nt
real, dimension(ntr,1-nbdy:idm+nbdy,1-nbdy:jdm+nbdy) ::
Expand All @@ -82,10 +83,12 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c
real intp1d
external intp1d
c
A_cgs2mks=1./(L_mks2cgs**2)
c
c --- Due to conservation, the ratio of ice and snow density must be
c --- equal to the ratio of ice and snow heat of fusion
if (abs(fuss/fusi-rhosnw/rhoice).gt.epsil) then
if (abs(fuss/fusi-rhosnw/rhoice).gt.epsilt) then
if (mnproc.eq.1) then
write (lp,*)
. 'thermf: check consistency between snow/ice densities'
Expand All @@ -97,7 +100,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c
c --- Set various constants
dt=baclin ! Time step
cpsw=spcifh*1.e3 ! Specific heat of seawater
cpsw=spcifh*M_mks2cgs ! Specific heat of seawater
rnf_fac=baclin/real(nrfets*86400) ! Runoff reservoar detrainment rate
sag_fac=exp(-sagets*dt) ! Snow aging rate
c
Expand Down Expand Up @@ -326,7 +329,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c --- ----- Ice volume that has to freeze to balance the heat budget
volice=-(qsww+qnsw-q)*(1.-fice)*dt/fusi
c
if (volice.gt.epsil) then
if (volice.gt.epsilt) then
c
c --- ------- New ice in the lead is formed with a specified thickness.
c --- ------- Estimate the change in ice fraction
Expand All @@ -344,7 +347,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c --- ----- If the lead is warming, let the fraction (1 - fice) go to
c --- ----- warm the lead, and the fraction fice to melt ice laterally
fice=fice-(swfac*qsww+qnsw)*fice*dt
. /max(hice*fusi+hsnw*fuss,epsil)
. /max(hice*fusi+hsnw*fuss,epsilt)
if (fice.lt.0.) then
fice=0.
hice=0.
Expand Down Expand Up @@ -398,14 +401,14 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
fwflx=eva(i,j)+lip(i,j)+sop(i,j)+rnf(i,j)+rfi(i,j)+fmltfz(i,j)
c
c --- --- Salt flux [kg m-2 s-1] (positive downwards)
sfl(i,j)=-sice*dvi*rhoice/dt*1.e-3
sfl(i,j)=-sice*dvi*rhoice/dt*g2kg
c
c --- --- Salt flux due to brine rejection of freezing sea
c --- --- ice [kg m-2 m-1] (positive downwards)
brnflx(i,j)=max(0.,-sotl*fmltfz(i,j)*1.e-3+sfl(i,j))
brnflx(i,j)=max(0.,-sotl*fmltfz(i,j)*g2kg+sfl(i,j))
c
c --- --- Virtual salt flux [kg m-2 s-1] (positive downwards)
vrtsfl(i,j)=-sotl*fwflx*1.e-3
vrtsfl(i,j)=-sotl*fwflx*g2kg
c
c --- --- Store area weighted virtual salt flux and fresh water flux
util1(i,j)=vrtsfl(i,j)*scp2(i,j)
Expand All @@ -415,11 +418,11 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
hmltfz(i,j)=(dvi*fusi+dvs*fuss)/dt
c
c --- --- Total heat flux in BLOM units [W cm-2] (positive upwards)
surflx(i,j)=-(swa(i,j)+nsf(i,j)+hmltfz(i,j))*1.e-4
surflx(i,j)=-(swa(i,j)+nsf(i,j)+hmltfz(i,j))*A_cgs2mks
c
c --- --- Short-wave heat flux in BLOM units [W cm-2] (positive
c --- --- upwards)
sswflx(i,j)=-qsww*(1.-fice0)*1.e-4
sswflx(i,j)=-qsww*(1.-fice0)*A_cgs2mks
c
#ifdef TRC
c --- ------------------------------------------------------------------
Expand Down Expand Up @@ -452,7 +455,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
endif
# endif
# endif
trflx(nt,i,j)=-trc(i,j,k1n,nt)*fwflx*1.e-3
trflx(nt,i,j)=-trc(i,j,k1n,nt)*fwflx*g2kg
ttrsf(nt,i,j)=trflx(nt,i,j)*scp2(i,j)
ttrav(nt,i,j)=trc(i,j,k1n,nt)*scp2(i,j)
enddo
Expand All @@ -465,16 +468,16 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
surrlx(i,j)=0.
c
c --- --- If trxday>0 , apply relaxation towards observed sst
if (trxday.gt.epsil) then
if (trxday.gt.epsilt) then
sstc=intp1d(sstclm(i,j,l1mi),sstclm(i,j,l2mi),
. sstclm(i,j,l3mi),sstclm(i,j,l4mi),
. sstclm(i,j,l5mi),xmi)
rice=intp1d(ricclm(i,j,l1mi),ricclm(i,j,l2mi),
. ricclm(i,j,l3mi),ricclm(i,j,l4mi),
. ricclm(i,j,l5mi),xmi)
sstc=(1.-rice)*max(sstc,tice_f)+rice*tice_f
trxflx=spcifh*100.*min(hmxl,trxdpt)/(trxday*86400.)
. *min(trxlim,max(-trxlim,sstc-tmxl))
trxflx=spcifh*L_mks2cgs*min(hmxl,trxdpt)/(trxday*86400.)
. *min(trxlim,max(-trxlim,sstc-tmxl))/alpha0
surrlx(i,j)=-trxflx
else
trxflx=0.
Expand All @@ -496,12 +499,12 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
salrlx(i,j)=0.
c
c --- --- if srxday>0 , apply relaxation towards observed sss
if (srxday.gt.epsil) then
if (srxday.gt.epsilt) then
sssc=intp1d(sssclm(i,j,l1mi),sssclm(i,j,l2mi),
. sssclm(i,j,l3mi),sssclm(i,j,l4mi),
. sssclm(i,j,l5mi),xmi)
srxflx=100.*min(hmxl,srxdpt)/(srxday*86400.)
. *min(srxlim,max(-srxlim,sssc-smxl))
srxflx=L_mks2cgs*min(hmxl,srxdpt)/(srxday*86400.)
. *min(srxlim,max(-srxlim,sssc-smxl))/alpha0
salrlx(i,j)=-srxflx
util3(i,j)=max(0.,salrlx(i,j))*scp2(i,j)
util4(i,j)=min(0.,salrlx(i,j))*scp2(i,j)
Expand Down Expand Up @@ -538,7 +541,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c --- -------------------------------------------------------------------
c
ustar(i,j)=(min(ustari(i,j),.8e-2)*fice0
. +ustarw(i,j)*(1.-fice0))*1.e2
. +ustarw(i,j)*(1.-fice0))*L_mks2cgs
c
enddo
enddo
Expand All @@ -556,7 +559,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
call xcsum(totwfl,util2,ips)
c
c --- Correction for the virtual salt flux [kg m-2 s-1]
sflxc=(-sref*totwfl*1.e-3-totsfl)/area
sflxc=(-sref*totwfl*g2kg-totsfl)/area
if (mnproc.eq.1) then
write (lp,*) 'thermf: totsfl/area,sflxc',totsfl/area,sflxc
endif
Expand All @@ -567,8 +570,10 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
do j=1,jj
do l=1,isp(j)
do i=max(1,ifp(j,l)),min(ii,ilp(j,l))
salflx(i,j)=-(vrtsfl(i,j)+sflxc+sfl(i,j))*1.e2
brnflx(i,j)=-brnflx(i,j)*1.e2
salflx(i,j)=-(vrtsfl(i,j)+sflxc+sfl(i,j))
. *(kg2g*(M_mks2cgs/L_mks2cgs**2))
brnflx(i,j)=-brnflx(i,j)
. *(kg2g*(M_mks2cgs/L_mks2cgs**2))
enddo
enddo
enddo
Expand All @@ -577,7 +582,7 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c --- if srxday>0 and srxbal=.true. , balance the sss relaxation flux
c --- so the net input of salt in grid cells connected to the world
c --- ocean is zero
if (srxday.gt.epsil.and.srxbal) then
if (srxday.gt.epsilt.and.srxbal) then
call xcsum(totsrp,util3,ipwocn)
call xcsum(totsrn,util4,ipwocn)
if (abs(totsrp).gt.abs(totsrn)) then
Expand Down Expand Up @@ -632,14 +637,14 @@ subroutine thermf_ben02(m,n,mm,nn,k1m,k1n)
c tottrav=tottrav/area
cc
c trflxc=(-tottrsf)/area
c trflxc=(-tottrav*totwfl*1.e-3-tottrsf)/area
c trflxc=(-tottrav*totwfl*g2kg-tottrsf)/area
trflxc=0.
c
c$OMP PARALLEL DO PRIVATE(l,i)
do j=1,jj
do l=1,isp(j)
do i=max(1,ifp(j,l)),min(ii,ilp(j,l))
trflx(nt,i,j)=-(trflx(nt,i,j)+trflxc)*1.e2
trflx(nt,i,j)=-(trflx(nt,i,j)+trflxc)*L_mks2cgs
enddo
enddo
enddo
Expand Down
7 changes: 4 additions & 3 deletions cesm/sfcstr_cesm.F
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
! ------------------------------------------------------------------------------
! Copyright (C) 2015-2020 Mats Bentsen
! Copyright (C) 2015-2022 Mats Bentsen, Mehmet Ilicak
!
! This file is part of BLOM.
!
Expand All @@ -24,6 +24,7 @@ subroutine sfcstr_cesm(m,n,mm,nn,k1m,k1n)
c --- ------------------------------------------------------------------
c
use mod_xc
use mod_constants, only: P_mks2cgs
use mod_forcing, only: ztx, mty, taux, tauy
use mod_checksum, only: csdiag, chksummsk
c
Expand All @@ -37,12 +38,12 @@ subroutine sfcstr_cesm(m,n,mm,nn,k1m,k1n)
do j=1,jj
do l=1,isu(j)
do i=max(1,ifu(j,l)),min(ii,ilu(j,l))
taux(i,j)=10.*ztx(i,j)
taux(i,j)=P_mks2cgs*ztx(i,j)
enddo
enddo
do l=1,isv(j)
do i=max(1,ifv(j,l)),min(ii,ilv(j,l))
tauy(i,j)=10.*mty(i,j)
tauy(i,j)=P_mks2cgs*mty(i,j)
enddo
enddo
enddo
Expand Down
Loading

0 comments on commit c0ae6c2

Please sign in to comment.