Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • nemo/nemo
  • sparonuz/nemo
  • hatfield/nemo
  • extdevs/nemo
4 results
Show changes
Showing
with 933 additions and 532 deletions
......@@ -36,6 +36,8 @@ MODULE p4zsed
REAL(wp), SAVE :: r1_rday
REAL(wp), SAVE :: sedsilfrac, sedcalfrac
LOGICAL :: l_dia_sdenit, l_dia_nfix, l_dia_sed
!! * Substitutions
# include "do_loop_substitute.h90"
# include "domzgr_substitute.h90"
......@@ -66,44 +68,41 @@ CONTAINS
REAL(wp) :: zsiloss, zcaloss, zws3, zws4, zwsc, zdep
REAL(wp) :: zwstpoc, zwstpon, zwstpop
REAL(wp) :: ztrfer, ztrpo4s, ztrdp, zwdust, zmudia, ztemp
REAL(wp) :: zsoufer, zlight, ztrpo4, ztrdop
REAL(wp) :: xdiano3, xdianh4
!
CHARACTER (len=25) :: charout
REAL(wp), DIMENSION(jpi,jpj ) :: zdenit2d, zbureff, zwork
REAL(wp), DIMENSION(jpi,jpj ) :: zwsbio3, zwsbio4
REAL(wp), DIMENSION(jpi,jpj ) :: zsedcal, zsedsi, zsedc
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsoufer, zlight
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrpo4, ztrdop, zirondep, zpdep
REAL(wp), DIMENSION(A2D(0)) :: zdenit2d, zbureff
REAL(wp), DIMENSION(A2D(0)) :: zwsbio3, zwsbio4
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zsedcal, zsedsi, zsedc
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p4z_sed')
!
! Allocate temporary workspace
ALLOCATE( ztrpo4(jpi,jpj,jpk) )
IF( ln_p5z ) ALLOCATE( ztrdop(jpi,jpj,jpk) )
zdenit2d(:,:) = 0.e0
zbureff (:,:) = 0.e0
zwork (:,:) = 0.e0
zsedsi (:,:) = 0.e0
zsedcal (:,:) = 0.e0
zsedc (:,:) = 0.e0
IF( kt == nittrc000 ) THEN
l_dia_nfix = iom_use( "Nfix" )
l_dia_sdenit = iom_use( "Sdenit" )
l_dia_sed = .NOT.lk_sed .AND. ( iom_use( "SedC" ) .OR. iom_use( "SedCal" ) .OR. iom_use( "SedSi" ) )
ENDIF
! OA: Warning, the following part is necessary to avoid CFL problems above the sediments
! --------------------------------------------------------------------
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = e3t(ji,jj,ikt,Kmm) / xstep
zwsbio4(ji,jj) = MIN( 0.99 * zdep, wsbio4(ji,jj,ikt) )
zwsbio3(ji,jj) = MIN( 0.99 * zdep, wsbio3(ji,jj,ikt) )
END_2D
!
zdenit2d(:,:) = 0.e0
zbureff (:,:) = 0.e0
!
IF( .NOT.lk_sed ) THEN
! Computation of the sediment denitrification proportion: The metamodel from midlleburg (2006) is being used
! Computation of the fraction of organic matter that is permanently buried from Dunne's model
! -------------------------------------------------------
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF( tmask(ji,jj,1) == 1 ) THEN
ikt = mbkt(ji,jj)
zflx = ( tr(ji,jj,ikt,jpgoc,Kbb) * zwsbio4(ji,jj) &
......@@ -129,7 +128,7 @@ CONTAINS
! ------------------------------------------------------
IF( .NOT.lk_sed ) zrivsil = 1._wp - sedsilfrac
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = xstep / e3t(ji,jj,ikt,Kmm)
zwsc = zwsbio4(ji,jj) * zdep
......@@ -141,7 +140,7 @@ CONTAINS
END_2D
!
IF( .NOT.lk_sed ) THEN
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = xstep / e3t(ji,jj,ikt,Kmm)
zwsc = zwsbio4(ji,jj) * zdep
......@@ -154,12 +153,11 @@ CONTAINS
zrivalk = sedcalfrac * zfactcal
tr(ji,jj,ikt,jptal,Krhs) = tr(ji,jj,ikt,jptal,Krhs) + zcaloss * zrivalk * 2.0
tr(ji,jj,ikt,jpdic,Krhs) = tr(ji,jj,ikt,jpdic,Krhs) + zcaloss * zrivalk
zsedcal(ji,jj) = (1.0 - zrivalk) * zcaloss * e3t(ji,jj,ikt,Kmm)
zsedsi (ji,jj) = (1.0 - zrivsil) * zsiloss * e3t(ji,jj,ikt,Kmm)
END_2D
ENDIF
!
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
!
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = xstep / e3t(ji,jj,ikt,Kmm)
zws4 = zwsbio4(ji,jj) * zdep
......@@ -171,7 +169,7 @@ CONTAINS
END_2D
!
IF( ln_p5z ) THEN
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = xstep / e3t(ji,jj,ikt,Kmm)
zws4 = zwsbio4(ji,jj) * zdep
......@@ -183,10 +181,10 @@ CONTAINS
END_2D
ENDIF
! The 0.5 factor in zpdenit is to avoid negative NO3 concentration after
! denitrification in the sediments. Not very clever, but simpliest option.
IF( .NOT.lk_sed ) THEN
! The 0.5 factor in zpdenit is to avoid negative NO3 concentration after
! denitrification in the sediments. Not very clever, but simpliest option.
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = xstep / e3t(ji,jj,ikt,Kmm)
zws4 = zwsbio4(ji,jj) * zdep
......@@ -204,27 +202,37 @@ CONTAINS
tr(ji,jj,ikt,jptal,Krhs) = tr(ji,jj,ikt,jptal,Krhs) + rno3 * (zolimit + (1.+rdenit) * zpdenit )
tr(ji,jj,ikt,jpdic,Krhs) = tr(ji,jj,ikt,jpdic,Krhs) + zpdenit + zolimit
sdenit(ji,jj) = rdenit * zpdenit * e3t(ji,jj,ikt,Kmm)
zsedc(ji,jj) = (1. - zrivno3) * zwstpoc * e3t(ji,jj,ikt,Kmm)
IF( ln_p5z ) THEN
zwstpop = tr(ji,jj,ikt,jpgop,Kbb) * zws4 + tr(ji,jj,ikt,jppop,Kbb) * zws3
zwstpon = tr(ji,jj,ikt,jpgon,Kbb) * zws4 + tr(ji,jj,ikt,jppon,Kbb) * zws3
END_2D
IF( ln_p5z ) THEN
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zdep = xstep / e3t(ji,jj,ikt,Kmm)
zws4 = zwsbio4(ji,jj) * zdep
zws3 = zwsbio3(ji,jj) * zdep
zrivno3 = 1. - zbureff(ji,jj)
zwstpoc = tr(ji,jj,ikt,jpgoc,Kbb) * zws4 + tr(ji,jj,ikt,jppoc,Kbb) * zws3
zpdenit = MIN( 0.5 * ( tr(ji,jj,ikt,jpno3,Kbb) - rtrn ) / rdenit, zdenit2d(ji,jj) * zwstpoc * zrivno3 )
z1pdenit = zwstpoc * zrivno3 - zpdenit
zolimit = MIN( ( tr(ji,jj,ikt,jpoxy,Kbb) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) )
zwstpop = tr(ji,jj,ikt,jpgop,Kbb) * zws4 + tr(ji,jj,ikt,jppop,Kbb) * zws3
zwstpon = tr(ji,jj,ikt,jpgon,Kbb) * zws4 + tr(ji,jj,ikt,jppon,Kbb) * zws3
tr(ji,jj,ikt,jpdon,Krhs) = tr(ji,jj,ikt,jpdon,Krhs) + ( z1pdenit - zolimit ) * zwstpon / (zwstpoc + rtrn)
tr(ji,jj,ikt,jpdop,Krhs) = tr(ji,jj,ikt,jpdop,Krhs) + ( z1pdenit - zolimit ) * zwstpop / (zwstpoc + rtrn)
ENDIF
END_2D
END_2D
ENDIF
ENDIF
! Nitrogen fixation process
! Small source iron from particulate inorganic iron
!-----------------------------------
DO jk = 1, jpkm1
zlight (:,:,jk) = ( 1.- EXP( -etot_ndcy(:,:,jk) / diazolight ) ) * ( 1. - fr_i(:,:) )
zsoufer(:,:,jk) = zlight(:,:,jk) * 2E-11 / ( 2E-11 + biron(:,:,jk) )
ENDDO
!
IF( ln_p4z ) THEN
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
! ! Potential nitrogen fixation dependant on temperature and iron
zlight = ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) ) * ( 1. - fr_i(ji,jj) )
zsoufer = zlight * 2E-11 / ( 2E-11 + biron(ji,jj,jk) )
!
ztemp = ts(ji,jj,jk,jp_tem,Kmm)
zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) / rno3
! Potential nitrogen fixation dependant on temperature and iron
......@@ -234,33 +242,11 @@ CONTAINS
IF( zlim <= 0.1 ) zlim = 0.01
zfact = zlim * rfact2
ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) )
ztrpo4(ji,jj,jk) = tr(ji,jj,jk,jppo4,Kbb) / ( 1E-6 + tr(ji,jj,jk,jppo4,Kbb) )
ztrdp = ztrpo4(ji,jj,jk)
nitrpot(ji,jj,jk) = zmudia * r1_rday * zfact * MIN( ztrfer, ztrdp ) * zlight(ji,jj,jk)
END_3D
ELSE ! p5z
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
! ! Potential nitrogen fixation dependant on temperature and iron
ztemp = ts(ji,jj,jk,jp_tem,Kmm)
zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) * 7.625
! Potential nitrogen fixation dependant on temperature and iron
xdianh4 = tr(ji,jj,jk,jpnh4,Kbb) / ( concnnh4 + tr(ji,jj,jk,jpnh4,Kbb) )
xdiano3 = tr(ji,jj,jk,jpno3,Kbb) / ( concnno3 + tr(ji,jj,jk,jpno3,Kbb) ) * (1. - xdianh4)
zlim = ( 1.- xdiano3 - xdianh4 )
IF( zlim <= 0.1 ) zlim = 0.01
zfact = zlim * rfact2
ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) )
ztrpo4(ji,jj,jk) = tr(ji,jj,jk,jppo4,Kbb) / ( 1E-6 + tr(ji,jj,jk,jppo4,Kbb) )
ztrdop(ji,jj,jk) = tr(ji,jj,jk,jpdop,Kbb) / ( 1E-6 + tr(ji,jj,jk,jpdop,Kbb) ) * (1. - ztrpo4(ji,jj,jk))
ztrdp = ztrpo4(ji,jj,jk) + ztrdop(ji,jj,jk)
nitrpot(ji,jj,jk) = zmudia * r1_rday * zfact * MIN( ztrfer, ztrdp ) * zlight(ji,jj,jk)
END_3D
ENDIF
! Nitrogen change due to nitrogen fixation
! ----------------------------------------
IF( ln_p4z ) THEN
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
ztrpo4 = tr(ji,jj,jk,jppo4,Kbb) / ( 1E-6 + tr(ji,jj,jk,jppo4,Kbb) )
nitrpot(ji,jj,jk) = zmudia * r1_rday * zfact * MIN( ztrfer, ztrpo4 ) * zlight
!
! Nitrogen change due to nitrogen fixation
! ----------------------------------------
zfact = nitrpot(ji,jj,jk) * nitrfix
tr(ji,jj,jk,jpnh4,Krhs) = tr(ji,jj,jk,jpnh4,Krhs) + zfact / 3.0
tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) + rno3 * zfact / 3.0
......@@ -273,23 +259,41 @@ CONTAINS
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - 30E-6 * zfact * 1.0 / 3.0
tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0
tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.002 * 4E-10 * zsoufer * rfact2 / rday
tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) + concdnh4 / ( concdnh4 + tr(ji,jj,jk,jppo4,Kbb) ) &
& * 0.001 * tr(ji,jj,jk,jpdoc,Kbb) * xstep
END_3D
ELSE ! p5z
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
ELSE ! p5z
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zlight = ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) ) * ( 1. - fr_i(ji,jj) )
zsoufer = zlight * 2E-11 / ( 2E-11 + biron(ji,jj,jk) )
!
! ! Potential nitrogen fixation dependant on temperature and iron
ztemp = ts(ji,jj,jk,jp_tem,Kmm)
zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) * 7.625
! Potential nitrogen fixation dependant on temperature and iron
xdianh4 = tr(ji,jj,jk,jpnh4,Kbb) / ( concnnh4 + tr(ji,jj,jk,jpnh4,Kbb) )
xdiano3 = tr(ji,jj,jk,jpno3,Kbb) / ( concnno3 + tr(ji,jj,jk,jpno3,Kbb) ) * (1. - xdianh4)
zlim = ( 1.- xdiano3 - xdianh4 )
IF( zlim <= 0.1 ) zlim = 0.01
zfact = zlim * rfact2
ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) )
ztrpo4 = tr(ji,jj,jk,jppo4,Kbb) / ( 1E-6 + tr(ji,jj,jk,jppo4,Kbb) )
ztrdop = tr(ji,jj,jk,jpdop,Kbb) / ( 1E-6 + tr(ji,jj,jk,jpdop,Kbb) ) * (1. - ztrpo4)
nitrpot(ji,jj,jk) = zmudia * r1_rday * zfact * MIN( ztrfer, ztrpo4 + ztrdop ) * zlight
! Nitrogen change due to nitrogen fixation
! ----------------------------------------
zfact = nitrpot(ji,jj,jk) * nitrfix
tr(ji,jj,jk,jpnh4,Krhs) = tr(ji,jj,jk,jpnh4,Krhs) + zfact / 3.0
tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) + rno3 * zfact / 3.0
tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) - zfact * 2.0 / 3.0
tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) - 16.0 / 46.0 * zfact * ( 1.0 - 1.0 / 3.0 ) &
& * ztrpo4(ji,jj,jk) / (ztrpo4(ji,jj,jk) + ztrdop(ji,jj,jk) + rtrn)
& * ztrpo4 / (ztrpo4 + ztrdop + rtrn)
tr(ji,jj,jk,jpdon,Krhs) = tr(ji,jj,jk,jpdon,Krhs) + zfact * 1.0 / 3.0
tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) + zfact * 1.0 / 3.0
tr(ji,jj,jk,jpdop,Krhs) = tr(ji,jj,jk,jpdop,Krhs) + 16.0 / 46.0 * zfact / 3.0 &
& - 16.0 / 46.0 * zfact * ztrdop(ji,jj,jk) &
& / (ztrpo4(ji,jj,jk) + ztrdop(ji,jj,jk) + rtrn)
& - 16.0 / 46.0 * zfact * ztrdop / ( ztrpo4 + ztrdop + rtrn)
tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) + zfact * 1.0 / 3.0 * 2.0 / 3.0
tr(ji,jj,jk,jppon,Krhs) = tr(ji,jj,jk,jppon,Krhs) + zfact * 1.0 / 3.0 * 2.0 /3.0
tr(ji,jj,jk,jppop,Krhs) = tr(ji,jj,jk,jppop,Krhs) + 16.0 / 46.0 * zfact * 1.0 / 3.0 * 2.0 /3.0
......@@ -300,18 +304,46 @@ CONTAINS
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - 30E-6 * zfact * 1.0 / 3.0
tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0
tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.002 * 4E-10 * zsoufer * rfact2 / rday
END_3D
!
ENDIF
IF( lk_iomput .AND. knt == nrdttrc ) THEN
zfact = 1.e+3 * rfact2r ! conversion from molC/l/kt to molN/m3/s
CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * rno3 * zfact * tmask(:,:,:) ) ! nitrogen fixation
CALL iom_put( "SedCal", zsedcal(:,:) * zfact )
CALL iom_put( "SedSi" , zsedsi (:,:) * zfact )
CALL iom_put( "SedC" , zsedc (:,:) * zfact )
CALL iom_put( "Sdenit", sdenit (:,:) * zfact * rno3 )
!
IF( l_dia_nfix ) THEN ! nitrogen fixation
zfact = rno3 * 1.e+3 * rfact2r ! conversion from molC/l/kt to molN/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = nitrpot(A2D(0),1:jpkm1) * nitrfix * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "Nfix", zw3d )
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_sed ) THEN
ALLOCATE( zsedcal(A2D(0)), zsedsi(A2D(0) ), zsedc(A2D(0) ) )
zfact = 1.e+3 * rfact2r ! conversion from molC/l/kt to molC/m3/s
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
zsedsi(ji,jj) = (1.0 - zrivsil) * tr(ji,jj,ikt,jpgsi,Kbb) * zwsbio4(ji,jj) * xstep
!
zfactcal = MAX(-0.1, MIN( excess(ji,jj,ikt), 0.2 ) )
zfactcal = 0.3 + 0.7 * MIN( 1., (0.1 + zfactcal) / ( 0.5 - zfactcal ) )
zrivalk = sedcalfrac * zfactcal
zsedcal(ji,jj) = (1.0 - zrivalk) * tr(ji,jj,ikt,jpcal,Kbb) * zwsbio4(ji,jj) * xstep
!
zrivno3 = 1. - zbureff(ji,jj)
zsedc(ji,jj) = (1. - zrivno3) * &
& ( tr(ji,jj,ikt,jpgoc,Kbb) * zwsbio4(ji,jj) + tr(ji,jj,ikt,jppoc,Kbb) * zwsbio3(ji,jj) ) * xstep
END_2D
CALL iom_put( "SedCal", zsedcal(:,:) * zfact )
CALL iom_put( "SedSi" , zsedsi (:,:) * zfact )
CALL iom_put( "SedC" , zsedc (:,:) * zfact )
DEALLOCATE( zsedcal, zsedsi, zsedc )
ENDIF
IF( l_dia_sdenit ) THEN
zfact = rno3 * 1.e+3 * rfact2r ! conversion from molC/l/kt to molN/m3/s
CALL iom_put( "Sdenit", sdenit(:,:) * zfact )
ENDIF
!
ENDIF
!
IF(sn_cfctl%l_prttrc) THEN ! print mean trneds (USEd for debugging)
......@@ -320,8 +352,6 @@ CONTAINS
CALL prt_ctl(tab4d_1=tr(:,:,:,:,Krhs), mask1=tmask, clinfo=ctrcnm)
ENDIF
!
IF( ln_p5z ) DEALLOCATE( ztrpo4, ztrdop )
!
IF( ln_timing ) CALL timing_stop('p4z_sed')
!
END SUBROUTINE p4z_sed
......@@ -367,7 +397,7 @@ CONTAINS
!
lk_sed = ln_sediment .AND. ln_sed_2way
!
nitrpot(:,:,jpk) = 0._wp ! define last level for iom_put
! nitrpot(:,:,jpk) = 0._wp ! define last level for iom_put
!
END SUBROUTINE p4z_sed_init
......@@ -375,7 +405,7 @@ CONTAINS
!!----------------------------------------------------------------------
!! *** ROUTINE p4z_sed_alloc ***
!!----------------------------------------------------------------------
ALLOCATE( nitrpot(jpi,jpj,jpk), sdenit(jpi,jpj), STAT=p4z_sed_alloc )
ALLOCATE( nitrpot(A2D(0),jpk), sdenit(A2D(0)), STAT=p4z_sed_alloc )
!
IF( p4z_sed_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p4z_sed_alloc: failed to allocate arrays' )
!
......
......@@ -40,6 +40,7 @@ MODULE p4zsink
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sinkfer2 !: Big iron sinking fluxes
INTEGER :: ik100
LOGICAL :: l_dia_sink2d, l_dia_sink3d, l_dia_tcexp
!! * Substitutions
# include "do_loop_substitute.h90"
......@@ -68,10 +69,20 @@ CONTAINS
INTEGER :: ji, jj, jk
CHARACTER (len=25) :: charout
REAL(wp) :: zmax, zfact
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
REAL(wp), ALLOCATABLE, DIMENSION(:,: ) :: zw2d
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p4z_sink')
IF( kt == nittrc000 ) THEN
l_dia_sink2d = iom_use( "EPC100" ) .OR. iom_use( "EPFE100" ) &
& .OR. iom_use( "EPCAL100" ) .OR. iom_use( "EPSI100" )
l_dia_sink3d = iom_use( "EXPC" ) .OR. iom_use( "EXPFE" ) &
& .OR. iom_use( "EXPCAL" ) .OR. iom_use( "EXPSI" )
l_dia_tcexp = iom_use( "tcexp" )
ENDIF
! Initialization of some global variables
! ---------------------------------------
prodpoc(:,:,:) = 0.
......@@ -86,7 +97,7 @@ CONTAINS
! CaCO3 and bSi are supposed to sink at the big particles speed
! due to their high density
! ---------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zmax = MAX( heup_01(ji,jj), hmld(ji,jj) )
zfact = MAX( 0., gdepw(ji,jj,jk+1,Kmm) - zmax ) / wsbio2scale
wsbio4(ji,jj,jk) = wsbio2 + MAX(0., ( wsbio2max - wsbio2 )) * zfact
......@@ -129,28 +140,61 @@ CONTAINS
ENDIF
! Total carbon export per year
IF( iom_use( "tcexp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) &
& t_oce_co2_exp = glob_sum( 'p4zsink', ( sinking(:,:,ik100) + sinking2(:,:,ik100) ) * e1e2t(:,:) * tmask(:,:,1) )
IF( l_dia_tcexp .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) THEN
ALLOCATE( zw2d(A2D(0)) )
zw2d(A2D(0)) = ( sinking(A2D(0),ik100) + sinking2(A2D(0),ik100) ) * e1e2t(A2D(0)) * tmask(A2D(0),1)
t_oce_co2_exp = glob_sum( 'p4zsink', zw2d(:,:) )
DEALLOCATE( zw2d )
ENDIF
!
IF( lk_iomput .AND. knt == nrdttrc ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
!
CALL iom_put( "EPC100" , ( sinking(:,:,ik100) + sinking2(:,:,ik100) ) * zfact * tmask(:,:,1) ) ! Export of carbon at 100m
CALL iom_put( "EPFE100" , ( sinkfer(:,:,ik100) + sinkfer2(:,:,ik100) ) * zfact * tmask(:,:,1) ) ! Export of iron at 100m
CALL iom_put( "EPCAL100", sinkcal(:,:,ik100) * zfact * tmask(:,:,1) ) ! Export of calcite at 100m
CALL iom_put( "EPSI100" , sinksil(:,:,ik100) * zfact * tmask(:,:,1) ) ! Export of bigenic silica at 100m
CALL iom_put( "EXPC" , ( sinking(:,:,:) + sinking2(:,:,:) ) * zfact * tmask(:,:,:) ) ! Export of carbon in the water column
CALL iom_put( "EXPFE" , ( sinkfer(:,:,:) + sinkfer2(:,:,:) ) * zfact * tmask(:,:,:) ) ! Export of iron
CALL iom_put( "EXPCAL" , sinkcal(:,:,:) * zfact * tmask(:,:,:) ) ! Export of calcite
CALL iom_put( "EXPSI" , sinksil(:,:,:) * zfact * tmask(:,:,:) ) ! Export of bigenic silica
CALL iom_put( "tcexp" , t_oce_co2_exp * zfact ) ! molC/s
IF( l_dia_sink2d ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw2d(A2D(0)) )
! Export of carbon at 100m
zw2d(A2D(0)) = ( sinking(A2D(0),ik100) + sinking2(A2D(0),ik100) ) * zfact * tmask(A2D(0),1)
CALL iom_put( "EPC100", zw2d )
! Export of iron at 100m
zw2d(A2D(0)) = ( sinkfer(A2D(0),ik100) + sinkfer2(A2D(0),ik100) ) * zfact * tmask(A2D(0),1)
CALL iom_put( "EPFE100", zw2d )
! Export of calcite at 100m
zw2d(A2D(0)) = sinkcal(A2D(0),ik100) * zfact * tmask(A2D(0),1)
CALL iom_put( "EPCAL100", zw2d )
! Export of bigenic silica at 100m
zw2d(A2D(0)) = sinksil(A2D(0),ik100) * zfact * tmask(A2D(0),1)
CALL iom_put( "EPSI100", zw2d )
!
DEALLOCATE( zw2d )
ENDIF
!
IF( l_dia_sink3d ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
! Export of carbon in the water column
zw3d(A2D(0),1:jpkm1) = ( sinking(A2D(0),1:jpkm1) + sinking2(A2D(0),1:jpkm1) ) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "EXPC", zw3d )
! Export of iron
zw3d(A2D(0),1:jpkm1) = ( sinkfer(A2D(0),1:jpkm1) + sinkfer2(A2D(0),1:jpkm1) ) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "EXPFE", zw3d )
! Export of calcite
zw3d(A2D(0),1:jpkm1) = sinkcal(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "EXPCAL", zw3d )
! Export of bigenic silica
zw3d(A2D(0),1:jpkm1) = sinksil(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "EXPSI", zw3d )
!
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_tcexp ) CALL iom_put( "tcexp", t_oce_co2_exp * 1.e+3 * rfact2r )
!
ENDIF
!
IF(sn_cfctl%l_prttrc) THEN ! print mean trends (used for debugging)
WRITE(charout, FMT="('sink')")
CALL prt_ctl_info( charout, cdcomp = 'top' )
CALL prt_ctl(tab4d_1=tr(:,:,:,:,Krhs), mask1=tmask, clinfo=ctrcnm)
CALL prt_ctl(tab4d_1=tr(:,:,:,:,Kbb), mask1=tmask, clinfo=ctrcnm)
ENDIF
!
IF( ln_timing ) CALL timing_stop('p4z_sink')
......@@ -192,13 +236,13 @@ CONTAINS
!
ierr(:) = 0
!
ALLOCATE( sinking(jpi,jpj,jpk) , sinking2(jpi,jpj,jpk) , &
& sinkcal(jpi,jpj,jpk) , sinksil (jpi,jpj,jpk) , &
& sinkfer2(jpi,jpj,jpk) , &
& sinkfer(jpi,jpj,jpk) , STAT=ierr(1) )
ALLOCATE( sinking(A2D(0),jpk) , sinking2(A2D(0),jpk) , &
& sinkcal(A2D(0),jpk) , sinksil (A2D(0),jpk) , &
& sinkfer2(A2D(0),jpk) , &
& sinkfer(A2D(0),jpk) , STAT=ierr(1) )
!
IF( ln_p5z ) ALLOCATE( sinkingn(jpi,jpj,jpk), sinking2n(jpi,jpj,jpk) , &
& sinkingp(jpi,jpj,jpk), sinking2p(jpi,jpj,jpk) , STAT=ierr(2) )
IF( ln_p5z ) ALLOCATE( sinkingn(A2D(0),jpk), sinking2n(A2D(0),jpk) , &
& sinkingp(A2D(0),jpk), sinking2p(A2D(0),jpk) , STAT=ierr(2) )
!
p4z_sink_alloc = MAXVAL( ierr )
IF( p4z_sink_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p4z_sink_alloc : failed to allocate arrays.' )
......
......@@ -140,7 +140,7 @@ CONTAINS
! ------------------------------------------------------------------
xnegtr(:,:,:) = 1.e0
DO jn = jp_pcs0, jp_pcs1
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk)
DO_3D( 0, 0, 0, 0, 1, jpk)
IF( ( tr(ji,jj,jk,jn,Kbb) + tr(ji,jj,jk,jn,Krhs) ) < 0.e0 ) THEN
ztra = ABS( tr(ji,jj,jk,jn,Kbb) ) / ( ABS( tr(ji,jj,jk,jn,Krhs) ) + rtrn )
xnegtr(ji,jj,jk) = MIN( xnegtr(ji,jj,jk), ztra )
......@@ -157,45 +157,56 @@ CONTAINS
IF( iom_use( 'INTdtAlk' ) .OR. iom_use( 'INTdtDIC' ) .OR. iom_use( 'INTdtFer' ) .OR. &
& iom_use( 'INTdtDIN' ) .OR. iom_use( 'INTdtDIP' ) .OR. iom_use( 'INTdtSil' ) ) THEN
!
ALLOCATE( zw3d(jpi,jpj,jpk), zw2d(jpi,jpj) )
zw3d(:,:,jpk) = 0.
DO jk = 1, jpkm1
zw3d(:,:,jk) = xnegtr(:,:,jk) * xfact * e3t(:,:,jk,Kmm) * tmask(:,:,jk)
ENDDO
ALLOCATE( zw3d(A2D(0),jpk), zw2d(A2D(0)) )
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zw3d(ji,jj,jk) = xnegtr(ji,jj,jk) * xfact * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
END_3D
!
zw2d(:,:) = 0.
DO jk = 1, jpkm1
zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jptal,Krhs)
DO_2D( 0, 0, 0, 0 )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) * tr(ji,jj,jk,jptal,Krhs)
END_2D
ENDDO
CALL iom_put( 'INTdtAlk', zw2d )
!
zw2d(:,:) = 0.
DO jk = 1, jpkm1
zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jpdic,Krhs)
DO_2D( 0, 0, 0, 0 )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) * tr(ji,jj,jk,jpdic,Krhs)
END_2D
ENDDO
CALL iom_put( 'INTdtDIC', zw2d )
!
zw2d(:,:) = 0.
DO jk = 1, jpkm1
zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * rno3 * ( tr(:,:,jk,jpno3,Krhs) + tr(:,:,jk,jpnh4,Krhs) )
DO_2D( 0, 0, 0, 0 )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) * rno3 * ( tr(ji,jj,jk,jpno3,Krhs) + tr(ji,jj,jk,jpnh4,Krhs) )
END_2D
ENDDO
CALL iom_put( 'INTdtDIN', zw2d )
!
zw2d(:,:) = 0.
DO jk = 1, jpkm1
zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * po4r * tr(:,:,jk,jppo4,Krhs)
DO_2D( 0, 0, 0, 0 )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) * po4r * tr(ji,jj,jk,jppo4,Krhs)
END_2D
ENDDO
CALL iom_put( 'INTdtDIP', zw2d )
!
zw2d(:,:) = 0.
DO jk = 1, jpkm1
zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jpfer,Krhs)
DO_2D( 0, 0, 0, 0 )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) * tr(ji,jj,jk,jpfer,Krhs)
END_2D
ENDDO
CALL iom_put( 'INTdtFer', zw2d )
!
zw2d(:,:) = 0.
DO jk = 1, jpkm1
zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jpsil,Krhs)
DO_2D( 0, 0, 0, 0 )
zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) * tr(ji,jj,jk,jpsil,Krhs)
END_2D
ENDDO
CALL iom_put( 'INTdtSil', zw2d )
!
......@@ -522,8 +533,9 @@ CONTAINS
INTEGER, INTENT( in ) :: Kmm ! time level indices
REAL(wp) :: zrdenittot, zsdenittot, znitrpottot
CHARACTER(LEN=100) :: cltxt
INTEGER :: jk
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwork
INTEGER :: ji, jj, jk
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
REAL(wp), ALLOCATABLE, DIMENSION(:,: ) :: zw2d
!!----------------------------------------------------------------------
!
IF( kt == nittrc000 ) THEN
......@@ -542,82 +554,113 @@ CONTAINS
! Compute the budget of NO3
IF( iom_use( "pno3tot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
IF( ln_p4z ) THEN
zwork(:,:,:) = tr(:,:,:,jpno3,Kmm) + tr(:,:,:,jpnh4,Kmm) &
& + tr(:,:,:,jpphy,Kmm) + tr(:,:,:,jpdia,Kmm) &
& + tr(:,:,:,jppoc,Kmm) + tr(:,:,:,jpgoc,Kmm) + tr(:,:,:,jpdoc,Kmm) &
& + tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm)
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jpno3,Kmm) + tr(ji,jj,jk,jpnh4,Kmm) &
& + tr(ji,jj,jk,jpphy,Kmm) + tr(ji,jj,jk,jpdia,Kmm) &
& + tr(ji,jj,jk,jppoc,Kmm) + tr(ji,jj,jk,jpgoc,Kmm) + tr(ji,jj,jk,jpdoc,Kmm) &
& + tr(ji,jj,jk,jpzoo,Kmm) + tr(ji,jj,jk,jpmes,Kmm) ) * cvol(ji,jj,jk)
END_3D
ELSE
zwork(:,:,:) = tr(:,:,:,jpno3,Kmm) + tr(:,:,:,jpnh4,Kmm) + tr(:,:,:,jpnph,Kmm) &
& + tr(:,:,:,jpndi,Kmm) + tr(:,:,:,jpnpi,Kmm) &
& + tr(:,:,:,jppon,Kmm) + tr(:,:,:,jpgon,Kmm) + tr(:,:,:,jpdon,Kmm) &
& + ( tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) ) * no3rat3
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jpno3,Kmm) + tr(ji,jj,jk,jpnh4,Kmm) + tr(ji,jj,jk,jpnph,Kmm) &
& + tr(ji,jj,jk,jpndi,Kmm) + tr(ji,jj,jk,jpnpi,Kmm) &
& + tr(ji,jj,jk,jppon,Kmm) + tr(ji,jj,jk,jpgon,Kmm) + tr(ji,jj,jk,jpdon,Kmm) &
& + ( tr(ji,jj,jk,jpzoo,Kmm) + tr(ji,jj,jk,jpmes,Kmm) ) * no3rat3 ) * cvol(ji,jj,jk)
END_3D
ENDIF
!
no3budget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:) )
no3budget = glob_sum( 'p4zsms', zw3d(:,:,:) )
no3budget = no3budget / areatot
CALL iom_put( "pno3tot", no3budget )
DEALLOCATE( zw3d )
ENDIF
!
! Compute the budget of PO4
IF( iom_use( "ppo4tot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
IF( ln_p4z ) THEN
zwork(:,:,:) = tr(:,:,:,jppo4,Kmm) &
& + tr(:,:,:,jpphy,Kmm) + tr(:,:,:,jpdia,Kmm) &
& + tr(:,:,:,jppoc,Kmm) + tr(:,:,:,jpgoc,Kmm) + tr(:,:,:,jpdoc,Kmm) &
& + tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm)
ELSE
zwork(:,:,:) = tr(:,:,:,jppo4,Kmm) + tr(:,:,:,jppph,Kmm) &
& + tr(:,:,:,jppdi,Kmm) + tr(:,:,:,jpppi,Kmm) &
& + tr(:,:,:,jppop,Kmm) + tr(:,:,:,jpgop,Kmm) + tr(:,:,:,jpdop,Kmm) &
& + ( tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) ) * po4rat3
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jppo4,Kmm) &
& + tr(ji,jj,jk,jpphy,Kmm) + tr(ji,jj,jk,jpdia,Kmm) &
& + tr(ji,jj,jk,jppoc,Kmm) + tr(ji,jj,jk,jpgoc,Kmm) + tr(ji,jj,jk,jpdoc,Kmm) &
& + tr(ji,jj,jk,jpzoo,Kmm) + tr(ji,jj,jk,jpmes,Kmm) ) * cvol(ji,jj,jk)
END_3D
ELSE
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jppo4,Kmm) + tr(ji,jj,jk,jppph,Kmm) &
& + tr(ji,jj,jk,jppdi,Kmm) + tr(ji,jj,jk,jpppi,Kmm) &
& + tr(ji,jj,jk,jppop,Kmm) + tr(ji,jj,jk,jpgop,Kmm) + tr(ji,jj,jk,jpdop,Kmm) &
& + ( tr(ji,jj,jk,jpzoo,Kmm) + tr(ji,jj,jk,jpmes,Kmm) ) * po4rat3 ) * cvol(ji,jj,jk)
END_3D
ENDIF
!
po4budget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:) )
po4budget = glob_sum( 'p4zsms', zw3d(:,:,:) )
po4budget = po4budget / areatot
CALL iom_put( "ppo4tot", po4budget )
DEALLOCATE( zw3d )
ENDIF
!
! Compute the budget of SiO3
IF( iom_use( "psiltot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
zwork(:,:,:) = tr(:,:,:,jpsil,Kmm) + tr(:,:,:,jpgsi,Kmm) + tr(:,:,:,jpdsi,Kmm)
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jpsil,Kmm) + tr(ji,jj,jk,jpgsi,Kmm) + tr(ji,jj,jk,jpdsi,Kmm) ) * cvol(ji,jj,jk)
END_3D
!
silbudget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:) )
silbudget = glob_sum( 'p4zsms', zw3d(:,:,:) )
silbudget = silbudget / areatot
CALL iom_put( "psiltot", silbudget )
DEALLOCATE( zw3d )
ENDIF
!
IF( iom_use( "palktot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
zwork(:,:,:) = tr(:,:,:,jpno3,Kmm) * rno3 + tr(:,:,:,jptal,Kmm) + tr(:,:,:,jpcal,Kmm) * 2.
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jpno3,Kmm) * rno3 + tr(ji,jj,jk,jptal,Kmm) + tr(ji,jj,jk,jpcal,Kmm) * 2. ) * cvol(ji,jj,jk)
END_3D
!
alkbudget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:) ) !
alkbudget = glob_sum( 'p4zsms', zw3d(:,:,:) ) !
alkbudget = alkbudget / areatot
CALL iom_put( "palktot", alkbudget )
DEALLOCATE( zw3d )
ENDIF
!
! Compute the budget of Iron
IF( iom_use( "pfertot" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
zwork(:,:,:) = tr(:,:,:,jpfer,Kmm) + tr(:,:,:,jpnfe,Kmm) + tr(:,:,:,jpdfe,Kmm) &
& + tr(:,:,:,jpbfe,Kmm) + tr(:,:,:,jpsfe,Kmm) &
& + ( tr(:,:,:,jpzoo,Kmm) * feratz + tr(:,:,:,jpmes,Kmm) ) * feratm
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpk)
zw3d(ji,jj,jk) = ( tr(ji,jj,jk,jpfer,Kmm) + tr(ji,jj,jk,jpnfe,Kmm) + tr(ji,jj,jk,jpdfe,Kmm) &
& + tr(ji,jj,jk,jpbfe,Kmm) + tr(ji,jj,jk,jpsfe,Kmm) &
& + tr(ji,jj,jk,jpzoo,Kmm) * feratz + tr(ji,jj,jk,jpmes,Kmm) * feratm ) * cvol(ji,jj,jk)
END_3D
!
ferbudget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:) )
ferbudget = glob_sum( 'p4zsms', zw3d(:,:,:) )
ferbudget = ferbudget / areatot
CALL iom_put( "pfertot", ferbudget )
DEALLOCATE( zw3d )
ENDIF
!
! Global budget of N SMS : denitrification in the water column and in the sediment
! nitrogen fixation by the diazotrophs
! --------------------------------------------------------------------------------
IF( iom_use( "tnfix" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
znitrpottot = glob_sum ( 'p4zsms', nitrpot(:,:,:) * nitrfix * cvol(:,:,:) )
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = nitrpot(A2D(0),1:jpkm1) * nitrfix * cvol(A2D(0),1:jpkm1)
znitrpottot = glob_sum ( 'p4zsms', zw3d)
CALL iom_put( "tnfix" , znitrpottot * xfact3 ) ! Global nitrogen fixation molC/l to molN/m3
DEALLOCATE( zw3d )
ENDIF
!
IF( iom_use( "tdenit" ) .OR. ( ln_check_mass .AND. kt == nitend ) ) THEN
zrdenittot = glob_sum ( 'p4zsms', denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) )
zsdenittot = glob_sum ( 'p4zsms', sdenit(:,:) * e1e2t(:,:) * tmask(:,:,1) )
ALLOCATE( zw3d(A2D(0),jpk), zw2d(A2D(0)) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = denitr(A2D(0),1:jpkm1) * rdenit * xnegtr(A2D(0),1:jpkm1) * cvol(A2D(0),1:jpkm1)
zw2d(A2D(0)) = sdenit(A2D(0)) * e1e2t(A2D(0)) * tmask(A2D(0),1)
zrdenittot = glob_sum ( 'p4zsms', zw3d )
zsdenittot = glob_sum ( 'p4zsms', Zw2d )
CALL iom_put( "tdenit" , ( zrdenittot + zsdenittot ) * xfact3 ) ! Total denitrification molC/l to molN/m3
DEALLOCATE( zw3d, zw2d )
ENDIF
!
IF( ln_check_mass .AND. kt == nitend ) THEN ! Compute the budget of NO3, ALK, Si, Fer
......
......@@ -94,6 +94,9 @@ MODULE p5zlim
REAL(wp) :: xcoef1 = 0.00167 / 55.85
REAL(wp) :: xcoef2 = 1.21E-5 * 14. / 55.85 / 7.625 * 0.5 * 1.5
REAL(wp) :: xcoef3 = 1.15E-4 * 14. / 55.85 / 7.625 * 0.5
LOGICAL :: l_dia_nut_lim, l_dia_iron_lim, l_dia_fracal
LOGICAL :: l_dia_size_lim, l_dia_size_pro
!! * Substitutions
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
......@@ -136,15 +139,23 @@ CONTAINS
REAL(wp) :: zfvn, zfvp, zfvf, zsizen, zsizep, zsized, znanochl, zpicochl, zdiatchl
REAL(wp) :: zqfemn, zqfemp, zqfemd, zbactno3, zbactnh4, zbiron
REAL(wp) :: znutlimtot, zlimno3, zlimnh4, zlim1f, zsizetmp
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrassn, zrassp, zrassd
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p5z_lim')
IF( kt == nittrc000 ) THEN
l_dia_nut_lim = iom_use( "LNnut" ) .OR. iom_use( "LDnut" ) .OR. iom_use( "LPnut" )
l_dia_iron_lim = iom_use( "LNFe" ) .OR. iom_use( "LDFe" ) .OR. iom_use( "LPFe" )
l_dia_size_lim = iom_use( "SIZEN" ) .OR. iom_use( "SIZED" ) .OR. iom_use( "SIZEP" )
l_dia_size_pro = iom_use( "RASSN" ) .OR. iom_use( "RASSP" ) .OR. iom_use( "RASSP" )
l_dia_fracal = iom_use( "xfracal" )
ENDIF
!
zratchl = 6.0
sizena(:,:,:) = 0.0 ; sizepa(:,:,:) = 0.0 ; sizeda(:,:,:) = 0.0
!
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
! Computation of the Chl/C ratio of each phytoplankton group
! -------------------------------------------------------
z1_trnphy = 1. / ( tr(ji,jj,jk,jpphy,Kbb) + rtrn )
......@@ -407,7 +418,7 @@ CONTAINS
! nutrient uptake pool and assembly machinery. DNA is assumed to represent 1% of the dry mass of
! phytoplankton (see Daines et al., 2013).
! --------------------------------------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
! Size estimation of nanophytoplankton based on total biomass
! Assumes that larger biomass implies addition of larger cells
! ------------------------------------------------------------
......@@ -419,7 +430,6 @@ CONTAINS
! Computed from Inomura et al. (2020) using Pavlova Lutheri
zrpho = 11.55 * tr(ji,jj,jk,jpnch,Kbb) / ( tr(ji,jj,jk,jpphy,Kbb) * 12. + rtrn )
zrass = MAX(0.62/4., ( 1. - zrpho - zfuptk ) * xlimnpn(ji,jj,jk) )
zrassn(ji,jj,jk) = zrass
xqpnmin(ji,jj,jk) = ( 0.0 + 0.0078 + 0.62/4. * 0.0783 ) * 16.
xqpnmax(ji,jj,jk) = ( zrpho * 0.0089 + zrass * 0.0783 ) * 16.
xqpnmax(ji,jj,jk) = xqpnmax(ji,jj,jk) + (0.033 + 0.0078 ) * 16.
......@@ -437,7 +447,6 @@ CONTAINS
! Computed from Inomura et al. (2020) using a synechococcus
zrpho = 13.4 * tr(ji,jj,jk,jppch,Kbb) / ( tr(ji,jj,jk,jppic,Kbb) * 12. + rtrn )
zrass = MAX(0.4/4., ( 1. - zrpho - zfuptk ) * xlimnpp(ji,jj,jk) )
zrassp(ji,jj,jk) = zrass
xqppmin(ji,jj,jk) = ( (0.0 + 0.0078 ) + 0.4/4. * 0.0517 ) * 16.
xqppmax(ji,jj,jk) = ( zrpho * 0.0076 + zrass * 0.0517 ) * 16.
xqppmax(ji,jj,jk) = xqppmax(ji,jj,jk) + (0.033 + 0.0078 ) * 16
......@@ -454,7 +463,6 @@ CONTAINS
! Computed from Inomura et al. (2020) using a synechococcus
zrpho = 8.08 * tr(ji,jj,jk,jpdch,Kbb) / ( tr(ji,jj,jk,jpndi,Kbb) * 12. + rtrn )
zrass = MAX(0.66/4., ( 1. - zrpho - zfuptk ) * xlimnpd(ji,jj,jk) )
zrassd(ji,jj,jk)=zrass
xqpdmin(ji,jj,jk) = ( ( 0.0 + 0.0078 ) + 0.66/4. * 0.0783 ) * 16.
xqpdmax(ji,jj,jk) = ( zrpho * 0.0135 + zrass * 0.0783 ) * 16.
xqpdmax(ji,jj,jk) = xqpdmax(ji,jj,jk) + ( 0.0078 + 0.033 ) * 16.
......@@ -465,7 +473,7 @@ CONTAINS
! This is a purely adhoc formulation described in Aumont et al. (2015)
! This fraction depends on nutrient limitation, light, temperature
! --------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zlim1 = tr(ji,jj,jk,jpnh4,Kbb) / ( tr(ji,jj,jk,jpnh4,Kbb) + concnnh4 ) + tr(ji,jj,jk,jpno3,Kbb) &
& / ( tr(ji,jj,jk,jpno3,Kbb) + concnno3 ) * ( 1.0 - tr(ji,jj,jk,jpnh4,Kbb) &
& / ( tr(ji,jj,jk,jpnh4,Kbb) + concnnh4 ) )
......@@ -482,7 +490,7 @@ CONTAINS
xfracal(ji,jj,jk) = MAX( 0.02, MIN( 0.8 , xfracal(ji,jj,jk) ) )
END_3D
!
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
! denitrification factor computed from O2 levels
nitrfac(ji,jj,jk) = MAX( 0.e0, 0.4 * ( 6.e-6 - tr(ji,jj,jk,jpoxy,Kbb) ) &
& / ( oxymin + tr(ji,jj,jk,jpoxy,Kbb) ) )
......@@ -490,19 +498,70 @@ CONTAINS
END_3D
!
IF( lk_iomput .AND. knt == nrdttrc ) THEN ! save output diagnostics
CALL iom_put( "xfracal", xfracal(:,:,:) * tmask(:,:,:) ) ! euphotic layer deptht
CALL iom_put( "LNnut" , xlimphy(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term
CALL iom_put( "LPnut" , xlimpic(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term
CALL iom_put( "LDnut" , xlimdia(:,:,:) * tmask(:,:,:) ) ! Nutrient limitation term
CALL iom_put( "LNFe" , xlimnfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "LPFe" , xlimpfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "LDFe" , xlimdfe(:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "SIZEN" , sizen (:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "SIZEP" , sizep (:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "SIZED" , sized (:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "RASSN" , zrassn (:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "RASSP" , zrassp (:,:,:) * tmask(:,:,:) ) ! Iron limitation term
CALL iom_put( "RASSD" , zrassd (:,:,:) * tmask(:,:,:) ) ! Iron limitation term
!
IF( l_dia_fracal ) THEN ! fraction of calcifiers
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = xfracal(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "xfracal", zw3d)
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_nut_lim ) THEN ! Nutrient limitation term
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = xlimphy(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LNnut", zw3d)
zw3d(A2D(0),1:jpkm1) = xlimdia(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LDnut", zw3d)
zw3d(A2D(0),1:jpkm1) = xlimpic(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LPnut", zw3d)
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_iron_lim ) THEN ! Iron limitation term
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = xlimnfe(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LNFe", zw3d)
zw3d(A2D(0),1:jpkm1) = xlimdfe(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LDFe", zw3d)
zw3d(A2D(0),1:jpkm1) = xlimpfe(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LPFe", zw3d)
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_size_lim ) THEN ! Size limitation term
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = sizen(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "SIZEN", zw3d)
zw3d(A2D(0),1:jpkm1) = sized(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "SIZED", zw3d)
zw3d(A2D(0),1:jpkm1) = sizep(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "SIZEP", zw3d)
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_size_pro ) THEN ! Size of the protein machinery
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zfuptk = 0.2 + 0.12 / ( 3.0 * sizen(ji,jj,jk) + rtrn )
zrpho = 11.55 * tr(ji,jj,jk,jpnch,Kbb) / ( tr(ji,jj,jk,jpphy,Kbb) * 12. + rtrn )
zw3d(ji,jj,jk) = MAX(0.62/4., ( 1. - zrpho - zfuptk ) * xlimnpn(ji,jj,jk) ) * tmask(ji,jj,jk)
END_3D
CALL iom_put( "RASSN", zw3d)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zfuptk = 0.2 + 0.12 / ( 3.0 * sizep(ji,jj,jk) + rtrn )
zrpho = 11.55 * tr(ji,jj,jk,jppch,Kbb) / ( tr(ji,jj,jk,jppic,Kbb) * 12. + rtrn )
zw3d(ji,jj,jk) = MAX(0.62/4., ( 1. - zrpho - zfuptk ) * xlimnpp(ji,jj,jk) ) * tmask(ji,jj,jk)
END_3D
CALL iom_put( "RASSP", zw3d)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zfuptk = 0.2 + 0.12 / ( 3.0 * sized(ji,jj,jk) + rtrn )
zrpho = 11.55 * tr(ji,jj,jk,jpdch,Kbb) / ( tr(ji,jj,jk,jpndi,Kbb) * 12. + rtrn )
zw3d(ji,jj,jk) = MAX(0.62/4., ( 1. - zrpho - zfuptk ) * xlimnpd(ji,jj,jk) ) * tmask(ji,jj,jk)
END_3D
CALL iom_put( "RASSD", zw3d)
DEALLOCATE( zw3d )
ENDIF
!
ENDIF
!
IF( ln_timing ) CALL timing_stop('p5z_lim')
......@@ -635,24 +694,24 @@ CONTAINS
ierr(:) = 0
!
!* Biological arrays for phytoplankton growth
ALLOCATE( xpicono3(jpi,jpj,jpk), xpiconh4(jpi,jpj,jpk), &
& xpicopo4(jpi,jpj,jpk), xpicodop(jpi,jpj,jpk), &
& xnanodop(jpi,jpj,jpk), xdiatdop(jpi,jpj,jpk), &
& xpicofer(jpi,jpj,jpk), xlimpfe (jpi,jpj,jpk), &
& fvnuptk (jpi,jpj,jpk), fvduptk (jpi,jpj,jpk), &
& xlimphys(jpi,jpj,jpk), xlimdias(jpi,jpj,jpk), &
& xlimnpp (jpi,jpj,jpk), xlimnpn (jpi,jpj,jpk), &
& xlimnpd (jpi,jpj,jpk), &
& xlimpics(jpi,jpj,jpk), xqfuncfecp(jpi,jpj,jpk), &
& fvpuptk (jpi,jpj,jpk), xlimpic (jpi,jpj,jpk), STAT=ierr(1) )
ALLOCATE( xpicono3(A2D(0),jpk), xpiconh4(A2D(0),jpk), &
& xpicopo4(A2D(0),jpk), xpicodop(A2D(0),jpk), &
& xnanodop(A2D(0),jpk), xdiatdop(A2D(0),jpk), &
& xpicofer(A2D(0),jpk), xlimpfe (A2D(0),jpk), &
& fvnuptk (A2D(0),jpk), fvduptk (A2D(0),jpk), &
& xlimphys(A2D(0),jpk), xlimdias(A2D(0),jpk), &
& xlimnpp (A2D(0),jpk), xlimnpn (A2D(0),jpk), &
& xlimnpd (A2D(0),jpk), &
& xlimpics(A2D(0),jpk), xqfuncfecp(A2D(0),jpk), &
& fvpuptk (A2D(0),jpk), xlimpic (A2D(0),jpk), STAT=ierr(1) )
!
!* Minimum/maximum quotas of phytoplankton
ALLOCATE( xqnnmin (jpi,jpj,jpk), xqnnmax(jpi,jpj,jpk), &
& xqpnmin (jpi,jpj,jpk), xqpnmax(jpi,jpj,jpk), &
& xqnpmin (jpi,jpj,jpk), xqnpmax(jpi,jpj,jpk), &
& xqppmin (jpi,jpj,jpk), xqppmax(jpi,jpj,jpk), &
& xqndmin (jpi,jpj,jpk), xqndmax(jpi,jpj,jpk), &
& xqpdmin (jpi,jpj,jpk), xqpdmax(jpi,jpj,jpk), STAT=ierr(2) )
ALLOCATE( xqnnmin (A2D(0),jpk), xqnnmax(A2D(0),jpk), &
& xqpnmin (A2D(0),jpk), xqpnmax(A2D(0),jpk), &
& xqnpmin (A2D(0),jpk), xqnpmax(A2D(0),jpk), &
& xqppmin (A2D(0),jpk), xqppmax(A2D(0),jpk), &
& xqndmin (A2D(0),jpk), xqndmax(A2D(0),jpk), &
& xqpdmin (A2D(0),jpk), xqpdmax(A2D(0),jpk), STAT=ierr(2) )
!
p5z_lim_alloc = MAXVAL( ierr )
!
......
......@@ -58,6 +58,7 @@ MODULE p5zmeso
REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: depmig !: DVM of mesozooplankton : migration depth
INTEGER , ALLOCATABLE, SAVE, DIMENSION(:,:) :: kmig !: Vertical indice of the the migration depth
LOGICAL :: l_dia_fezoo2, l_dia_graz2, l_dia_lprodz2
!! * Substitutions
# include "do_loop_substitute.h90"
# include "domzgr_substitute.h90"
......@@ -103,21 +104,40 @@ CONTAINS
REAL(wp) :: zmigreltime, zrum, zcodel, zargu, zval, zmigthick
CHARACTER (len=25) :: charout
REAL(wp) :: zrfact2, zmetexcess, zsigma, zdiffdn
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo2
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrarem, zgraref, zgrapoc, zgrapof
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrarep, zgraren, zgrapon, zgrapop
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgradoc, zgradon, zgradop
REAL(wp), DIMENSION(A2D(0),jpk) :: zgrarem, zgraref, zgrapoc, zgrapof
REAL(wp), DIMENSION(A2D(0),jpk) :: zgrarep, zgraren, zgrapon, zgrapop
REAL(wp), DIMENSION(A2D(0),jpk) :: zgradoc, zgradon, zgradop
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zgramigrem, zgramigref, zgramigpoc, zgramigpof
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zgramigrep, zgramigren, zgramigpop, zgramigpon
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zgramigdoc, zgramigdop, zgramigdon
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zgrazing2, zfezoo2, zzligprod2, zw3d
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p5z_meso')
!
IF( kt == nittrc000 ) THEN
l_dia_graz2 = iom_use( "GRAZ2" )
l_dia_fezoo2 = iom_use( "FEZOO2" )
l_dia_lprodz2 = ln_ligand .AND. iom_use( "LPRODZ2" )
ENDIF
IF( l_dia_lprodz2 ) THEN
ALLOCATE( zzligprod2(A2D(0),jpk) )
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zzligprod2(ji,jj,jk) = tr(ji,jj,jk,jplgw,Krhs)
END_3D
ENDIF
IF( l_dia_fezoo2 ) THEN
ALLOCATE( zfezoo2(A2D(0),jpk) )
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zfezoo2(ji,jj,jk) = tr(ji,jj,jk,jpfer,Krhs)
END_3D
ENDIF
IF( l_dia_graz2 ) THEN
ALLOCATE( zgrazing2(A2D(0),jpk) )
ENDIF
! Initialization of local arrays
zgrazing(:,:,:) = 0._wp ; zfezoo2(:,:,:) = 0._wp
zgrarem (:,:,:) = 0._wp ; zgraren(:,:,:) = 0._wp
zgrarep (:,:,:) = 0._wp ; zgraref(:,:,:) = 0._wp
zgrapoc (:,:,:) = 0._wp ; zgrapon(:,:,:) = 0._wp
......@@ -136,7 +156,7 @@ CONTAINS
zmetexcess = 0.0
IF ( bmetexc2 ) zmetexcess = 1.0
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zcompam = MAX( ( tr(ji,jj,jk,jpmes,Kbb) - 1.e-9 ), 0.e0 )
zfact = xstep * tgfunc2(ji,jj,jk) * zcompam
......@@ -284,7 +304,7 @@ CONTAINS
& + zgrazffpp + zgrazffpg
! Total grazing ( grazing by microzoo is already computed in p5zmicro )
zgrazing(ji,jj,jk) = zgraztotc
IF( l_dia_graz2 ) zgrazing2(ji,jj,jk) = zgraztotc
! Stoichiometruc ratios of the food ingested by zooplanton
! --------------------------------------------------------
......@@ -427,7 +447,7 @@ CONTAINS
! This fraction is sumed over the euphotic zone and is removed from
! the fluxes driven by mesozooplankton in the euphotic zone.
! --------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zmigreltime = (1. - strn(ji,jj))
IF( gdept(ji,jj,jk,Kmm) <= heup(ji,jj) ) THEN
zmigthick = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) * ( 1. - zmigreltime )
......@@ -460,7 +480,7 @@ CONTAINS
! The inorganic and organic fluxes induced by migrating organisms are added at the
! the migration depth (corresponding indice is set by kmig)
! --------------------------------------------------------------------------------
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF( tmask(ji,jj,1) == 1. ) THEN
jkt = kmig(ji,jj)
zdep = 1. / e3t(ji,jj,jkt,Kmm)
......@@ -490,7 +510,7 @@ CONTAINS
! This only concerns the variables which are affected by DVM (inorganic
! nutrients, DOC agands, and particulate organic carbon).
! ---------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) + zgrarep(ji,jj,jk)
tr(ji,jj,jk,jpnh4,Krhs) = tr(ji,jj,jk,jpnh4,Krhs) + zgraren(ji,jj,jk)
tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) + zgradoc(ji,jj,jk)
......@@ -512,12 +532,47 @@ CONTAINS
tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + zgrapof(ji,jj,jk)
END_3D
!
! Write the output
IF( lk_iomput .AND. knt == nrdttrc ) THEN
CALL iom_put( "PCAL" , prodcal(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ) ! Calcite production
CALL iom_put( "GRAZ2" , zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) ) ! Total grazing of phyto by zoo
CALL iom_put( "FEZOO2", zfezoo2(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) )
IF( ln_ligand ) &
& CALL iom_put( "LPRODZ2", zgradoc(:,:,:) * ldocz * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) )
!
IF( iom_use ( "PCAL" ) ) THEN ! Calcite production
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zw3d(ji,jj,jk) = prodcal(ji,jj,jk) * 1.e+3 * rfact2r * tmask(ji,jj,jk)
END_3D
CALL iom_put( "PCAL", zw3d )
DEALLOCATE( zw3d )
ENDIF
!
IF( l_dia_graz2 ) THEN ! Total grazing of phyto by zooplankton
zgrazing2(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zgrazing2(ji,jj,jk) = zgrazing2(ji,jj,jk) * 1.e+3 * rfact2r * tmask(ji,jj,jk) ! conversion in mol/m2/s
END_3D
CALL iom_put( "GRAZ2" , zgrazing2 )
DEALLOCATE( zgrazing2 )
ENDIF
!
IF( l_dia_fezoo2 ) THEN
zfezoo2(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zfezoo2(ji,jj,jk) = ( tr(ji,jj,jk,jpfer,Krhs) - zfezoo2(ji,jj,jk) ) &
& * 1e9 * 1.e+3 * rfact2r * tmask(ji,jj,jk) ! conversion in nmol/m2/s
END_3D
CALL iom_put( "FEZOO2", zfezoo2 )
DEALLOCATE( zfezoo2 )
ENDIF
!
IF( l_dia_lprodz2 ) THEN
zzligprod2(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zzligprod2(ji,jj,jk) = ( tr(ji,jj,jk,jplgw,Krhs) - zzligprod2(ji,jj,jk) ) &
& * 1e9 * 1.e+3 * rfact2r * tmask(ji,jj,jk) ! conversion in nmol/m2/s
END_3D
CALL iom_put( "LPRODZ2", zzligprod2 )
DEALLOCATE( zzligprod2 )
ENDIF
!
ENDIF
!
IF(sn_cfctl%l_prttrc) THEN ! print mean trends (used for debugging)
......@@ -626,7 +681,7 @@ CONTAINS
! Compute the averaged values of oxygen, temperature over the domain
! 150m to 500 m depth.
! ------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
DO_3D( 0, 0, 0, 0, 1, jpk )
IF( tmask(ji,jj,jk) == 1.) THEN
IF( gdept(ji,jj,jk,Kmm) >= 150. .AND. gdept(ji,jj,jk,kmm) <= 500.) THEN
oxymoy(ji,jj) = oxymoy(ji,jj) + tr(ji,jj,jk,jpoxy,Kbb) * 1E6 * e3t(ji,jj,jk,Kmm)
......@@ -639,7 +694,7 @@ CONTAINS
! Compute the difference between surface values and the mean values in the mesopelagic
! domain
! ------------------------------------------------------------------------------------
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
z1dep = 1. / ( zdepmoy(ji,jj) + rtrn )
oxymoy(ji,jj) = tr(ji,jj,1,jpoxy,Kbb) * 1E6 - oxymoy(ji,jj) * z1dep
tempmoy(ji,jj) = ts(ji,jj,1,jp_tem,Kmm) - tempmoy(ji,jj) * z1dep
......@@ -648,7 +703,7 @@ CONTAINS
! Computation of the migration depth based on the parameterization of
! Bianchi et al. (2013)
! -------------------------------------------------------------------
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF( tmask(ji,jj,1) == 1. ) THEN
ztotchl = ( tr(ji,jj,1,jppch,Kbb) + tr(ji,jj,1,jpnch,Kbb) + tr(ji,jj,1,jpdch,Kbb) ) * 1E6
depmig(ji,jj) = 398. - 0.56 * oxymoy(ji,jj) -115. * log10(ztotchl) + 0.36 * hmld(ji,jj) -2.4 * tempmoy(ji,jj)
......@@ -657,7 +712,7 @@ CONTAINS
! Computation of the corresponding jk indice
! ------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
IF( depmig(ji,jj) >= gdepw(ji,jj,jk,Kmm) .AND. depmig(ji,jj) < gdepw(ji,jj,jk+1,Kmm) ) THEN
kmig(ji,jj) = jk
ENDIF
......@@ -669,7 +724,7 @@ CONTAINS
! to 0. Thus, to avoid that problem, the migration depth is adjusted so
! that it falls above the OMZ
! -----------------------------------------------------------------------
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF( tr(ji,jj,kmig(ji,jj),jpoxy,Kbb) < 5E-6 ) THEN
DO jk = kmig(ji,jj),1,-1
IF( tr(ji,jj,jk,jpoxy,Kbb) >= 5E-6 .AND. tr(ji,jj,jk+1,jpoxy,Kbb) < 5E-6) THEN
......@@ -689,7 +744,7 @@ CONTAINS
!! *** ROUTINE p5z_meso_alloc ***
!!----------------------------------------------------------------------
!
ALLOCATE( depmig(jpi,jpj), kmig(jpi,jpj), STAT= p5z_meso_alloc )
ALLOCATE( depmig(A2D(0)), kmig(A2D(0)), STAT= p5z_meso_alloc )
!
IF( p5z_meso_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p5z_meso_alloc : failed to allocate arrays.' )
!
......
......@@ -53,6 +53,8 @@ MODULE p5zmicro
REAL(wp), PUBLIC :: xsigmadel !: Maximum additional width of the grazing window at low food density
LOGICAL, PUBLIC :: bmetexc !: Use of excess carbon for respiration
LOGICAL :: l_dia_fezoo, l_dia_graz1, l_dia_lprodz
!! * Substitutions
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
......@@ -88,18 +90,39 @@ CONTAINS
REAL(wp) :: zgraznc, zgraznn, zgraznp, zgrazpoc, zgrazpon, zgrazpop, zgrazpof
REAL(wp) :: zgrazdc, zgrazdn, zgrazdp, zgrazdf, zgraznf, zgrazz
REAL(wp) :: zgrazpc, zgrazpn, zgrazpp, zgrazpf, zbeta, zrfact2, zmetexcess
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo, zzligprod
REAL(wp) :: zsigma, zdiffdn, zdiffpn, zdiffdp, zproport, zproport2
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zgrazing, zfezoo, zzligprod, zw3d
CHARACTER (len=25) :: charout
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p5z_micro')
!
IF( kt == nittrc000 ) THEN
l_dia_graz1 = iom_use( "GRAZ1" )
l_dia_fezoo = iom_use( "FEZOO" )
l_dia_lprodz = ln_ligand .AND. iom_use( "LPRODZ" )
ENDIF
IF( l_dia_lprodz ) THEN
ALLOCATE( zzligprod(A2D(0),jpk) )
DO_3D( 0, 0, 0, 0, 1, jpk)
zzligprod(ji,jj,jk) = tr(ji,jj,jk,jplgw,Krhs)
END_3D
ENDIF
IF( l_dia_fezoo ) THEN
ALLOCATE( zfezoo(A2D(0),jpk) )
DO_3D( 0, 0, 0, 0, 1, jpk)
zfezoo(ji,jj,jk) = tr(ji,jj,jk,jpfer,Krhs)
END_3D
ENDIF
IF( l_dia_graz1 ) THEN
ALLOCATE( zgrazing(A2D(0),jpk) )
ENDIF
! Use of excess carbon for metabolism
zmetexcess = 0.0
IF ( bmetexc ) zmetexcess = 1.0
!
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zcompaz = MAX( ( tr(ji,jj,jk,jpzoo,Kbb) - 1.e-9 ), 0.e0 )
zfact = xstep * tgfunc2(ji,jj,jk) * zcompaz
! Proportion of nano and diatoms that are within the size range
......@@ -207,14 +230,13 @@ CONTAINS
zgrazdf = zgrazdc * tr(ji,jj,jk,jpdfe,Kbb) / (tr(ji,jj,jk,jpdia,Kbb) + rtrn)
!
! Total ingestion rates in C, P, Fe, N
zgraztotc = zgraznc + zgrazpoc + zgrazdc + zgrazz + zgrazpc
zgraztotc = zgraznc + zgrazpoc + zgrazdc + zgrazz + zgrazpc ! Grazing by microzooplankton
IF( l_dia_graz1 ) zgrazing(ji,jj,jk) = zgraztotc
zgraztotn = zgraznn + zgrazpn + zgrazpon + zgrazdn + zgrazz * no3rat3
zgraztotp = zgraznp + zgrazpp + zgrazpop + zgrazdp + zgrazz * po4rat3
zgraztotf = zgraznf + zgrazpf + zgrazpof + zgrazdf + zgrazz * feratz
!
! Grazing by microzooplankton
zgrazing(ji,jj,jk) = zgraztotc
! Stoichiometruc ratios of the food ingested by zooplanton
! --------------------------------------------------------
zgrasratf = (zgraztotf + rtrn) / ( zgraztotc + rtrn )
......@@ -298,14 +320,12 @@ CONTAINS
!
IF( ln_ligand ) THEN
tr(ji,jj,jk,jplgw,Krhs) = tr(ji,jj,jk,jplgw,Krhs) + zgradoc * ldocz
zzligprod(ji,jj,jk) = zgradoc * ldocz
ENDIF
!
tr(ji,jj,jk,jpdon,Krhs) = tr(ji,jj,jk,jpdon,Krhs) + zgradon
tr(ji,jj,jk,jpdop,Krhs) = tr(ji,jj,jk,jpdop,Krhs) + zgradop
tr(ji,jj,jk,jpoxy,Krhs) = tr(ji,jj,jk,jpoxy,Krhs) - o2ut * zgrarem
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + zgraref
zfezoo(ji,jj,jk) = zgraref
tr(ji,jj,jk,jpzoo,Krhs) = tr(ji,jj,jk,jpzoo,Krhs) + zepsherv * zgraztotc - zrespirc - ztortz - zgrazz
tr(ji,jj,jk,jpphy,Krhs) = tr(ji,jj,jk,jpphy,Krhs) - zgraznc
tr(ji,jj,jk,jpnph,Krhs) = tr(ji,jj,jk,jpnph,Krhs) - zgraznn
......@@ -342,15 +362,36 @@ CONTAINS
END_3D
!
IF( lk_iomput .AND. knt == nrdttrc ) THEN
IF( iom_use("GRAZ1") ) THEN ! Total grazing of phyto by zooplankton
zgrazing(:,:,jpk) = 0._wp ; CALL iom_put( "GRAZ1" , zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:) )
ENDIF
IF( iom_use("FEZOO") ) THEN
zfezoo (:,:,jpk) = 0._wp ; CALL iom_put( "FEZOO" , zfezoo(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) )
ENDIF
IF( ln_ligand ) THEN
zzligprod(:,:,jpk) = 0._wp ; CALL iom_put( "LPRODZ", zzligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:))
ENDIF
!
IF( l_dia_graz1 ) THEN ! Total grazing of phyto by zooplankton
zgrazing(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zgrazing(ji,jj,jk) = zgrazing(ji,jj,jk) * 1.e+3 * rfact2r * tmask(ji,jj,jk) ! conversion in mol/m2/s
END_3D
CALL iom_put( "GRAZ1" , zgrazing )
DEALLOCATE( zgrazing )
ENDIF
!
IF( l_dia_fezoo ) THEN
zfezoo(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zfezoo(ji,jj,jk) = ( tr(ji,jj,jk,jpfer,Krhs) - zfezoo(ji,jj,jk) ) &
& * 1e9 * 1.e+3 * rfact2r * tmask(ji,jj,jk) ! conversion in nmol/m2/s
END_3D
CALL iom_put( "FEZOO", zfezoo )
DEALLOCATE( zfezoo )
ENDIF
!
IF( l_dia_lprodz ) THEN
zzligprod(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zzligprod(ji,jj,jk) = ( tr(ji,jj,jk,jplgw,Krhs) - zzligprod(ji,jj,jk) ) &
& * 1e9 * 1.e+3 * rfact2r * tmask(ji,jj,jk) ! conversion in nmol/m2/s
END_3D
CALL iom_put( "LPRODZ", zzligprod )
DEALLOCATE( zzligprod )
ENDIF
!
ENDIF
!
IF(sn_cfctl%l_prttrc) THEN ! print mean trends (used for debugging)
......
......@@ -80,7 +80,7 @@ CONTAINS
IF( ln_timing ) CALL timing_start('p5z_mort_nano')
!
prodcal(:,:,:) = 0. !: calcite production variable set to zero
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zcompaph = MAX( ( tr(ji,jj,jk,jpphy,Kbb) - 1e-9 ), 0.e0 )
! Quadratic mortality of nano due to aggregation during
......@@ -151,7 +151,7 @@ CONTAINS
!
IF( ln_timing ) CALL timing_start('p5z_mort_pico')
!
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zcompaph = MAX( ( tr(ji,jj,jk,jppic,Kbb) - 1e-9 ), 0.e0 )
! Quadratic mortality of pico due to aggregation during
......@@ -215,7 +215,7 @@ CONTAINS
IF( ln_timing ) CALL timing_start('p5z_mort_diat')
!
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zcompadi = MAX( ( tr(ji,jj,jk,jpdia,Kbb) - 1E-9), 0. )
......
......@@ -26,7 +26,6 @@ MODULE p5zprod
PUBLIC p5z_prod ! called in p5zbio.F90
PUBLIC p5z_prod_init ! called in trcsms_pisces.F90
PUBLIC p5z_prod_alloc
!! * Shared module variables
REAL(wp), PUBLIC :: pislopen !: P-I slope of nanophytoplankton
......@@ -43,12 +42,16 @@ MODULE p5zprod
REAL(wp), PUBLIC :: chlcmin !: Minimum Chl/C ratio of phytoplankton
REAL(wp), PUBLIC :: grosip !: Mean Si/C ratio of diatoms
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: zdaylen ! day length
REAL(wp) :: r1_rday !: 1 / rday
REAL(wp) :: texcretn !: 1 - excretn
REAL(wp) :: texcretp !: 1 - excretp
REAL(wp) :: texcretd !: 1 - excretd
REAL(wp) :: r1_rday !: 1 / rday
REAL(wp) :: texcretn !: 1 - excretn
REAL(wp) :: texcretp !: 1 - excretp
REAL(wp) :: texcretd !: 1 - excretd
REAL(wp) :: xq10_n !: q10 coef for nano = 1. + xpsino3 * qnnmax
REAL(wp) :: xq10_p !: q10 coef for pico = 1. + xpsino3 * qnpmax
REAL(wp) :: xq10_d !: q10 coef for diat = 1. + xpsino3 * qndmax
LOGICAL :: l_dia_ppphy, l_dia_ppnew, l_dia_ppbfe, l_dia_ppbsi
LOGICAL :: l_dia_mu, l_dia_light, l_dia_lprod
!! * Substitutions
# include "do_loop_substitute.h90"
......@@ -76,51 +79,60 @@ CONTAINS
INTEGER :: ji, jj, jk
REAL(wp) :: zsilfac, znanotot, zpicotot, zdiattot, zconctemp, zconctemp2
REAL(wp) :: zration, zratiop, zratiof, zmax, ztn, zadap
REAL(wp) :: zpronmax, zpropmax, zprofmax, zratio
REAL(wp) :: zprofmax, zratio
REAL(wp) :: zpronewn, zpronewp, zpronewd
REAL(wp) :: zproregn, zproregp, zproregd
REAL(wp) :: zpropo4n, zpropo4p, zpropo4d
REAL(wp) :: zprodopn, zprodopp, zprodopd
REAL(wp) :: zlim, zsilfac2, zsiborn, zprod, zprontot, zproptot, zprodtot
REAL(wp) :: zproddoc, zproddon, zproddop, zprodsil, zprodfer, zprodlig, zresptot
REAL(wp) :: zprnutmax, zprochln, zprochld, zprochlp
REAL(wp) :: zpislopen, zpislopep, zpisloped
REAL(wp) :: zval, zpptot, zpnewtot, zpregtot
REAL(wp) :: zmxl_chl, zmxl_fac
REAL(wp) :: zqfpmax, zqfnmax, zqfdmax
REAL(wp) :: zfact, zrfact2, zmaxsi, zratiosi, zsizetmp, zlimfac, zsilim
CHARACTER (len=25) :: charout
REAL(wp), DIMENSION(jpi,jpj ) :: zmixnano, zmixpico, zmixdiat
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadp, zpislopeadd
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprnut, zprmaxp, zprmaxn, zprmaxd
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprbio, zprpic, zprdia, zysopt
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprchln, zprchlp, zprchld
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprorcan, zprorcap, zprorcad
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprofed, zprofep, zprofen
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronewn, zpronewp, zpronewd
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zproregn, zproregp, zproregd
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpropo4n, zpropo4p, zpropo4d
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprodopn, zprodopp, zprodopd
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrespn, zrespp, zrespd
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zpligprod
REAL(wp), DIMENSION(A2D(0),jpk) :: zprorcan, zprorcap, zprorcad
REAL(wp), DIMENSION(A2D(0),jpk) :: zpislopeadn, zpislopeadp, zpislopeadd
REAL(wp), DIMENSION(A2D(0),jpk) :: zprnut, zprbio, zprpic, zprdia, zysopt
REAL(wp), DIMENSION(A2D(0),jpk) :: zprchln, zprchlp, zprchld
REAL(wp), DIMENSION(A2D(0),jpk) :: zprofed, zprofep, zprofen
REAL(wp), DIMENSION(A2D(0),jpk) :: zpronmaxn, zpronmaxp,zpronmaxd
REAL(wp), DIMENSION(A2D(0),jpk) :: zpropmaxn, zpropmaxp,zpropmaxd
! REAL(wp), DIMENSION(A2D(0),jpk) :: zrespn, zrespp, zrespd
REAL(wp), DIMENSION(A2D(0),jpk) :: zprmaxn, zprmaxd, zprmaxp, zmxl
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p5z_prod')
!
IF( kt == nittrc000 ) THEN
l_dia_ppphy = iom_use( "PPPHYN" ) .OR. iom_use( "PPPHYD" ) .OR. iom_use( "PPPHYP" ) .OR. iom_use( "TPP" )
l_dia_ppnew = iom_use( "PPNEWN" ) .OR. iom_use( "PPNEWD" ) .OR. iom_use( "PPNEWP" ) .OR. iom_use( "TPNEW")
l_dia_ppbfe = iom_use( "PFeN" ) .OR. iom_use( "PFeD" ) .OR. iom_use( "PFeP" ) .OR. iom_use( "TPBFE")
l_dia_ppbsi = iom_use( "PBSi" )
l_dia_mu = iom_use( "Mumax" ) .OR. iom_use( "MuN" ) .OR. iom_use( "MuD") .OR. iom_use( "MuP")
l_dia_light = iom_use( "LNlight") .OR. iom_use( "LDlight") .OR. iom_use( "LPlight")
l_dia_lprod = ln_ligand .AND. ( iom_use( "LPRODP") .OR. iom_use( "LDETP") )
ENDIF
! Initialize the local arrays
zprorcan(:,:,:) = 0._wp ; zprorcap(:,:,:) = 0._wp ; zprorcad(:,:,:) = 0._wp
zprofed (:,:,:) = 0._wp ; zprofep (:,:,:) = 0._wp ; zprofen (:,:,:) = 0._wp
zpronewn(:,:,:) = 0._wp ; zpronewp(:,:,:) = 0._wp ; zpronewd(:,:,:) = 0._wp
zproregn(:,:,:) = 0._wp ; zproregp(:,:,:) = 0._wp ; zproregd(:,:,:) = 0._wp
zpropo4n(:,:,:) = 0._wp ; zpropo4p(:,:,:) = 0._wp ; zpropo4d(:,:,:) = 0._wp
zprdia (:,:,:) = 0._wp ; zprpic (:,:,:) = 0._wp ; zprbio (:,:,:) = 0._wp
zprodopn(:,:,:) = 0._wp ; zprodopp(:,:,:) = 0._wp ; zprodopd(:,:,:) = 0._wp
zysopt (:,:,:) = 0._wp
zrespn (:,:,:) = 0._wp ; zrespp (:,:,:) = 0._wp ; zrespd (:,:,:) = 0._wp
zmxl_fac(:,:,:) = 0._wp ; zmxl_chl(:,:,:) = 0._wp
zprorcan (:,:,:) = 0._wp ; zprorcap (:,:,:) = 0._wp ; zprorcad (:,:,:) = 0._wp
zprofen (:,:,:) = 0._wp ; zprofep (:,:,:) = 0._wp ; zprofed (:,:,:) = 0._wp
zprbio (:,:,:) = 0._wp ; zprpic (:,:,:) = 0._wp ; zprdia (:,:,:) = 0._wp
zpronmaxn(:,:,:) = 0._wp ; zpronmaxp(:,:,:) = 0._wp ; zpronmaxd(:,:,:) = 0._wp
zpropmaxn(:,:,:) = 0._wp ; zpropmaxp(:,:,:) = 0._wp ; zpropmaxd(:,:,:) = 0._wp
zmxl (:,:,:) = 0._wp ; zysopt (:,:,:) = 0._wp
! zrespn (:,:,:) = 0._wp ; zrespp (:,:,:) = 0._wp ; zrespd (:,:,:) = 0._wp
! Computation of the optimal production rates and nutrient uptake
! rates. Based on a Q10 description of the thermal dependency.
zprnut (:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)
zprmaxn(:,:,:) = 0.8_wp * (1. + xpsino3 * qnnmax ) * r1_rday * tgfunc(:,:,:)
zprmaxd(:,:,:) = 0.8_wp * (1. + xpsino3 * qndmax ) * r1_rday * tgfunc(:,:,:)
zprmaxp(:,:,:) = 0.6_wp * (1. + xpsino3 * qnpmax ) * r1_rday * tgfunc(:,:,:)
zprnut (:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)
zprmaxn(:,:,:) = 0.8_wp * xq10_n * r1_rday * tgfunc(:,:,:)
zprmaxd(:,:,:) = 0.8_wp * xq10_d * r1_rday * tgfunc(:,:,:)
zprmaxp(:,:,:) = 0.6_wp * xq10_p * r1_rday * tgfunc(:,:,:)
! Impact of the day duration and light intermittency on phytoplankton growth
! Intermittency is supposed to have a similar effect on production as
......@@ -131,42 +143,39 @@ CONTAINS
! -------------------------------------------------------------------------
IF ( ln_p4z_dcyc ) THEN ! Diurnal cycle in PISCES
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN
zval = MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn ))
ENDIF
zmxl_chl(ji,jj,jk) = zval / 24.
zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval )
zmxl(ji,jj,jk) = zval
ENDIF
END_3D
ELSE ! No diurnal cycle in PISCES
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
zval = MAX( 1., strn(ji,jj) )
IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN
zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn ))
ENDIF
zmxl_chl(ji,jj,jk) = zval / 24.
zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval )
zmxl(ji,jj,jk) = zval
ENDIF
END_3D
ENDIF
zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:)
zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:)
zprpic(:,:,:) = zprmaxp(:,:,:) * zmxl_fac(:,:,:)
! Maximum light intensity
zdaylen(:,:) = MAX(1., strn(:,:)) / 24.
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
zmxl_fac = 1.0 - EXP( -0.26 * zmxl(ji,jj,jk) )
zprbio(ji,jj,jk) = zprmaxn(ji,jj,jk) * zmxl_fac
zprdia(ji,jj,jk) = zprmaxd(ji,jj,jk) * zmxl_fac
zprpic(ji,jj,jk) = zprmaxp(ji,jj,jk) * zmxl_fac
ENDIF
END_3D
! Computation of the P-I slope for nanos, picos and diatoms
! The formulation proposed by Geider et al. (1997) has been used.
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
! Computation of the P-I slope for nanos and diatoms
ztn = MAX( 0., ts(ji,jj,jk,jp_tem,Kmm) - 15. )
......@@ -191,26 +200,31 @@ CONTAINS
! Computation of production function for Carbon
! Actual light levels are used here
! ---------------------------------------------
zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) / zmxl_fac(ji,jj,jk) ) )
zprpic(ji,jj,jk) = zprpic(ji,jj,jk) * ( 1.- EXP( -zpislopep * epico(ji,jj,jk) / zmxl_fac(ji,jj,jk) ) )
zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) / zmxl_fac(ji,jj,jk) ) )
zmxl_fac = 1.0 - EXP( -0.26 * zmxl(ji,jj,jk) )
zmxl_chl = zmxl(ji,jj,jk) / 24.
!
zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) / zmxl_fac ) )
zprpic(ji,jj,jk) = zprpic(ji,jj,jk) * ( 1.- EXP( -zpislopep * epico(ji,jj,jk) / zmxl_fac ) )
zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) / zmxl_fac ) )
! Computation of production function for Chlorophyll
! Mean light level in the mixed layer (when appropriate)
! is used here (acclimation is in general slower than
! the characteristic time scales of vertical mixing)
! ------------------------------------------------------
zpislopen = zpislopen * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
zpislopep = zpislopep * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
zpislopen = zpislopen * zmxl_fac / ( zmxl_chl + rtrn )
zpisloped = zpisloped * zmxl_fac / ( zmxl_chl + rtrn )
zpislopep = zpislopep * zmxl_fac / ( zmxl_chl + rtrn )
!
zprchln(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) )
zprchlp(ji,jj,jk) = zprmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep * epicom(ji,jj,jk) ) )
zprchld(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) )
ENDIF
END_3D
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
!
! Si/C of diatoms
! ------------------------
! Si/C increases with iron stress and silicate availability (zsilfac)
......@@ -242,7 +256,7 @@ CONTAINS
! Sea-ice effect on production
! No production is assumed below sea ice
! --------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
zprpic(ji,jj,jk) = zprpic(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
......@@ -256,7 +270,7 @@ CONTAINS
! quota, uptake is downregulated according to a sigmoidal function
! (power 2), as proposed by Flynn (2003)
! ---------------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
! production terms for nanophyto.
zprorcan(ji,jj,jk) = zprbio(ji,jj,jk) * xlimphy(ji,jj,jk) * tr(ji,jj,jk,jpphy,Kbb) * rfact2
......@@ -278,17 +292,13 @@ CONTAINS
! Uptake of nitrogen
zratio = 1.0 - MIN( 1., zration / (xqnnmax(ji,jj,jk) + rtrn) )
zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpnmin(ji,jj,jk) ) &
zpronmaxn(ji,jj,jk) = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpnmin(ji,jj,jk) ) &
& / ( xqpnmax(ji,jj,jk) - xqpnmin(ji,jj,jk) + rtrn ), xlimnfe(ji,jj,jk) ) )
zpronmax = zpronmax * xqnnmin(ji,jj,jk) / qnnmin
zpronewn(ji,jj,jk) = zpronmax * xnanono3(ji,jj,jk)
zproregn(ji,jj,jk) = zpronmax * xnanonh4(ji,jj,jk)
zpronmaxn(ji,jj,jk) = zpronmaxn(ji,jj,jk) * xqnnmin(ji,jj,jk) / qnnmin
! Uptake of phosphorus and DOP
zratio = 1.0 - MIN( 1., zratiop / (xqpnmax(ji,jj,jk) + rtrn) )
zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
zpropmax = zprnutmax * zmax * xlimnfe(ji,jj,jk)
zpropo4n(ji,jj,jk) = zpropmax * xnanopo4(ji,jj,jk)
zprodopn(ji,jj,jk) = zpropmax * xnanodop(ji,jj,jk)
zpropmaxn(ji,jj,jk) = zprnutmax * zmax * xlimnfe(ji,jj,jk)
! Uptake of iron
zqfnmax = xqfuncfecn(ji,jj,jk) + ( qfnmax - xqfuncfecn(ji,jj,jk) ) * xlimnpn(ji,jj,jk)
zratio = 1.0 - MIN( 1., zratiof / zqfnmax )
......@@ -307,8 +317,9 @@ CONTAINS
! quota, uptake is downregulated according to a sigmoidal function
! (power 2), as proposed by Flynn (2003)
! ---------------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
!
! production terms for picophyto.
zprorcap(ji,jj,jk) = zprpic(ji,jj,jk) * xlimpic(ji,jj,jk) * tr(ji,jj,jk,jppic,Kbb) * rfact2
! Size computation
......@@ -328,17 +339,13 @@ CONTAINS
! Uptake of nitrogen
zratio = 1.0 - MIN( 1., zration / (xqnpmax(ji,jj,jk) + rtrn) )
zmax = MAX(0., MIN(1., zratio**2/ (0.05**2 + zratio**2) ) )
zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqppmin(ji,jj,jk) ) &
zpronmaxp(ji,jj,jk) = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqppmin(ji,jj,jk) ) &
& / ( xqppmax(ji,jj,jk) - xqppmin(ji,jj,jk) + rtrn ), xlimpfe(ji,jj,jk) ) )
zpronmax = zpronmax * xqnpmin(ji,jj,jk) / qnnmin
zpronewp(ji,jj,jk) = zpronmax * xpicono3(ji,jj,jk)
zproregp(ji,jj,jk) = zpronmax * xpiconh4(ji,jj,jk)
zpronmaxp(ji,jj,jk) = zpronmaxp(ji,jj,jk) * xqnpmin(ji,jj,jk) / qnnmin
! Uptake of phosphorus
zratio = 1.0 - MIN( 1., zratiop / (xqppmax(ji,jj,jk) + rtrn) )
zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
zpropmax = zprnutmax * zmax * xlimpfe(ji,jj,jk)
zpropo4p(ji,jj,jk) = zpropmax * xpicopo4(ji,jj,jk)
zprodopp(ji,jj,jk) = zpropmax * xpicodop(ji,jj,jk)
zpropmaxp(ji,jj,jk) = zprnutmax * zmax * xlimpfe(ji,jj,jk)
! Uptake of iron
zqfpmax = xqfuncfecp(ji,jj,jk) + ( qfpmax - xqfuncfecp(ji,jj,jk) ) * xlimnpp(ji,jj,jk)
zratio = 1.0 - MIN( 1., zratiof / zqfpmax )
......@@ -357,8 +364,9 @@ CONTAINS
! quota, uptake is downregulated according to a sigmoidal function
! (power 2), as proposed by Flynn (2003)
! ---------------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
!
! production terms for diatomees
zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * tr(ji,jj,jk,jpdia,Kbb) * rfact2
! Size computation
......@@ -378,17 +386,13 @@ CONTAINS
! Uptake of nitrogen
zratio = 1.0 - MIN( 1., zration / (xqndmax(ji,jj,jk) + rtrn) )
zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpdmin(ji,jj,jk) ) &
zpronmaxd(ji,jj,jk) = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpdmin(ji,jj,jk) ) &
& / ( xqpdmax(ji,jj,jk) - xqpdmin(ji,jj,jk) + rtrn ), xlimdfe(ji,jj,jk) ) )
zpronmax = zpronmax * xqndmin(ji,jj,jk) / qnnmin
zpronewd(ji,jj,jk) = zpronmax * xdiatno3(ji,jj,jk)
zproregd(ji,jj,jk) = zpronmax * xdiatnh4(ji,jj,jk)
zpronmaxd(ji,jj,jk) = zpronmaxd(ji,jj,jk) * xqndmin(ji,jj,jk) / qnnmin
! Uptake of phosphorus
zratio = 1.0 - MIN( 1., zratiop / (xqpdmax(ji,jj,jk) + rtrn) )
zmax = MAX(0., MIN(1., zratio**2/ (0.05**2 + zratio**2) ) )
zpropmax = zprnutmax * zmax * xlimdfe(ji,jj,jk)
zpropo4d(ji,jj,jk) = zpropmax * xdiatpo4(ji,jj,jk)
zprodopd(ji,jj,jk) = zpropmax * xdiatdop(ji,jj,jk)
zpropmaxd(ji,jj,jk) = zprnutmax * zmax * xlimdfe(ji,jj,jk)
! Uptake of iron
zqfdmax = xqfuncfecd(ji,jj,jk) + ( qfdmax - xqfuncfecd(ji,jj,jk) ) * xlimnpd(ji,jj,jk)
zratio = 1.0 - MIN( 1., zratiof / zqfdmax )
......@@ -403,21 +407,28 @@ CONTAINS
! Production of Chlorophyll. The formulation proposed by Geider et al.
! is adopted here.
! --------------------------------------------------------------------
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
! production terms for nanophyto. ( chlorophyll )
znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
zprod = rday * (zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk)
zmxl_chl = zmxl(ji,jj,jk) / 24.
! production terms for nanophyto. ( chlorophyll )
zpronewn = zpronmaxn(ji,jj,jk) * xnanono3(ji,jj,jk)
zproregn = zpronmaxn(ji,jj,jk) * xnanonh4(ji,jj,jk)
znanotot = enanom(ji,jj,jk) / ( zmxl_chl + rtrn )
zprod = rday * (zpronewn + zproregn) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk)
zprochln = thetannm * zprod / ( zpislopeadn(ji,jj,jk) * znanotot + rtrn )
zprochln = MAX(zprochln, chlcmin * 12. * zprorcan (ji,jj,jk) )
! production terms for picophyto. ( chlorophyll )
zpicotot = epicom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
zprod = rday * (zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk)
! production terms for picophyto. ( chlorophyll )
zpronewp = zpronmaxp(ji,jj,jk) * xpicono3(ji,jj,jk)
zproregp = zpronmaxp(ji,jj,jk) * xpiconh4(ji,jj,jk)
zpicotot = epicom(ji,jj,jk) / ( zmxl_chl + rtrn )
zprod = rday * (zpronewp + zproregp) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk)
zprochlp = thetanpm * zprod / ( zpislopeadp(ji,jj,jk) * zpicotot + rtrn )
zprochlp = MAX(zprochlp, chlcmin * 12. * zprorcap(ji,jj,jk) )
! production terms for diatoms ( chlorophyll )
zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
zprod = rday * (zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk)
zpronewd = zpronmaxd(ji,jj,jk) * xdiatno3(ji,jj,jk)
zproregd = zpronmaxd(ji,jj,jk) * xdiatnh4(ji,jj,jk)
zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl + rtrn )
zprod = rday * (zpronewd + zproregd) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk)
zprochld = thetandm * zprod / ( zpislopeadd(ji,jj,jk) * zdiattot + rtrn )
zprochld = MAX(zprochld, chlcmin * 12. * zprorcad(ji,jj,jk) )
! Update the arrays TRA which contain the Chla sources and sinks
......@@ -428,83 +439,104 @@ CONTAINS
END_3D
! Update the arrays TRA which contain the biological sources and sinks
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
zpptot = zpropo4n(ji,jj,jk) + zpropo4d(ji,jj,jk) + zpropo4p(ji,jj,jk)
zpnewtot = zpronewn(ji,jj,jk) + zpronewd(ji,jj,jk) + zpronewp(ji,jj,jk)
zpregtot = zproregn(ji,jj,jk) + zproregd(ji,jj,jk) + zproregp(ji,jj,jk)
zprontot = zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)
zproptot = zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)
zprodtot = zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)
!
zproddoc = excretd * zprorcad(ji,jj,jk) &
& + excretn * zprorcan(ji,jj,jk) &
& + excretp * zprorcap(ji,jj,jk)
!
zproddop = excretd * zpropo4d(ji,jj,jk) - texcretd * zprodopd(ji,jj,jk) &
& + excretn * zpropo4n(ji,jj,jk) - texcretn * zprodopn(ji,jj,jk) &
& + excretp * zpropo4p(ji,jj,jk) - texcretp * zprodopp(ji,jj,jk)
zproddon = excretd * zprodtot + excretn * zprontot + excretp * zproptot
zprodfer = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk)
zresptot = zrespn(ji,jj,jk) + zrespp(ji,jj,jk) + zrespd(ji,jj,jk)
!
tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) - zpptot
tr(ji,jj,jk,jpno3,Krhs) = tr(ji,jj,jk,jpno3,Krhs) - zpnewtot
tr(ji,jj,jk,jpnh4,Krhs) = tr(ji,jj,jk,jpnh4,Krhs) - zpregtot
!
tr(ji,jj,jk,jpphy,Krhs) = tr(ji,jj,jk,jpphy,Krhs) &
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zpronewn = zpronmaxn(ji,jj,jk) * xnanono3(ji,jj,jk)
zpronewp = zpronmaxp(ji,jj,jk) * xpicono3(ji,jj,jk)
zpronewd = zpronmaxd(ji,jj,jk) * xdiatno3(ji,jj,jk)
!
zproregn = zpronmaxn(ji,jj,jk) * xnanonh4(ji,jj,jk)
zproregp = zpronmaxp(ji,jj,jk) * xpiconh4(ji,jj,jk)
zproregd = zpronmaxd(ji,jj,jk) * xdiatnh4(ji,jj,jk)
!
zpropo4n = zpropmaxn(ji,jj,jk) * xnanopo4(ji,jj,jk)
zpropo4p = zpropmaxp(ji,jj,jk) * xpicopo4(ji,jj,jk)
zpropo4d = zpropmaxd(ji,jj,jk) * xdiatpo4(ji,jj,jk)
!
zprodopn = zpropmaxn(ji,jj,jk) * xnanodop(ji,jj,jk)
zprodopp = zpropmaxp(ji,jj,jk) * xpicodop(ji,jj,jk)
zprodopd = zpropmaxd(ji,jj,jk) * xdiatdop(ji,jj,jk)
!
zpptot = zpropo4n + zpropo4d + zpropo4p
zpnewtot = zpronewn + zpronewd + zpronewp
zpregtot = zproregn + zproregd + zproregp
zprontot = zpronewn + zproregn
zproptot = zpronewp + zproregp
zprodtot = zpronewd + zproregd
!
zproddoc = excretd * zprorcad(ji,jj,jk) &
& + excretn * zprorcan(ji,jj,jk) &
& + excretp * zprorcap(ji,jj,jk)
!
zproddop = excretd * zpropo4d - texcretd * zprodopd &
& + excretn * zpropo4n - texcretn * zprodopn &
& + excretp * zpropo4p - texcretp * zprodopp
zproddon = excretd * zprodtot + excretn * zprontot + excretp * zproptot
zprodfer = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk)
! CE : zrespn/d/p ????
! zresptot = zrespn(ji,jj,jk) + zrespp(ji,jj,jk) + zrespd(ji,jj,jk)
zresptot = 0._wp
!
tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) - zpptot
tr(ji,jj,jk,jpno3,Krhs) = tr(ji,jj,jk,jpno3,Krhs) - zpnewtot
tr(ji,jj,jk,jpnh4,Krhs) = tr(ji,jj,jk,jpnh4,Krhs) - zpregtot
!
tr(ji,jj,jk,jpphy,Krhs) = tr(ji,jj,jk,jpphy,Krhs) &
& + zprorcan(ji,jj,jk) * texcretn &
& - xpsino3 * zpronewn(ji,jj,jk) &
& - xpsinh4 * zproregn(ji,jj,jk) &
& - zrespn(ji,jj,jk)
tr(ji,jj,jk,jpnph,Krhs) = tr(ji,jj,jk,jpnph,Krhs) + zprontot * texcretn
tr(ji,jj,jk,jppph,Krhs) = tr(ji,jj,jk,jppph,Krhs) + ( zpropo4n(ji,jj,jk) + zprodopn(ji,jj,jk) ) * texcretn
tr(ji,jj,jk,jpnfe,Krhs) = tr(ji,jj,jk,jpnfe,Krhs) + zprofen(ji,jj,jk) * texcretn
!
tr(ji,jj,jk,jppic,Krhs) = tr(ji,jj,jk,jppic,Krhs) &
& + zprorcap(ji,jj,jk) * texcretp &
& - xpsino3 * zpronewp(ji,jj,jk) &
& - xpsinh4 * zproregp(ji,jj,jk) &
& - zrespp(ji,jj,jk)
tr(ji,jj,jk,jpnpi,Krhs) = tr(ji,jj,jk,jpnpi,Krhs) + zproptot * texcretp
tr(ji,jj,jk,jpppi,Krhs) = tr(ji,jj,jk,jpppi,Krhs) + ( zpropo4p(ji,jj,jk) + zprodopp(ji,jj,jk) ) * texcretp
tr(ji,jj,jk,jppfe,Krhs) = tr(ji,jj,jk,jppfe,Krhs) + zprofep(ji,jj,jk) * texcretp
!
tr(ji,jj,jk,jpdia,Krhs) = tr(ji,jj,jk,jpdia,Krhs) &
& + zprorcad(ji,jj,jk) * texcretd &
& - xpsino3 * zpronewd(ji,jj,jk) &
& - xpsinh4 * zproregd(ji,jj,jk) &
& - zrespd(ji,jj,jk)
tr(ji,jj,jk,jpndi,Krhs) = tr(ji,jj,jk,jpndi,Krhs) + zprodtot * texcretd
tr(ji,jj,jk,jppdi,Krhs) = tr(ji,jj,jk,jppdi,Krhs) + ( zpropo4d(ji,jj,jk) + zprodopd(ji,jj,jk) ) * texcretd
tr(ji,jj,jk,jpdfe,Krhs) = tr(ji,jj,jk,jpdfe,Krhs) + zprofed(ji,jj,jk) * texcretd
tr(ji,jj,jk,jpdsi,Krhs) = tr(ji,jj,jk,jpdsi,Krhs) + zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * tr(ji,jj,jk,jpdia,Kbb)
tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) + zproddoc
tr(ji,jj,jk,jpdon,Krhs) = tr(ji,jj,jk,jpdon,Krhs) + zproddon
tr(ji,jj,jk,jpdop,Krhs) = tr(ji,jj,jk,jpdop,Krhs) + zproddop
& - xpsino3 * zpronewn &
& - xpsinh4 * zproregn
! & - zrespn(ji,jj,jk)
tr(ji,jj,jk,jpnph,Krhs) = tr(ji,jj,jk,jpnph,Krhs) + zprontot * texcretn
tr(ji,jj,jk,jppph,Krhs) = tr(ji,jj,jk,jppph,Krhs) + ( zpropo4n + zprodopn ) * texcretn
tr(ji,jj,jk,jpnfe,Krhs) = tr(ji,jj,jk,jpnfe,Krhs) + zprofen(ji,jj,jk) * texcretn
!
tr(ji,jj,jk,jppic,Krhs) = tr(ji,jj,jk,jppic,Krhs) &
& + zprorcap(ji,jj,jk) * texcretp &
& - xpsino3 * zpronewp &
& - xpsinh4 * zproregp
! & - zrespp(ji,jj,jk)
tr(ji,jj,jk,jpnpi,Krhs) = tr(ji,jj,jk,jpnpi,Krhs) + zproptot * texcretp
tr(ji,jj,jk,jpppi,Krhs) = tr(ji,jj,jk,jpppi,Krhs) + ( zpropo4p + zprodopp ) * texcretp
tr(ji,jj,jk,jppfe,Krhs) = tr(ji,jj,jk,jppfe,Krhs) + zprofep(ji,jj,jk) * texcretp
!
tr(ji,jj,jk,jpdia,Krhs) = tr(ji,jj,jk,jpdia,Krhs) &
& + zprorcad(ji,jj,jk) * texcretd &
& - xpsino3 * zpronewd &
& - xpsinh4 * zproregd
! & - zrespd(ji,jj,jk)
!
zprodsil = zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * tr(ji,jj,jk,jpdia,Kbb)
!
tr(ji,jj,jk,jpndi,Krhs) = tr(ji,jj,jk,jpndi,Krhs) + zprodtot * texcretd
tr(ji,jj,jk,jppdi,Krhs) = tr(ji,jj,jk,jppdi,Krhs) + ( zpropo4d + zprodopd ) * texcretd
tr(ji,jj,jk,jpdfe,Krhs) = tr(ji,jj,jk,jpdfe,Krhs) + zprofed(ji,jj,jk) * texcretd
tr(ji,jj,jk,jpdsi,Krhs) = tr(ji,jj,jk,jpdsi,Krhs) + zprodsil
tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) + zproddoc
tr(ji,jj,jk,jpdon,Krhs) = tr(ji,jj,jk,jpdon,Krhs) + zproddon
tr(ji,jj,jk,jpdop,Krhs) = tr(ji,jj,jk,jpdop,Krhs) + zproddop
tr(ji,jj,jk,jpoxy,Krhs) = tr(ji,jj,jk,jpoxy,Krhs) &
tr(ji,jj,jk,jpoxy,Krhs) = tr(ji,jj,jk,jpoxy,Krhs) &
& + o2ut * zpregtot + ( o2ut + o2nit ) * zpnewtot - o2ut * zresptot
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - zprodfer
consfe3(ji,jj,jk) = zprodfer * 75.0 / ( rtrn + ( plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) ) &
& * tr(ji,jj,jk,jpfer,Kbb) ) / rfact2
tr(ji,jj,jk,jpsil,Krhs) = tr(ji,jj,jk,jpsil,Krhs) - zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * tr(ji,jj,jk,jpdia,Kbb)
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - zprodfer
consfe3(ji,jj,jk) = zprodfer * 75.0 / ( rtrn + ( plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) ) &
& * tr(ji,jj,jk,jpfer,Kbb) ) / rfact2
tr(ji,jj,jk,jpsil,Krhs) = tr(ji,jj,jk,jpsil,Krhs) - zprodsil
tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) - zpptot &
& + xpsino3 * zpronewn(ji,jj,jk) + xpsinh4 * zproregn(ji,jj,jk) &
& + xpsino3 * zpronewp(ji,jj,jk) + xpsinh4 * zproregp(ji,jj,jk) &
& + xpsino3 * zpronewd(ji,jj,jk) + xpsinh4 * zproregd(ji,jj,jk)
tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) - zpptot &
& + xpsino3 * zpronewn + xpsinh4 * zproregn &
& + xpsino3 * zpronewp + xpsinh4 * zproregp &
& + xpsino3 * zpronewd + xpsinh4 * zproregd
tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) + rno3 * ( zpnewtot - zpregtot )
!
tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) + rno3 * ( zpnewtot - zpregtot )
!
END_3D
! Production and uptake of ligands by phytoplankton. This part is activated
......@@ -513,7 +545,7 @@ CONTAINS
! Shaked and Lis (2012)
! -------------------------------------------------------------------------
IF( ln_ligand ) THEN
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1)
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zproddoc = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk) + excretp * zprorcap(ji,jj,jk)
zprodfer = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk)
zprodlig = plig(ji,jj,jk) / ( rtrn + plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) ) * lthet
......@@ -522,47 +554,141 @@ CONTAINS
END_3D
ENDIF
! Output of the diagnostics
! Total primary production per year
IF( iom_use( "tintpp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) &
& tpp = glob_sum( 'p5zprod', ( zprorcan(:,:,:) + zprorcad(:,:,:) + zprorcap(:,:,:) ) * cvol(:,:,:) )
IF( l_dia_ppphy .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) THEN
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1)
zw3d(ji,jj,jk) = ( zprorcan(ji,jj,jk) + zprorcad(ji,jj,jk) + zprorcap(ji,jj,jk) ) * cvol(ji,jj,jk)
END_3D
tpp = glob_sum( 'p5zprod', zw3d )
DEALLOCATE ( zw3d )
ENDIF
IF( lk_iomput .AND. knt == nrdttrc ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
!
CALL iom_put( "PPPHYP" , zprorcap(:,:,:) * zfact * tmask(:,:,:) ) ! primary production by picophyto
CALL iom_put( "PPPHYN" , zprorcan(:,:,:) * zfact * tmask(:,:,:) ) ! primary production by nanophyto
CALL iom_put( "PPPHYD" , zprorcad(:,:,:) * zfact * tmask(:,:,:) ) ! primary production by diatomes
CALL iom_put( "PPNEWN" , zpronewp(:,:,:) * zfact * tmask(:,:,:) ) ! new primary production by picophyto
CALL iom_put( "PPNEWN" , zpronewn(:,:,:) * zfact * tmask(:,:,:) ) ! new primary production by nanophyto
CALL iom_put( "PPNEWD" , zpronewd(:,:,:) * zfact * tmask(:,:,:) ) ! new primary production by diatomes
CALL iom_put( "PBSi" , zprmaxd (:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:) ) ! biogenic silica production
CALL iom_put( "PFeP" , zprofep (:,:,:) * zfact * tmask(:,:,:) ) ! biogenic iron production by picophyto
CALL iom_put( "PFeN" , zprofen(:,:,:) * zfact * tmask(:,:,:) ) ! biogenic iron production by nanophyto
CALL iom_put( "PFeD" , zprofed(:,:,:) * zfact * tmask(:,:,:) ) ! biogenic iron production by diatomes
IF( ln_ligand .AND. ( iom_use( "LPRODP" ) .OR. iom_use( "LDETP" ) ) ) THEN
ALLOCATE( zpligprod(jpi,jpj,jpk) )
zpligprod(:,:,:) = excretd * zprorcad(:,:,:) + excretn * zprorcan(:,:,:) + excretp * zprorcap(:,:,:)
CALL iom_put( "LPRODP" , zpligprod(:,:,:) * ldocp * 1e9 * zfact * tmask(:,:,:) )
!
zpligprod(:,:,:) = ( texcretn * zprofen(:,:,:) + texcretd * zprofed(:,:,:) + texcretp * zprofep(:,:,:) ) &
& * plig(:,:,:) / ( rtrn + plig(:,:,:) + 75.0 * (1.0 - plig(:,:,:) ) )
CALL iom_put( "LDETP" , zpligprod(:,:,:) * lthet * 1e9 * zfact * tmask(:,:,:) )
DEALLOCATE( zpligprod )
IF( l_dia_ppphy ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
! primary production by nanophyto
zw3d(A2D(0),1:jpkm1) = zprorcan(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PPPHYN", zw3d )
! primary production by diatomes
zw3d(A2D(0),1:jpkm1) = zprorcad(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PPPHYD", zw3d )
! primary production by pico
zw3d(A2D(0),1:jpkm1) = zprorcap(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PPPHYP", zw3d )
! total primary production
zw3d(A2D(0),1:jpkm1) = ( zprorcan(A2D(0),1:jpkm1) + zprorcad(A2D(0),1:jpkm1) + zprorcap(A2D(0),1:jpkm1) ) &
& * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "TPP", zw3d )
CALL iom_put( "tintpp" , tpp * zfact ) ! global total integrated primary production molC/s
DEALLOCATE ( zw3d )
ENDIF
!
IF( l_dia_ppnew ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
! new primary production by nano
zw3d(A2D(0),1:jpkm1) = zpronmaxn(A2D(0),1:jpkm1) * xnanono3(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PPNEWN", zw3d )
! new primary production by diatomes
zw3d(A2D(0),1:jpkm1) = zpronmaxd(A2D(0),1:jpkm1) * xdiatno3(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PPNEWD", zw3d )
! new primary production by pico
zw3d(A2D(0),1:jpkm1) = zpronmaxp(A2D(0),1:jpkm1) * xpicono3(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PPNEWP", zw3d )
! total new production
zw3d(A2D(0),1:jpkm1) = ( zpronmaxn(A2D(0),1:jpkm1) * xnanono3(A2D(0),1:jpkm1) + &
& zpronmaxd(A2D(0),1:jpkm1) * xdiatno3(A2D(0),1:jpkm1) + &
& zpronmaxp(A2D(0),1:jpkm1) * xpicono3(A2D(0),1:jpkm1) ) &
& * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "TPNEW", zw3d )
DEALLOCATE ( zw3d )
ENDIF
CALL iom_put( "Mumax" , zprmaxn(:,:,:) * tmask(:,:,:) ) ! Maximum growth rate
CALL iom_put( "MuP" , zprpic(:,:,:) * xlimpic(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for picophyto
CALL iom_put( "MuN" , zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for nanophyto
CALL iom_put( "MuD" , zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for diatoms
CALL iom_put( "LPlight" , zprpic(:,:,:) / (zprmaxp(:,:,:) + rtrn) * tmask(:,:,:) ) ! light limitation term
CALL iom_put( "LNlight" , zprbio(:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ) ! light limitation term
CALL iom_put( "LDlight" , zprdia(:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:) )
CALL iom_put( "MunetP" , ( tr(:,:,:,jppic,Krhs)/rfact2/(tr(:,:,:,jppic,Kbb)+ rtrn ) * tmask(:,:,:)) ) ! Realized growth rate for picophyto
CALL iom_put( "MunetN" , ( tr(:,:,:,jpphy,Krhs)/rfact2/(tr(:,:,:,jpphy,Kbb)+ rtrn ) * tmask(:,:,:)) ) ! Realized growth rate for picophyto
CALL iom_put( "MunetD" , ( tr(:,:,:,jpdia,Krhs)/rfact2/(tr(:,:,:,jpdia,Kbb)+ rtrn ) * tmask(:,:,:)) ) ! Realized growth rate for picophyto
CALL iom_put( "TPP" , ( zprorcap(:,:,:) + zprorcan(:,:,:) + zprorcad(:,:,:) ) * zfact * tmask(:,:,:) ) ! total primary production
CALL iom_put( "TPNEW" , ( zpronewp(:,:,:) + zpronewn(:,:,:) + zpronewd(:,:,:) ) * zfact * tmask(:,:,:) ) ! total new production
CALL iom_put( "TPBFE" , ( zprofep (:,:,:) + zprofen (:,:,:) + zprofed (:,:,:) ) * zfact * tmask(:,:,:) ) ! total biogenic iron production
CALL iom_put( "tintpp" , tpp * zfact ) ! global total integrated primary production molC/s
!
IF( l_dia_ppbsi ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
! biogenic silica production
zw3d(A2D(0),1:jpkm1) = zprmaxd(A2D(0),1:jpkm1) * zysopt(A2D(0),1:jpkm1) &
& * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PBSi", zw3d )
DEALLOCATE ( zw3d )
ENDIF
!
!
IF( l_dia_ppbfe ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
! biogenic iron production by nanophyto
zw3d(A2D(0),1:jpkm1) = zprofen(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PFeN", zw3d )
! biogenic iron production by diatomes
zw3d(A2D(0),1:jpkm1) = zprofed(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PFeD", zw3d )
! biogenic iron production by pico
zw3d(A2D(0),1:jpkm1) = zprofep(A2D(0),1:jpkm1) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "PFeP", zw3d )
! total biogenic iron production
zw3d(A2D(0),1:jpkm1) = ( zprofen(A2D(0),1:jpkm1) + zprofed(A2D(0),1:jpkm1) + zprofep(A2D(0),1:jpkm1) ) &
& * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "TPBFE", zw3d )
DEALLOCATE ( zw3d )
ENDIF
!
IF( l_dia_mu ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = zprmaxn(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "Mumax", zw3d )
! Realized growth rate for nanophyto
zw3d(A2D(0),1:jpkm1) = zprbio(A2D(0),1:jpkm1) * xlimphy(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "MuN", zw3d )
! Realized growth rate for diatoms
zw3d(A2D(0),1:jpkm1) = zprdia(A2D(0),1:jpkm1) * xlimdia(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "MuD", zw3d )
! Realized growth rate for pico
zw3d(A2D(0),1:jpkm1) = zprpic(A2D(0),1:jpkm1) * xlimpic(A2D(0),1:jpkm1) * tmask(A2D(0),1:jpkm1)
CALL iom_put( "MuP", zw3d )
DEALLOCATE ( zw3d )
ENDIF
!
!
IF( l_dia_light ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
! light limitation term for nano
zw3d(A2D(0),1:jpkm1) = zprbio(A2D(0),1:jpkm1) / (zprmaxn(A2D(0),1:jpkm1)+rtrn) &
& * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LNlight", zw3d )
! light limitation term for diatomes
zw3d(A2D(0),1:jpkm1) = zprdia(A2D(0),1:jpkm1) / (zprmaxd(A2D(0),1:jpkm1)+rtrn) &
& * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LDlight", zw3d )
! light limitation term for pico
zw3d(A2D(0),1:jpkm1) = zprpic(A2D(0),1:jpkm1) / (zprmaxp(A2D(0),1:jpkm1)+rtrn) &
& * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LPlight", zw3d )
DEALLOCATE ( zw3d )
ENDIF
!
IF( l_dia_lprod ) THEN
zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = ( excretd * zprorcad(A2D(0),1:jpkm1) + excretn * zprorcan(A2D(0),1:jpkm1) + &
& excretp * zprorcap(A2D(0),1:jpkm1) ) * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LPRODP" , zw3d * ldocp * 1e9 )
!
zw3d(A2D(0),1:jpkm1) = ( texcretn * zprofen(A2D(0),1:jpkm1) + texcretd * zprofed(A2D(0),1:jpkm1) + &
& texcretp * zprofep(A2D(0),1:jpkm1) ) * plig(A2D(0),1:jpkm1) &
& / ( rtrn + plig(A2D(0),1:jpkm1) + 75.0 * (1.0 - plig(A2D(0),1:jpkm1) ) ) &
& * zfact * tmask(A2D(0),1:jpkm1)
CALL iom_put( "LDETP" , zw3d * lthet * 1e9 )
DEALLOCATE ( zw3d )
ENDIF
!
ENDIF
IF(sn_cfctl%l_prttrc) THEN ! print mean trends (used for debugging)
......@@ -625,17 +751,11 @@ CONTAINS
texcretd = 1._wp - excretd
tpp = 0._wp
!
xq10_n = 1. + xpsino3 * qnnmax
xq10_d = 1. + xpsino3 * qndmax
xq10_p = 1. + xpsino3 * qnpmax
!
END SUBROUTINE p5z_prod_init
INTEGER FUNCTION p5z_prod_alloc()
!!----------------------------------------------------------------------
!! *** ROUTINE p5z_prod_alloc ***
!!----------------------------------------------------------------------
ALLOCATE( zdaylen(jpi,jpj), STAT = p5z_prod_alloc )
!
IF( p5z_prod_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p5z_prod_alloc : failed to allocate arrays.' )
!
END FUNCTION p5z_prod_alloc
!!======================================================================
END MODULE p5zprod
......@@ -138,7 +138,7 @@ CONTAINS
IF (ln_sediment_offline) THEN
CALL sed_chem_cst
ELSE
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
IF ( tmask(ji,jj,ikt) == 1 ) THEN
zchem_data(ji,jj,1) = ak13 (ji,jj,ikt)
......
......@@ -49,7 +49,7 @@ CONTAINS
CALL unpack_arr ( jpoce, trc_data(1:jpi,1:jpj,8), iarroce(1:jpoce), pwcp(1:jpoce,1,jwfe2) )
CALL unpack_arr ( jpoce, trc_data(1:jpi,1:jpj,9), iarroce(1:jpoce), pwcp(1:jpoce,1,jwlgw) )
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
IF ( tmask(ji,jj,ikt) == 1 ) THEN
tr(ji,jj,ikt,jptal,Kbb) = trc_data(ji,jj,1)
......
......@@ -92,7 +92,7 @@ CONTAINS
jl = n_trc_index(jn)
CALL trc_dta( kt, jl, ztrcdta ) ! read tracer data at nit000
!
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
ikt = mbkt(ji,jj)
tr(ji,jj,ikt,jn,Kbb) = ztrcdta(ji,jj,ikt) + ( tr(ji,jj,ikt,jn,Kbb) - ztrcdta(ji,jj,ikt) ) &
& * exp( -restosed(ji,jj,ikt) * dtsed )
......
......@@ -124,6 +124,8 @@ MODULE sms_pisces
LOGICAL, SAVE :: lk_sed
!! * Substitutions
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/TOP 4.0 , NEMO Consortium (2018)
!! $Id: sms_pisces.F90 15459 2021-10-29 08:19:18Z cetlod $
......@@ -140,52 +142,52 @@ CONTAINS
!!----------------------------------------------------------------------
ierr(:) = 0
!* Biological fluxes for light : shared variables for pisces & lobster
ALLOCATE( xksi(jpi,jpj), strn(jpi,jpj), STAT=ierr(1) )
ALLOCATE( xksi(A2D(0)), strn(A2D(0)), STAT=ierr(1) )
IF( ln_p4z .OR. ln_p5z ) THEN
!* Optics
ALLOCATE( enano(jpi,jpj,jpk) , ediat(jpi,jpj,jpk) , &
& enanom(jpi,jpj,jpk), ediatm(jpi,jpj,jpk), &
& emoy(jpi,jpj,jpk) , etotm(jpi,jpj,jpk), STAT=ierr(2) )
ALLOCATE( enano(A2D(0),jpk) , ediat(A2D(0),jpk) , &
& enanom(A2D(0),jpk), ediatm(A2D(0),jpk), &
& emoy(A2D(0),jpk) , etotm(A2D(0),jpk), STAT=ierr(2) )
!* Biological SMS
ALLOCATE( xksimax(jpi,jpj) , biron(jpi,jpj,jpk) , STAT=ierr(3) )
ALLOCATE( xksimax(A2D(0)) , biron(A2D(0),jpk) , STAT=ierr(3) )
! Biological SMS
ALLOCATE( xfracal (jpi,jpj,jpk), orem (jpi,jpj,jpk), &
& nitrfac (jpi,jpj,jpk), nitrfac2(jpi,jpj,jpk), &
& prodcal (jpi,jpj,jpk), xdiss (jpi,jpj,jpk), &
& prodpoc (jpi,jpj,jpk), conspoc (jpi,jpj,jpk), &
& prodgoc (jpi,jpj,jpk), consgoc (jpi,jpj,jpk), &
& blim (jpi,jpj,jpk), consfe3 (jpi,jpj,jpk), &
& xfecolagg(jpi,jpj,jpk), xcoagfe (jpi,jpj,jpk), STAT=ierr(4) )
ALLOCATE( xfracal (A2D(0),jpk), orem (A2D(0),jpk), &
& nitrfac (A2D(0),jpk), nitrfac2(A2D(0),jpk), &
& prodcal (A2D(0),jpk), xdiss (A2D(0),jpk), &
& prodpoc (A2D(0),jpk), conspoc (A2D(0),jpk), &
& prodgoc (A2D(0),jpk), consgoc (A2D(0),jpk), &
& blim (A2D(0),jpk), consfe3 (A2D(0),jpk), &
& xfecolagg(A2D(0),jpk), xcoagfe (A2D(0),jpk), STAT=ierr(4) )
!* Carbonate chemistry
ALLOCATE( ak13 (jpi,jpj,jpk) , &
& ak23(jpi,jpj,jpk) , aksp (jpi,jpj,jpk) , &
& hi (jpi,jpj,jpk) , excess(jpi,jpj,jpk) , &
& aphscale(jpi,jpj,jpk), STAT=ierr(5) )
ALLOCATE( ak13(A2D(0),jpk), &
& ak23(A2D(0),jpk), aksp (A2D(0),jpk) , &
& hi (A2D(0),jpk), excess(A2D(0),jpk) , &
& aphscale(A2D(0),jpk), STAT=ierr(5) )
!
!* Temperature dependency of SMS terms
ALLOCATE( tgfunc (jpi,jpj,jpk) , tgfunc2(jpi,jpj,jpk), STAT=ierr(6) )
ALLOCATE( tgfunc (A2D(0),jpk) , tgfunc2(A2D(0),jpk), STAT=ierr(6) )
!
!* Sinking speed
ALLOCATE( wsbio3 (jpi,jpj,jpk) , wsbio4 (jpi,jpj,jpk), STAT=ierr(7) )
ALLOCATE( wsbio3 (A2D(0),jpk) , wsbio4 (A2D(0),jpk), STAT=ierr(7) )
!* Size of phytoplankton cells
ALLOCATE( sizen (jpi,jpj,jpk), sized (jpi,jpj,jpk), &
& sizena(jpi,jpj,jpk), sizeda(jpi,jpj,jpk), STAT=ierr(8) )
ALLOCATE( sizen (A2D(0),jpk), sized (A2D(0),jpk), &
& sizena(A2D(0),jpk), sizeda(A2D(0),jpk), STAT=ierr(8) )
!
ALLOCATE( plig(jpi,jpj,jpk) , STAT=ierr(9) )
ALLOCATE( plig(A2D(0),jpk) , STAT=ierr(9) )
ENDIF
!
IF( ln_p5z ) THEN
! PISCES-QUOTA specific part
ALLOCATE( epico(jpi,jpj,jpk) , epicom(jpi,jpj,jpk) , STAT=ierr(10) )
ALLOCATE( epico(A2D(0),jpk) , epicom(A2D(0),jpk) , STAT=ierr(10) )
!* Size of phytoplankton cells
ALLOCATE( sizep(jpi,jpj,jpk), sizepa(jpi,jpj,jpk), STAT=ierr(11) )
ALLOCATE( sizep(A2D(0),jpk), sizepa(A2D(0),jpk), STAT=ierr(11) )
ENDIF
!
sms_pisces_alloc = MAXVAL( ierr )
......
......@@ -19,6 +19,8 @@ MODULE trcice_pisces
PUBLIC trc_ice_ini_pisces ! called by trcini.F90 module
!! * Substitutions
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/TOP 4.0 , NEMO Consortium (2018)
!! $Id: trcice_pisces.F90 10794 2019-03-22 09:25:28Z cetlod $
......@@ -283,15 +285,15 @@ CONTAINS
ENDIF
!
DO jn = jp_pcs0, jp_pcs1
IF( cn_trc_o(jn) == 'GL ' ) trc_o(:,:,jn) = zpisc(jn,1) ! Global case
IF( cn_trc_o(jn) == 'GL ' ) trc_o(A2D(0),jn) = zpisc(jn,1) ! Global case
IF( cn_trc_o(jn) == 'AA ' ) THEN
WHERE( gphit(:,:) >= 0._wp ) ; trc_o(:,:,jn) = zpisc(jn,2) ; END WHERE ! Arctic
WHERE( gphit(:,:) < 0._wp ) ; trc_o(:,:,jn) = zpisc(jn,3) ; END WHERE ! Antarctic
WHERE( gphit(A2D(0)) >= 0._wp ) ; trc_o(A2D(0),jn) = zpisc(jn,2) ; END WHERE ! Arctic
WHERE( gphit(A2D(0)) < 0._wp ) ; trc_o(A2D(0),jn) = zpisc(jn,3) ; END WHERE ! Antarctic
ENDIF
IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA" ) THEN ! Baltic Sea particular case for ORCA configurations
WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND. &
54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp )
trc_o(:,:,jn) = zpisc(jn,4)
WHERE( 14._wp <= glamt(A2D(0)) .AND. glamt(A2D(0)) <= 32._wp .AND. &
54._wp <= gphit(A2D(0)) .AND. gphit(A2D(0)) <= 66._wp )
trc_o(A2D(0),jn) = zpisc(jn,4)
END WHERE
ENDIF
ENDDO
......@@ -321,16 +323,16 @@ CONTAINS
DO jn = jp_pcs0, jp_pcs1
!-- Everywhere but in the Baltic
IF ( trc_ice_ratio(jn) >= -1._wp ) THEN ! no prescribed conc. ; typically everything but iron)
trc_i(:,:,jn) = zratio(jn,1) * trc_o(:,:,jn)
trc_i(A2D(0),jn) = zratio(jn,1) * trc_o(A2D(0),jn)
ELSE ! prescribed concentration
trc_i(:,:,jn) = trc_ice_prescr(jn)
trc_i(A2D(0),jn) = trc_ice_prescr(jn)
ENDIF
!-- Baltic
IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA" ) THEN
IF ( trc_ice_ratio(jn) >= - 1._wp ) THEN ! no prescribed conc. ; typically everything but iron)
WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND. &
54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp )
trc_i(:,:,jn) = zratio(jn,2) * trc_o(:,:,jn)
WHERE( 14._wp <= glamt(A2D(0)) .AND. glamt(A2D(0)) <= 32._wp .AND. &
54._wp <= gphit(A2D(0)) .AND. gphit(A2D(0)) <= 66._wp )
trc_i(A2D(0),jn) = zratio(jn,2) * trc_o(A2D(0),jn)
END WHERE
ENDIF
ENDIF
......
......@@ -123,7 +123,6 @@ CONTAINS
ELSE
! PISCES-QUOTA part
ierr = ierr + p5z_lim_alloc()
ierr = ierr + p5z_prod_alloc()
ierr = ierr + p5z_meso_alloc()
ENDIF
ierr = ierr + p4z_rem_alloc()
......
......@@ -38,7 +38,7 @@ CONTAINS
CHARACTER (len=20) :: cltra
REAL(wp) :: zfact
INTEGER :: ji, jj, jk, jn
REAL(wp), DIMENSION(jpi,jpj) :: zdic, zo2min, zdepo2min
REAL(wp), DIMENSION(A2D(0)) :: zdic, zo2min, zdepo2min
!!---------------------------------------------------------------------
! write the tracer concentrations in the file
......@@ -60,15 +60,19 @@ CONTAINS
IF( iom_use( "INTDIC" ) ) THEN ! DIC content in kg/m2
zdic(:,:) = 0.
DO jk = 1, jpkm1
zdic(:,:) = zdic(:,:) + tr(:,:,jk,jpdic,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) * 12.
DO_2D( 0, 0, 0, 0 )
zdic(ji,jj) = zdic(ji,jj) + tr(ji,jj,jk,jpdic,Kmm) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) * 12.
END_2D
ENDDO
CALL iom_put( 'INTDIC', zdic )
CALL iom_put( 'INTDIC', zdic )
ENDIF
!
IF( iom_use( "O2MIN" ) .OR. iom_use ( "ZO2MIN" ) ) THEN ! Oxygen minimum concentration and depth
zo2min (:,:) = tr(:,:,1,jpoxy,Kmm) * tmask(:,:,1)
zdepo2min(:,:) = gdepw(:,:,1,Kmm) * tmask(:,:,1)
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 2, jpkm1 )
IF( iom_use( "O2MIN" ) .OR. iom_use ( "ZO2MIN" ) ) THEN ! Oxygen minimum concentration and depth
DO_2D( 0, 0, 0, 0 )
zo2min (ji,jj) = tr(ji,jj,1,jpoxy,Kmm) * tmask(ji,jj,1)
zdepo2min(ji,jj) = gdepw(ji,jj,1,Kmm) * tmask(ji,jj,1)
END_2D
DO_3D( 0, 0, 0, 0, 2, jpkm1 )
IF( tmask(ji,jj,jk) == 1 ) then
IF( tr(ji,jj,jk,jpoxy,Kmm) < zo2min(ji,jj) ) then
zo2min (ji,jj) = tr(ji,jj,jk,jpoxy,Kmm)
......
......@@ -123,19 +123,19 @@ CONTAINS
!
IF( ln_wave .AND. ln_sdw ) THEN
DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 ) ! eulerian transport + Stokes Drift
zuu(ji,jj,jk) = e2u (ji,jj) * e3u(ji,jj,jk,Kmm) * ( zptu(ji,jj,jk) + usd(ji,jj,jk) )
zvv(ji,jj,jk) = e1v (ji,jj) * e3v(ji,jj,jk,Kmm) * ( zptv(ji,jj,jk) + vsd(ji,jj,jk) )
zuu(ji,jj,jk) = e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * ( zptu(ji,jj,jk) + usd(ji,jj,jk) )
zvv(ji,jj,jk) = e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * ( zptv(ji,jj,jk) + vsd(ji,jj,jk) )
END_3D
DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )
zww(ji,jj,jk) = e1e2t(ji,jj) * ( zptw(ji,jj,jk) + wsd(ji,jj,jk) )
zww(ji,jj,jk) = e1e2t(ji,jj) * ( zptw(ji,jj,jk) + wsd(ji,jj,jk) )
END_3D
ELSE
DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )
zuu(ji,jj,jk) = e2u (ji,jj) * e3u(ji,jj,jk,Kmm) * zptu(ji,jj,jk) ! eulerian transport
zvv(ji,jj,jk) = e1v (ji,jj) * e3v(ji,jj,jk,Kmm) * zptv(ji,jj,jk)
zuu(ji,jj,jk) = e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * zptu(ji,jj,jk) ! eulerian transport
zvv(ji,jj,jk) = e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * zptv(ji,jj,jk)
END_3D
DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )
zww(ji,jj,jk) = e1e2t(ji,jj) * zptw(ji,jj,jk)
zww(ji,jj,jk) = e1e2t(ji,jj) * zptw(ji,jj,jk)
END_3D
ENDIF
!
......
......@@ -57,7 +57,7 @@ CONTAINS
!!----------------------------------------------------------------------
!! *** ROUTINE trc_dmp_alloc ***
!!----------------------------------------------------------------------
ALLOCATE( restotr(jpi,jpj,jpk) , STAT=trc_dmp_alloc )
ALLOCATE( restotr(A2D(0),jpk) , STAT=trc_dmp_alloc )
!
IF( trc_dmp_alloc /= 0 ) CALL ctl_warn('trc_dmp_alloc: failed to allocate array')
!
......@@ -329,11 +329,11 @@ CONTAINS
! convert the position in local domain indices
! --------------------------------------------
DO jc = 1, npncts
nctsi1(jc) = mi0( nctsi1(jc) )
nctsj1(jc) = mj0( nctsj1(jc) )
nctsi1(jc) = mi0( nctsi1(jc), nn_hls )
nctsj1(jc) = mj0( nctsj1(jc), nn_hls )
!
nctsi2(jc) = mi1( nctsi2(jc) )
nctsj2(jc) = mj1( nctsj2(jc) )
nctsi2(jc) = mi1( nctsi2(jc), nn_hls )
nctsj2(jc) = mj1( nctsj2(jc), nn_hls )
END DO
!
ENDIF
......
......@@ -115,14 +115,14 @@ CONTAINS
CASE ( -1 ) ! ! No tracers in sea ice ( trc_i = 0 )
!
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
sbc_trc(ji,jj,jn) = 0._wp
END_2D
END DO
!
IF( ln_linssh ) THEN !* linear free surface
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
sbc_trc(ji,jj,jn) = sbc_trc(ji,jj,jn) + r1_rho0 * emp(ji,jj) * ptr(ji,jj,1,jn,Kmm) !==>> add concentration/dilution effect due to constant volume cell
END_2D
END DO
......@@ -131,14 +131,14 @@ CONTAINS
CASE ( 0 ) ! Same concentration in sea ice and in the ocean ( trc_i = ptr(...,Kmm) )
!
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
sbc_trc(ji,jj,jn) = - fmmflx(ji,jj) * r1_rho0 * ptr(ji,jj,1,jn,Kmm)
END_2D
END DO
!
IF( ln_linssh ) THEN !* linear free surface
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
sbc_trc(ji,jj,jn) = sbc_trc(ji,jj,jn) + r1_rho0 * emp(ji,jj) * ptr(ji,jj,1,jn,Kmm) !==>> add concentration/dilution effect due to constant volume cell
END_2D
END DO
......@@ -147,21 +147,21 @@ CONTAINS
CASE ( 1 ) ! Specific treatment of sea ice fluxes with an imposed concentration in sea ice
!
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
sbc_trc(ji,jj,jn) = - fmmflx(ji,jj) * r1_rho0 * trc_i(ji,jj,jn)
END_2D
END DO
!
IF( ln_linssh ) THEN !* linear free surface
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
sbc_trc(ji,jj,jn) = sbc_trc(ji,jj,jn) + r1_rho0 * emp(ji,jj) * ptr(ji,jj,1,jn,Kmm) !==>> add concentration/dilution effect due to constant volume cell
END_2D
END DO
ENDIF
!
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
zse3t = rDt_trc / e3t(ji,jj,1,Kmm)
zdtra = ptr(ji,jj,1,jn,Kmm) + sbc_trc(ji,jj,jn) * zse3t
IF( zdtra < 0. ) sbc_trc(ji,jj,jn) = MAX( zdtra, -ptr(ji,jj,1,jn,Kmm) / zse3t ) ! avoid negative concentration that can occurs if trc_i > ptr
......@@ -176,7 +176,7 @@ CONTAINS
!
IF( l_trdtrc ) ztrtrd(:,:,:) = ptr(:,:,:,jn,Krhs) ! save trends
!
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
zse3t = zfact / e3t(ji,jj,1,Kmm)
ptr(ji,jj,1,jn,Krhs) = ptr(ji,jj,1,jn,Krhs) + ( sbc_trc_b(ji,jj,jn) + sbc_trc(ji,jj,jn) ) * zse3t
END_2D
......@@ -295,7 +295,7 @@ CONTAINS
CASE ( 0 ) ! Same concentration in sea ice and in the ocean fmm contribution to concentration/dilution effect has to be removed
!
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
z1_rho0_e3t = r1_rho0 / e3t(ji,jj,1,Kmm)
ptr(ji,jj,1,jn,Krhs) = ptr(ji,jj,1,jn,Krhs) + ( emp(ji,jj) - fmmflx(ji,jj) ) * r1_rho0 * ptr(ji,jj,1,jn,Kmm)
END_2D
......@@ -331,7 +331,7 @@ CONTAINS
CASE ( 0 ) ! Same concentration in sea ice and in the ocean : correct concentration/dilution effect due to "freezing - melting"
!
DO jn = 1, jptra
DO_2D( 0, 0, 0, 1 )
DO_2D( 0, 0, 0, 0 )
z1_rho0_e3t = r1_rho0 / e3t(ji,jj,1,Kmm)
ptr(ji,jj,1,jn,Krhs) = ptr(ji,jj,1,jn,Krhs) - fmmflx(ji,jj) * r1_rho0 * ptr(ji,jj,1,jn,Kmm)
END_2D
......
......@@ -50,12 +50,12 @@ CONTAINS
INTEGER , INTENT(in) :: Kbb, Kmm
INTEGER , INTENT(in) :: jp_tra ! tracer index index
REAL(wp), INTENT(in) :: rsfact ! time step duration
REAL(wp), INTENT(in) , DIMENSION(jpi,jpj,jpk) :: pwsink
REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: psinkflx
REAL(wp), INTENT(in) , DIMENSION(A2D(0),jpk) :: pwsink
REAL(wp), INTENT(inout), DIMENSION(A2D(0),jpk) :: psinkflx
INTEGER :: ji, jj, jk
INTEGER, DIMENSION(jpi, jpj) :: iiter
INTEGER, DIMENSION(A2D(0)) :: iiter
REAL(wp) :: zfact, zwsmax, zmax
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwsink
REAL(wp), DIMENSION(A2D(0),jpk) :: zwsink
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('trc_sink')
......@@ -73,7 +73,7 @@ CONTAINS
IF( nitermax == 1 ) THEN
iiter(:,:) = 1
ELSE
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
iiter(ji,jj) = 1
DO jk = 1, jpkm1
IF( tmask(ji,jj,jk) == 1.0 ) THEN
......@@ -85,7 +85,7 @@ CONTAINS
iiter(:,:) = MIN( iiter(:,:), nitermax )
ENDIF
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
IF( tmask(ji,jj,jk) == 1.0 ) THEN
zwsmax = 0.5 * e3t(ji,jj,jk,Kmm) * rday / rsfact
zwsink(ji,jj,jk) = MIN( pwsink(ji,jj,jk), zwsmax * REAL( iiter(ji,jj), wp ) )
......@@ -121,23 +121,25 @@ CONTAINS
INTEGER, INTENT(in ) :: Kbb, Kmm ! time level indices
INTEGER, INTENT(in ) :: jp_tra ! tracer index index
REAL(wp), INTENT(in ) :: rsfact ! duration of time step
INTEGER, INTENT(in ), DIMENSION(jpi,jpj) :: kiter ! number of iterations for time-splitting
REAL(wp), INTENT(in ), DIMENSION(jpi,jpj,jpk) :: pwsink ! sinking speed
REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: psinkflx ! sinking fluxe
INTEGER, INTENT(in ), DIMENSION(A2D(0)) :: kiter ! number of iterations for time-splitting
REAL(wp), INTENT(in ), DIMENSION(A2D(0),jpk) :: pwsink ! sinking speed
REAL(wp), INTENT(inout), DIMENSION(A2D(0),jpk) :: psinkflx ! sinking fluxe
!
INTEGER :: ji, jj, jk, jn, jt
REAL(wp) :: zigma,zew,zign, zflx, zstep
REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztraz, zakz, zwsink2, ztrb, psinking
REAL(wp), DIMENSION(A2D(0),jpk) :: ztraz, zakz, zwsink2, ztrb, psinking
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('trc_sink2')
!
DO jk = 1, jpkm1
zwsink2(:,:,jk+1) = -pwsink(:,:,jk) / rday * tmask(:,:,jk+1)
END DO
zwsink2(:,:,1) = 0.e0
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zwsink2(ji,jj,jk+1) = -pwsink(ji,jj,jk) / rday * tmask(ji,jj,jk+1)
END_3D
DO_2D( 0, 0, 0, 0 )
zwsink2(ji,jj,1) = 0.e0
END_2D
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
! Vertical advective flux
zstep = rsfact / REAL( kiter(ji,jj), wp ) / 2.
DO jt = 1, kiter(ji,jj)
......
......@@ -158,19 +158,19 @@ CONTAINS
!!-------------------------------------------------------------------
ierr(:) = 0
!
ALLOCATE( tr(jpi,jpj,jpk,jptra,jpt) , &
& trc_i(jpi,jpj,jptra) , trc_o(jpi,jpj,jptra) , &
& gtru (jpi,jpj,jptra) , gtrv (jpi,jpj,jptra) , &
& gtrui(jpi,jpj,jptra) , gtrvi(jpi,jpj,jptra) , &
& trc_ice_ratio(jptra) , trc_ice_prescr(jptra) , cn_trc_o(jptra) , &
& neln(jpi,jpj) , heup(jpi,jpj) , heup_01(jpi,jpj) , &
& etot(jpi,jpj,jpk) , etot_ndcy(jpi,jpj,jpk) , &
& sbc_trc_b(jpi,jpj,jptra), sbc_trc(jpi,jpj,jptra) , &
& cvol(jpi,jpj,jpk) , trai(jptra) , &
& ctrcnm(jptra) , ctrcln(jptra) , ctrcun(jptra) , &
& ln_trc_ini(jptra) , &
& ln_trc_sbc(jptra) , ln_trc_cbc(jptra) , ln_trc_obc(jptra) , &
& ln_trc_ais(jptra) , &
ALLOCATE( tr(jpi,jpj,jpk,jptra,jpt) , &
& gtru (jpi,jpj,jptra) , gtrv (jpi,jpj,jptra) , &
& gtrui(jpi,jpj,jptra) , gtrvi(jpi,jpj,jptra) , &
& trc_i(A2D(0),jptra) , trc_o(A2D(0),jptra) , &
& trc_ice_ratio(jptra) , trc_ice_prescr(jptra) , cn_trc_o(jptra) , &
& neln(A2D(0)) , heup(A2D(0)) , heup_01(A2D(0)) , &
& etot(A2D(0),jpk) , etot_ndcy(A2D(0),jpk) , &
& sbc_trc_b(A2D(0),jptra), sbc_trc(A2D(0),jptra) , &
& cvol(jpi,jpj,jpk) , trai(jptra) , &
& ctrcnm(jptra) , ctrcln(jptra) , ctrcun(jptra) , &
& ln_trc_ini(jptra) , &
& ln_trc_sbc(jptra) , ln_trc_cbc(jptra) , ln_trc_obc(jptra) , &
& ln_trc_ais(jptra) , &
& STAT = ierr(1) )
!
IF( ln_bdy ) ALLOCATE( trcdta_bdy(jptra, jp_bdy) , STAT = ierr(2) )
......