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 @ 9304

Last change on this file since 9304 was 9304, checked in by jpalmier, 6 years ago
JPALM -- Update trc_bio_check print statement - does not check surface values only
File size: 43.5 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
[9298]	184	REAL(wp) :: fq3,fq4
	185	REAL(wp) :: this_y, this_d, this_s, fyear
[8074]	186	!!
[9257]	187	!! T and S check temporary variable
	188	REAL(wp) :: sumtsn, tsnavg
	189	INTEGER :: summask
	190	CHARACTER(40) :: charout, charout2, charout3, charout4, charout5
	191	!!
[8441]	192	!!------------------------------------------------------------------
[5726]	193
[5841]	194	# if defined key_debug_medusa
[8074]	195	IF (lwp) write (numout,*) 'trc_bio_medusa: variables defined'
	196	CALL flush(numout)
[5841]	197	# endif
	198
[8074]	199	!! AXY (20/11/14): alter this to report on first MEDUSA call
	200	!! IF( kt == nit000 ) THEN
	201	IF( kt == nittrc000 ) THEN
	202	IF(lwp) WRITE(numout,*)
	203	IF(lwp) WRITE(numout,*) ' trc_bio: MEDUSA bio-model'
	204	IF(lwp) WRITE(numout,*) ' ~~~~~~~'
	205	IF(lwp) WRITE(numout,*) ' kt =',kt
	206	ENDIF
[5726]	207
[8074]	208	!! AXY (13/01/12): is benthic model properly interactive? 0 = no, 1 = yes
	209	ibenthic = 1
[5726]	210
[8441]	211	!!------------------------------------------------------------------
[8074]	212	!! b0 is present for debugging purposes; using b0 = 0 sets the tendency
	213	!! terms of all biological equations to 0.
[8441]	214	!!------------------------------------------------------------------
[8074]	215	!!
	216	!! AXY (03/09/14): probably not the smartest move ever, but it'll fit
	217	!! the bill for now; another item on the things-to-sort-
	218	!! out-in-the-future list ...
[5726]	219	# if defined key_kill_medusa
[8074]	220	b0 = 0.
[5726]	221	# else
[8074]	222	b0 = 1.
[5726]	223	# endif
[8441]	224	!!------------------------------------------------------------------
[8074]	225	!! fast detritus ballast scheme (0 = no; 1 = yes)
	226	!! alternative to ballast scheme is same scheme but with no ballast
	227	!! protection (not dissimilar to Martin et al., 1987)
[8441]	228	!!------------------------------------------------------------------
[8074]	229	!!
	230	iball = 1
[5726]	231
[8441]	232	!!------------------------------------------------------------------
[8074]	233	!! full flux diagnostics (0 = no; 1 = yes); appear in ocean.output
	234	!! these should only be used in 1D since they give comprehensive
	235	!! output for ecological functions in the model; primarily used in
	236	!! debugging
[8441]	237	!!------------------------------------------------------------------
[8074]	238	!!
	239	idf = 0
	240	!!
	241	!! timer mechanism
	242	if (kt/120*120.eq.kt) then
	243	idfval = 1
	244	else
	245	idfval = 0
	246	endif
[5726]	247
[8441]	248	!!------------------------------------------------------------------
	249	!! Initialise arrays to zero and set up arrays for diagnostics
	250	!!------------------------------------------------------------------
	251	CALL bio_medusa_init( kt )
[5726]	252
	253	# if defined key_axy_nancheck
[8441]	254	DO jn = jp_msa0,jp_msa1
[8074]	255	!! fq0 = MINVAL(trn(:,:,:,jn))
	256	!! fq1 = MAXVAL(trn(:,:,:,jn))
	257	fq2 = SUM(trn(:,:,:,jn))
[8441]	258	!! if (lwp) write (numout,'(a,2i6,3(1x,1pe15.5))') 'NAN-CHECK', &
	259	!! kt, jn, fq0, fq1, fq2
	260	!! AXY (30/01/14): much to our surprise, the next line doesn't
	261	!! work on HECTOR and has been replaced here with
	262	!! a specialist routine
[8074]	263	!! if (fq2 /= fq2 ) then
	264	if ( ieee_is_nan( fq2 ) ) then
	265	!! there's a NaN here
[8441]	266	if (lwp) write(numout,*) 'NAN detected in field', jn, &
	267	'at time', kt, 'at position:'
[8074]	268	DO jk = 1,jpk
	269	DO jj = 1,jpj
	270	DO ji = 1,jpi
	271	!! AXY (30/01/14): "isnan" problem on HECTOR
	272	!! if (trn(ji,jj,jk,jn) /= trn(ji,jj,jk,jn)) then
	273	if ( ieee_is_nan( trn(ji,jj,jk,jn) ) ) then
[8441]	274	if (lwp) write (numout,'(a,1pe12.2,4i6)') &
	275	'NAN-CHECK', tmask(ji,jj,jk), ji, jj, jk, jn
[8074]	276	endif
	277	enddo
[7224]	278	enddo
	279	enddo
[8074]	280	CALL ctl_stop( 'trcbio_medusa, NAN in incoming tracer field' )
	281	endif
	282	ENDDO
	283	CALL flush(numout)
[5726]	284	# endif
	285
[5841]	286	# if defined key_debug_medusa
[8441]	287	IF (lwp) write (numout,*) &
	288	'trc_bio_medusa: variables initialised and checked'
[8074]	289	CALL flush(numout)
[5841]	290	# endif
	291
[5726]	292	# if defined key_roam
[8441]	293	!!------------------------------------------------------------------
[8074]	294	!! calculate atmospheric pCO2
[8441]	295	!!------------------------------------------------------------------
[8074]	296	!!
	297	!! what's atmospheric pCO2 doing? (data start in 1859)
[9298]	298	!! AXY - JPALM new interpolation scheme usinf nyear_len
	299	this_y = real(nyear)
	300	this_d = real(nday_year)
	301	this_s = real(nsec_day)
	302	!!
	303	fyear = this_y + ((this_d - 1) + (this_s / (60. * 60. * 24.))) / real(nyear_len(1))
	304	!!
	305	IF ( ( kt == nittrc000 ) .AND. lwp ) THEN
	306	WRITE(numout,*) ' MEDUSA year =', this_y
	307	WRITE(numout,*) ' Year length =', real(nyear_len(1))
	308	WRITE(numout,*) ' MEDUSA nday_year =', this_d
	309	WRITE(numout,*) ' MEDUSA nsec_day =', this_s
	310	ENDIF
	311	!!
	312	!! different case test
	313	IF (fyear .LE. co2_yinit) THEN
	314	!! before first record -- pre-industrial value
	315	f_xco2a(:,:) = hist_pco2(1)
	316	ELSEIF (fyear .GE. co2_yend) THEN
	317	!! after last record - continue to use the last value
	318	f_xco2a(:,:) = hist_pco2(int(co2_yend - co2_yinit + 1.) )
	319	ELSE
[8074]	320	!! just right
[9298]	321	iyr1 = int(fyear - co2_yinit) + 1
	322	iyr2 = iyr1 + 1
	323	fq3 = fyear - real(iyr1) - co2_yinit + 1.
	324	fq4 = ((1 - fq3) * hist_pco2(iyr1)) + (fq3 * hist_pco2(iyr2))
[8441]	325	f_xco2a(:,:) = fq4
[9298]	326	!!
	327	IF ( ( kt == nittrc000 ) .AND. lwp ) THEN
	328	WRITE(numout,*) ' MEDUSA year1 =', iyr1
	329	WRITE(numout,*) ' MEDUSA year2 =', iyr2
	330	WRITE(numout,*) ' xCO2 year1 =', hist_pco2(iyr1)
	331	WRITE(numout,*) ' xCO2 year2 =', hist_pco2(iyr2)
	332	WRITE(numout,*) ' Year2 weight =', fq3
	333	ENDIF
	334	ENDIF
[8074]	335	# if defined key_axy_pi_co2
[8441]	336	!! OCMIP pre-industrial pCO2
	337	!! f_xco2a(:,:) = 284.725 !! CMIP5 pre-industrial pCO2
[9285]	338	f_xco2a(:,:) = 284.317 !! CMIP6 pre-industrial pCO2
[8074]	339	# endif
	340	!! AXY (29/01/14): remove surplus diagnostics
[9298]	341	IF ( ( kt == nittrc000 ) .AND. lwp ) &
	342	WRITE(numout,*) ' final atm xCO2 =', f_xco2a(1,1)
[5726]	343	# endif
	344
[5841]	345	# if defined key_debug_medusa
[8074]	346	IF (lwp) write (numout,*) 'trc_bio_medusa: ready for carbonate chemistry'
	347	IF (lwp) write (numout,*) 'trc_bio_medusa: kt = ', kt
	348	IF (lwp) write (numout,*) 'trc_bio_medusa: nittrc000 = ', nittrc000
	349	CALL flush(numout)
[5841]	350	# endif
	351
[5726]	352	# if defined key_roam
[8074]	353	!! AXY (20/11/14): alter to call on first MEDUSA timestep and then every
	354	!! month (this is hardwired as 960 timesteps but should
	355	!! be calculated and done properly
	356	!! IF( kt == nit000 .or. mod(kt,1920) == 0 ) THEN
	357	!! IF( kt == nittrc000 .or. mod(kt,960) == 0 ) THEN
	358	!!=============================
	359	!! Jpalm -- 07-10-2016 -- need to change carb-chem frequency call :
[8441]	360	!! we don't want to call on the first time-step of all run
	361	!! submission, but only on the very first time-step, and
	362	!! then every month. So we call on nittrc000 if not
	363	!! restarted run, else if one month after last call.
	364	!! Assume one month is 30d --> 36002430 : 2592000s
	365	!! try to call carb-chem at 1st month's tm-stp :
	366	!! x * 30d + 1*rdt(i.e: mod = rdt)
[8074]	367	!! ++ need to pass carb-chem output var through restarts
[8224]	368	If ( (kt == nittrc000 .AND. .NOT.ln_rsttr) .OR. &
	369	((86400mod(nn_date0,100) + mod(ktrdt,2592000.)) == rdt) ) THEN
[8441]	370	!!---------------------------------------------------------------
[8074]	371	!! Calculate the carbonate chemistry for the whole ocean on the first
	372	!! simulation timestep and every month subsequently; the resulting 3D
	373	!! field of omega calcite is used to determine the depth of the CCD
[8441]	374	!!---------------------------------------------------------------
	375	CALL carb_chem( kt )
	376
[8074]	377	ENDIF
[5726]	378	# endif
	379
[5841]	380	# if defined key_debug_medusa
[8074]	381	IF (lwp) write (numout,*) 'trc_bio_medusa: ready for full domain calculations'
	382	CALL flush(numout)
[5841]	383	# endif
	384
[8441]	385	!!------------------------------------------------------------------
[8074]	386	!! MEDUSA has unified equation through the water column
	387	!! (Diff. from LOBSTER which has two sets: bio- and non-bio layers)
	388	!! Statement below in LOBSTER is different: DO jk = 1, jpkbm1
[8441]	389	!!------------------------------------------------------------------
[8074]	390	!!
	391	!! NOTE: the ordering of the loops below differs from that of some other
	392	!! models; looping over the vertical dimension is the outermost loop and
	393	!! this complicates some calculations (e.g. storage of vertical fluxes
	394	!! that can otherwise be done via a singular variable require 2D fields
	395	!! here); however, these issues are relatively easily resolved, but the
	396	!! loops CANNOT be reordered without potentially causing code efficiency
	397	!! problems (e.g. array indexing means that reordering the loops would
	398	!! require skipping between widely-spaced memory location; potentially
	399	!! outside those immediately cached)
	400	!!
	401	!! OPEN vertical loop
	402	DO jk = 1,jpk
	403	!! OPEN horizontal loops
	404	DO jj = 2,jpjm1
[8441]	405	DO ji = 2,jpim1
	406	!! OPEN wet point IF..THEN loop
	407	if (tmask(ji,jj,jk) == 1) then
	408	!!======================================================
	409	!! SETUP LOCAL GRID CELL
	410	!!======================================================
[8074]	411	!!
[8441]	412	!!------------------------------------------------------
	413	!! Some notes on grid vertical structure
	414	!! - fsdepw(ji,jj,jk) is the depth of the upper surface of
	415	!! level jk
	416	!! - fsde3w(ji,jj,jk) is approximately the midpoint of
	417	!! level jk
	418	!! - fse3t(ji,jj,jk) is the thickness of level jk
	419	!!------------------------------------------------------
[8074]	420	!!
[8441]	421	!! AXY (01/03/10): set up level depth (bottom of level)
	422	fdep1(ji,jj) = fsdepw(ji,jj,jk) + fse3t(ji,jj,jk)
[8074]	423	!!
[8441]	424	!! set up model tracers
	425	!! negative values of state variables are not allowed to
	426	!! contribute to the calculated fluxes
	427	!! non-diatom chlorophyll
	428	zchn(ji,jj) = max(0.,trn(ji,jj,jk,jpchn))
	429	!! diatom chlorophyll
	430	zchd(ji,jj) = max(0.,trn(ji,jj,jk,jpchd))
	431	!! non-diatoms
	432	zphn(ji,jj) = max(0.,trn(ji,jj,jk,jpphn))
	433	!! diatoms
	434	zphd(ji,jj) = max(0.,trn(ji,jj,jk,jpphd))
	435	!! diatom silicon
	436	zpds(ji,jj) = max(0.,trn(ji,jj,jk,jppds))
	437	!! AXY (28/01/10): probably need to take account of
	438	!! chl/biomass connection
	439	if (zchn(ji,jj).eq.0.) zphn(ji,jj) = 0.
	440	if (zchd(ji,jj).eq.0.) zphd(ji,jj) = 0.
	441	if (zphn(ji,jj).eq.0.) zchn(ji,jj) = 0.
	442	if (zphd(ji,jj).eq.0.) zchd(ji,jj) = 0.
	443	!! AXY (23/01/14): duh - why did I forget diatom silicon?
	444	if (zpds(ji,jj).eq.0.) zphd(ji,jj) = 0.
	445	if (zphd(ji,jj).eq.0.) zpds(ji,jj) = 0.
	446	ENDIF
	447	ENDDO
	448	ENDDO
[7224]	449
[8441]	450	DO jj = 2,jpjm1
	451	DO ji = 2,jpim1
	452	if (tmask(ji,jj,jk) == 1) then
	453	!! microzooplankton
	454	zzmi(ji,jj) = max(0.,trn(ji,jj,jk,jpzmi))
	455	!! mesozooplankton
	456	zzme(ji,jj) = max(0.,trn(ji,jj,jk,jpzme))
	457	!! detrital nitrogen
	458	zdet(ji,jj) = max(0.,trn(ji,jj,jk,jpdet))
	459	!! dissolved inorganic nitrogen
	460	zdin(ji,jj) = max(0.,trn(ji,jj,jk,jpdin))
	461	!! dissolved silicic acid
	462	zsil(ji,jj) = max(0.,trn(ji,jj,jk,jpsil))
	463	!! dissolved "iron"
	464	zfer(ji,jj) = max(0.,trn(ji,jj,jk,jpfer))
[8074]	465	ENDIF
[8441]	466	ENDDO
	467	ENDDO
[7224]	468
[8441]	469	# if defined key_roam
[9257]	470	!! extra MEDUSA-2 tracers
[8441]	471	DO jj = 2,jpjm1
	472	DO ji = 2,jpim1
	473	if (tmask(ji,jj,jk) == 1) then
	474	!! detrital carbon
	475	zdtc(ji,jj) = max(0.,trn(ji,jj,jk,jpdtc))
	476	!! dissolved inorganic carbon
[9257]	477	zdic(ji,jj) = trn(ji,jj,jk,jpdic)
[8074]	478	!! alkalinity
[9257]	479	zalk(ji,jj) = trn(ji,jj,jk,jpalk)
[8441]	480	!! oxygen
	481	zoxy(ji,jj) = max(0.,trn(ji,jj,jk,jpoxy))
	482	# if defined key_axy_carbchem && defined key_mocsy
	483	!! phosphate via DIN and Redfield
	484	zpho(ji,jj) = max(0.,trn(ji,jj,jk,jpdin)) / 16.0
[8074]	485	# endif
	486	!!
[8441]	487	!! also need physical parameters for gas exchange
	488	!! calculations
	489	ztmp(ji,jj) = tsn(ji,jj,jk,jp_tem)
	490	zsal(ji,jj) = tsn(ji,jj,jk,jp_sal)
	491	ENDIF
	492	ENDDO
	493	ENDDO
[5726]	494	# else
[9257]	495	!! diagnostic MEDUSA-1 detrital carbon tracer
[8441]	496	DO jj = 2,jpjm1
	497	DO ji = 2,jpim1
[9257]	498	IF (tmask(ji,jj,jk) == 1) THEN
[8441]	499	!! implicit detrital carbon
	500	zdtc(ji,jj) = zdet(ji,jj) * xthetad
	501	ENDIF
	502	ENDDO
	503	ENDDO
[5726]	504	# endif
[9257]	505
	506	# if defined key_roam
	507	!! ---------------------------------------------
	508	!! JPALM -- 14-12-2017 -- Here, before any exeptionnal crazy value is
	509	!! removed, we want to tell to the Master Processor that this
	510	!! Exceptionnal value did exist
	511	!!
[9304]	512	Call trc_bio_check(kt, jk)
[9257]	513
	514	!!================================================================
	515	!! AXY (03/11/17): check input fields
	516	!! tracer values that exceed thresholds can cause carbonate system
	517	!! failures when passed to MOCSY; temporary temperature excursions
	518	!! in recent UKESM0.8 runs trigger such failures but are too short
	519	!! to have physical consequences; this section checks for such
	520	!! values and amends them using neighbouring values
	521	!!
	522	!! the check on temperature here is also carried out at the end of
	523	!! each model time step and anomalies are reported in the master
	524	!! ocean.output file; the error reporting below is strictly local
	525	!! to the relevant ocean.output_XXXX file so will not be visible
	526	!! unless all processors are reporting output
	527	!!================================================================
	528	!!
	529	DO jj = 2,jpjm1
	530	DO ji = 2,jpim1
	531	if (tmask(ji,jj,jk) == 1) then
	532	!! the thresholds for the four tracers are ...
	533	IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0 ) .OR. &
	534	(zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT. 50.0 ) .OR. &
	535	(zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3 ) .OR. &
	536	(zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3 ) ) THEN
	537	!!
	538	!! all tracer values are reported in the event of any excursion
	539	IF (lwp) THEN
	540	WRITE(charout,*) ' Tmp = ', ztmp(ji,jj)
	541	WRITE(charout2,*) ' Sal = ', zsal(ji,jj)
	542	WRITE(charout3,*) ' DIC = ', zdic(ji,jj)
	543	WRITE(charout4,*) ' Alk = ', zalk(ji,jj)
	544	WRITE(charout5,*) mig(ji), mjg(jj), jk, kt
	545	CALL ctl_warn( 'trc_bio_medusa: carbonate chemistry WARNING:', &
	546	TRIM(charout),TRIM(charout2),TRIM(charout3),TRIM(charout4), &
	547	'at i, j, k, kt:', TRIM(charout5) )
	548	ENDIF
	549	!!
	550	!! Detect, report and correct tracer excursions
	551	IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0) ) &
	552	CALL trc_bio_exceptionnal_fix( &
	553	tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_tem), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	554	'Tmp', -3.0, 40.0, ztmp(ji,jj) )
	555	!!
	556	IF ( (zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT. 50.0) ) &
	557	CALL trc_bio_exceptionnal_fix( &
	558	tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_sal), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	559	'Sal', 1.0, 50.0, zsal(ji,jj) )
	560	!!
	561	IF ( (zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3) ) &
	562	CALL trc_bio_exceptionnal_fix( &
	563	trn(ji-1:ji+1,jj-1:jj+1,jk,jpdic), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	564	'DIC', 100.0, 4.0E3, zdic(ji,jj) )
	565	!!
	566	IF ( (zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3) ) &
	567	CALL trc_bio_exceptionnal_fix( &
	568	trn(ji-1:ji+1,jj-1:jj+1,jk,jpalk), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
	569	'Alk', 100.0, 4.0E3, zalk(ji,jj) )
	570	ENDIF
	571	ENDIF
	572	ENDDO
	573	ENDDO
	574	# endif
	575
[5726]	576	# if defined key_debug_medusa
[8441]	577	DO jj = 2,jpjm1
	578	DO ji = 2,jpim1
	579	if (tmask(ji,jj,jk) == 1) then
	580	if (idf.eq.1) then
	581	!! AXY (15/01/10)
	582	if (trn(ji,jj,jk,jpdin).lt.0.) then
	583	IF (lwp) write (numout,*) &
	584	'------------------------------'
	585	IF (lwp) write (numout,*) 'NEGATIVE DIN ERROR =', &
	586	trn(ji,jj,jk,jpdin)
	587	IF (lwp) write (numout,*) 'NEGATIVE DIN ERROR @', &
	588	ji, jj, jk, kt
[8074]	589	endif
[8441]	590	if (trn(ji,jj,jk,jpsil).lt.0.) then
	591	IF (lwp) write (numout,*) &
	592	'------------------------------'
	593	IF (lwp) write (numout,*) 'NEGATIVE SIL ERROR =', &
	594	trn(ji,jj,jk,jpsil)
	595	IF (lwp) write (numout,*) 'NEGATIVE SIL ERROR @', &
	596	ji, jj, jk, kt
[8074]	597	endif
[8441]	598	# if defined key_roam
	599	if (trn(ji,jj,jk,jpdic).lt.0.) then
	600	IF (lwp) write (numout,*) &
	601	'------------------------------'
	602	IF (lwp) write (numout,*) 'NEGATIVE DIC ERROR =', &
	603	trn(ji,jj,jk,jpdic)
	604	IF (lwp) write (numout,*) 'NEGATIVE DIC ERROR @', &
	605	ji, jj, jk, kt
[8074]	606	endif
[8441]	607	if (trn(ji,jj,jk,jpalk).lt.0.) then
	608	IF (lwp) write (numout,*) &
	609	'------------------------------'
	610	IF (lwp) write (numout,*) 'NEGATIVE ALK ERROR =', &
	611	trn(ji,jj,jk,jpalk)
	612	IF (lwp) write (numout,*) 'NEGATIVE ALK ERROR @', &
	613	ji, jj, jk, kt
[8074]	614	endif
[8441]	615	if (trn(ji,jj,jk,jpoxy).lt.0.) then
	616	IF (lwp) write (numout,*) &
	617	'------------------------------'
	618	IF (lwp) write (numout,*) 'NEGATIVE OXY ERROR =', &
	619	trn(ji,jj,jk,jpoxy)
	620	IF (lwp) write (numout,*) 'NEGATIVE OXY ERROR @', &
	621	ji, jj, jk, kt
[8074]	622	endif
[5726]	623	# endif
[8074]	624	endif
[8441]	625	ENDIF
	626	ENDDO
	627	ENDDO
[5726]	628	# endif
[8074]	629	# if defined key_debug_medusa
[8441]	630	! I'M NOT SURE THIS IS USEFUL NOW THAT I'VE SPLIT THE DO LOOP - marc 8/5/17
	631	! if (idf.eq.1.AND.idfval.eq.1) then
	632	! DO jj = 2,jpjm1
	633	! DO ji = 2,jpim1
	634	! if (tmask(ji,jj,jk) == 1) then
	635	! !! report state variable values
	636	! IF (lwp) write (numout,*) &
	637	! '------------------------------'
	638	! IF (lwp) write (numout,*) 'fthk(',jk,') = ', &
	639	! fse3t(ji,jj,jk)
	640	! IF (lwp) write (numout,*) 'zphn(',jk,') = ', zphn(ji,jj)
	641	! IF (lwp) write (numout,*) 'zphd(',jk,') = ', zphd(ji,jj)
	642	! IF (lwp) write (numout,*) 'zpds(',jk,') = ', zpds(ji,jj)
	643	! IF (lwp) write (numout,*) 'zzmi(',jk,') = ', zzmi(ji,jj)
	644	! IF (lwp) write (numout,*) 'zzme(',jk,') = ', zzme(ji,jj)
	645	! IF (lwp) write (numout,*) 'zdet(',jk,') = ', zdet(ji,jj)
	646	! IF (lwp) write (numout,*) 'zdin(',jk,') = ', zdin(ji,jj)
	647	! IF (lwp) write (numout,*) 'zsil(',jk,') = ', zsil(ji,jj)
	648	! IF (lwp) write (numout,*) 'zfer(',jk,') = ', zfer(ji,jj)
[8074]	649	# if defined key_roam
[8441]	650	! IF (lwp) write (numout,*) 'zdtc(',jk,') = ', zdtc(ji,jj)
	651	! IF (lwp) write (numout,*) 'zdic(',jk,') = ', zdic(ji,jj)
	652	! IF (lwp) write (numout,*) 'zalk(',jk,') = ', zalk(ji,jj)
	653	! IF (lwp) write (numout,*) 'zoxy(',jk,') = ', zoxy(ji,jj)
[8074]	654	# endif
[8441]	655	! ENDIF
	656	! ENDDO
	657	! ENDDO
	658	! endif
[8074]	659	# endif
	660
[8441]	661	# if defined key_debug_medusa
	662	! I'M NOT SURE THIS IS USEFUL NOW THAT I'VE SPLIT THE DO LOOP - marc 8/5/17
	663	! if (idf.eq.1.AND.idfval.eq.1.AND.jk.eq.1) then
	664	! DO jj = 2,jpjm1
	665	! DO ji = 2,jpim1
	666	! if (tmask(ji,jj,jk) == 1) then
	667	! IF (lwp) write (numout,*) &
	668	! '------------------------------'
	669	! IF (lwp) write (numout,*) 'dust = ', dust(ji,jj)
	670	! ENDIF
	671	! ENDDO
	672	! ENDDO
	673	! endif
[5726]	674	# endif
	675
[8441]	676	!!---------------------------------------------------------------
	677	!! Calculate air-sea gas exchange and river inputs
	678	!!---------------------------------------------------------------
	679	IF ( jk == 1 ) THEN
	680	CALL air_sea( kt )
	681	ENDIF
[5726]	682
[8441]	683	!!---------------------------------------------------------------
	684	!! Phytoplankton growth, zooplankton grazing and miscellaneous
	685	!! plankton losses.
	686	!!---------------------------------------------------------------
	687	CALL plankton( jk )
[7224]	688
[8441]	689	!!---------------------------------------------------------------
	690	!! Iron chemistry and scavenging
	691	!!---------------------------------------------------------------
	692	CALL iron_chem_scav( jk )
[6022]	693
[8441]	694	!!---------------------------------------------------------------
	695	!! Detritus processes
	696	!!---------------------------------------------------------------
	697	CALL detritus( jk, iball )
[5726]	698
[8441]	699	!!---------------------------------------------------------------
	700	!! Updating tracers
	701	!!---------------------------------------------------------------
	702	CALL bio_medusa_update( kt, jk )
[5726]	703
[8441]	704	!!---------------------------------------------------------------
	705	!! Diagnostics
	706	!!---------------------------------------------------------------
	707	CALL bio_medusa_diag( jk )
[8074]	708
[8441]	709	!!-------------------------------------------------------
	710	!! 2d specific k level diags
	711	!!-------------------------------------------------------
[8442]	712	IF( lk_iomput ) THEN
[8441]	713	CALL bio_medusa_diag_slice( jk )
	714	ENDIF
[5726]	715
[8074]	716	!! CLOSE vertical loop
	717	ENDDO
	718
[8441]	719	!!------------------------------------------------------------------
	720	!! Final calculations for diagnostics
	721	!!------------------------------------------------------------------
	722	CALL bio_medusa_fin( kt )
[8074]	723
[5726]	724	# if defined key_debug_medusa
[8074]	725	IF(lwp) WRITE(numout,*) ' MEDUSA exiting trc_bio_medusa at kt =', kt
	726	CALL flush(numout)
[5726]	727	# endif
	728
[8074]	729	END SUBROUTINE trc_bio_medusa
[5726]	730
[9257]	731
	732
	733	SUBROUTINE trc_bio_exceptionnal_fix(tiny_var, tiny_mask, var_nm, mini, maxi, varout)
	734	!! JPALM (27/10/17): This function is called only when abnormal values that
	735	!! could break the model's carbonate system are fed to MEDUSA
	736	!!
	737	!! The function calculates an average tracer value based on the 3x3 cell
	738	!! neighbourhood around the abnormal cell, and reports this back
	739	!!
	740	!! Tracer array values are not modified, but MEDUSA uses "corrected" values
	741	!! in its calculations
	742	!!
	743	!! temporary variables
	744	REAL(wp), INTENT( in ), DIMENSION(3,3) :: tiny_var, tiny_mask
	745	CHARACTER(25), INTENT( in ) :: var_nm
	746	REAL(wp), INTENT( in ) :: mini, maxi
	747	REAL(wp), INTENT( out ) :: varout
	748	REAL(wp) :: sumtsn, tsnavg
	749	INTEGER :: summask
	750	CHARACTER(25) :: charout1, charout2
	751	CHARACTER(60) :: charout3, charout4
	752	INTEGER :: ii,ij
	753
	754	!! point to the center of the 3*3 zoom-grid, to check around
	755	ii = 2
	756	ij = 2
	757	!! Print surounding values to check if isolated Crazy value or
	758	!! General error
	759	IF(lwp) THEN
	760	WRITE(numout,*) &
	761	'----------------------------------------------------------------------'
	762	WRITE(numout,*) &
	763	'trc_bio_medusa: 3x3 neighbourhood surrounding abnormal ', TRIM(var_nm)
	764	WRITE(numout,9100) &
	765	3, tiny_var(ii-1,ij+1), tiny_var(ii ,ij+1), tiny_var(ii+1,ij+1)
	766	WRITE(numout,9100) &
	767	2, tiny_var(ii-1,ij ), tiny_var(ii ,ij ), tiny_var(ii+1,ij )
	768	WRITE(numout,9100) &
	769	1, tiny_var(ii-1,ij-1), tiny_var(ii ,ij-1), tiny_var(ii+1,ij-1)
	770	WRITE(numout,*) &
	771	'trc_bio_medusa: 3x3 land-sea neighbourhood, tmask'
	772	WRITE(numout,9100) &
	773	3, tiny_mask(ii-1,ij+1), tiny_mask(ii ,ij+1), tiny_mask(ii+1,ij+1)
	774	WRITE(numout,9100) &
	775	2, tiny_mask(ii-1,ij ), tiny_mask(ii ,ij ), tiny_mask(ii+1,ij )
	776	WRITE(numout,9100) &
	777	1, tiny_mask(ii-1,ij-1), tiny_mask(ii ,ij-1), tiny_mask(ii+1,ij-1)
	778	ENDIF
	779	!! Correct out of range values
	780	sumtsn = ( 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	( tiny_mask(ii-1,ij ) * tiny_var(ii-1,ij ) ) + &
	784	( tiny_mask(ii+1,ij ) * tiny_var(ii+1,ij ) ) + &
	785	( tiny_mask(ii-1,ij-1) * tiny_var(ii-1,ij-1) ) + &
	786	( tiny_mask(ii ,ij-1) * tiny_var(ii ,ij-1) ) + &
	787	( tiny_mask(ii+1,ij-1) * tiny_var(ii+1,ij-1) )
	788	!!
	789	summask = tiny_mask(ii-1,ij+1) + tiny_mask(ii ,ij+1) + &
	790	tiny_mask(ii+1,ij+1) + tiny_mask(ii-1,ij ) + &
	791	tiny_mask(ii+1,ij ) + tiny_mask(ii-1,ij-1) + &
	792	tiny_mask(ii ,ij-1) + tiny_mask(ii+1,ij-1)
	793	!!
	794	IF ( summask .GT. 0 ) THEN
	795	tsnavg = ( sumtsn / summask )
	796	varout = MAX( MIN( maxi, tsnavg), mini )
	797	ELSE
	798	IF (ztmp(ii,ij) .LT. mini ) varout = mini
	799	IF (ztmp(ii,ij) .GT. maxi ) varout = maxi
	800	ENDIF
	801	!!
	802	IF (lwp) THEN
	803	WRITE(charout1,9200) tiny_var(ii,ij)
	804	WRITE(charout2,9200) varout
	805	WRITE(charout3,*) ' ', charout1, ' -> ', charout2
	806	WRITE(charout4,*) ' Tracer: ', trim(var_nm)
	807	!!
	808	WRITE(numout,) 'trc_bio_medusa: * EXCEPTIONAL VALUE SWITCHING **'
	809	WRITE(numout,*) charout4
	810	WRITE(numout,*) charout3
	811	WRITE(numout,*) '----------------------------------------------------------------------'
	812	WRITE(numout,*) ' '
	813	ENDIF
	814
	815	9100 FORMAT('Row:', i1, ' ', e16.6, ' ', e16.6, ' ', e16.6)
	816	9200 FORMAT(e16.6)
	817
	818	END SUBROUTINE trc_bio_exceptionnal_fix
	819
[9304]	820	SUBROUTINE trc_bio_check(kt, jk)
[9257]	821	!!-----------------------------------
	822	!! JPALM -- 14-12-2017 -- Still dealing with this micro-boil/carb failure
	823	!! problem. The model is now able to correct a local
	824	!! crazy value. but does it silently.
	825	!! We need to spread the word to the master processor. we
	826	!! don't want the model to correct values without telling us
	827	!! This module will tell at least when crazy DIC or
	828	!! ALK appears.
	829	!!-------------------------------------
	830	REAL(wp) :: zmax, zmin ! temporary scalars
	831	INTEGER :: ji,jj ! dummy loop
	832	INTEGER :: ii,ij ! temporary scalars
	833	INTEGER, DIMENSION(2) :: ilocs !
[9304]	834	INTEGER, INTENT( in ) :: kt, jk
[9257]	835	!!
	836	!!==========================
	837	!! DIC Check
	838	zmax = -5.0 ! arbitrary low maximum value
	839	zmin = 4.0E4 ! arbitrary high minimum value
	840	DO jj = 2, jpjm1
	841	DO ji = 2,jpim1
	842	IF( tmask(ji,jj,1) == 1) THEN
	843	zmax = MAX(zmax,zdic(ji,jj)) ! find local maximum
	844	zmin = MIN(zmin,zdic(ji,jj)) ! find local minimum
	845	ENDIF
	846	END DO
	847	END DO
	848	IF( lk_mpp ) CALL mpp_max( zmax ) ! max over the global domain
	849	IF( lk_mpp ) CALL mpp_min( zmin ) ! min over the global domain
	850	!
	851	IF( zmax .GT. 4.0E3) THEN ! we've got a problem
	852	IF (lk_mpp) THEN
	853	CALL mpp_maxloc ( zdic(:,:),tmask(:,:,1), zmax, ii,ij )
	854	ELSE
	855	ilocs = MAXLOC( zdic(:,:), mask = tmask(:,:,1) == 1. )
	856	ii = ilocs(1) + nimpp - 1
	857	ij = ilocs(2) + njmpp - 1
	858	ENDIF
	859	!
	860	IF(lwp) THEN
	861	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
[9304]	862	WRITE(numout,*) 'trc_bio:tracer anomaly: DIC concentration > 4000 '
	863	WRITE(numout,9600) kt, zmax, ii, ij, jk
[9257]	864	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	865	ENDIF
	866	ENDIF
	867	!
	868	IF( zmin .LE. 0.0) THEN ! we've got a problem
	869	IF (lk_mpp) THEN
	870	CALL mpp_minloc ( zdic(:,:),tmask(:,:,1), zmin, ii,ij )
	871	ELSE
	872	ilocs = MINLOC( zdic(:,:), mask = tmask(:,:,1) == 1. )
	873	ii = ilocs(1) + nimpp - 1
	874	ij = ilocs(2) + njmpp - 1
	875	ENDIF
	876	!
	877	IF(lwp) THEN
	878	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
[9304]	879	WRITE(numout,*) 'trc_bio:tracer anomaly: DIC concentration <= 0 '
	880	WRITE(numout,9700) kt, zmin, ii, ij, jk
[9257]	881	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	882	ENDIF
	883	ENDIF
	884	!!
	885	!!==========================
	886	!! ALKALINITY Check
	887	zmax = -5.0 ! arbitrary low maximum value
	888	zmin = 4.0E4 ! arbitrary high minimum value
	889	DO jj = 2, jpjm1
	890	DO ji = 2,jpim1
	891	IF( tmask(ji,jj,1) == 1) THEN
	892	zmax = MAX(zmax,zalk(ji,jj)) ! find local maximum
	893	zmin = MIN(zmin,zalk(ji,jj)) ! find local minimum
	894	ENDIF
	895	END DO
	896	END DO
	897	IF( lk_mpp ) CALL mpp_max( zmax ) ! max over the global domain
	898	IF( lk_mpp ) CALL mpp_min( zmin ) ! min over the global domain
	899	!
	900	IF( zmax .GT. 4.0E3) THEN ! we've got a problem
	901	IF (lk_mpp) THEN
	902	CALL mpp_maxloc ( zalk(:,:),tmask(:,:,1), zmax, ii,ij )
	903	ELSE
	904	ilocs = MAXLOC( zalk(:,:), mask = tmask(:,:,1) == 1. )
	905	ii = ilocs(1) + nimpp - 1
	906	ij = ilocs(2) + njmpp - 1
	907	ENDIF
	908	!
	909	IF(lwp) THEN
	910	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
[9304]	911	WRITE(numout,*) 'trc_bio:tracer anomaly: ALK concentration > 4000 '
	912	WRITE(numout,9800) kt, zmax, ii, ij, jk
[9257]	913	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	914	ENDIF
	915	ENDIF
	916	!
	917	IF( zmin .LE. 0.0) THEN ! we've got a problem
	918	IF (lk_mpp) THEN
	919	CALL mpp_minloc ( zalk(:,:),tmask(:,:,1), zmin, ii,ij )
	920	ELSE
	921	ilocs = MINLOC( zalk(:,:), mask = tmask(:,:,1) == 1. )
	922	ii = ilocs(1) + nimpp - 1
	923	ij = ilocs(2) + njmpp - 1
	924	ENDIF
	925	!
	926	IF(lwp) THEN
	927	WRITE(numout,) 'trc_bio:tracer anomaly: ** WARNING ***'
[9304]	928	WRITE(numout,*) 'trc_bio:tracer anomaly: ALK concentration <= 0 '
	929	WRITE(numout,9900) kt, zmin, ii, ij, jk
[9257]	930	WRITE(numout,) 'trc_bio:tracer anomaly: ** END OF WARNING ***'
	931	ENDIF
	932	ENDIF
	933
	934
[9304]	935	9600 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max DIC: ',f16.10,', i j k: ',3i5)
	936	9700 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min DIC: ',f16.10,', i j k: ',3i5)
	937	9800 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max ALK: ',f16.10,', i j k: ',3i5)
	938	9900 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min ALK: ',f16.10,', i j k: ',3i5)
[9257]	939
	940	END SUBROUTINE trc_bio_check
	941
	942
[5726]	943	#else
[8441]	944	!!=====================================================================
[5726]	945	!! Dummy module : No MEDUSA bio-model
[8441]	946	!!=====================================================================
[5726]	947	CONTAINS
	948	SUBROUTINE trc_bio_medusa( kt ) ! Empty routine
	949	INTEGER, INTENT( in ) :: kt
	950	WRITE(,) 'trc_bio_medusa: You should not have seen this print! error?', kt
	951	END SUBROUTINE trc_bio_medusa
	952	#endif
	953
[8441]	954	!!=====================================================================
[5726]	955	END MODULE trcbio_medusa

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

Download in other formats: