New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

p4zsed.F90 in branches/2011/dev_NEMO_MERGE_2011/NEMOGCM/NEMO/TOP_SRC/PISCES – NEMO

Context Navigation

source: branches/2011/dev_NEMO_MERGE_2011/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsed.F90 @ 3124

Last change on this file since 3124 was 3124, checked in by cetlod, 13 years ago
dev_NEMO_MERGE_2011/NEMOGCM:minor modifications on the use of nittrc000 + style corrections
Property svn:keywords set to `Id`
File size: 29.8 KB

Rev	Line
[935]	1	MODULE p4zsed
	2	!!======================================================================
	3	!! * MODULE p4sed *
	4	!! TOP : PISCES Compute loss of organic matter in the sediments
	5	!!======================================================================
	6	!! History : 1.0 ! 2004-03 (O. Aumont) Original code
	7	!! 2.0 ! 2007-12 (C. Ethe, G. Madec) F90
[2977]	8	!! 3.4 ! 2011-06 (O. Aumont, C. Ethe) USE of fldread
[935]	9	!!----------------------------------------------------------------------
	10	#if defined key_pisces
	11	!!----------------------------------------------------------------------
	12	!! 'key_pisces' PISCES bio-model
	13	!!----------------------------------------------------------------------
	14	!! p4z_sed : Compute loss of organic matter in the sediments
	15	!! p4z_sbc : Read and interpolate time-varying nutrients fluxes
	16	!! p4z_sed_init : Initialization of p4z_sed
	17	!!----------------------------------------------------------------------
[2977]	18	USE oce_trc ! shared variables between ocean and passive tracers
	19	USE trc ! passive tracers common variables
	20	USE sms_pisces ! PISCES Source Minus Sink variables
	21	USE p4zsink ! vertical flux of particulate matter due to sinking
	22	USE p4zopt ! optical model
	23	USE p4zlim ! Co-limitations of differents nutrients
	24	USE p4zrem ! Remineralisation of organic matter
	25	USE p4zint ! interpolation and computation of various fields
	26	USE iom ! I/O manager
	27	USE fldread ! time interpolation
	28	USE prtctl_trc ! print control for debugging
[935]	29
	30	IMPLICIT NONE
	31	PRIVATE
	32
[1073]	33	PUBLIC p4z_sed
[2528]	34	PUBLIC p4z_sed_init
[2715]	35	PUBLIC p4z_sed_alloc
[935]	36
	37	!! * Shared module variables
[2977]	38	LOGICAL :: ln_dust = .FALSE. !: boolean for dust input from the atmosphere
	39	LOGICAL :: ln_river = .FALSE. !: boolean for river input of nutrients
	40	LOGICAL :: ln_ndepo = .FALSE. !: boolean for atmospheric deposition of N
	41	LOGICAL :: ln_ironsed = .FALSE. !: boolean for Fe input from sediments
[935]	42
[2977]	43	REAL(wp) :: sedfeinput = 1.E-9_wp !: Coastal release of Iron
	44	REAL(wp) :: dustsolub = 0.014_wp !: Solubility of the dust
	45	REAL(wp) :: wdust = 2.0_wp !: Sinking speed of the dust
	46	REAL(wp) :: nitrfix = 1E-7_wp !: Nitrogen fixation rate
	47	REAL(wp) :: diazolight = 50._wp !: Nitrogen fixation sensitivty to light
	48	REAL(wp) :: concfediaz = 1.E-10_wp !: Fe half-saturation Cste for diazotrophs
[935]	49
[2977]	50
[935]	51	!! * Module variables
[2715]	52	REAL(wp) :: ryyss !: number of seconds per year
[2977]	53	REAL(wp) :: r1_ryyss !: inverse of ryyss
[2715]	54	REAL(wp) :: rmtss !: number of seconds per month
[2977]	55	REAL(wp) :: r1_rday !: inverse of rday
	56	LOGICAL :: ll_sbc
[1735]	57
[2977]	58	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_dust ! structure of input dust
	59	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_riverdic ! structure of input riverdic
	60	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_riverdoc ! structure of input riverdoc
	61	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_ndepo ! structure of input nitrogen deposition
	62	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_ironsed ! structure of input iron from sediment
[2715]	63
[2977]	64	INTEGER , PARAMETER :: nbtimes = 365 !: maximum number of times record in a file
	65	INTEGER :: ntimes_dust, ntimes_riv, ntimes_ndep ! number of time steps in a file
[2715]	66
	67	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: dust !: dust fields
	68	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: rivinp, cotdep !: river input fields
	69	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: nitdep !: atmospheric N deposition
	70	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ironsed !: Coastal supply of iron
	71
[935]	72	REAL(wp) :: sumdepsi, rivalkinput, rivpo4input, nitdepinput
	73
	74	!!* Substitution
[1503]	75	# include "top_substitute.h90"
[935]	76	!!----------------------------------------------------------------------
[2528]	77	!! NEMO/TOP 3.3 , NEMO Consortium (2010)
[1180]	78	!! $Header:$
[2528]	79	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[935]	80	!!----------------------------------------------------------------------
	81
	82	CONTAINS
	83
[2715]	84
[2528]	85	SUBROUTINE p4z_sed( kt, jnt )
[935]	86	!!---------------------------------------------------------------------
	87	!! * ROUTINE p4z_sed *
	88	!!
	89	!! ** Purpose : Compute loss of organic matter in the sediments. This
	90	!! is by no way a sediment model. The loss is simply
	91	!! computed to balance the inout from rivers and dust
	92	!!
	93	!! ** Method : - ???
	94	!!---------------------------------------------------------------------
[2977]	95	USE wrk_nemo, ONLY: wrk_in_USE, wrk_not_released
	96	USE wrk_nemo, ONLY: zsidep => wrk_2d_11
	97	USE wrk_nemo, ONLY: zwork1 => wrk_2d_12, zwork2 => wrk_2d_13, zwork3 => wrk_2d_14
[2715]	98	USE wrk_nemo, ONLY: znitrpot => wrk_3d_2, zirondep => wrk_3d_3
	99	!
[935]	100	INTEGER, INTENT(in) :: kt, jnt ! ocean time step
[2528]	101	INTEGER :: ji, jj, jk, ikt
[1180]	102	#if ! defined key_sed
[935]	103	REAL(wp) :: zsumsedsi, zsumsedpo4, zsumsedcal
[2528]	104	REAL(wp) :: zrivalk, zrivsil, zrivpo4
[1180]	105	#endif
[2977]	106	REAL(wp) :: zdenitot, znitrpottot, zlim, zfact, zfactcal
	107	REAL(wp) :: zsiloss, zcaloss, zwsbio3, zwsbio4, zwscal, zdep
[935]	108	CHARACTER (len=25) :: charout
	109	!!---------------------------------------------------------------------
	110
[2977]	111	IF( ( wrk_in_USE(2, 11,12,13,14) ) .OR. ( wrk_in_USE(3, 2,3) ) ) THEN
[2715]	112	CALL ctl_stop('p4z_sed: requested workspace arrays unavailable') ; RETURN
	113	END IF
	114
[2977]	115	IF( jnt == 1 .AND. ll_sbc ) CALL p4z_sbc( kt )
[935]	116
[2977]	117	zirondep(:,:,:) = 0.e0 ! Initialisation of variables USEd to compute deposition
	118	zsidep (:,:) = 0.e0
	119
[935]	120	! Iron and Si deposition at the surface
	121	! -------------------------------------
	122	DO jj = 1, jpj
	123	DO ji = 1, jpi
[2977]	124	zdep = rfact2 / fse3t(ji,jj,1)
	125	zirondep(ji,jj,1) = ( dustsolub * dust(ji,jj) / ( 55.85 * rmtss ) + 3.e-10 * r1_ryyss ) * zdep
	126	zsidep (ji,jj) = 8.8 * 0.075 * dust(ji,jj) * zdep / ( 28.1 * rmtss )
[935]	127	END DO
	128	END DO
	129
	130	! Iron solubilization of particles in the water column
	131	! ----------------------------------------------------
	132	DO jk = 2, jpkm1
[2977]	133	zirondep(:,:,jk) = dust(:,:) / ( wdust * 55.85 * rmtss ) * rfact2 * 1.e-4 * EXP( -fsdept(:,:,jk) / 1000. )
[935]	134	END DO
	135
	136	! Add the external input of nutrients, carbon and alkalinity
	137	! ----------------------------------------------------------
	138	trn(:,:,1,jppo4) = trn(:,:,1,jppo4) + rivinp(:,:) * rfact2
	139	trn(:,:,1,jpno3) = trn(:,:,1,jpno3) + (rivinp(:,:) + nitdep(:,:)) * rfact2
	140	trn(:,:,1,jpfer) = trn(:,:,1,jpfer) + rivinp(:,:) * 3.e-5 * rfact2
	141	trn(:,:,1,jpsil) = trn(:,:,1,jpsil) + zsidep (:,:) + cotdep(:,:) * rfact2 / 6.
	142	trn(:,:,1,jpdic) = trn(:,:,1,jpdic) + rivinp(:,:) * 2.631 * rfact2
	143	trn(:,:,1,jptal) = trn(:,:,1,jptal) + (cotdep(:,:) - rno3(rivinp(:,:) + nitdep(:,:) ) ) rfact2
	144
	145
	146	! Add the external input of iron which is 3D distributed
	147	! (dust, river and sediment mobilization)
	148	! ------------------------------------------------------
	149	DO jk = 1, jpkm1
[1457]	150	trn(:,:,jk,jpfer) = trn(:,:,jk,jpfer) + zirondep(:,:,jk) + ironsed(:,:,jk) * rfact2
[935]	151	END DO
	152
[1180]	153	#if ! defined key_sed
[935]	154	! Loss of biogenic silicon, Caco3 organic carbon in the sediments.
	155	! First, the total loss is computed.
	156	! The factor for calcite comes from the alkalinity effect
	157	! -------------------------------------------------------------
	158	DO jj = 1, jpj
	159	DO ji = 1, jpi
[2528]	160	ikt = mbkt(ji,jj)
[935]	161	# if defined key_kriest
[2977]	162	zwork1(ji,jj) = trn(ji,jj,ikt,jpdsi) * wscal (ji,jj,ikt)
	163	zwork2(ji,jj) = trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt)
[935]	164	# else
[2977]	165	zwork1(ji,jj) = trn(ji,jj,ikt,jpdsi) * wsbio4(ji,jj,ikt)
	166	zwork2(ji,jj) = trn(ji,jj,ikt,jpgoc) * wsbio4(ji,jj,ikt) + trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt)
[935]	167	# endif
[2977]	168	! For calcite, burial efficiency is made a function of saturation
	169	zfactcal = MIN( excess(ji,jj,ikt), 0.2 )
	170	zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) )
	171	zwork3(ji,jj) = trn(ji,jj,ikt,jpcal) * wscal (ji,jj,ikt) * 2.e0 * zfactcal
[935]	172	END DO
	173	END DO
[2977]	174	zsumsedsi = glob_sum( zwork1(:,:) * e1e2t(:,:) ) * r1_rday
	175	zsumsedpo4 = glob_sum( zwork2(:,:) * e1e2t(:,:) ) * r1_rday
	176	zsumsedcal = glob_sum( zwork3(:,:) * e1e2t(:,:) ) * r1_rday
[1180]	177	#endif
	178
[2977]	179	! THEN this loss is scaled at each bottom grid cell for
[935]	180	! equilibrating the total budget of silica in the ocean.
	181	! Thus, the amount of silica lost in the sediments equal
	182	! the supply at the surface (dust+rivers)
	183	! ------------------------------------------------------
[2977]	184	#if ! defined key_sed
	185	zrivsil = 1._wp - ( sumdepsi + rivalkinput * r1_ryyss / 6. ) / zsumsedsi
	186	zrivpo4 = 1._wp - ( rivpo4input * r1_ryyss ) / zsumsedpo4
	187	#endif
[935]	188
	189	DO jj = 1, jpj
	190	DO ji = 1, jpi
[2977]	191	ikt = mbkt(ji,jj)
	192	zdep = xstep / fse3t(ji,jj,ikt)
	193	zwsbio4 = wsbio4(ji,jj,ikt) * zdep
	194	zwscal = wscal (ji,jj,ikt) * zdep
[2528]	195	# if defined key_kriest
[2977]	196	zsiloss = trn(ji,jj,ikt,jpdsi) * zwsbio4
[935]	197	# else
[2977]	198	zsiloss = trn(ji,jj,ikt,jpdsi) * zwscal
[935]	199	# endif
[2977]	200	zcaloss = trn(ji,jj,ikt,jpcal) * zwscal
	201	!
	202	trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) - zsiloss
	203	trn(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) - zcaloss
	204	#if ! defined key_sed
	205	trn(ji,jj,ikt,jpsil) = trn(ji,jj,ikt,jpsil) + zsiloss * zrivsil
	206	zfactcal = MIN( excess(ji,jj,ikt), 0.2 )
	207	zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) )
	208	zrivalk = 1._wp - ( rivalkinput * r1_ryyss ) * zfactcal / zsumsedcal
	209	trn(ji,jj,ikt,jptal) = trn(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0
	210	trn(ji,jj,ikt,jpdic) = trn(ji,jj,ikt,jpdic) + zcaloss * zrivalk
	211	#endif
[935]	212	END DO
	213	END DO
	214
	215	DO jj = 1, jpj
	216	DO ji = 1, jpi
[2977]	217	ikt = mbkt(ji,jj)
	218	zdep = xstep / fse3t(ji,jj,ikt)
	219	zwsbio4 = wsbio4(ji,jj,ikt) * zdep
	220	zwsbio3 = wsbio3(ji,jj,ikt) * zdep
	221	# if ! defined key_kriest
	222	trn(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) - trn(ji,jj,ikt,jpgoc) * zwsbio4
	223	trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - trn(ji,jj,ikt,jppoc) * zwsbio3
	224	trn(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) - trn(ji,jj,ikt,jpbfe) * zwsbio4
	225	trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * zwsbio3
	226	#if ! defined key_sed
	227	trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc) &
	228	& + ( trn(ji,jj,ikt,jpgoc) * zwsbio4 + trn(ji,jj,ikt,jppoc) * zwsbio3 ) * zrivpo4
	229	#endif
	230
[935]	231	# else
[2977]	232	trn(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum) - trn(ji,jj,ikt,jpnum) * zwsbio4
	233	trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - trn(ji,jj,ikt,jppoc) * zwsbio3
	234	trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * zwsbio3
	235	#if ! defined key_sed
	236	trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc) &
	237	& + ( trn(ji,jj,ikt,jpnum) * zwsbio4 + trn(ji,jj,ikt,jppoc) * zwsbio3 ) * zrivpo4
	238	#endif
	239
[935]	240	# endif
	241	END DO
	242	END DO
	243
[2977]	244
[935]	245	! Nitrogen fixation (simple parameterization). The total gain
	246	! from nitrogen fixation is scaled to balance the loss by
	247	! denitrification
	248	! -------------------------------------------------------------
	249
[2977]	250	zdenitot = glob_sum( ( denitr(:,:,:) * rdenit + denitnh4(:,:,:) * rdenita ) * cvol(:,:,:) )
[935]	251
[1678]	252	! Potential nitrogen fixation dependant on temperature and iron
[935]	253	! -------------------------------------------------------------
	254
	255	!CDIR NOVERRCHK
	256	DO jk = 1, jpk
	257	!CDIR NOVERRCHK
	258	DO jj = 1, jpj
	259	!CDIR NOVERRCHK
	260	DO ji = 1, jpi
	261	zlim = ( 1.- xnanono3(ji,jj,jk) - xnanonh4(ji,jj,jk) )
	262	IF( zlim <= 0.2 ) zlim = 0.01
[2977]	263	#if defined key_degrad
	264	zfact = zlim * rfact2 * facvol(ji,jj,jk)
	265	#else
	266	zfact = zlim * rfact2
	267	#endif
	268	znitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday ) &
	269	& * zfact * trn(ji,jj,jk,jpfer) / ( concfediaz + trn(ji,jj,jk,jpfer) ) &
	270	& * ( 1.- EXP( -etot(ji,jj,jk) / diazolight ) )
[935]	271	END DO
	272	END DO
	273	END DO
	274
[2528]	275	znitrpottot = glob_sum( znitrpot(:,:,:) * cvol(:,:,:) )
[935]	276
	277	! Nitrogen change due to nitrogen fixation
	278	! ----------------------------------------
	279	DO jk = 1, jpk
	280	DO jj = 1, jpj
	281	DO ji = 1, jpi
[2977]	282	zfact = znitrpot(ji,jj,jk) * nitrfix
[935]	283	trn(ji,jj,jk,jpnh4) = trn(ji,jj,jk,jpnh4) + zfact
[2977]	284	trn(ji,jj,jk,jptal) = trn(ji,jj,jk,jptal) + rno3 * zfact
[935]	285	trn(ji,jj,jk,jpoxy) = trn(ji,jj,jk,jpoxy) + zfact * o2nit
[2977]	286	trn(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) + 30. / 46. * zfact
	287	! trn(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) + zfact
	288	END DO
	289	END DO
[935]	290	END DO
[2977]	291	!
	292	IF( ln_diatrc ) THEN
	293	zfact = 1.e+3 * rfact2r
	294	IF( lk_iomput ) THEN
	295	zwork1(:,:) = ( zirondep(:,:,1) + ironsed(:,:,1) * rfact2 ) * zfact * fse3t(:,:,1) * tmask(:,:,1)
	296	zwork2(:,:) = znitrpot(:,:,1) * nitrfix * zfact * fse3t(:,:,1) * tmask(:,:,1)
	297	IF( jnt == nrdttrc ) THEN
	298	CALL iom_put( "Irondep", zwork1 ) ! surface downward net flux of iron
	299	CALL iom_put( "Nfix" , zwork2 ) ! nitrogen fixation at surface
	300	ENDIF
	301	ELSE
	302	trc2d(:,:,jp_pcs0_2d + 11) = zirondep(:,:,1) * zfact * fse3t(:,:,1) * tmask(:,:,1)
	303	trc2d(:,:,jp_pcs0_2d + 12) = znitrpot(:,:,1) * nitrfix * zfact * fse3t(:,:,1) * tmask(:,:,1)
	304	ENDIF
[2528]	305	ENDIF
[935]	306	!
[2977]	307	IF(ln_ctl) THEN ! print mean trends (USEd for debugging)
	308	WRITE(charout, fmt="('sed ')")
[935]	309	CALL prt_ctl_trc_info(charout)
	310	CALL prt_ctl_trc(tab4d=trn, mask=tmask, clinfo=ctrcnm)
[2977]	311	ENDIF
[935]	312
[2977]	313	IF( ( wrk_not_released(2, 11,12,13,14) ) .OR. ( wrk_not_released(3, 2,3) ) ) &
[2715]	314	& CALL ctl_stop('p4z_sed: failed to release workspace arrays')
	315
[935]	316	END SUBROUTINE p4z_sed
	317
[2528]	318	SUBROUTINE p4z_sbc( kt )
[935]	319	!!----------------------------------------------------------------------
[2977]	320	!! * routine p4z_sbc *
[935]	321	!!
[2977]	322	!! ** purpose : read and interpolate the external sources of
[935]	323	!! nutrients
	324	!!
[2977]	325	!! ** method : read the files and interpolate the appropriate variables
[935]	326	!!
	327	!! ** input : external netcdf files
	328	!!
	329	!!----------------------------------------------------------------------
	330	!! * arguments
	331	INTEGER, INTENT( in ) :: kt ! ocean time step
	332
[2977]	333	!! * local declarations
	334	INTEGER :: ji,jj
	335	REAL(wp) :: zcoef
[935]	336
	337	!!---------------------------------------------------------------------
	338
[3124]	339	! Compute dust at nittrc000 or only if there is more than 1 time record in dust file
[2977]	340	IF( ln_dust ) THEN
[3124]	341	IF( kt == nittrc000 .OR. ( kt /= nittrc000 .AND. ntimes_dust > 1 ) ) THEN
[2977]	342	CALL fld_read( kt, 1, sf_dust )
	343	dust(:,:) = sf_dust(1)%fnow(:,:,1)
	344	ENDIF
	345	ENDIF
[935]	346
[2977]	347	! N/P and Si releases due to coastal rivers
[3124]	348	! Compute river at nittrc000 or only if there is more than 1 time record in river file
[2977]	349	! -----------------------------------------
	350	IF( ln_river ) THEN
[3124]	351	IF( kt == nittrc000 .OR. ( kt /= nittrc000 .AND. ntimes_riv > 1 ) ) THEN
[2977]	352	CALL fld_read( kt, 1, sf_riverdic )
	353	CALL fld_read( kt, 1, sf_riverdoc )
	354	DO jj = 1, jpj
	355	DO ji = 1, jpi
	356	zcoef = ryyss * cvol(ji,jj,1)
	357	cotdep(ji,jj) = sf_riverdic(1)%fnow(ji,jj,1) * 1E9 / ( 12. * zcoef + rtrn )
	358	rivinp(ji,jj) = ( sf_riverdic(1)%fnow(ji,jj,1) + sf_riverdoc(1)%fnow(ji,jj,1) ) * 1E9 / ( 31.6* zcoef + rtrn )
	359	END DO
	360	END DO
	361	ENDIF
	362	ENDIF
[935]	363
[3124]	364	! Compute N deposition at nittrc000 or only if there is more than 1 time record in N deposition file
[2977]	365	IF( ln_ndepo ) THEN
[3124]	366	IF( kt == nittrc000 .OR. ( kt /= nittrc000 .AND. ntimes_ndep > 1 ) ) THEN
[2977]	367	CALL fld_read( kt, 1, sf_ndepo )
	368	DO jj = 1, jpj
	369	DO ji = 1, jpi
	370	nitdep(ji,jj) = 7.6 * sf_ndepo(1)%fnow(ji,jj,1) / ( 14E6 * ryyss * fse3t(ji,jj,1) + rtrn )
	371	END DO
	372	END DO
	373	ENDIF
[935]	374	ENDIF
[2977]	375	!
[935]	376	END SUBROUTINE p4z_sbc
	377
	378	SUBROUTINE p4z_sed_init
	379
	380	!!----------------------------------------------------------------------
[2977]	381	!! * routine p4z_sed_init *
[935]	382	!!
[2977]	383	!! ** purpose : initialization of the external sources of nutrients
[935]	384	!!
[2977]	385	!! ** method : read the files and compute the budget
[3124]	386	!! called at the first timestep (nittrc000)
[935]	387	!!
	388	!! ** input : external netcdf files
	389	!!
	390	!!----------------------------------------------------------------------
[2774]	391	!
[2977]	392	INTEGER :: ji, jj, jk, jm
	393	INTEGER :: numdust, numriv, numiron, numdepo
	394	INTEGER :: ierr, ierr1, ierr2, ierr3
	395	REAL(wp) :: zexpide, zdenitide, zmaskt
	396	REAL(wp), DIMENSION(nbtimes) :: zsteps ! times records
	397	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zdust, zndepo, zriverdic, zriverdoc, zcmask
[2774]	398	!
[2977]	399	CHARACTER(len=100) :: cn_dir ! Root directory for location of ssr files
	400	TYPE(FLD_N) :: sn_dust, sn_riverdoc, sn_riverdic, sn_ndepo, sn_ironsed ! informations about the fields to be read
	401	NAMELIST/nampissed/cn_dir, sn_dust, sn_riverdic, sn_riverdoc, sn_ndepo, sn_ironsed, &
	402	& ln_dust, ln_river, ln_ndepo, ln_ironsed, &
	403	& sedfeinput, dustsolub, wdust, nitrfix, diazolight, concfediaz
[2774]	404	!!----------------------------------------------------------------------
[2977]	405	! ! number of seconds per year and per month
	406	ryyss = nyear_len(1) * rday
	407	rmtss = ryyss / raamo
	408	r1_rday = 1. / rday
	409	r1_ryyss = 1. / ryyss
	410	! !* set file information
	411	cn_dir = './' ! directory in which the model is executed
	412	! ... default values (NB: frequency positive => hours, negative => months)
	413	! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation !
	414	! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs !
	415	sn_dust = FLD_N( 'dust' , -1 , 'dust' , .true. , .true. , 'yearly' , '' , '' )
	416	sn_riverdic = FLD_N( 'river' , -12 , 'riverdic' , .false. , .true. , 'yearly' , '' , '' )
	417	sn_riverdoc = FLD_N( 'river' , -12 , 'riverdoc' , .false. , .true. , 'yearly' , '' , '' )
	418	sn_ndepo = FLD_N( 'ndeposition', -12 , 'ndep' , .false. , .true. , 'yearly' , '' , '' )
	419	sn_ironsed = FLD_N( 'ironsed' , -12 , 'bathy' , .false. , .true. , 'yearly' , '' , '' )
[935]	420
[2977]	421	REWIND( numnatp ) ! read numnatp
	422	READ ( numnatp, nampissed )
[935]	423
	424	IF(lwp) THEN
	425	WRITE(numout,*) ' '
[2977]	426	WRITE(numout,*) ' namelist : nampissed '
[935]	427	WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~ '
[2977]	428	WRITE(numout,*) ' dust input from the atmosphere ln_dust = ', ln_dust
	429	WRITE(numout,*) ' river input of nutrients ln_river = ', ln_river
	430	WRITE(numout,*) ' atmospheric deposition of n ln_ndepo = ', ln_ndepo
	431	WRITE(numout,*) ' fe input from sediments ln_sedinput = ', ln_ironsed
	432	WRITE(numout,*) ' coastal release of iron sedfeinput = ', sedfeinput
	433	WRITE(numout,*) ' solubility of the dust dustsolub = ', dustsolub
	434	WRITE(numout,*) ' sinking speed of the dust wdust = ', wdust
	435	WRITE(numout,*) ' nitrogen fixation rate nitrfix = ', nitrfix
	436	WRITE(numout,*) ' nitrogen fixation sensitivty to light diazolight = ', diazolight
	437	WRITE(numout,*) ' fe half-saturation cste for diazotrophs concfediaz = ', concfediaz
	438	END IF
	439
	440	IF( ln_dust .OR. ln_river .OR. ln_ndepo ) THEN
	441	ll_sbc = .TRUE.
	442	ELSE
	443	ll_sbc = .FALSE.
[935]	444	ENDIF
	445
[2977]	446	! dust input from the atmosphere
[935]	447	! ------------------------------
[2977]	448	IF( ln_dust ) THEN
	449	IF(lwp) WRITE(numout,*) ' initialize dust input from atmosphere '
[935]	450	IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
[2977]	451	!
	452	ALLOCATE( sf_dust(1), STAT=ierr ) !* allocate and fill sf_sst (forcing structure) with sn_sst
	453	IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'p4z_sed_init: unable to allocate sf_apr structure' )
	454	!
	455	CALL fld_fill( sf_dust, (/ sn_dust /), cn_dir, 'p4z_sed_init', 'Iron from sediment ', 'nampissed' )
	456	ALLOCATE( sf_dust(1)%fnow(jpi,jpj,1) )
	457	IF( sn_dust%ln_tint ) ALLOCATE( sf_dust(1)%fdta(jpi,jpj,1,2) )
	458	!
	459	! Get total input dust ; need to compute total atmospheric supply of Si in a year
	460	CALL iom_open ( TRIM( sn_dust%clname ) , numdust )
	461	CALL iom_gettime( numdust, zsteps, kntime=ntimes_dust) ! get number of record in file
	462	ALLOCATE( zdust(jpi,jpj,ntimes_dust) )
	463	DO jm = 1, ntimes_dust
	464	CALL iom_get( numdust, jpdom_data, TRIM( sn_dust%clvar ), zdust(:,:,jm), jm )
[935]	465	END DO
	466	CALL iom_close( numdust )
[2977]	467	sumdepsi = 0.e0
	468	DO jm = 1, ntimes_dust
	469	sumdepsi = sumdepsi + glob_sum( zdust(:,:,jm) * e1e2t(:,:) * tmask(:,:,1) )
	470	ENDDO
	471	sumdepsi = sumdepsi * r1_ryyss * 8.8 * 0.075 / 28.1
	472	DEALLOCATE( zdust)
[935]	473	ELSE
[2977]	474	dust(:,:) = 0._wp
	475	sumdepsi = 0._wp
	476	END IF
[935]	477
[2977]	478	! nutrient input from rivers
[935]	479	! --------------------------
[1511]	480	IF( ln_river ) THEN
[2977]	481	ALLOCATE( sf_riverdic(1), STAT=ierr1 ) !* allocate and fill sf_sst (forcing structure) with sn_sst
	482	ALLOCATE( sf_riverdoc(1), STAT=ierr2 ) !* allocate and fill sf_sst (forcing structure) with sn_sst
	483	IF( ierr1 + ierr2 > 0 ) CALL ctl_stop( 'STOP', 'p4z_sed_init: unable to allocate sf_apr structure' )
	484	!
	485	CALL fld_fill( sf_riverdic, (/ sn_riverdic /), cn_dir, 'p4z_sed_init', 'Input DOC from river ', 'nampissed' )
	486	CALL fld_fill( sf_riverdoc, (/ sn_riverdoc /), cn_dir, 'p4z_sed_init', 'Input DOC from river ', 'nampissed' )
	487	ALLOCATE( sf_riverdic(1)%fnow(jpi,jpj,1) )
	488	ALLOCATE( sf_riverdoc(1)%fnow(jpi,jpj,1) )
	489	IF( sn_riverdic%ln_tint ) ALLOCATE( sf_riverdic(1)%fdta(jpi,jpj,1,2) )
	490	IF( sn_riverdoc%ln_tint ) ALLOCATE( sf_riverdoc(1)%fdta(jpi,jpj,1,2) )
	491	! Get total input rivers ; need to compute total river supply in a year
	492	CALL iom_open ( TRIM( sn_riverdic%clname ), numriv )
	493	CALL iom_gettime( numriv, zsteps, kntime=ntimes_riv)
	494	ALLOCATE( zriverdic(jpi,jpj,ntimes_riv) ) ; ALLOCATE( zriverdoc(jpi,jpj,ntimes_riv) )
	495	DO jm = 1, ntimes_riv
	496	CALL iom_get( numriv, jpdom_data, TRIM( sn_riverdic%clvar ), zriverdic(:,:,jm), jm )
	497	CALL iom_get( numriv, jpdom_data, TRIM( sn_riverdoc%clvar ), zriverdoc(:,:,jm), jm )
	498	END DO
[935]	499	CALL iom_close( numriv )
[2977]	500	! N/P and Si releases due to coastal rivers
	501	! -----------------------------------------
	502	rivpo4input = 0._wp
	503	rivalkinput = 0._wp
	504	DO jm = 1, ntimes_riv
	505	rivpo4input = rivpo4input + glob_sum( ( zriverdic(:,:,jm) + zriverdoc(:,:,jm) ) * tmask(:,:,1) )
	506	rivalkinput = rivalkinput + glob_sum( zriverdic(:,:,jm) * tmask(:,:,1) )
	507	END DO
	508	rivpo4input = rivpo4input * 1E9 / 31.6_wp
	509	rivalkinput = rivalkinput * 1E9 / 12._wp
	510	DEALLOCATE( zriverdic) ; DEALLOCATE( zriverdoc)
[935]	511	ELSE
[2977]	512	rivinp(:,:) = 0._wp
	513	cotdep(:,:) = 0._wp
	514	rivpo4input = 0._wp
	515	rivalkinput = 0._wp
	516	END IF
[935]	517
[2977]	518	! nutrient input from dust
[935]	519	! ------------------------
[1511]	520	IF( ln_ndepo ) THEN
[2977]	521	IF(lwp) WRITE(numout,*) ' initialize the nutrient input by dust from ndeposition.orca.nc'
[935]	522	IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
[2977]	523	ALLOCATE( sf_ndepo(1), STAT=ierr3 ) !* allocate and fill sf_sst (forcing structure) with sn_sst
	524	IF( ierr3 > 0 ) CALL ctl_stop( 'STOP', 'p4z_sed_init: unable to allocate sf_apr structure' )
	525	!
	526	CALL fld_fill( sf_ndepo, (/ sn_ndepo /), cn_dir, 'p4z_sed_init', 'Iron from sediment ', 'nampissed' )
	527	ALLOCATE( sf_ndepo(1)%fnow(jpi,jpj,1) )
	528	IF( sn_ndepo%ln_tint ) ALLOCATE( sf_ndepo(1)%fdta(jpi,jpj,1,2) )
	529	!
	530	! Get total input dust ; need to compute total atmospheric supply of N in a year
	531	CALL iom_open ( TRIM( sn_ndepo%clname ), numdepo )
	532	CALL iom_gettime( numdepo, zsteps, kntime=ntimes_ndep)
	533	ALLOCATE( zndepo(jpi,jpj,ntimes_ndep) )
	534	DO jm = 1, ntimes_ndep
	535	CALL iom_get( numdepo, jpdom_data, TRIM( sn_ndepo%clvar ), zndepo(:,:,jm), jm )
	536	END DO
	537	CALL iom_close( numdepo )
	538	nitdepinput = 0._wp
	539	DO jm = 1, ntimes_ndep
	540	nitdepinput = nitdepinput + glob_sum( zndepo(:,:,jm) * e1e2t(:,:) * tmask(:,:,1) )
	541	ENDDO
	542	nitdepinput = nitdepinput * 7.6 / 14E6
	543	DEALLOCATE( zndepo)
[935]	544	ELSE
[2977]	545	nitdep(:,:) = 0._wp
	546	nitdepinput = 0._wp
[935]	547	ENDIF
	548
[2977]	549	! coastal and island masks
[935]	550	! ------------------------
[2977]	551	IF( ln_ironsed ) THEN
	552	IF(lwp) WRITE(numout,*) ' computation of an island mask to enhance coastal supply of iron'
[935]	553	IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
[2977]	554	CALL iom_open ( TRIM( sn_ironsed%clname ), numiron )
	555	ALLOCATE( zcmask(jpi,jpj,jpk) )
	556	CALL iom_get ( numiron, jpdom_data, TRIM( sn_ironsed%clvar ), zcmask(:,:,:), 1 )
	557	CALL iom_close( numiron )
[935]	558	!
	559	DO jk = 1, 5
	560	DO jj = 2, jpjm1
[1503]	561	DO ji = fs_2, fs_jpim1
[935]	562	IF( tmask(ji,jj,jk) /= 0. ) THEN
	563	zmaskt = tmask(ji+1,jj,jk) * tmask(ji-1,jj,jk) * tmask(ji,jj+1,jk) &
	564	& * tmask(ji,jj-1,jk) * tmask(ji,jj,jk+1)
[2774]	565	IF( zmaskt == 0. ) zcmask(ji,jj,jk ) = MAX( 0.1, zcmask(ji,jj,jk) )
[2977]	566	END IF
[935]	567	END DO
	568	END DO
	569	END DO
[2977]	570	CALL lbc_lnk( zcmask , 'T', 1. ) ! lateral boundary conditions on cmask (sign unchanged)
[935]	571	DO jk = 1, jpk
	572	DO jj = 1, jpj
	573	DO ji = 1, jpi
[2977]	574	zexpide = MIN( 8.,( fsdept(ji,jj,jk) / 500. )**(-1.5) )
	575	zdenitide = -0.9543 + 0.7662 * LOG( zexpide ) - 0.235 * LOG( zexpide )**2
	576	zcmask(ji,jj,jk) = zcmask(ji,jj,jk) * MIN( 1., EXP( zdenitide ) / 0.5 )
[935]	577	END DO
	578	END DO
	579	END DO
[2977]	580	! Coastal supply of iron
	581	! -------------------------
	582	ironsed(:,:,jpk) = 0._wp
	583	DO jk = 1, jpkm1
	584	ironsed(:,:,jk) = sedfeinput * zcmask(:,:,jk) / ( fse3t(:,:,jk) * rday )
	585	END DO
	586	DEALLOCATE( zcmask)
[935]	587	ELSE
[2977]	588	ironsed(:,:,:) = 0._wp
[935]	589	ENDIF
[2977]	590	!
	591	IF(lwp) THEN
	592	WRITE(numout,*)
	593	WRITE(numout,*) ' Total input of elements from river supply'
	594	WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
	595	WRITE(numout,) ' N Supply : ', rivpo4input/7.61E3/1E12*14.,' TgN/yr'
	596	WRITE(numout,) ' Si Supply : ', rivalkinput/6.1E3/1E12*32.,' TgSi/yr'
	597	WRITE(numout,) ' Alk Supply : ', rivalkinput1E3/1E12,' Teq/yr'
	598	WRITE(numout,) ' DIC Supply : ', rivpo4input2.6311E312./1E12,'TgC/yr'
	599	WRITE(numout,*)
	600	WRITE(numout,*) ' Total input of elements from atmospheric supply'
	601	WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
	602	WRITE(numout,) ' N Supply : ', nitdepinput/7.61E3/1E12*14.,' TgN/yr'
	603	WRITE(numout,*)
	604	ENDIF
	605	!
[935]	606	END SUBROUTINE p4z_sed_init
	607
[2715]	608	INTEGER FUNCTION p4z_sed_alloc()
	609	!!----------------------------------------------------------------------
	610	!! * ROUTINE p4z_sed_alloc *
	611	!!----------------------------------------------------------------------
	612
[2977]	613	ALLOCATE( dust (jpi,jpj), rivinp(jpi,jpj) , cotdep(jpi,jpj), &
	614	& nitdep(jpi,jpj), ironsed(jpi,jpj,jpk), STAT=p4z_sed_alloc )
[2715]	615
	616	IF( p4z_sed_alloc /= 0 ) CALL ctl_warn('p4z_sed_alloc : failed to allocate arrays.')
	617
	618	END FUNCTION p4z_sed_alloc
[935]	619	#else
	620	!!======================================================================
	621	!! Dummy module : No PISCES bio-model
	622	!!======================================================================
	623	CONTAINS
	624	SUBROUTINE p4z_sed ! Empty routine
	625	END SUBROUTINE p4z_sed
	626	#endif
	627
	628	!!======================================================================
	629	END MODULE p4zsed

Note: See TracBrowser for help on using the repository browser.

Download in other formats: