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.

p4zprod.F90 in trunk/NEMO/TOP_SRC/PISCES – NEMO

Context Navigation

source: trunk/NEMO/TOP_SRC/PISCES/p4zprod.F90 @ 1073

Last change on this file since 1073 was 1073, checked in by cetlod, 16 years ago
update PISCES model, see ticket:190
Property svn:executable set to ``*
File size: 18.1 KB

Line
1	MODULE p4zprod
2	!!======================================================================
3	!! * MODULE p4zprod *
4	!! TOP : PISCES
5	!!======================================================================
6	!! History : 1.0 ! 2004 (O. Aumont) Original code
7	!! 2.0 ! 2007-12 (C. Ethe, G. Madec) F90
8	!!----------------------------------------------------------------------
9	#if defined key_pisces
10	!!----------------------------------------------------------------------
11	!! 'key_pisces' PISCES bio-model
12	!!----------------------------------------------------------------------
13	!! p4z_prod :
14	!!----------------------------------------------------------------------
15	USE trc
16	USE oce_trc !
17	USE trp_trc !
18	USE sms_pisces !
19	USE prtctl_trc
20	USE p4zint
21	USE p4zlim
22	USE p4zopt
23
24	IMPLICIT NONE
25	PRIVATE
26
27	PUBLIC p4z_prod ! called in p4zbio.F90
28
29	!! * Shared module variables
30	REAL(wp), PUBLIC :: &
31	pislope = 3.0_wp , & !:
32	pislope2 = 3.0_wp , & !:
33	excret = 10.e-5_wp , & !:
34	excret2 = 0.05_wp , & !:
35	chlcnm = 0.033_wp , & !:
36	chlcdm = 0.05_wp , & !:
37	fecnm = 10.E-6_wp , & !:
38	fecdm = 15.E-6_wp , & !:
39	grosip = 0.151_wp
40
41	REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: &
42	& prmax
43
44	REAL(wp) :: &
45	texcret , & !: 1 - excret
46	texcret2 , & !: 1 - excret2
47	rpis180 , & !: rpi / 180
48	tpp = 0. !: Total primary production
49
50	!!* Substitution
51	# include "domzgr_substitute.h90"
52	!!----------------------------------------------------------------------
53	!! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
54	!! $Header:$
55	!! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
56	!!----------------------------------------------------------------------
57
58	CONTAINS
59
60	SUBROUTINE p4z_prod( kt , jnt )
61	!!---------------------------------------------------------------------
62	!! * ROUTINE p4z_prod *
63	!!
64	!! ** Purpose : Compute the phytoplankton production depending on
65	!! light, temperature and nutrient availability
66	!!
67	!! ** Method : - ???
68	!!---------------------------------------------------------------------
69	INTEGER, INTENT(in) :: kt, jnt
70	INTEGER :: ji, jj, jk, nspyr
71	REAL(wp) :: zsilfac, zfact, zrfact2
72	REAL(wp) :: zprdiachl, zprbiochl, zsilim, ztn, zadap, zadap2
73	REAL(wp) :: zlim, zsilfac2, zsiborn, zprod, zetot2, zmax, zproreg, zproreg2
74	REAL(wp) :: zmxltst, zmxlday, zlim1
75	REAL(wp) :: zpislopen , zpislope2n
76	REAL(wp) :: zrum, zcodel, zargu
77	REAL(wp), DIMENSION(jpi,jpj) :: zmixnano , zmixdiat, zstrn
78	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopead , zpislopead2
79	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprdia , zprbio, zysopt
80	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprorca , zprorcad, zprofed
81	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprofen , zprochln, zprochld
82	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronew , zpronewd
83	CHARACTER (len=25) :: charout
84	!!---------------------------------------------------------------------
85
86
87	IF( ( kt * jnt ) == nittrc000 ) CALL p4z_prod_init ! Initialization (first time-step only)
88
89
90	zprorca (:,:,:) = 0.0
91	zprorcad(:,:,:) = 0.0
92	zprofed(:,:,:) = 0.0
93	zprofen(:,:,:) = 0.0
94	zprochln(:,:,:) = 0.0
95	zprochld(:,:,:) = 0.0
96	zpronew (:,:,:) = 0.0
97	zpronewd(:,:,:) = 0.0
98	zprdia (:,:,:) = 0.0
99	zprbio (:,:,:) = 0.0
100	zysopt (:,:,:) = 0.0
101
102	nspyr = INT( raass / rdt )
103
104
105	! Computation of the optimal production
106	! -------------------------------------
107
108	# if defined key_off_degrad
109	prmax(:,:,:) = 0.6 / rjjss * tgfunc(:,:,:) * facvol(:,:,:)
110	# else
111	prmax(:,:,:) = 0.6 / rjjss * tgfunc(:,:,:)
112	# endif
113
114	! compute the day length depending on latitude and the day
115	!--------------------------------------------------------
116	IF(lwp) write(numout,*)
117	IF(lwp) write(numout,*) 'p4zday : - Julian day ', nday_year
118	IF(lwp) write(numout,*) '~~~~~~'
119
120	IF( nleapy == 1 .AND. MOD( nyear, 4 ) == 0 ) THEN
121	zrum = FLOAT( nday_year - 80 ) / 366.
122	ELSE
123	zrum = FLOAT( nday_year - 80 ) / 365.
124	ENDIF
125	zcodel = ASIN( SIN( zrum * rpi * 2. ) * SIN( rpis180 * 23.5 ) )
126
127	! day length in hours
128	zstrn(:,:) = 0.
129	DO jj = 1, jpj
130	DO ji = 1, jpi
131	zargu = TAN( zcodel ) * TAN( gphit(ji,jj) * rpis180 )
132	zargu = MAX( -1., MIN( 1., zargu ) )
133	zstrn(ji,jj) = MAX( 0.0, 24. - 2. * ACOS( zargu ) / rpis180 / 15. )
134	END DO
135	END DO
136
137
138	!CDIR NOVERRCHK
139	DO jk = 1, jpkm1
140	!CDIR NOVERRCHK
141	DO jj = 1, jpj
142	!CDIR NOVERRCHK
143	DO ji = 1, jpi
144
145	! Computation of the P-I slope for nanos and diatoms
146	! --------------------------------------------------
147	IF( etot(ji,jj,jk) > 1.E-3 ) THEN
148	ztn = MAX( 0., tn(ji,jj,jk) - 15. )
149	zadap = 0.+ 1.* ztn / ( 2.+ ztn )
150	zadap2 = 0.e0
151
152	zfact = EXP( -0.21 * emoy(ji,jj,jk) )
153
154	zpislopead (ji,jj,jk) = pislope * ( 1.+ zadap * zfact )
155	zpislopead2(ji,jj,jk) = pislope2 * ( 1.+ zadap2 * zfact )
156
157	zpislopen = zpislopead(ji,jj,jk) * trn(ji,jj,jk,jpnch) &
158	& / ( trn(ji,jj,jk,jpphy) * 12. + rtrn ) &
159	& / ( prmax(ji,jj,jk) * rjjss * xlimphy(ji,jj,jk) + rtrn )
160
161	zpislope2n = zpislopead2(ji,jj,jk) * trn(ji,jj,jk,jpdch) &
162	& / ( trn(ji,jj,jk,jpdia) * 12. + rtrn ) &
163	& / ( prmax(ji,jj,jk) * rjjss * xlimdia(ji,jj,jk) + rtrn )
164
165	! Computation of production function
166	! ----------------------------------
167
168	zprbio(ji,jj,jk) = prmax(ji,jj,jk) * &
169	& ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) )
170	zprdia(ji,jj,jk) = prmax(ji,jj,jk) * &
171	& ( 1.- EXP( -zpislope2n * ediat(ji,jj,jk) ) )
172	ENDIF
173	END DO
174	END DO
175	END DO
176
177
178	DO jk = 1, jpkm1
179	DO jj = 1, jpj
180	DO ji = 1, jpi
181
182	IF( etot(ji,jj,jk) > 1.E-3 ) THEN
183	! Si/C of diatoms
184	! ------------------------
185	! Si/C increases with iron stress and silicate availability
186	! Si/C is arbitrariliy increased for very high Si concentrations
187	! to mimic the very high ratios observed in the Southern Ocean (silpot2)
188
189	zlim1 = trn(ji,jj,jk,jpsil) / ( trn(ji,jj,jk,jpsil) + xksi1 )
190	zlim = xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk)
191
192	zsilim = MIN( zprdia(ji,jj,jk) / ( rtrn + prmax(ji,jj,jk) ), &
193	& trn(ji,jj,jk,jpfer) / ( concdfe(ji,jj,jk) + trn(ji,jj,jk,jpfer) ), &
194	& trn(ji,jj,jk,jppo4) / ( concdnh4 + trn(ji,jj,jk,jppo4) ), &
195	& zlim )
196	zsilfac = 5.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim1 - 0.5 ) ) ) + 1.e0
197	zsiborn = MAX( 0.e0, ( trn(ji,jj,jk,jpsil) - 15.e-6 ) )
198	zsilfac2 = 1.+ 3.* zsiborn / ( zsiborn + xksi2 )
199	zsilfac = MIN( 6.4,zsilfac * zsilfac2)
200	zysopt(ji,jj,jk) = grosip * zlim1 * zsilfac
201
202	ENDIF
203	END DO
204	END DO
205	END DO
206
207	! Computation of the limitation term due to
208	! A mixed layer deeper than the euphotic depth
209	! --------------------------------------------
210
211	DO jj = 1, jpj
212	DO ji = 1, jpi
213	zmxltst = MAX( 0.e0, hmld(ji,jj) - heup(ji,jj) )
214	zmxlday = zmxltst**2 / rjjss
215	zmixnano(ji,jj) = 1.- zmxlday / ( 1.+ zmxlday )
216	zmixdiat(ji,jj) = 1.- zmxlday / ( 3.+ zmxlday )
217	END DO
218	END DO
219
220	DO jk = 1, jpkm1
221	DO jj = 1, jpj
222	DO ji = 1, jpi
223	IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
224
225	! Mixed-layer effect on production
226	! --------------------------------
227	zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * zmixnano(ji,jj)
228	zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * zmixdiat(ji,jj)
229	ENDIF
230	END DO
231	END DO
232	END DO
233
234	! Computation of the fractionnal day length
235	! -----------------------------------------
236	! zstrn(:,:) = strn(:,:)
237
238	WHERE( zstrn(:,:) < 1.e0 ) zstrn(:,:) = 24.
239	zstrn(:,:) = 24. / zstrn(:,:)
240	! DO jj = 1, jpj
241	! DO ji = 1, jpi
242	!
243	! Computation of the maximum light intensity
244	! ------------------------------------------
245	! IF( zstrn(ji,jj) < 1.e0 ) zstrn(ji,jj) = 24.
246	! zstrn(ji,jj) = 24. / zstrn(ji,jj)
247	! END DO
248	! END DO
249
250	!CDIR NOVERRCHK
251	DO jk = 1, jpkm1
252	!CDIR NOVERRCHK
253	DO jj = 1, jpj
254	!CDIR NOVERRCHK
255	DO ji = 1, jpi
256
257	IF( etot(ji,jj,jk) > 1.E-3 ) THEN
258	! Computation of the various production terms for nanophyto.
259	! ----------------------------------------------------------
260	zetot2 = enano(ji,jj,jk) * zstrn(ji,jj)
261	zmax = MAX( 0.1, xlimphy(ji,jj,jk) )
262	zpislopen = zpislopead(ji,jj,jk) &
263	& * trn(ji,jj,jk,jpnch) / ( rtrn + trn(ji,jj,jk,jpphy) * 12.) &
264	& / ( prmax(ji,jj,jk) * rjjss * zmax + rtrn )
265
266	zprbiochl = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * zetot2 ) )
267
268	zprorca(ji,jj,jk) = zprbio(ji,jj,jk) * xlimphy(ji,jj,jk) * trn(ji,jj,jk,jpphy) * rfact2
269
270	zpronew(ji,jj,jk) = zprorca(ji,jj,jk) * xnanono3(ji,jj,jk) &
271	& / ( xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) + rtrn )
272	zprod = rjjss * zprorca(ji,jj,jk) * zprbiochl * trn(ji,jj,jk,jpphy) *zmax
273
274	zprofen(ji,jj,jk) = (fecnm)*2 zprod / chlcnm &
275	& / ( zpislopead(ji,jj,jk) * zetot2 * trn(ji,jj,jk,jpnfe) + rtrn )
276
277	zprochln(ji,jj,jk) = chlcnm * 144. * zprod &
278	& / ( zpislopead(ji,jj,jk) * zetot2 * trn(ji,jj,jk,jpnch) + rtrn )
279	ENDIF
280	END DO
281	END DO
282	END DO
283
284	!CDIR NOVERRCHK
285	DO jk = 1, jpkm1
286	!CDIR NOVERRCHK
287	DO jj = 1, jpj
288	!CDIR NOVERRCHK
289	DO ji = 1, jpi
290	IF( etot(ji,jj,jk) > 1.E-3 ) THEN
291	! Computation of the various production terms for diatoms
292	! -------------------------------------------------------
293	zetot2 = ediat(ji,jj,jk) * zstrn(ji,jj)
294	zmax = MAX( 0.1, xlimdia(ji,jj,jk) )
295	zpislope2n = zpislopead2(ji,jj,jk) * trn(ji,jj,jk,jpdch) &
296	& / ( rtrn + trn(ji,jj,jk,jpdia) * 12.) &
297	& / ( prmax(ji,jj,jk) * rjjss * zmax + rtrn )
298
299	zprdiachl = prmax(ji,jj,jk) * ( 1.- EXP( -zetot2 * zpislope2n ) )
300
301	zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * trn(ji,jj,jk,jpdia) * rfact2
302
303	zpronewd(ji,jj,jk) = zprorcad(ji,jj,jk) * xdiatno3(ji,jj,jk) &
304	& / ( xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) + rtrn )
305
306	zprod = rjjss * zprorcad(ji,jj,jk) * zprdiachl * trn(ji,jj,jk,jpdia) * zmax
307
308	zprofed(ji,jj,jk) = (fecdm)*2 zprod / chlcdm &
309	& / ( zpislopead2(ji,jj,jk) * zetot2 * trn(ji,jj,jk,jpdfe) + rtrn )
310
311	zprochld(ji,jj,jk) = chlcdm * 144. * zprod &
312	& / ( zpislopead2(ji,jj,jk) * zetot2 * trn(ji,jj,jk,jpdch) + rtrn )
313
314	ENDIF
315	END DO
316	END DO
317	END DO
318	!
319
320	!
321	! Update the arrays TRA which contain the biological sources and sinks
322	! --------------------------------------------------------------------
323
324	DO jk = 1, jpkm1
325	DO jj = 1, jpj
326	DO ji =1 ,jpi
327	zproreg = zprorca(ji,jj,jk) - zpronew(ji,jj,jk)
328	zproreg2 = zprorcad(ji,jj,jk) - zpronewd(ji,jj,jk)
329	tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zprorca(ji,jj,jk) - zprorcad(ji,jj,jk)
330	tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - zpronew(ji,jj,jk) - zpronewd(ji,jj,jk)
331	tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) - zproreg - zproreg2
332	tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) + zprorca(ji,jj,jk) * texcret
333	tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) + zprochln(ji,jj,jk) * texcret
334	tra(ji,jj,jk,jpnfe) = tra(ji,jj,jk,jpnfe) + zprofen(ji,jj,jk) * texcret
335	tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) + zprorcad(ji,jj,jk) * texcret2
336	tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) + zprochld(ji,jj,jk) * texcret2
337	tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) + zprofed(ji,jj,jk) * texcret2
338	tra(ji,jj,jk,jpbsi) = tra(ji,jj,jk,jpbsi) + zprorcad(ji,jj,jk) * zysopt(ji,jj,jk) * texcret2
339	tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + &
340	& excret2 * zprorcad(ji,jj,jk) + excret * zprorca(ji,jj,jk)
341	tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2ut * ( zproreg + zproreg2) &
342	& + ( o2ut + o2nit ) * ( zpronew(ji,jj,jk) + zpronewd(ji,jj,jk) )
343	tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) &
344	& - texcret * zprofen(ji,jj,jk) - texcret2 * zprofed(ji,jj,jk)
345	tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) &
346	& - texcret2 * zprorcad(ji,jj,jk) * zysopt(ji,jj,jk)
347	tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) - zprorca(ji,jj,jk) - zprorcad(ji,jj,jk)
348	tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) &
349	& + rno3 * ( zpronew(ji,jj,jk) + zpronewd(ji,jj,jk) )
350	END DO
351	END DO
352	END DO
353
354	! Total primary production per year
355	DO jk = 1, jpkm1
356	DO jj = 1, jpj
357	DO ji = 1, jpi
358	tpp = tpp + ( zprorca(ji,jj,jk) + zprorcad(ji,jj,jk) ) &
359	#if defined key_off_degrad
360	& * facvol(ji,jj,jk) &
361	#endif
362	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) * tmask_i(ji,jj)
363	END DO
364	END DO
365	END DO
366
367	IF( MOD( kt, nspyr ) == 0 ) THEN
368	WRITE(numout,*) 'Total PP :'
369	WRITE(numout,) '-------------------- : ',tpp 12. / 1.E12
370	WRITE(numout,*) '(GtC/an)'
371	tpp = 0.
372	ENDIF
373
374	#if defined key_trc_dia3d
375	zrfact2 = 1.e3 * rfact2r
376	! Supplementary diagnostics
377	! -------------------------
378	trc3d(:,:,:,4) = etot(:,:,:)
379	trc3d(:,:,:,5) = zprorca(:,:,:) * zrfact2
380	trc3d(:,:,:,6) = zprorcad(:,:,:) * zrfact2
381	trc3d(:,:,:,7) = zpronew(:,:,:) * zrfact2
382	trc3d(:,:,:,8) = zpronewd(:,:,:) * zrfact2
383	trc3d(:,:,:,9) = zprorcad(:,:,:) * zysopt(:,:,:) * zrfact2
384	trc3d(:,:,:,10) = zprofed(:,:,:) * zrfact2
385	#if ! defined key_kriest
386	trc3d(:,:,:,11) = zprofen(:,:,:) * zrfact2
387	#endif
388	#endif
389
390	IF(ln_ctl) THEN ! print mean trends (used for debugging)
391	WRITE(charout, FMT="('prod')")
392	CALL prt_ctl_trc_info(charout)
393	CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
394	ENDIF
395
396	END SUBROUTINE p4z_prod
397
398	SUBROUTINE p4z_prod_init
399
400	!!----------------------------------------------------------------------
401	!! * ROUTINE p4z_prod_init *
402	!!
403	!! ** Purpose : Initialization of phytoplankton production parameters
404	!!
405	!! ** Method : Read the natprod namelist and check the parameters
406	!! called at the first timestep (nittrc000)
407	!!
408	!! ** input : Namelist natprod
409	!!
410	!!----------------------------------------------------------------------
411
412	NAMELIST/natprod/ pislope, pislope2, excret, excret2, chlcnm, chlcdm, &
413	& fecnm, fecdm, grosip
414
415	REWIND( numnat ) ! read numnat
416	READ ( numnat, natprod )
417
418	IF(lwp) THEN ! control print
419	WRITE(numout,*) ' '
420	WRITE(numout,*) ' Namelist parameters for phytoplankton growth, natprod'
421	WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
422	WRITE(numout,*) ' mean Si/C ratio grosip =', grosip
423	WRITE(numout,*) ' P-I slope pislope =', pislope
424	WRITE(numout,*) ' excretion ratio of nanophytoplankton excret =', excret
425	WRITE(numout,*) ' excretion ratio of diatoms excret2 =', excret2
426	WRITE(numout,*) ' P-I slope for diatoms pislope2 =', pislope2
427	WRITE(numout,*) ' Minimum Chl/C in nanophytoplankton chlcnm =', chlcnm
428	WRITE(numout,*) ' Minimum Chl/C in diatoms chlcdm =', chlcdm
429	WRITE(numout,*) ' Maximum Fe/C in nanophytoplankton fecnm =', fecnm
430	WRITE(numout,*) ' Minimum Fe/C in diatoms fecdm =', fecdm
431	ENDIF
432
433	rpis180 = rpi / 180.
434	texcret = 1.0 - excret
435	texcret2 = 1.0 - excret2
436
437	END SUBROUTINE p4z_prod_init
438
439
440
441	#else
442	!!======================================================================
443	!! Dummy module : No PISCES bio-model
444	!!======================================================================
445	CONTAINS
446	SUBROUTINE p4z_prod ! Empty routine
447	END SUBROUTINE p4z_prod
448	#endif
449
450	!!======================================================================
451	END MODULE p4zprod

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

Download in other formats: