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.

trcbio_medusa.F90 in branches/NERC/dev_r5518_GO6_CO2_cmip/NEMOGCM/NEMO/TOP_SRC/MEDUSA – NEMO

Context Navigation

source: branches/NERC/dev_r5518_GO6_CO2_cmip/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90 @ 9285

Last change on this file since 9285 was 9285, checked in by jpalmier, 6 years ago
JPALM --26-01-2018 -- change MEDUSA atm co2 input method. now reads a file a variable length
File size: 43.3 KB

Rev	Line
[5726]	1	MODULE trcbio_medusa
	2	!!======================================================================
	3	!! * MODULE trcbio *
	4	!! TOP : MEDUSA
	5	!!======================================================================
	6	!! History :
	7	!! - ! 1999-07 (M. Levy) original code
[8441]	8	!! - ! 2000-12 (E. Kestenare) assign parameters to name
	9	!! individual tracers
[5726]	10	!! - ! 2001-03 (M. Levy) LNO3 + dia2d
	11	!! 2.0 ! 2007-12 (C. Deltel, G. Madec) F90
	12	!! - ! 2008-08 (K. Popova) adaptation for MEDUSA
	13	!! - ! 2008-11 (A. Yool) continuing adaptation for MEDUSA
	14	!! - ! 2010-03 (A. Yool) updated for branch inclusion
	15	!! - ! 2011-08 (A. Yool) updated for ROAM (see below)
	16	!! - ! 2013-03 (A. Yool) updated for iMARNET
	17	!! - ! 2013-05 (A. Yool) updated for v3.5
	18	!! - ! 2014-08 (A. Yool, J. Palm) Add DMS module for UKESM1 model
[5841]	19	!! - ! 2015-06 (A. Yool) Update to include MOCSY
	20	!! - ! 2015-07 (A. Yool) Update for rolling averages
[8441]	21	!! - ! 2015-10 (J. Palm) Update for diag outputs through
	22	!! iom_use
[7224]	23	!! - ! 2016-11 (A. Yool) Updated diags for CMIP6
[8131]	24	!! - ! 2017-05 (A. Yool) Added extra DMS calculation
[9257]	25	!! - ! 2017-11 (J. Palm, A. Yool) Diagnose tracer excursions
[5726]	26	!!----------------------------------------------------------------------
	27	!!
	28	#if defined key_roam
	29	!!----------------------------------------------------------------------
	30	!! Updates for the ROAM project include:
	31	!! - addition of DIC, alkalinity, detrital carbon and oxygen tracers
[8074]	32	!! - addition of air-sea fluxes of CO2 and oxygen
[5726]	33	!! - periodic (monthly) calculation of full 3D carbonate chemistry
	34	!! - detrital C:N ratio now free to evolve dynamically
	35	!! - benthic storage pools
	36	!!
	37	!! Opportunity also taken to add functionality:
	38	!! - switch for Liebig Law (= most-limiting) nutrient uptake
	39	!! - switch for accelerated seafloor detritus remineralisation
	40	!! - switch for fast -> slow detritus transfer at seafloor
	41	!! - switch for ballast vs. Martin vs. Henson fast detritus remin.
	42	!! - per GMD referee remarks, xfdfrac3 introduced for grazed PDS
[8074]	43	!!----------------------------------------------------------------------
	44	#endif
[5726]	45	!!
[8074]	46	#if defined key_mocsy
	47	!!----------------------------------------------------------------------
[5841]	48	!! Updates with the addition of MOCSY include:
	49	!! - option to use PML or MOCSY carbonate chemistry (the latter is
	50	!! preferred)
	51	!! - central calculation of gas transfer velocity, f_kw660; previously
	52	!! this was done separately for CO2 and O2 with predictable results
	53	!! - distribution of f_kw660 to both PML and MOCSY CO2 air-sea flux
	54	!! calculations and to those for O2 air-sea flux
	55	!! - extra diagnostics included for MOCSY
	56	!!----------------------------------------------------------------------
	57	#endif
	58	!!
[5726]	59	#if defined key_medusa
	60	!!----------------------------------------------------------------------
	61	!! MEDUSA bio-model
	62	!!----------------------------------------------------------------------
	63	!! trc_bio_medusa :
	64	!!----------------------------------------------------------------------
	65	!! AXY (30/01/14): necessary to find NaNs on HECTOR
	66	USE, INTRINSIC :: ieee_arithmetic
	67
[8441]	68	USE bio_medusa_mod, ONLY: b0, fdep1, &
	69	ibenthic, idf, idfval, &
	70	# if defined key_roam
	71	f_xco2a, &
	72	zalk, zdic, zoxy, zsal, ztmp, &
	73	# endif
	74	# if defined key_mocsy
	75	zpho, &
	76	# endif
	77	zchd, zchn, zdet, zdin, zdtc, &
	78	zfer, zpds, zphd, zphn, zsil, &
	79	zzme, zzmi
	80	USE dom_oce, ONLY: e3t_0, e3t_n, &
	81	gdept_0, gdept_n, &
	82	gdepw_0, gdepw_n, &
[9285]	83	nday_year, nsec_day, &
	84	nyear, nyear_len, &
[9257]	85	rdt, tmask, mig, mjg, nimpp, &
	86	njmpp
[8441]	87	USE in_out_manager, ONLY: lwp, numout, nn_date0
	88	# if defined key_iomput
	89	USE iom, ONLY: lk_iomput
	90	# endif
	91	USE lbclnk, ONLY: lbc_lnk
[9257]	92	USE lib_mpp, ONLY: mpp_max, mpp_maxloc, &
	93	mpp_min, mpp_minloc, &
	94	ctl_stop, ctl_warn, lk_mpp
[8441]	95	USE oce, ONLY: tsb, tsn
	96	USE par_kind, ONLY: wp
	97	USE par_medusa, ONLY: jpalk, jpchd, jpchn, jpdet, &
	98	jpdic, jpdin, jpdtc, jpfer, &
	99	jpoxy, jppds, jpphd, jpphn, &
	100	jpsil, jpzme, jpzmi
	101	USE par_oce, ONLY: jp_sal, jp_tem, jpi, jpim1, &
	102	jpj, jpjm1, jpk
[6744]	103	!! JPALM (27-06-2016): add lk_oasis for CO2 and DMS coupling with atm
[8441]	104	USE sbc_oce, ONLY: lk_oasis
[9285]	105	USE sms_medusa, ONLY: hist_pco2, co2_yinit, co2_yend
[8442]	106	USE trc, ONLY: ln_rsttr, nittrc000, trn
[8441]	107	USE bio_medusa_init_mod, ONLY: bio_medusa_init
	108	USE carb_chem_mod, ONLY: carb_chem
	109	USE air_sea_mod, ONLY: air_sea
	110	USE plankton_mod, ONLY: plankton
	111	USE iron_chem_scav_mod, ONLY: iron_chem_scav
	112	USE detritus_mod, ONLY: detritus
	113	USE bio_medusa_update_mod, ONLY: bio_medusa_update
	114	USE bio_medusa_diag_mod, ONLY: bio_medusa_diag
	115	USE bio_medusa_diag_slice_mod, ONLY: bio_medusa_diag_slice
	116	USE bio_medusa_fin_mod, ONLY: bio_medusa_fin
[6744]	117
[5726]	118	IMPLICIT NONE
	119	PRIVATE
	120
[8441]	121	PUBLIC trc_bio_medusa ! called in trcsms_medusa.F90
[9257]	122	PUBLIC trc_bio_exceptionnal_fix ! here
[5726]	123
	124	!!* Substitution
	125	# include "domzgr_substitute.h90"
	126	!!----------------------------------------------------------------------
	127	!! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
	128	!! $Id$
	129	!! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
	130	!!----------------------------------------------------------------------
	131
	132	CONTAINS
	133
[8074]	134	SUBROUTINE trc_bio_medusa( kt )
[8441]	135	!!------------------------------------------------------------------
[8074]	136	!! * ROUTINE trc_bio *
	137	!!
[8441]	138	!! ** Purpose : compute the now trend due to biogeochemical processes
	139	!! and add it to the general trend of passive tracers
	140	!! equations
[8074]	141	!!
[8441]	142	!! ** Method : each now biological flux is calculated in function of
	143	!! now concentrations of tracers.
	144	!! depending on the tracer, these fluxes are sources or
	145	!! sinks.
	146	!! The total of the sources and sinks for each tracer
[8074]	147	!! is added to the general trend.
	148	!!
	149	!! tra = tra + zf...tra - zftra...
	150	!! \| \|
	151	!! \| \|
	152	!! source sink
	153	!!
	154	!! IF 'key_trc_diabio' defined , the biogeochemical trends
	155	!! for passive tracers are saved for futher diagnostics.
[8441]	156	!!------------------------------------------------------------------
[8074]	157	!!
	158	!!
[8441]	159	!!------------------------------------------------------------------
[8074]	160	!! Variable conventions
[8441]	161	!!------------------------------------------------------------------
[8074]	162	!!
	163	!! names: z*** - state variable
[8441]	164	!! f*** - function (or temporary variable used in part of
	165	!! a function)
[8074]	166	!! x*** - parameter
	167	!! b*** - right-hand part (sources and sinks)
	168	!! i*** - integer variable (usually used in yes/no flags)
	169	!!
	170	!! time (integer timestep)
	171	INTEGER, INTENT( in ) :: kt
	172	!!
	173	!! spatial array indices
	174	INTEGER :: ji,jj,jk,jn
	175	!!
	176	INTEGER :: iball
[5726]	177	# if defined key_roam
[8074]	178	!!
	179	INTEGER :: iyr1, iyr2
	180	!!
[5726]	181	# endif
[8074]	182	!!
[8441]	183	!! temporary variables
	184	REAL(wp) :: fq0,fq1,fq2,fq3,fq4
[8074]	185	!!
[9257]	186	!! T and S check temporary variable
	187	REAL(wp) :: sumtsn, tsnavg
	188	INTEGER :: summask
	189	CHARACTER(40) :: charout, charout2, charout3, charout4, charout5
	190	!!
[8441]	191	!!------------------------------------------------------------------
[5726]	192
[5841]	193	# if defined key_debug_medusa
[8074]	194	IF (lwp) write (numout,*) 'trc_bio_medusa: variables defined'
	195	CALL flush(numout)
[5841]	196	# endif
	197
[8074]	198	!! AXY (20/11/14): alter this to report on first MEDUSA call
	199	!! IF( kt == nit000 ) THEN
	200	IF( kt == nittrc000 ) THEN
	201	IF(lwp) WRITE(numout,*)
	202	IF(lwp) WRITE(numout,*) ' trc_bio: MEDUSA bio-model'
	203	IF(lwp) WRITE(numout,*) ' ~~~~~~~'
	204	IF(lwp) WRITE(numout,*) ' kt =',kt
	205	ENDIF
[5726]	206
[8074]	207	!! AXY (13/01/12): is benthic model properly interactive? 0 = no, 1 = yes
	208	ibenthic = 1
[5726]	209
[8441]	210	!!------------------------------------------------------------------
[8074]	211	!! b0 is present for debugging purposes; using b0 = 0 sets the tendency
	212	!! terms of all biological equations to 0.
[8441]	213	!!------------------------------------------------------------------
[8074]	214	!!
	215	!! AXY (03/09/14): probably not the smartest move ever, but it'll fit
	216	!! the bill for now; another item on the things-to-sort-
	217	!! out-in-the-future list ...
[5726]	218	# if defined key_kill_medusa
[8074]	219	b0 = 0.
[5726]	220	# else
[8074]	221	b0 = 1.
[5726]	222	# endif
[8441]	223	!!------------------------------------------------------------------
[8074]	224	!! fast detritus ballast scheme (0 = no; 1 = yes)
	225	!! alternative to ballast scheme is same scheme but with no ballast
	226	!! protection (not dissimilar to Martin et al., 1987)
[8441]	227	!!------------------------------------------------------------------
[8074]	228	!!
	229	iball = 1
[5726]	230
[8441]	231	!!------------------------------------------------------------------
[8074]	232	!! full flux diagnostics (0 = no; 1 = yes); appear in ocean.output
	233	!! these should only be used in 1D since they give comprehensive
	234	!! output for ecological functions in the model; primarily used in
	235	!! debugging
[8441]	236	!!------------------------------------------------------------------
[8074]	237	!!
	238	idf = 0
	239	!!
	240	!! timer mechanism
	241	if (kt/120*120.eq.kt) then
	242	idfval = 1
	243	else
	244	idfval = 0
	245	endif
[5726]	246
[8441]	247	!!------------------------------------------------------------------
	248	!! Initialise arrays to zero and set up arrays for diagnostics
	249	!!------------------------------------------------------------------
	250	CALL bio_medusa_init( kt )
[5726]	251
	252	# if defined key_axy_nancheck
[8441]	253	DO jn = jp_msa0,jp_msa1
[8074]	254	!! fq0 = MINVAL(trn(:,:,:,jn))
	255	!! fq1 = MAXVAL(trn(:,:,:,jn))
	256	fq2 = SUM(trn(:,:,:,jn))
[8441]	257	!! if (lwp) write (numout,'(a,2i6,3(1x,1pe15.5))') 'NAN-CHECK', &
	258	!! kt, jn, fq0, fq1, fq2
	259	!! AXY (30/01/14): much to our surprise, the next line doesn't
	260	!! work on HECTOR and has been replaced here with
	261	!! a specialist routine
[8074]	262	!! if (fq2 /= fq2 ) then
	263	if ( ieee_is_nan( fq2 ) ) then
	264	!! there's a NaN here
[8441]	265	if (lwp) write(numout,*) 'NAN detected in field', jn, &
	266	'at time', kt, 'at position:'
[8074]	267	DO jk = 1,jpk
	268	DO jj = 1,jpj
	269	DO ji = 1,jpi
	270	!! AXY (30/01/14): "isnan" problem on HECTOR
	271	!! if (trn(ji,jj,jk,jn) /= trn(ji,jj,jk,jn)) then
	272	if ( ieee_is_nan( trn(ji,jj,jk,jn) ) ) then
[8441]	273	if (lwp) write (numout,'(a,1pe12.2,4i6)') &
	274	'NAN-CHECK', tmask(ji,jj,jk), ji, jj, jk, jn
[8074]	275	endif
	276	enddo
[7224]	277	enddo
	278	enddo
[8074]	279	CALL ctl_stop( 'trcbio_medusa, NAN in incoming tracer field' )
	280	endif
	281	ENDDO
	282	CALL flush(numout)
[5726]	283	# endif
	284
[5841]	285	# if defined key_debug_medusa
[8441]	286	IF (lwp) write (numout,*) &
	287	'trc_bio_medusa: variables initialised and checked'
[8074]	288	CALL flush(numout)
[5841]	289	# endif
	290
[5726]	291	# if defined key_roam
[8441]	292	!!------------------------------------------------------------------
[8074]	293	!! calculate atmospheric pCO2
[8441]	294	!!------------------------------------------------------------------
[8074]	295	!!
	296	!! what's atmospheric pCO2 doing? (data start in 1859)
[9285]	297	iyr1 = nyear
[8074]	298	iyr2 = iyr1 + 1
[9285]	299	if (iyr1 .le. co2_yinit) then
	300	!! before 1850
	301	f_xco2a(:,:) = 284.317
	302	elseif (iyr2 .ge. co2_yend) then
	303	!! after 2010
	304	f_xco2a(:,:) = hist_pco2(co2_yend - co2_yinit +1)
[8074]	305	else
	306	!! just right
	307	fq0 = hist_pco2(iyr1)
	308	fq1 = hist_pco2(iyr2)
	309	fq2 = real(nsec_day) / (60.0 * 60.0 * 24.0)
[9285]	310	fq3 = (real(nday_year) - 1.0 + fq2) / real(nyear_len(1))
[8074]	311	!! AXY (14/06/12): tweaked to make more sense (and be correct)
[9285]	312	!!# if defined key_bs_axy_yrlen
	313	!! !! bugfix: for 360d year with HadGEM2-ES forcing
	314	!! fq3 = (real(nday_year) - 1.0 + fq2) / 360.0
	315	!!# else
	316	!! !! original use of 365 days (not accounting for leap year or
	317	!! !! 360d year)
	318	!! fq3 = (real(nday_year) - 1.0 + fq2) / 365.0
	319	!!# endif
[8074]	320	fq4 = (fq0 * (1.0 - fq3)) + (fq1 * fq3)
[8441]	321	f_xco2a(:,:) = fq4
[8074]	322	endif
	323	# if defined key_axy_pi_co2
[8441]	324	!! OCMIP pre-industrial pCO2
	325	!! f_xco2a(:,:) = 284.725 !! CMIP5 pre-industrial pCO2
[9285]	326	f_xco2a(:,:) = 284.317 !! CMIP6 pre-industrial pCO2
[8074]	327	# endif
[9285]	328	IF(lwp) WRITE(numout,*) ' MEDUSA nyear =', nyear
	329	IF(lwp) WRITE(numout,*) ' MEDUSA nyear =', real(nyear_len(1))
	330	IF(lwp) WRITE(numout,*) ' MEDUSA nsec_day =', real(nsec_day)
	331	IF(lwp) WRITE(numout,*) ' MEDUSA nday_year =', real(nday_year)
[8074]	332	!! AXY (29/01/14): remove surplus diagnostics
[9285]	333	IF(lwp) WRITE(numout,*) ' MEDUSA fq0 =', fq0
	334	IF(lwp) WRITE(numout,*) ' MEDUSA fq1 =', fq1
	335	IF(lwp) WRITE(numout,*) ' MEDUSA fq2 =', fq2
	336	IF(lwp) WRITE(numout,*) ' MEDUSA fq3 =', fq3
	337	IF(lwp) WRITE(numout,*) ' MEDUSA atm pCO2 =', f_xco2a(1,1)
[5726]	338	# endif
	339
[5841]	340	# if defined key_debug_medusa
[8074]	341	IF (lwp) write (numout,*) 'trc_bio_medusa: ready for carbonate chemistry'
	342	IF (lwp) write (numout,*) 'trc_bio_medusa: kt = ', kt
	343	IF (lwp) write (numout,*) 'trc_bio_medusa: nittrc000 = ', nittrc000
	344	CALL flush(numout)
[5841]	345	# endif
	346
[5726]	347	# if defined key_roam
[8074]	348	!! AXY (20/11/14): alter to call on first MEDUSA timestep and then every
	349	!! month (this is hardwired as 960 timesteps but should
	350	!! be calculated and done properly
	351	!! IF( kt == nit000 .or. mod(kt,1920) == 0 ) THEN
	352	!! IF( kt == nittrc000 .or. mod(kt,960) == 0 ) THEN
	353	!!=============================
	354	!! Jpalm -- 07-10-2016 -- need to change carb-chem frequency call :
[8441]	355	!! we don't want to call on the first time-step of all run
	356	!! submission, but only on the very first time-step, and
	357	!! then every month. So we call on nittrc000 if not
	358	!! restarted run, else if one month after last call.
	359	!! Assume one month is 30d --> 36002430 : 2592000s
	360	!! try to call carb-chem at 1st month's tm-stp :
	361	!! x * 30d + 1*rdt(i.e: mod = rdt)
[8074]	362	!! ++ need to pass carb-chem output var through restarts
[8224]	363	If ( (kt == nittrc000 .AND. .NOT.ln_rsttr) .OR. &
	364	((86400mod(nn_date0,100) + mod(ktrdt,2592000.)) == rdt) ) THEN
[8441]	365	!!---------------------------------------------------------------
[8074]	366	!! Calculate the carbonate chemistry for the whole ocean on the first
	367	!! simulation timestep and every month subsequently; the resulting 3D
	368	!! field of omega calcite is used to determine the depth of the CCD
[8441]	369	!!---------------------------------------------------------------
	370	CALL carb_chem( kt )
	371
[8074]	372	ENDIF
[5726]	373	# endif
	374
[5841]	375	# if defined key_debug_medusa
[8074]	376	IF (lwp) write (numout,*) 'trc_bio_medusa: ready for full domain calculations'
	377	CALL flush(numout)
[5841]	378	# endif
	379
[8441]	380	!!------------------------------------------------------------------
[8074]	381	!! MEDUSA has unified equation through the water column
	382	!! (Diff. from LOBSTER which has two sets: bio- and non-bio layers)
	383	!! Statement below in LOBSTER is different: DO jk = 1, jpkbm1
[8441]	384	!!------------------------------------------------------------------
[8074]	385	!!
	386	!! NOTE: the ordering of the loops below differs from that of some other
	387	!! models; looping over the vertical dimension is the outermost loop and
	388	!! this complicates some calculations (e.g. storage of vertical fluxes
	389	!! that can otherwise be done via a singular variable require 2D fields
	390	!! here); however, these issues are relatively easily resolved, but the
	391	!! loops CANNOT be reordered without potentially causing code efficiency
	392	!! problems (e.g. array indexing means that reordering the loops would
	393	!! require skipping between widely-spaced memory location; potentially
	394	!! outside those immediately cached)
	395	!!
	396	!! OPEN vertical loop
	397	DO jk = 1,jpk
	398	!! OPEN horizontal loops
	399	DO jj = 2,jpjm1
[8441]	400	DO ji = 2,jpim1
	401	!! OPEN wet point IF..THEN loop
	402	if (tmask(ji,jj,jk) == 1) then
	403	!!======================================================
	404	!! SETUP LOCAL GRID CELL
	405	!!======================================================
[8074]	406	!!
[8441]	407	!!------------------------------------------------------
	408	!! Some notes on grid vertical structure
	409	!! - fsdepw(ji,jj,jk) is the depth of the upper surface of
	410	!! level jk
	411	!! - fsde3w(ji,jj,jk) is approximately the midpoint of
	412	!! level jk
	413	!! - fse3t(ji,jj,jk) is the thickness of level jk
	414	!!------------------------------------------------------
[8074]	415	!!
[8441]	416	!! AXY (01/03/10): set up level depth (bottom of level)
	417	fdep1(ji,jj) = fsdepw(ji,jj,jk) + fse3t(ji,jj,jk)
[8074]	418	!!
[8441]	419	!! set up model tracers
	420	!! negative values of state variables are not allowed to
	421	!! contribute to the calculated fluxes
	422	!! non-diatom chlorophyll
	423	zchn(ji,jj) = max(0.,trn(ji,jj,jk,jpchn))
	424	!! diatom chlorophyll
	425	zchd(ji,jj) = max(0.,trn(ji,jj,jk,jpchd))
	426	!! non-diatoms
	427	zphn(ji,jj) = max(0.,trn(ji,jj,jk,jpphn))
	428	!! diatoms
	429	zphd(ji,jj) = max(0.,trn(ji,jj,jk,jpphd))
	430	!! diatom silicon
	431	zpds(ji,jj) = max(0.,trn(ji,jj,jk,jppds))
	432	!! AXY (28/01/10): probably need to take account of
	433	!! chl/biomass connection
	434	if (zchn(ji,jj).eq.0.) zphn(ji,jj) = 0.
	435	if (zchd(ji,jj).eq.0.) zphd(ji,jj) = 0.
	436	if (zphn(ji,jj).eq.0.) zchn(ji,jj) = 0.
	437	if (zphd(ji,jj).eq.0.) zchd(ji,jj) = 0.
	438	!! AXY (23/01/14): duh - why did I forget diatom silicon?
	439	if (zpds(ji,jj).eq.0.) zphd(ji,jj) = 0.
	440	if (zphd(ji,jj).eq.0.) zpds(ji,jj) = 0.
	441	ENDIF
	442	ENDDO
	443	ENDDO
[7224]	444
[8441]	445	DO jj = 2,jpjm1
	446	DO ji = 2,jpim1
	447	if (tmask(ji,jj,jk) == 1) then
	448	!! microzooplankton
	449	zzmi(ji,jj) = max(0.,trn(ji,jj,jk,jpzmi))
	450	!! mesozooplankton
	451	zzme(ji,jj) = max(0.,trn(ji,jj,jk,jpzme))
	452	!! detrital nitrogen
	453	zdet(ji,jj) = max(0.,trn(ji,jj,jk,jpdet))
	454	!! dissolved inorganic nitrogen
	455	zdin(ji,jj) = max(0.,trn(ji,jj,jk,jpdin))
	456	!! dissolved silicic acid
	457	zsil(ji,jj) = max(0.,trn(ji,jj,jk,jpsil))
	458	!! dissolved "iron"
	459	zfer(ji,jj) = max(0.,trn(ji,jj,jk,jpfer))
[8074]	460	ENDIF
[8441]	461	ENDDO
	462	ENDDO
[7224]	463
[8441]	464	# if defined key_roam
[9257]	465	!! extra MEDUSA-2 tracers
[8441]	466	DO jj = 2,jpjm1
	467	DO ji = 2,jpim1
	468	if (tmask(ji,jj,jk) == 1) then
	469	!! detrital carbon
	470	zdtc(ji,jj) = max(0.,trn(ji,jj,jk,jpdtc))
	471	!! dissolved inorganic carbon
[9257]	472	zdic(ji,jj) = trn(ji,jj,jk,jpdic)
[8074]	473	!! alkalinity
[9257]	474	zalk(ji,jj) = trn(ji,jj,jk,jpalk)
[8441]	475	!! oxygen
	476	zoxy(ji,jj) = max(0.,trn(ji,jj,jk,jpoxy))
	477	# if defined key_axy_carbchem && defined key_mocsy
	478	!! phosphate via DIN and Redfield
	479	zpho(ji,jj) = max(0.,trn(ji,jj,jk,jpdin)) / 16.0
[8074]	480	# endif
	481	!!
[8441]	482	!! also need physical parameters for gas exchange
	483	!! calculations
	484	ztmp(ji,jj) = tsn(ji,jj,jk,jp_tem)
	485	zsal(ji,jj) = tsn(ji,jj,jk,jp_sal)
	486	ENDIF
	487	ENDDO
	488	ENDDO
[5726]	489	# else
[9257]	490	!! diagnostic MEDUSA-1 detrital carbon tracer
[8441]	491	DO jj = 2,jpjm1
	492	DO ji = 2,jpim1
[9257]	493	IF (tmask(ji,jj,jk) == 1) THEN
[8441]	494	!! implicit detrital carbon
	495	zdtc(ji,jj) = zdet(ji,jj) * xthetad
	496	ENDIF
	497	ENDDO
	498	ENDDO
[5726]	499	# endif
[9257]	500
	501	# if defined key_roam
	502	!! ---------------------------------------------
	503	!! JPALM -- 14-12-2017 -- Here, before any exeptionnal crazy value is
	504	!! removed, we want to tell to the Master Processor that this
	505	!! Exceptionnal value did exist
	506	!!
	507	Call trc_bio_check(kt)
	508
	509	!!================================================================
	510	!! AXY (03/11/17): check input fields
	511	!! tracer values that exceed thresholds can cause carbonate system
	512	!! failures when passed to MOCSY; temporary temperature excursions
	513	!! in recent UKESM0.8 runs trigger such failures but are too short
	514	!! to have physical consequences; this section checks for such
	515	!! values and amends them using neighbouring values
	516	!!
	517	!! the check on temperature here is also carried out at the end of
	518	!! each model time step and anomalies are reported in the master
	519	!! ocean.output file; the error reporting below is strictly local
	520	!! to the relevant ocean.output_XXXX file so will not be visible
	521	!! unless all processors are reporting output
	522	!!================================================================
	523	!!
	524	DO jj = 2,jpjm1
	525	DO ji = 2,jpim1
	526	if (tmask(ji,jj,jk) == 1) then
	527	!! the thresholds for the four tracers are ...
	528	IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0 ) .OR. &
	529	(zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT. 50.0 ) .OR. &
	530	(zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3 ) .OR. &
	531	(zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3 ) ) THEN
	532	!!
	533	!! all tracer values are reported in the event of any excursion
	534	IF (lwp) THEN
	535	WRITE(charout,*) ' Tmp = ', ztmp(ji,jj)
	536	WRITE(charout2,*) ' Sal = ', zsal(ji,jj)
	537	WRITE(charout3,*) ' DIC = ', zdic(ji,jj)
	538	WRITE(charout4,*) ' Alk = ', zalk(ji,jj)
	539	WRITE(charout5,*) mig(ji), mjg(jj), jk, kt
	540	CALL ctl_warn( 'trc_bio_medusa: carbonate chemistry WARNING:', &
	541	TRIM(charout),TRIM(charout2),TRIM(charout3),TRIM(charout4), &
	542	'at i, j, k, kt:', TRIM(charout5) )
	543	ENDIF
	544	!!
	545	!! Detect, report and correct tracer excursions
	546	IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0) ) &
	547	CALL trc_bio_exceptionnal_fix( &
	548	tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_tem), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	549	'Tmp', -3.0, 40.0, ztmp(ji,jj) )
	550	!!
	551	IF ( (zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT. 50.0) ) &
	552	CALL trc_bio_exceptionnal_fix( &
	553	tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_sal), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	554	'Sal', 1.0, 50.0, zsal(ji,jj) )
	555	!!
	556	IF ( (zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3) ) &
	557	CALL trc_bio_exceptionnal_fix( &
	558	trn(ji-1:ji+1,jj-1:jj+1,jk,jpdic), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	559	'DIC', 100.0, 4.0E3, zdic(ji,jj) )
	560	!!
	561	IF ( (zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3) ) &
	562	CALL trc_bio_exceptionnal_fix( &
	563	trn(ji-1:ji+1,jj-1:jj+1,jk,jpalk), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	564	'Alk', 100.0, 4.0E3, zalk(ji,jj) )
	565	ENDIF
	566	ENDIF
	567	ENDDO
	568	ENDDO
	569	# endif
	570
[5726]	571	# if defined key_debug_medusa
[8441]	572	DO jj = 2,jpjm1
	573	DO ji = 2,jpim1
	574	if (tmask(ji,jj,jk) == 1) then
	575	if (idf.eq.1) then
	576	!! AXY (15/01/10)
	577	if (trn(ji,jj,jk,jpdin).lt.0.) then
	578	IF (lwp) write (numout,*) &
	579	'------------------------------'
	580	IF (lwp) write (numout,*) 'NEGATIVE DIN ERROR =', &
	581	trn(ji,jj,jk,jpdin)
	582	IF (lwp) write (numout,*) 'NEGATIVE DIN ERROR @', &
	583	ji, jj, jk, kt
[8074]	584	endif
[8441]	585	if (trn(ji,jj,jk,jpsil).lt.0.) then
	586	IF (lwp) write (numout,*) &
	587	'------------------------------'
	588	IF (lwp) write (numout,*) 'NEGATIVE SIL ERROR =', &
	589	trn(ji,jj,jk,jpsil)
	590	IF (lwp) write (numout,*) 'NEGATIVE SIL ERROR @', &
	591	ji, jj, jk, kt
[8074]	592	endif
[8441]	593	# if defined key_roam
	594	if (trn(ji,jj,jk,jpdic).lt.0.) then
	595	IF (lwp) write (numout,*) &
	596	'------------------------------'
	597	IF (lwp) write (numout,*) 'NEGATIVE DIC ERROR =', &
	598	trn(ji,jj,jk,jpdic)
	599	IF (lwp) write (numout,*) 'NEGATIVE DIC ERROR @', &
	600	ji, jj, jk, kt
[8074]	601	endif
[8441]	602	if (trn(ji,jj,jk,jpalk).lt.0.) then
	603	IF (lwp) write (numout,*) &
	604	'------------------------------'
	605	IF (lwp) write (numout,*) 'NEGATIVE ALK ERROR =', &
	606	trn(ji,jj,jk,jpalk)
	607	IF (lwp) write (numout,*) 'NEGATIVE ALK ERROR @', &
	608	ji, jj, jk, kt
[8074]	609	endif
[8441]	610	if (trn(ji,jj,jk,jpoxy).lt.0.) then
	611	IF (lwp) write (numout,*) &
	612	'------------------------------'
	613	IF (lwp) write (numout,*) 'NEGATIVE OXY ERROR =', &
	614	trn(ji,jj,jk,jpoxy)
	615	IF (lwp) write (numout,*) 'NEGATIVE OXY ERROR @', &
	616	ji, jj, jk, kt
[8074]	617	endif
[5726]	618	# endif
[8074]	619	endif
[8441]	620	ENDIF
	621	ENDDO
	622	ENDDO
[5726]	623	# endif
[8074]	624	# if defined key_debug_medusa
[8441]	625	! I'M NOT SURE THIS IS USEFUL NOW THAT I'VE SPLIT THE DO LOOP - marc 8/5/17
	626	! if (idf.eq.1.AND.idfval.eq.1) then
	627	! DO jj = 2,jpjm1
	628	! DO ji = 2,jpim1
	629	! if (tmask(ji,jj,jk) == 1) then
	630	! !! report state variable values
	631	! IF (lwp) write (numout,*) &
	632	! '------------------------------'
	633	! IF (lwp) write (numout,*) 'fthk(',jk,') = ', &
	634	! fse3t(ji,jj,jk)
	635	! IF (lwp) write (numout,*) 'zphn(',jk,') = ', zphn(ji,jj)
	636	! IF (lwp) write (numout,*) 'zphd(',jk,') = ', zphd(ji,jj)
	637	! IF (lwp) write (numout,*) 'zpds(',jk,') = ', zpds(ji,jj)
	638	! IF (lwp) write (numout,*) 'zzmi(',jk,') = ', zzmi(ji,jj)
	639	! IF (lwp) write (numout,*) 'zzme(',jk,') = ', zzme(ji,jj)
	640	! IF (lwp) write (numout,*) 'zdet(',jk,') = ', zdet(ji,jj)
	641	! IF (lwp) write (numout,*) 'zdin(',jk,') = ', zdin(ji,jj)
	642	! IF (lwp) write (numout,*) 'zsil(',jk,') = ', zsil(ji,jj)
	643	! IF (lwp) write (numout,*) 'zfer(',jk,') = ', zfer(ji,jj)
[8074]	644	# if defined key_roam
[8441]	645	! IF (lwp) write (numout,*) 'zdtc(',jk,') = ', zdtc(ji,jj)
	646	! IF (lwp) write (numout,*) 'zdic(',jk,') = ', zdic(ji,jj)
	647	! IF (lwp) write (numout,*) 'zalk(',jk,') = ', zalk(ji,jj)
	648	! IF (lwp) write (numout,*) 'zoxy(',jk,') = ', zoxy(ji,jj)
[8074]	649	# endif
[8441]	650	! ENDIF
	651	! ENDDO
	652	! ENDDO
	653	! endif
[8074]	654	# endif
	655
[8441]	656	# if defined key_debug_medusa
	657	! I'M NOT SURE THIS IS USEFUL NOW THAT I'VE SPLIT THE DO LOOP - marc 8/5/17
	658	! if (idf.eq.1.AND.idfval.eq.1.AND.jk.eq.1) then
	659	! DO jj = 2,jpjm1
	660	! DO ji = 2,jpim1
	661	! if (tmask(ji,jj,jk) == 1) then
	662	! IF (lwp) write (numout,*) &
	663	! '------------------------------'
	664	! IF (lwp) write (numout,*) 'dust = ', dust(ji,jj)
	665	! ENDIF
	666	! ENDDO
	667	! ENDDO
	668	! endif
[5726]	669	# endif
	670
[8441]	671	!!---------------------------------------------------------------
	672	!! Calculate air-sea gas exchange and river inputs
	673	!!---------------------------------------------------------------
	674	IF ( jk == 1 ) THEN
	675	CALL air_sea( kt )
	676	ENDIF
[5726]	677
[8441]	678	!!---------------------------------------------------------------
	679	!! Phytoplankton growth, zooplankton grazing and miscellaneous
	680	!! plankton losses.
	681	!!---------------------------------------------------------------
	682	CALL plankton( jk )
[7224]	683
[8441]	684	!!---------------------------------------------------------------
	685	!! Iron chemistry and scavenging
	686	!!---------------------------------------------------------------
	687	CALL iron_chem_scav( jk )
[6022]	688
[8441]	689	!!---------------------------------------------------------------
	690	!! Detritus processes
	691	!!---------------------------------------------------------------
	692	CALL detritus( jk, iball )
[5726]	693
[8441]	694	!!---------------------------------------------------------------
	695	!! Updating tracers
	696	!!---------------------------------------------------------------
	697	CALL bio_medusa_update( kt, jk )
[5726]	698
[8441]	699	!!---------------------------------------------------------------
	700	!! Diagnostics
	701	!!---------------------------------------------------------------
	702	CALL bio_medusa_diag( jk )
[8074]	703
[8441]	704	!!-------------------------------------------------------
	705	!! 2d specific k level diags
	706	!!-------------------------------------------------------
[8442]	707	IF( lk_iomput ) THEN
[8441]	708	CALL bio_medusa_diag_slice( jk )
	709	ENDIF
[5726]	710
[8074]	711	!! CLOSE vertical loop
	712	ENDDO
	713
[8441]	714	!!------------------------------------------------------------------
	715	!! Final calculations for diagnostics
	716	!!------------------------------------------------------------------
	717	CALL bio_medusa_fin( kt )
[8074]	718
[5726]	719	# if defined key_debug_medusa
[8074]	720	IF(lwp) WRITE(numout,*) ' MEDUSA exiting trc_bio_medusa at kt =', kt
	721	CALL flush(numout)
[5726]	722	# endif
	723
[8074]	724	END SUBROUTINE trc_bio_medusa
[5726]	725
[9257]	726
	727
	728	SUBROUTINE trc_bio_exceptionnal_fix(tiny_var, tiny_mask, var_nm, mini, maxi, varout)
	729	!! JPALM (27/10/17): This function is called only when abnormal values that
	730	!! could break the model's carbonate system are fed to MEDUSA
	731	!!
	732	!! The function calculates an average tracer value based on the 3x3 cell
	733	!! neighbourhood around the abnormal cell, and reports this back
	734	!!
	735	!! Tracer array values are not modified, but MEDUSA uses "corrected" values
	736	!! in its calculations
	737	!!
	738	!! temporary variables
	739	REAL(wp), INTENT( in ), DIMENSION(3,3) :: tiny_var, tiny_mask
	740	CHARACTER(25), INTENT( in ) :: var_nm
	741	REAL(wp), INTENT( in ) :: mini, maxi
	742	REAL(wp), INTENT( out ) :: varout
	743	REAL(wp) :: sumtsn, tsnavg
	744	INTEGER :: summask
	745	CHARACTER(25) :: charout1, charout2
	746	CHARACTER(60) :: charout3, charout4
	747	INTEGER :: ii,ij
	748
	749	!! point to the center of the 3*3 zoom-grid, to check around
	750	ii = 2
	751	ij = 2
	752	!! Print surounding values to check if isolated Crazy value or
	753	!! General error
	754	IF(lwp) THEN
	755	WRITE(numout,*) &
	756	'----------------------------------------------------------------------'
	757	WRITE(numout,*) &
	758	'trc_bio_medusa: 3x3 neighbourhood surrounding abnormal ', TRIM(var_nm)
	759	WRITE(numout,9100) &
	760	3, tiny_var(ii-1,ij+1), tiny_var(ii ,ij+1), tiny_var(ii+1,ij+1)
	761	WRITE(numout,9100) &
	762	2, tiny_var(ii-1,ij ), tiny_var(ii ,ij ), tiny_var(ii+1,ij )
	763	WRITE(numout,9100) &
	764	1, tiny_var(ii-1,ij-1), tiny_var(ii ,ij-1), tiny_var(ii+1,ij-1)
	765	WRITE(numout,*) &
	766	'trc_bio_medusa: 3x3 land-sea neighbourhood, tmask'
	767	WRITE(numout,9100) &
	768	3, tiny_mask(ii-1,ij+1), tiny_mask(ii ,ij+1), tiny_mask(ii+1,ij+1)
	769	WRITE(numout,9100) &
	770	2, tiny_mask(ii-1,ij ), tiny_mask(ii ,ij ), tiny_mask(ii+1,ij )
	771	WRITE(numout,9100) &
	772	1, tiny_mask(ii-1,ij-1), tiny_mask(ii ,ij-1), tiny_mask(ii+1,ij-1)
	773	ENDIF
	774	!! Correct out of range values
	775	sumtsn = ( tiny_mask(ii-1,ij+1) * tiny_var(ii-1,ij+1) ) + &
	776	( tiny_mask(ii ,ij+1) * tiny_var(ii ,ij+1) ) + &
	777	( tiny_mask(ii+1,ij+1) * tiny_var(ii+1,ij+1) ) + &
	778	( tiny_mask(ii-1,ij ) * tiny_var(ii-1,ij ) ) + &
	779	( tiny_mask(ii+1,ij ) * tiny_var(ii+1,ij ) ) + &
	780	( tiny_mask(ii-1,ij-1) * tiny_var(ii-1,ij-1) ) + &
	781	( tiny_mask(ii ,ij-1) * tiny_var(ii ,ij-1) ) + &
	782	( tiny_mask(ii+1,ij-1) * tiny_var(ii+1,ij-1) )
	783	!!
	784	summask = tiny_mask(ii-1,ij+1) + tiny_mask(ii ,ij+1) + &
	785	tiny_mask(ii+1,ij+1) + tiny_mask(ii-1,ij ) + &
	786	tiny_mask(ii+1,ij ) + tiny_mask(ii-1,ij-1) + &
	787	tiny_mask(ii ,ij-1) + tiny_mask(ii+1,ij-1)
	788	!!
	789	IF ( summask .GT. 0 ) THEN
	790	tsnavg = ( sumtsn / summask )
	791	varout = MAX( MIN( maxi, tsnavg), mini )
	792	ELSE
	793	IF (ztmp(ii,ij) .LT. mini ) varout = mini
	794	IF (ztmp(ii,ij) .GT. maxi ) varout = maxi
	795	ENDIF
	796	!!
	797	IF (lwp) THEN
	798	WRITE(charout1,9200) tiny_var(ii,ij)
	799	WRITE(charout2,9200) varout
	800	WRITE(charout3,*) ' ', charout1, ' -> ', charout2
	801	WRITE(charout4,*) ' Tracer: ', trim(var_nm)
	802	!!
	803	WRITE(numout,) 'trc_bio_medusa: * EXCEPTIONAL VALUE SWITCHING **'
	804	WRITE(numout,*) charout4
	805	WRITE(numout,*) charout3
	806	WRITE(numout,*) '----------------------------------------------------------------------'
	807	WRITE(numout,*) ' '
	808	ENDIF
	809
	810	9100 FORMAT('Row:', i1, ' ', e16.6, ' ', e16.6, ' ', e16.6)
	811	9200 FORMAT(e16.6)
	812
	813	END SUBROUTINE trc_bio_exceptionnal_fix
	814
	815	SUBROUTINE trc_bio_check(kt)
	816	!!-----------------------------------
	817	!! JPALM -- 14-12-2017 -- Still dealing with this micro-boil/carb failure
	818	!! problem. The model is now able to correct a local
	819	!! crazy value. but does it silently.
	820	!! We need to spread the word to the master processor. we
	821	!! don't want the model to correct values without telling us
	822	!! This module will tell at least when crazy DIC or
	823	!! ALK appears.
	824	!!-------------------------------------
	825	REAL(wp) :: zmax, zmin ! temporary scalars
	826	INTEGER :: ji,jj ! dummy loop
	827	INTEGER :: ii,ij ! temporary scalars
	828	INTEGER, DIMENSION(2) :: ilocs !
	829	INTEGER, INTENT( in ) :: kt
	830	!!
	831	!!==========================
	832	!! DIC Check
	833	zmax = -5.0 ! arbitrary low maximum value
	834	zmin = 4.0E4 ! arbitrary high minimum value
	835	DO jj = 2, jpjm1
	836	DO ji = 2,jpim1
	837	IF( tmask(ji,jj,1) == 1) THEN
	838	zmax = MAX(zmax,zdic(ji,jj)) ! find local maximum
	839	zmin = MIN(zmin,zdic(ji,jj)) ! find local minimum
	840	ENDIF
	841	END DO
	842	END DO
	843	IF( lk_mpp ) CALL mpp_max( zmax ) ! max over the global domain
	844	IF( lk_mpp ) CALL mpp_min( zmin ) ! min over the global domain
	845	!
	846	IF( zmax .GT. 4.0E3) THEN ! we've got a problem
	847	IF (lk_mpp) THEN
	848	CALL mpp_maxloc ( zdic(:,:),tmask(:,:,1), zmax, ii,ij )
	849	ELSE
	850	ilocs = MAXLOC( zdic(:,:), mask = tmask(:,:,1) == 1. )
	851	ii = ilocs(1) + nimpp - 1
	852	ij = ilocs(2) + njmpp - 1
	853	ENDIF
	854	!
	855	IF(lwp) THEN
	856	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
	857	WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC > 4000 '
	858	WRITE(numout,9600) kt, zmax, ii, ij
	859	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	860	ENDIF
	861	ENDIF
	862	!
	863	IF( zmin .LE. 0.0) THEN ! we've got a problem
	864	IF (lk_mpp) THEN
	865	CALL mpp_minloc ( zdic(:,:),tmask(:,:,1), zmin, ii,ij )
	866	ELSE
	867	ilocs = MINLOC( zdic(:,:), mask = tmask(:,:,1) == 1. )
	868	ii = ilocs(1) + nimpp - 1
	869	ij = ilocs(2) + njmpp - 1
	870	ENDIF
	871	!
	872	IF(lwp) THEN
	873	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
	874	WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC <= 0 '
	875	WRITE(numout,9700) kt, zmin, ii, ij
	876	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	877	ENDIF
	878	ENDIF
	879	!!
	880	!!==========================
	881	!! ALKALINITY Check
	882	zmax = -5.0 ! arbitrary low maximum value
	883	zmin = 4.0E4 ! arbitrary high minimum value
	884	DO jj = 2, jpjm1
	885	DO ji = 2,jpim1
	886	IF( tmask(ji,jj,1) == 1) THEN
	887	zmax = MAX(zmax,zalk(ji,jj)) ! find local maximum
	888	zmin = MIN(zmin,zalk(ji,jj)) ! find local minimum
	889	ENDIF
	890	END DO
	891	END DO
	892	IF( lk_mpp ) CALL mpp_max( zmax ) ! max over the global domain
	893	IF( lk_mpp ) CALL mpp_min( zmin ) ! min over the global domain
	894	!
	895	IF( zmax .GT. 4.0E3) THEN ! we've got a problem
	896	IF (lk_mpp) THEN
	897	CALL mpp_maxloc ( zalk(:,:),tmask(:,:,1), zmax, ii,ij )
	898	ELSE
	899	ilocs = MAXLOC( zalk(:,:), mask = tmask(:,:,1) == 1. )
	900	ii = ilocs(1) + nimpp - 1
	901	ij = ilocs(2) + njmpp - 1
	902	ENDIF
	903	!
	904	IF(lwp) THEN
	905	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
	906	WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity > 4000 '
	907	WRITE(numout,9800) kt, zmax, ii, ij
	908	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	909	ENDIF
	910	ENDIF
	911	!
	912	IF( zmin .LE. 0.0) THEN ! we've got a problem
	913	IF (lk_mpp) THEN
	914	CALL mpp_minloc ( zalk(:,:),tmask(:,:,1), zmin, ii,ij )
	915	ELSE
	916	ilocs = MINLOC( zalk(:,:), mask = tmask(:,:,1) == 1. )
	917	ii = ilocs(1) + nimpp - 1
	918	ij = ilocs(2) + njmpp - 1
	919	ENDIF
	920	!
	921	IF(lwp) THEN
	922	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
	923	WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity <= 0 '
	924	WRITE(numout,9900) kt, zmin, ii, ij
	925	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	926	ENDIF
	927	ENDIF
	928
	929
	930	9600 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max DIC: ',f16.10,', i j: ',2i5)
	931	9700 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min DIC: ',f16.10,', i j: ',2i5)
	932	9800 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max ALK: ',f16.10,', i j: ',2i5)
	933	9900 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min ALK: ',f16.10,', i j: ',2i5)
	934
	935	END SUBROUTINE trc_bio_check
	936
	937
[5726]	938	#else
[8441]	939	!!=====================================================================
[5726]	940	!! Dummy module : No MEDUSA bio-model
[8441]	941	!!=====================================================================
[5726]	942	CONTAINS
	943	SUBROUTINE trc_bio_medusa( kt ) ! Empty routine
	944	INTEGER, INTENT( in ) :: kt
	945	WRITE(,) 'trc_bio_medusa: You should not have seen this print! error?', kt
	946	END SUBROUTINE trc_bio_medusa
	947	#endif
	948
[8441]	949	!!=====================================================================
[5726]	950	END MODULE trcbio_medusa

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

Download in other formats: