Newer
Older
MODULE p4zbc
!!======================================================================
!! *** MODULE p4sbc ***
!! TOP : PISCES surface boundary conditions of external inputs of nutrients
!!======================================================================
!! History : 3.5 ! 2012-07 (O. Aumont, C. Ethe) Original code
!!----------------------------------------------------------------------
!! p4z_bc : Read and interpolate time-varying nutrients fluxes
!! p4z_bc_init : Initialization of p4z_bc
!!----------------------------------------------------------------------
USE oce_trc ! shared variables between ocean and passive tracers
USE trc ! passive tracers common variables
USE sms_pisces ! PISCES Source Minus Sink variables
USE iom ! I/O manager
USE fldread ! time interpolation
USE prtctl ! print control for debugging
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
USE trcbc
IMPLICIT NONE
PRIVATE
PUBLIC p4z_bc
PUBLIC p4z_bc_init
LOGICAL , PUBLIC :: ln_ironsed !: boolean for Fe input from sediments
LOGICAL , PUBLIC :: ln_hydrofe !: boolean for Fe input from hydrothermal vents
REAL(wp), PUBLIC :: sedfeinput !: Coastal release of Iron
REAL(wp), PUBLIC :: icefeinput !: Iron concentration in sea ice
REAL(wp), PUBLIC :: wdust !: Sinking speed of the dust
REAL(wp), PUBLIC :: mfrac !: Mineral Content of the dust
REAL(wp) :: hratio !: Fe:3He ratio assumed for vent iron supply
REAL(wp) :: distcoast !: Distance off the coast for Iron from sediments
REAL(wp), PUBLIC :: lgw_rath !: Weak ligand ratio from hydro sources
LOGICAL , PUBLIC :: ll_bc
LOGICAL , PUBLIC :: ll_dust !: boolean for dust input from the atmosphere
LOGICAL , PUBLIC :: ll_river !: boolean for river input of nutrients
LOGICAL , PUBLIC :: ll_ndepo !: boolean for atmospheric deposition of N
TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_dust ! structure of input dust
TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_hydrofe ! structure of input iron from sediment
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: dust !: dust fields
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ironsed !: Coastal supply of iron
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hydrofe !: Hydrothermal vent supply of iron
REAL(wp), PUBLIC :: sedsilfrac, sedcalfrac
!! * Substitutions
# include "do_loop_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/TOP 4.0 , NEMO Consortium (2018)
!! $Id: p4zbc.F90 10869 2019-04-15 10:34:03Z cetlod $
!! Software governed by the CeCILL license (see ./LICENSE)
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE p4z_bc( kt, Kbb, Kmm, Krhs )
!!----------------------------------------------------------------------
!! *** routine p4z_bc ***
!!
!! ** purpose : read and interpolate the external sources of nutrients
!!
!! ** method : read the files and interpolate the appropriate variables
!!
!! ** input : external netcdf files
!!
!!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt ! ocean time step
INTEGER, INTENT(in) :: Kbb, Kmm, Krhs ! time level index
!
INTEGER :: ji, jj, jk, jl
REAL(wp) :: zdep, zwflux, zironice
REAL(wp) :: zcoef, zwdust, zrivdin, zdustdep, zndep
!
CHARACTER (len=25) :: charout
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
REAL(wp), ALLOCATABLE, DIMENSION(:,: ) :: zw2d
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('p4z_bc')
! Add the external input of nutrients from dust deposition in the water column
! The inputs at surface have already been added
! ----------------------------------------------------------
IF( ll_dust ) THEN
!
CALL fld_read( kt, 1, sf_dust )
DO_2D( 0, 0, 0, 0 )
dust(ji,jj) = MAX( rtrn, sf_dust(1)%fnow(ji,jj,1) )
END_2D
!
! Iron solubilization of particles in the water column
! dust in kg/m2/s ---> 1/55.85 to put in mol/Fe ; wdust in m/d
! Dust are supposed to sink at wdust sinking speed. 3% of the iron
! in dust is hypothesized to be soluble at a dissolution rate set to
! 1/(250 days). The vertical distribution of iron in dust is computed
! from a steady state assumption. Parameters are very uncertain and
! are estimated from the literature quoted in Raiswell et al. (2011)
! -------------------------------------------------------------------
zwdust = 0.03 / ( wdust / rday ) / ( 250. * rday )
! Atmospheric input of Iron dissolves in the water column
IF ( ln_trc_sbc(jpfer) ) THEN
zdustdep = dust(ji,jj) * zwdust * rfact * EXP( -gdept(ji,jj,jk,Kmm) /( 250. * wdust ) )
!
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + zdustdep * mfrac / mMass_Fe
END_3D
IF( lk_iomput ) THEN
! surface downward dust depo of iron
zw2d(A2D(0)) = rf_trsfac(jl) * ( sf_trcsbc(jl)%fnow(A2D(0),1) / rn_sbc_time ) * 1.e+3 * tmask(A2D(0),1)
CALL iom_put( "Irondep", zw2d )
DEALLOCATE( zw2d )
ENDIF
ENDIF
! Atmospheric input of PO4 dissolves in the water column
IF ( ln_trc_sbc(jppo4) ) THEN
zdustdep = dust(ji,jj) * zwdust * rfact * EXP( -gdept(ji,jj,jk,Kmm) /( 250. * wdust ) )
!
tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) + zdustdep * 1.e-3 / mMass_P
END_3D
ENDIF
! Atmospheric input of Si dissolves in the water column
IF ( ln_trc_sbc(jpsil) ) THEN
zdustdep = dust(ji,jj) * zwdust * rfact * EXP( -gdept(ji,jj,jk,Kmm) /( 250. * wdust ) )
!
tr(ji,jj,jk,jpsil,Krhs) = tr(ji,jj,jk,jpsil,Krhs) + zdustdep * 0.269 / mMass_Si
END_3D
ENDIF
!
IF( lk_iomput ) THEN
! dust concentration at surface
ALLOCATE( zw2d(A2D(0)) )
zw2d(A2D(0)) = dust(A2D(0)) / ( wdust / rday ) * tmask(A2D(0),1)
CALL iom_put( "pdust", zw2d )
DEALLOCATE( zw2d )
ENDIF
ENDIF
! -----------------------------------------
! Add the external input of nutrients from river
! ----------------------------------------------------------
IF( ll_river ) THEN
jl = n_trc_indcbc(jpno3)
DO jk = 1, nk_rnf(ji,jj)
zcoef = rn_rfact / ( e1e2t(ji,jj) * h_rnf(ji,jj) * rn_cbc_time ) * tmask(ji,jj,1)
zrivdin = rf_trcfac(jl) * sf_trccbc(jl)%fnow(ji,jj,1) * zcoef
tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) - rno3 * zrivdin * rfact
ENDDO
END_2D
ENDIF
! Add the external input of nutrients from nitrogen deposition
! ----------------------------------------------------------
IF( ll_ndepo ) THEN
IF( ln_trc_sbc(jpno3) ) THEN
jl = n_trc_indsbc(jpno3)
zndep = rf_trsfac(jl) * sf_trcsbc(jl)%fnow(ji,jj,1) / e3t(ji,jj,1,Kmm) / rn_sbc_time
tr(ji,jj,1,jptal,Krhs) = tr(ji,jj,1,jptal,Krhs) - rno3 * zndep * rfact
END_2D
ENDIF
IF( ln_trc_sbc(jpnh4) ) THEN
jl = n_trc_indsbc(jpnh4)
zndep = rf_trsfac(jl) * sf_trcsbc(jl)%fnow(ji,jj,1) / e3t(ji,jj,1,Kmm) / rn_sbc_time
tr(ji,jj,1,jptal,Krhs) = tr(ji,jj,1,jptal,Krhs) + rno3 * zndep * rfact
END_2D
ENDIF
ENDIF
!
! Iron input/uptake due to sea ice : Crude parameterization based on
! Lancelot et al. Iron concentration in sea-ice is constant and set
! in the namelist_pisces (icefeinput). ln_ironice is forced to false
! when nn_ice_tr = 1
! ----------------------------------------------------
IF( ln_ironice ) THEN
!
! Compute the iron flux between sea ice and sea water
! Simple parameterization assuming a fixed constant concentration in
! sea-ice (icefeinput)
! ------------------------------------------------------------------
zdep = rfact / e3t(ji,jj,1,Kmm)
zwflux = fmmflx(ji,jj) / 1000._wp
zironice = MAX( -0.99 * tr(ji,jj,1,jpfer,Kbb), -zwflux * icefeinput * zdep )
tr(ji,jj,1,jpfer,Krhs) = tr(ji,jj,1,jpfer,Krhs) + zironice
END_2D
!
IF( lk_iomput ) THEN
! iron flux from ice
ALLOCATE( zw2d(A2D(0)) )
zw2d(A2D(0)) = MAX( -0.99 * tr(A2D(0),1,jpfer,Kbb), (-1.*fmmflx(A2D(0))/1000._wp )*icefeinput*1.e+3*tmask(A2D(0),1))
CALL iom_put( "Ironice", zw2d )
DEALLOCATE( zw2d )
ENDIF
!
ENDIF
! Add the external input of iron from sediment mobilization
! ------------------------------------------------------
IF( ln_ironsed .AND. .NOT.lk_sed ) THEN
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + ironsed(ji,jj,jk) * rfact
END_3D
!
IF( lk_iomput ) THEN
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = ironsed(A2D(0),1:jpkm1) * 1.e+3 * tmask(A2D(0),1:jpkm1)
CALL iom_put( "Ironsed", zw3d )
DEALLOCATE( zw3d )
ENDIF
ENDIF
! Add the external input of iron from hydrothermal vents
! ------------------------------------------------------
IF( ln_hydrofe ) THEN
CALL fld_read( kt, 1, sf_hydrofe )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
hydrofe(ji,jj,jk) = ( MAX( rtrn, sf_hydrofe(1)%fnow(ji,jj,jk) ) * hratio ) &
& / ( e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * ryyss + rtrn ) / 1000._wp &
& * tmask(ji,jj,jk)
END_3D
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + hydrofe(ji,jj,jk) * rfact
END_3D
IF( ln_ligand ) THEN
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
tr(ji,jj,jk,jplgw,Krhs) = tr(ji,jj,jk,jplgw,Krhs) + ( hydrofe(ji,jj,jk) * lgw_rath ) * rfact
END_3D
ENDIF
IF( lk_iomput ) THEN
! hydrothermal iron input
ALLOCATE( zw3d(A2D(0),jpk) ) ; zw3d(A2D(0),jpk) = 0._wp
zw3d(A2D(0),1:jpkm1) = hydrofe(A2D(0),1:jpkm1) * 1.e+3 * tmask(A2D(0),1:jpkm1)
CALL iom_put( "HYDR", zw3d )
DEALLOCATE( zw3d )
ENDIF
!
IF(sn_cfctl%l_prttrc) THEN ! print mean trends (used for debugging)
WRITE(charout, FMT="('bc')")
CALL prt_ctl_info( charout, cdcomp = 'top' )
CALL prt_ctl(tab4d_1=tr(:,:,:,:,Krhs), mask1=tmask, clinfo=ctrcnm)
ENDIF
!
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
IF( ln_timing ) CALL timing_stop('p4z_bc')
!
END SUBROUTINE p4z_bc
SUBROUTINE p4z_bc_init( Kmm )
!!----------------------------------------------------------------------
!! *** routine p4z_bc_init ***
!!
!! ** purpose : initialization of the external sources of nutrients
!!
!! ** method : read the files and compute the budget
!! called at the first timestep (nittrc000)
!!
!! ** input : external netcdf files
!!
!!----------------------------------------------------------------------
INTEGER, INTENT( in ) :: Kmm ! time level index
INTEGER :: ji, jj, jk, jm
INTEGER :: ii0, ii1, ij0, ij1
INTEGER :: numiron
INTEGER :: ierr, ierr1, ierr2, ierr3
INTEGER :: ios ! Local integer output status for namelist read
INTEGER :: ik50 ! last level where depth less than 50 m
REAL(wp) :: zexpide, zdenitide, zmaskt, zsurfc, zsurfp,ze3t, ze3t2, zcslp
REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zriver, zcmask
!
CHARACTER(len=100) :: cn_dir ! Root directory for location of ssr files
TYPE(FLD_N) :: sn_dust, sn_ironsed, sn_hydrofe ! informations about the fields to be read
!!
NAMELIST/nampisbc/cn_dir, sn_dust, sn_ironsed, sn_hydrofe, &
& ln_ironsed, ln_ironice, ln_hydrofe, &
& sedfeinput, distcoast, icefeinput, wdust, mfrac, &
& hratio, lgw_rath
!!----------------------------------------------------------------------
!
IF(lwp) THEN
WRITE(numout,*)
WRITE(numout,*) 'p4z_bc_init : initialization of the external sources of nutrients '
WRITE(numout,*) '~~~~~~~~~~~~ '
ENDIF
! !* set file information
READ ( numnatp_ref, nampisbc, IOSTAT = ios, ERR = 901)
901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisbc in reference namelist' )
READ ( numnatp_cfg, nampisbc, IOSTAT = ios, ERR = 902 )
902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nampisbc in configuration namelist' )
IF(lwm) WRITE ( numonp, nampisbc )
IF(lwp) THEN
WRITE(numout,*) ' Namelist : nampissbc '
WRITE(numout,*) ' Fe input from sediments ln_ironsed = ', ln_ironsed
WRITE(numout,*) ' Fe input from seaice ln_ironice = ', ln_ironice
WRITE(numout,*) ' fe input from hydrothermal vents ln_hydrofe = ', ln_hydrofe
IF( ln_ironsed ) THEN
WRITE(numout,*) ' coastal release of iron sedfeinput = ', sedfeinput
WRITE(numout,*) ' distance off the coast distcoast = ', distcoast
ENDIF
IF( ln_ligand ) THEN
WRITE(numout,*) ' Weak ligand ratio from sed hydro sources lgw_rath = ', lgw_rath
ENDIF
IF( ln_ironice ) THEN
WRITE(numout,*) ' Iron concentration in sea ice icefeinput = ', icefeinput
ENDIF
IF( ln_trc_sbc(jpfer) ) THEN
WRITE(numout,*) ' Mineral Fe content of the dust mfrac = ', mfrac
WRITE(numout,*) ' sinking speed of the dust wdust = ', wdust
ENDIF
IF( ln_hydrofe ) THEN
WRITE(numout,*) ' Fe to 3He ratio assumed for vent iron supply hratio = ', hratio
ENDIF
END IF
ll_bc = ( ln_trcbc .AND. lltrcbc ) .OR. ln_hydrofe .OR. ln_ironsed .OR. ln_ironice
ll_dust = ln_trc_sbc(jpfer) .OR. ln_trc_sbc(jppo4) .OR. ln_trc_sbc(jpsil) .OR. ln_sediment
ll_ndepo = ln_trc_sbc(jpno3) .OR. ln_trc_sbc(jpnh4)
ll_river = ln_trc_cbc(jpno3)
! dust input from the atmosphere
! ------------------------------
IF( ll_dust ) THEN
!
IF(lwp) WRITE(numout,*) ' initialize dust input from atmosphere '
IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
!
!
ALLOCATE( sf_dust(1), STAT=ierr ) !* allocate and fill sf_sst (forcing structure) with sn_sst
IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'p4z_bc_init: unable to allocate sf_dust structure' )
!
CALL fld_fill( sf_dust, (/ sn_dust /), cn_dir, 'p4z_bc_init', 'Atmospheric dust deposition', 'nampisbc' )
ALLOCATE( sf_dust(1)%fnow(jpi,jpj,1) )
IF( sn_dust%ln_tint ) ALLOCATE( sf_dust(1)%fdta(jpi,jpj,1,2) )
!
END IF
! coastal and island masks
! ------------------------
IF( ln_ironsed ) THEN
!
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' ==>>> ln_ironsed=T , computation of an island mask to enhance coastal supply of iron'
!
ALLOCATE( ironsed(A2D(0),jpk) ) ! allocation
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
!
CALL iom_open ( TRIM( sn_ironsed%clname ), numiron )
ALLOCATE( zcmask(jpi,jpj,jpk) )
CALL iom_get ( numiron, jpdom_global, TRIM( sn_ironsed%clvar ), zcmask(:,:,:), 1 )
CALL iom_close( numiron )
!
ik50 = 5 ! last level where depth less than 50 m
DO jk = jpkm1, 1, -1
IF( gdept_1d(jk) > 50. ) ik50 = jk - 1
END DO
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' Level corresponding to 50m depth ', ik50,' ', gdept_1d(ik50+1)
DO_3D( 0, 0, 0, 0, 1, ik50 )
ze3t = e3t_0(ji,jj,jk)
zsurfc = e1u(ji,jj) * ( 1. - umask(ji ,jj ,jk) ) &
+ e1u(ji,jj) * ( 1. - umask(ji-1,jj ,jk) ) &
+ e2v(ji,jj) * ( 1. - vmask(ji ,jj ,jk) ) &
+ e2v(ji,jj) * ( 1. - vmask(ji ,jj-1,jk) )
zsurfp = zsurfc * ze3t / e1e2t(ji,jj)
! estimation of the coastal slope : 5 km off the coast
ze3t2 = ze3t * ze3t
zcslp = SQRT( ( distcoast*distcoast + ze3t2 ) / ze3t2 )
!
zcmask(ji,jj,jk) = zcmask(ji,jj,jk) + zcslp * zsurfp
END_3D
!
CALL lbc_lnk( 'p4zbc', zcmask , 'T', 1.0_wp ) ! lateral boundary conditions on cmask (sign unchanged)
!
zexpide = MIN( 8.,( gdept(ji,jj,jk,Kmm) / 500. )**(-1.5) )
zdenitide = -0.9543 + 0.7662 * LOG( zexpide ) - 0.235 * LOG( zexpide )**2
zcmask(ji,jj,jk) = zcmask(ji,jj,jk) * MIN( 1., EXP( zdenitide ) / 0.5 )
END_3D
! Coastal supply of iron
! -------------------------
ironsed(:,:,jpk) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
ironsed(ji,jj,jk) = sedfeinput * zcmask(ji,jj,jk) / ( e3t_0(ji,jj,jk) * rday )
END_3D
DEALLOCATE( zcmask)
ENDIF
!
! Iron from Hydrothermal vents
! ------------------------
IF( ln_hydrofe ) THEN
!
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' ==>>> ln_hydrofe=T , Input of iron from hydrothermal vents'
!
ALLOCATE( hydrofe(A2D(0),jpk) ) ! allocation
!
ALLOCATE( sf_hydrofe(1), STAT=ierr ) !* allocate and fill sf_sst (forcing structure) with sn_sst
IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'p4z_bc_init: unable to allocate sf_hydro structure' )
!
CALL fld_fill( sf_hydrofe, (/ sn_hydrofe /), cn_dir, 'p4z_bc_init', 'Input of iron from hydrothermal vents', 'nampisbc' )
ALLOCATE( sf_hydrofe(1)%fnow(jpi,jpj,jpk) )
IF( sn_hydrofe%ln_tint ) ALLOCATE( sf_hydrofe(1)%fdta(jpi,jpj,jpk,2) )
!
ENDIF
!
END SUBROUTINE p4z_bc_init
!!======================================================================
END MODULE p4zbc