Index: /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml
===================================================================
--- /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml (revision 9256)
+++ /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml (revision 9257)
@@ -1645,34 +1645,53 @@
-
+
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-
-
-
-
-
-
-
-
-
-
-
-
@@ -1685,11 +1704,14 @@
-
+
+
+
+
@@ -1697,11 +1719,8 @@
-
-
-
Index: /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/stpctl.F90
===================================================================
--- /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/stpctl.F90 (revision 9256)
+++ /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/stpctl.F90 (revision 9257)
@@ -58,4 +58,6 @@
INTEGER :: ii, ij, ik ! temporary integers
REAL(wp) :: zumax, zsmin, zssh2 ! temporary scalars
+ REAL(wp) :: ztmax, ztmin ! Scalar to get temperature extreme
+ ! values and warn if they're out of Range
INTEGER, DIMENSION(3) :: ilocu !
INTEGER, DIMENSION(2) :: ilocs !
@@ -148,4 +150,59 @@
9500 FORMAT (' kt=',i6,' min SSS: ',1pg11.4,', i j: ',2i5)
+! ====================================================================================================
+! ====================================================================================================
+ ! !AXY (25/10/17)
+ ! !* Test max/min limits of temperature
+ ! ! ----------------------------------
+ ztmax = -5.e0 ! arbitrary low maximum value
+ ztmin = 100.e0 ! arbitrary high minimum value
+ DO jj = 2, jpjm1
+ DO ji = 2, jpim1
+ IF( tmask(ji,jj,1) == 1) THEN
+ ztmax = MAX(ztmax,tsn(ji,jj,1,jp_tem)) ! find local maximum
+ ztmin = MIN(ztmin,tsn(ji,jj,1,jp_tem)) ! find local minimum
+ ENDIF
+ END DO
+ END DO
+ IF( lk_mpp ) CALL mpp_max( ztmax ) ! max over the global domain
+ IF( lk_mpp ) CALL mpp_min( ztmin ) ! min over the global domain
+ !
+ IF( ztmax > 40.) THEN ! we've got a problem
+ IF (lk_mpp) THEN
+ CALL mpp_maxloc ( tsn(:,:,1,jp_tem),tmask(:,:,1), ztmax, ii,ij )
+ ELSE
+ ilocs = MAXLOC( tsn(:,:,1,jp_tem), mask = tmask(:,:,1) == 1.e0 )
+ ii = ilocs(1) + nimpp - 1
+ ij = ilocs(2) + njmpp - 1
+ ENDIF
+ !
+ IF(lwp) THEN
+ WRITE(numout,*) 'stp_ctl:tracer anomaly: ***** WARNING *****'
+ WRITE(numout,*) 'stp_ctl:tracer anomaly: sea surface temperature > 40C'
+ WRITE(numout,9600) kt, ztmax, ii, ij
+ WRITE(numout,*) 'stp_ctl:tracer anomaly: ***** END OF WARNING *****'
+ ENDIF
+ ENDIF
+ !
+ IF( ztmin < -3.) THEN ! we've got a problem
+ IF (lk_mpp) THEN
+ CALL mpp_minloc ( tsn(:,:,1,jp_tem),tmask(:,:,1), ztmin, ii,ij )
+ ELSE
+ ilocs = MINLOC( tsn(:,:,1,jp_tem), mask = tmask(:,:,1) == 1.e0 )
+ ii = ilocs(1) + nimpp - 1
+ ij = ilocs(2) + njmpp - 1
+ ENDIF
+ !
+ IF(lwp) THEN
+ WRITE(numout,*) 'stp_ctl:tracer anomaly: ***** WARNING *****'
+ WRITE(numout,*) 'stp_ctl:tracer anomaly: sea surface temperature < -3C'
+ WRITE(numout,9700) kt, ztmin, ii, ij
+ WRITE(numout,*) 'stp_ctl:tracer anomaly: ***** END OF WARNING *****'
+ ENDIF
+ ENDIF
+9600 FORMAT ('stp_ctl:tracer anomaly: kt=',i6,' max SST: ',f16.10,', i j: ',2i5)
+9700 FORMAT ('stp_ctl:tracer anomaly: kt=',i6,' min SST: ',f16.10,', i j: ',2i5)
+! ====================================================================================================
+! ====================================================================================================
IF( lk_c1d ) RETURN ! No log file in case of 1D vertical configuration
Index: /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90
===================================================================
--- /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90 (revision 9256)
+++ /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90 (revision 9257)
@@ -62,5 +62,5 @@
# endif
zchd, zchn, zdin, zsil
- USE dom_oce, ONLY: e3t_0, e3t_n, gphit, tmask
+ USE dom_oce, ONLY: e3t_0, e3t_n, gphit, tmask, mig, mjg
# if defined key_iomput
USE iom, ONLY: lk_iomput
@@ -91,4 +91,6 @@
USE trcoxy_medusa, ONLY: trc_oxy_medusa
# endif
+ USE lib_mpp, ONLY: ctl_stop
+ USE trcstat, ONLY: trc_rst_dia_stat
!!* Substitution
@@ -122,4 +124,19 @@
# if defined key_roam
+ !! init
+ f_fco2w(:,:) = 0.0
+ f_fco2atm(:,:) = 0.0
+ f_schmidtco2(:,:) = 0.0
+ f_kwco2(:,:) = 0.0
+ f_co2starair(:,:) = 0.0
+ f_dpco2(:,:) = 0.0
+ f_rhosw(:,:) = 0.0
+ f_K0(:,:) = 0.0
+ !! air pressure (atm); ultimately this will use air
+ !! pressure at the base of the UKESM1 atmosphere
+ !!
+ f_pp0(:,:) = 1.0
+
+
!!-----------------------------------------------------------
!! Air-sea gas exchange
@@ -134,8 +151,17 @@
DO ji = 2,jpim1
!! OPEN wet point IF..THEN loop
- if (tmask(ji,jj,1) == 1) then
+ IF (tmask(ji,jj,1) == 1) then
IF (lk_oasis) THEN
!! use 2D atm xCO2 from atm coupling
f_xco2a(ji,jj) = PCO2a_in_cpl(ji,jj)
+ !!!
+ !!! Jpalm test on atm xCO2
+ IF ( (f_xco2a(ji,jj) > 1500.0 ).OR.(f_xco2a(ji,jj) < 100.0 ) ) THEN
+ IF(lwp) THEN
+ WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj), &
+ ' -- ji =', mig(ji),' jj = ', mjg(jj)
+ ENDIF
+ CALL ctl_stop( 'MEDUSA - Air-Sea :', 'unrealistic atm xCO2 ' )
+ ENDIF
ENDIF
!!
@@ -162,19 +188,49 @@
'air-sea: carb-chem kt = ', kt
CALL flush(numout)
+ !! JPALM add carb print:
+ call trc_rst_dia_stat(f_xco2a(:,:), 'f_xco2a')
+ call trc_rst_dia_stat(wndm(:,:), 'wndm')
+ call trc_rst_dia_stat(f_kw660(:,:), 'f_kw660')
+ call trc_rst_dia_stat(ztmp(:,:), 'ztmp')
+ call trc_rst_dia_stat(zsal(:,:), 'zsal')
+ call trc_rst_dia_stat(zalk(:,:), 'zalk')
+ call trc_rst_dia_stat(zdic(:,:), 'zdic')
+ call trc_rst_dia_stat(zsil(:,:), 'zsil')
+ call trc_rst_dia_stat(zpho(:,:), 'zpho')
# endif
DO jj = 2,jpjm1
DO ji = 2,jpim1
if (tmask(ji,jj,1) == 1) then
- !! air pressure (atm); ultimately this will use air
- !! pressure at the base of the UKESM1 atmosphere
- !!
- f_pp0(ji,jj) = 1.0
- !!
- !! IF(lwp) WRITE(numout,*) ' MEDUSA ztmp =', ztmp(ji,jj)
- !! IF(lwp) WRITE(numout,*) ' MEDUSA wndm =', wndm(ji,jj)
- !! IF(lwp) WRITE(numout,*) ' MEDUSA fr_i =', fr_i(ji,jj)
!!
# if defined key_axy_carbchem
# if defined key_mocsy
+ !! Jpalm -- 12-09-2017 -- add extra check after reccurent
+ !! carbonate failure in the coupled run.
+ !! must be associated to air-sea flux or air xCO2...
+ !! Check MOCSY inputs
+ IF ( (zsal(ji,jj) > 75.0 ).OR.(zsal(ji,jj) < 0.0 ) .OR. &
+ (ztmp(ji,jj) > 50.0 ).OR.(ztmp(ji,jj) < -20.0 ) .OR. &
+ (zalk(ji,jj) > 35.0E2 ).OR.(zalk(ji,jj) <= 0.0 ) .OR. &
+ (zdic(ji,jj) > 35.0E2 ).OR.(zdic(ji,jj) <= 0.0 ) .OR. &
+ (f_kw660(ji,jj) > 1.0E-2 ).OR.(f_kw660(ji,jj) < 0.0 ) ) THEN
+ IF(lwp) THEN
+ WRITE(numout,*) ' surface T = ',ztmp(ji,jj)
+ WRITE(numout,*) ' surface S = ',zsal(ji,jj)
+ WRITE(numout,*) ' surface ALK = ',zalk(ji,jj)
+ WRITE(numout,*) ' surface DIC = ',zdic(ji,jj)
+ WRITE(numout,*) ' KW660 = ',f_kw660(ji,jj)
+ WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj)
+ WRITE(numout,*) ' surface pco2w = ',f_pco2w(ji,jj)
+ WRITE(numout,*) ' surface fco2w = ',f_fco2w(ji,jj)
+ WRITE(numout,*) ' surface fco2a = ',f_fco2atm(ji,jj)
+ WRITE(numout,*) ' surface co2flx = ',f_co2flux(ji,jj)
+ WRITE(numout,*) ' surface dpco2 = ',f_dpco2(ji,jj)
+ WRITE(numout,*) ' MOCSY input: ji =', mig(ji),' jj = ', mjg(jj), &
+ ' kt = ', kt
+ WRITE(numout,*) 'MEDUSA - Air-Sea INPUT: unrealistic surface Carb. Chemistry'
+ ENDIF
+ CALL ctl_stop( 'MEDUSA - Air-Sea INPUT: ', &
+ 'unrealistic surface Carb. Chemistry -- INPUTS' )
+ ENDIF
!!
!! AXY (22/06/15): use Orr & Epitalon (2015) MOCSY-2 carbonate
@@ -201,8 +257,47 @@
f_TALK(ji,jj) = (zalk(ji,jj) / f_rhosw(ji,jj)) * 1000.
f_dcf(ji,jj) = f_rhosw(ji,jj)
+ !! Jpalm -- 12-09-2017 -- add extra check after reccurent
+ !! carbonate failure in the coupled run.
+ !! must be associated to air-sea flux or air xCO2...
+ !! Check MOCSY outputs
+ IF ( (f_pco2w(ji,jj) > 1.E4 ).OR.(f_pco2w(ji,jj) < 0.0 ) .OR. &
+ (f_fco2w(ji,jj) > 1.E4 ).OR.(f_fco2w(ji,jj) < 0.0 ) .OR. &
+ (f_fco2atm(ji,jj) > 1.E4 ).OR.(f_fco2atm(ji,jj) < 0.0 ) .OR. &
+ (f_co2flux(ji,jj) > 1.E-2 ).OR.(f_co2flux(ji,jj) < -1.E-2 ) .OR. &
+ (f_dpco2(ji,jj) > 1.E4 ).OR.(f_dpco2(ji,jj) < -1.E4 ) ) THEN
+ IF(lwp) THEN
+ WRITE(numout,*) ' surface T = ',ztmp(ji,jj)
+ WRITE(numout,*) ' surface S = ',zsal(ji,jj)
+ WRITE(numout,*) ' surface ALK = ',zalk(ji,jj)
+ WRITE(numout,*) ' surface DIC = ',zdic(ji,jj)
+ WRITE(numout,*) ' KW660 = ',f_kw660(ji,jj)
+ WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj)
+ WRITE(numout,*) ' surface pco2w = ',f_pco2w(ji,jj)
+ WRITE(numout,*) ' surface fco2w = ',f_fco2w(ji,jj)
+ WRITE(numout,*) ' surface fco2a = ',f_fco2atm(ji,jj)
+ WRITE(numout,*) ' surface co2flx = ',f_co2flux(ji,jj)
+ WRITE(numout,*) ' surface dpco2 = ',f_dpco2(ji,jj)
+ WRITE(numout,*) ' MOCSY output: ji =', mig(ji),' jj = ', mjg(jj), &
+ ' kt = ', kt
+ WRITE(numout,*) 'MEDUSA - Air-Sea OUTPUT: unrealistic surface Carb. Chemistry'
+ ENDIF
+ CALL ctl_stop( 'MEDUSA - Air-Sea OUTPUT: ', &
+ 'unrealistic surface Carb. Chemistry -- OUTPUTS' )
+ ENDIF
ENDIF
ENDDO
ENDDO
+# if defined key_debug_medusa
+ !! JPALM add carb print:
+ call trc_rst_dia_stat(f_pco2w(:,:), 'f_pco2w')
+ call trc_rst_dia_stat(f_fco2w(:,:), 'f_fco2w')
+ call trc_rst_dia_stat(f_fco2atm(:,:), 'f_fco2atm')
+ call trc_rst_dia_stat(f_schmidtco2(:,:), 'f_schmidtco2')
+ call trc_rst_dia_stat(f_kwco2(:,:), 'f_kwco2')
+ call trc_rst_dia_stat(f_co2starair(:,:), 'f_co2starair')
+ call trc_rst_dia_stat(f_co2flux(:,:), 'f_co2flux')
+ call trc_rst_dia_stat(f_dpco2(:,:), 'f_dpco2')
+# endif
# else
Index: /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90
===================================================================
--- /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90 (revision 9256)
+++ /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90 (revision 9257)
@@ -167,4 +167,23 @@
fslownflux(:,:) = 0.0
fslowcflux(:,:) = 0.0
+ !!
+ !! JPALM -- 21-09-2017 -- needed to debug air-sea carb
+ f_xco2a(:,:) = 0.0
+ f_pco2w(:,:) = 0.0
+ f_ph(:,:) = 0.0
+ f_kw660(:,:) = 0.0
+ ztmp(:,:) = 0.0
+ zsal(:,:) = 0.0
+ zalk(:,:) = 0.0
+ zdic(:,:) = 0.0
+ zsil(:,:) = 0.0
+ zpho(:,:) = 0.0
+ f_co2flux(:,:) = 0.0
+ f_pco2atm(:,:) = 0.0
+ f_h2co3(:,:) = 0.0
+ f_hco3(:,:) = 0.0
+ f_co3(:,:) = 0.0
+ f_omarg(:,:) = 0.0
+ f_omcal(:,:) = 0.0
!!
!! AXY (08/08/17): zero slow detritus fluxes
Index: /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90
===================================================================
--- /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90 (revision 9256)
+++ /branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90 (revision 9257)
@@ -23,4 +23,5 @@
!! - ! 2016-11 (A. Yool) Updated diags for CMIP6
!! - ! 2017-05 (A. Yool) Added extra DMS calculation
+ !! - ! 2017-11 (J. Palm, A. Yool) Diagnose tracer excursions
!!----------------------------------------------------------------------
!!
@@ -81,5 +82,6 @@
gdepw_0, gdepw_n, &
nday_year, nsec_day, nyear, &
- rdt, tmask
+ rdt, tmask, mig, mjg, nimpp, &
+ njmpp
USE in_out_manager, ONLY: lwp, numout, nn_date0
# if defined key_iomput
@@ -87,5 +89,7 @@
# endif
USE lbclnk, ONLY: lbc_lnk
- USE lib_mpp, ONLY: ctl_stop
+ USE lib_mpp, ONLY: mpp_max, mpp_maxloc, &
+ mpp_min, mpp_minloc, &
+ ctl_stop, ctl_warn, lk_mpp
USE oce, ONLY: tsb, tsn
USE par_kind, ONLY: wp
@@ -115,4 +119,5 @@
PUBLIC trc_bio_medusa ! called in trcsms_medusa.F90
+ PUBLIC trc_bio_exceptionnal_fix ! here
!!* Substitution
@@ -177,4 +182,9 @@
!! temporary variables
REAL(wp) :: fq0,fq1,fq2,fq3,fq4
+ !!
+ !! T and S check temporary variable
+ REAL(wp) :: sumtsn, tsnavg
+ INTEGER :: summask
+ CHARACTER(40) :: charout, charout2, charout3, charout4, charout5
!!
!!------------------------------------------------------------------
@@ -450,4 +460,5 @@
# if defined key_roam
+ !! extra MEDUSA-2 tracers
DO jj = 2,jpjm1
DO ji = 2,jpim1
@@ -456,7 +467,7 @@
zdtc(ji,jj) = max(0.,trn(ji,jj,jk,jpdtc))
!! dissolved inorganic carbon
- zdic(ji,jj) = max(0.,trn(ji,jj,jk,jpdic))
+ zdic(ji,jj) = trn(ji,jj,jk,jpdic)
!! alkalinity
- zalk(ji,jj) = max(0.,trn(ji,jj,jk,jpalk))
+ zalk(ji,jj) = trn(ji,jj,jk,jpalk)
!! oxygen
zoxy(ji,jj) = max(0.,trn(ji,jj,jk,jpoxy))
@@ -470,30 +481,12 @@
ztmp(ji,jj) = tsn(ji,jj,jk,jp_tem)
zsal(ji,jj) = tsn(ji,jj,jk,jp_sal)
- !!
- !! AXY (28/02/14): check input fields
- if (ztmp(ji,jj) .lt. -3.0 .or. ztmp(ji,jj) .gt. 40.0 ) then
- IF(lwp) WRITE(numout,*) &
- ' trc_bio_medusa: T WARNING 2D, ', &
- tsb(ji,jj,jk,jp_tem), tsn(ji,jj,jk,jp_tem), &
- ' at (', ji, ',', jj, ',', jk, ') at time', kt
- IF(lwp) WRITE(numout,*) &
- ' trc_bio_medusa: T SWITCHING 2D, ', &
- tsn(ji,jj,jk,jp_tem), ' -> ', tsb(ji,jj,jk,jp_tem)
- !! temperatur
- ztmp(ji,jj) = tsb(ji,jj,jk,jp_tem)
- endif
- if (zsal(ji,jj) .lt. 0.0 .or. zsal(ji,jj) .gt. 45.0 ) then
- IF(lwp) WRITE(numout,*) &
- ' trc_bio_medusa: S WARNING 2D, ', &
- tsb(ji,jj,jk,jp_sal), tsn(ji,jj,jk,jp_sal), &
- ' at (', ji, ',', jj, ',', jk, ') at time', kt
- endif
ENDIF
ENDDO
ENDDO
# else
+ !! diagnostic MEDUSA-1 detrital carbon tracer
DO jj = 2,jpjm1
DO ji = 2,jpim1
- if (tmask(ji,jj,jk) == 1) then
+ IF (tmask(ji,jj,jk) == 1) THEN
!! implicit detrital carbon
zdtc(ji,jj) = zdet(ji,jj) * xthetad
@@ -502,4 +495,75 @@
ENDDO
# endif
+
+# if defined key_roam
+ !! ---------------------------------------------
+ !! JPALM -- 14-12-2017 -- Here, before any exeptionnal crazy value is
+ !! removed, we want to tell to the Master Processor that this
+ !! Exceptionnal value did exist
+ !!
+ Call trc_bio_check(kt)
+
+ !!================================================================
+ !! AXY (03/11/17): check input fields
+ !! tracer values that exceed thresholds can cause carbonate system
+ !! failures when passed to MOCSY; temporary temperature excursions
+ !! in recent UKESM0.8 runs trigger such failures but are too short
+ !! to have physical consequences; this section checks for such
+ !! values and amends them using neighbouring values
+ !!
+ !! the check on temperature here is also carried out at the end of
+ !! each model time step and anomalies are reported in the master
+ !! ocean.output file; the error reporting below is strictly local
+ !! to the relevant ocean.output_XXXX file so will not be visible
+ !! unless all processors are reporting output
+ !!================================================================
+ !!
+ DO jj = 2,jpjm1
+ DO ji = 2,jpim1
+ if (tmask(ji,jj,jk) == 1) then
+ !! the thresholds for the four tracers are ...
+ IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0 ) .OR. &
+ (zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT. 50.0 ) .OR. &
+ (zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3 ) .OR. &
+ (zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3 ) ) THEN
+ !!
+ !! all tracer values are reported in the event of any excursion
+ IF (lwp) THEN
+ WRITE(charout,*) ' Tmp = ', ztmp(ji,jj)
+ WRITE(charout2,*) ' Sal = ', zsal(ji,jj)
+ WRITE(charout3,*) ' DIC = ', zdic(ji,jj)
+ WRITE(charout4,*) ' Alk = ', zalk(ji,jj)
+ WRITE(charout5,*) mig(ji), mjg(jj), jk, kt
+ CALL ctl_warn( 'trc_bio_medusa: carbonate chemistry WARNING:', &
+ TRIM(charout),TRIM(charout2),TRIM(charout3),TRIM(charout4), &
+ 'at i, j, k, kt:', TRIM(charout5) )
+ ENDIF
+ !!
+ !! Detect, report and correct tracer excursions
+ IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0) ) &
+ CALL trc_bio_exceptionnal_fix( &
+ tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_tem), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
+ 'Tmp', -3.0, 40.0, ztmp(ji,jj) )
+ !!
+ IF ( (zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT. 50.0) ) &
+ CALL trc_bio_exceptionnal_fix( &
+ tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_sal), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
+ 'Sal', 1.0, 50.0, zsal(ji,jj) )
+ !!
+ IF ( (zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3) ) &
+ CALL trc_bio_exceptionnal_fix( &
+ trn(ji-1:ji+1,jj-1:jj+1,jk,jpdic), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
+ 'DIC', 100.0, 4.0E3, zdic(ji,jj) )
+ !!
+ IF ( (zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3) ) &
+ CALL trc_bio_exceptionnal_fix( &
+ trn(ji-1:ji+1,jj-1:jj+1,jk,jpalk), tmask(ji-1:ji+1,jj-1:jj+1,jk), &
+ 'Alk', 100.0, 4.0E3, zalk(ji,jj) )
+ ENDIF
+ ENDIF
+ ENDDO
+ ENDDO
+# endif
+
# if defined key_debug_medusa
DO jj = 2,jpjm1
@@ -657,4 +721,216 @@
END SUBROUTINE trc_bio_medusa
+
+
+ SUBROUTINE trc_bio_exceptionnal_fix(tiny_var, tiny_mask, var_nm, mini, maxi, varout)
+ !! JPALM (27/10/17): This function is called only when abnormal values that
+ !! could break the model's carbonate system are fed to MEDUSA
+ !!
+ !! The function calculates an average tracer value based on the 3x3 cell
+ !! neighbourhood around the abnormal cell, and reports this back
+ !!
+ !! Tracer array values are not modified, but MEDUSA uses "corrected" values
+ !! in its calculations
+ !!
+ !! temporary variables
+ REAL(wp), INTENT( in ), DIMENSION(3,3) :: tiny_var, tiny_mask
+ CHARACTER(25), INTENT( in ) :: var_nm
+ REAL(wp), INTENT( in ) :: mini, maxi
+ REAL(wp), INTENT( out ) :: varout
+ REAL(wp) :: sumtsn, tsnavg
+ INTEGER :: summask
+ CHARACTER(25) :: charout1, charout2
+ CHARACTER(60) :: charout3, charout4
+ INTEGER :: ii,ij
+
+ !! point to the center of the 3*3 zoom-grid, to check around
+ ii = 2
+ ij = 2
+ !! Print surounding values to check if isolated Crazy value or
+ !! General error
+ IF(lwp) THEN
+ WRITE(numout,*) &
+ '----------------------------------------------------------------------'
+ WRITE(numout,*) &
+ 'trc_bio_medusa: 3x3 neighbourhood surrounding abnormal ', TRIM(var_nm)
+ WRITE(numout,9100) &
+ 3, tiny_var(ii-1,ij+1), tiny_var(ii ,ij+1), tiny_var(ii+1,ij+1)
+ WRITE(numout,9100) &
+ 2, tiny_var(ii-1,ij ), tiny_var(ii ,ij ), tiny_var(ii+1,ij )
+ WRITE(numout,9100) &
+ 1, tiny_var(ii-1,ij-1), tiny_var(ii ,ij-1), tiny_var(ii+1,ij-1)
+ WRITE(numout,*) &
+ 'trc_bio_medusa: 3x3 land-sea neighbourhood, tmask'
+ WRITE(numout,9100) &
+ 3, tiny_mask(ii-1,ij+1), tiny_mask(ii ,ij+1), tiny_mask(ii+1,ij+1)
+ WRITE(numout,9100) &
+ 2, tiny_mask(ii-1,ij ), tiny_mask(ii ,ij ), tiny_mask(ii+1,ij )
+ WRITE(numout,9100) &
+ 1, tiny_mask(ii-1,ij-1), tiny_mask(ii ,ij-1), tiny_mask(ii+1,ij-1)
+ ENDIF
+ !! Correct out of range values
+ sumtsn = ( tiny_mask(ii-1,ij+1) * tiny_var(ii-1,ij+1) ) + &
+ ( tiny_mask(ii ,ij+1) * tiny_var(ii ,ij+1) ) + &
+ ( tiny_mask(ii+1,ij+1) * tiny_var(ii+1,ij+1) ) + &
+ ( tiny_mask(ii-1,ij ) * tiny_var(ii-1,ij ) ) + &
+ ( tiny_mask(ii+1,ij ) * tiny_var(ii+1,ij ) ) + &
+ ( tiny_mask(ii-1,ij-1) * tiny_var(ii-1,ij-1) ) + &
+ ( tiny_mask(ii ,ij-1) * tiny_var(ii ,ij-1) ) + &
+ ( tiny_mask(ii+1,ij-1) * tiny_var(ii+1,ij-1) )
+ !!
+ summask = tiny_mask(ii-1,ij+1) + tiny_mask(ii ,ij+1) + &
+ tiny_mask(ii+1,ij+1) + tiny_mask(ii-1,ij ) + &
+ tiny_mask(ii+1,ij ) + tiny_mask(ii-1,ij-1) + &
+ tiny_mask(ii ,ij-1) + tiny_mask(ii+1,ij-1)
+ !!
+ IF ( summask .GT. 0 ) THEN
+ tsnavg = ( sumtsn / summask )
+ varout = MAX( MIN( maxi, tsnavg), mini )
+ ELSE
+ IF (ztmp(ii,ij) .LT. mini ) varout = mini
+ IF (ztmp(ii,ij) .GT. maxi ) varout = maxi
+ ENDIF
+ !!
+ IF (lwp) THEN
+ WRITE(charout1,9200) tiny_var(ii,ij)
+ WRITE(charout2,9200) varout
+ WRITE(charout3,*) ' ', charout1, ' -> ', charout2
+ WRITE(charout4,*) ' Tracer: ', trim(var_nm)
+ !!
+ WRITE(numout,*) 'trc_bio_medusa: ** EXCEPTIONAL VALUE SWITCHING **'
+ WRITE(numout,*) charout4
+ WRITE(numout,*) charout3
+ WRITE(numout,*) '----------------------------------------------------------------------'
+ WRITE(numout,*) ' '
+ ENDIF
+
+9100 FORMAT('Row:', i1, ' ', e16.6, ' ', e16.6, ' ', e16.6)
+9200 FORMAT(e16.6)
+
+ END SUBROUTINE trc_bio_exceptionnal_fix
+
+ SUBROUTINE trc_bio_check(kt)
+ !!-----------------------------------
+ !! JPALM -- 14-12-2017 -- Still dealing with this micro-boil/carb failure
+ !! problem. The model is now able to correct a local
+ !! crazy value. but does it silently.
+ !! We need to spread the word to the master processor. we
+ !! don't want the model to correct values without telling us
+ !! This module will tell at least when crazy DIC or
+ !! ALK appears.
+ !!-------------------------------------
+ REAL(wp) :: zmax, zmin ! temporary scalars
+ INTEGER :: ji,jj ! dummy loop
+ INTEGER :: ii,ij ! temporary scalars
+ INTEGER, DIMENSION(2) :: ilocs !
+ INTEGER, INTENT( in ) :: kt
+ !!
+ !!==========================
+ !! DIC Check
+ zmax = -5.0 ! arbitrary low maximum value
+ zmin = 4.0E4 ! arbitrary high minimum value
+ DO jj = 2, jpjm1
+ DO ji = 2,jpim1
+ IF( tmask(ji,jj,1) == 1) THEN
+ zmax = MAX(zmax,zdic(ji,jj)) ! find local maximum
+ zmin = MIN(zmin,zdic(ji,jj)) ! find local minimum
+ ENDIF
+ END DO
+ END DO
+ IF( lk_mpp ) CALL mpp_max( zmax ) ! max over the global domain
+ IF( lk_mpp ) CALL mpp_min( zmin ) ! min over the global domain
+ !
+ IF( zmax .GT. 4.0E3) THEN ! we've got a problem
+ IF (lk_mpp) THEN
+ CALL mpp_maxloc ( zdic(:,:),tmask(:,:,1), zmax, ii,ij )
+ ELSE
+ ilocs = MAXLOC( zdic(:,:), mask = tmask(:,:,1) == 1. )
+ ii = ilocs(1) + nimpp - 1
+ ij = ilocs(2) + njmpp - 1
+ ENDIF
+ !
+ IF(lwp) THEN
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
+ WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC > 4000 '
+ WRITE(numout,9600) kt, zmax, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
+ ENDIF
+ ENDIF
+ !
+ IF( zmin .LE. 0.0) THEN ! we've got a problem
+ IF (lk_mpp) THEN
+ CALL mpp_minloc ( zdic(:,:),tmask(:,:,1), zmin, ii,ij )
+ ELSE
+ ilocs = MINLOC( zdic(:,:), mask = tmask(:,:,1) == 1. )
+ ii = ilocs(1) + nimpp - 1
+ ij = ilocs(2) + njmpp - 1
+ ENDIF
+ !
+ IF(lwp) THEN
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
+ WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC <= 0 '
+ WRITE(numout,9700) kt, zmin, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
+ ENDIF
+ ENDIF
+ !!
+ !!==========================
+ !! ALKALINITY Check
+ zmax = -5.0 ! arbitrary low maximum value
+ zmin = 4.0E4 ! arbitrary high minimum value
+ DO jj = 2, jpjm1
+ DO ji = 2,jpim1
+ IF( tmask(ji,jj,1) == 1) THEN
+ zmax = MAX(zmax,zalk(ji,jj)) ! find local maximum
+ zmin = MIN(zmin,zalk(ji,jj)) ! find local minimum
+ ENDIF
+ END DO
+ END DO
+ IF( lk_mpp ) CALL mpp_max( zmax ) ! max over the global domain
+ IF( lk_mpp ) CALL mpp_min( zmin ) ! min over the global domain
+ !
+ IF( zmax .GT. 4.0E3) THEN ! we've got a problem
+ IF (lk_mpp) THEN
+ CALL mpp_maxloc ( zalk(:,:),tmask(:,:,1), zmax, ii,ij )
+ ELSE
+ ilocs = MAXLOC( zalk(:,:), mask = tmask(:,:,1) == 1. )
+ ii = ilocs(1) + nimpp - 1
+ ij = ilocs(2) + njmpp - 1
+ ENDIF
+ !
+ IF(lwp) THEN
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
+ WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity > 4000 '
+ WRITE(numout,9800) kt, zmax, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
+ ENDIF
+ ENDIF
+ !
+ IF( zmin .LE. 0.0) THEN ! we've got a problem
+ IF (lk_mpp) THEN
+ CALL mpp_minloc ( zalk(:,:),tmask(:,:,1), zmin, ii,ij )
+ ELSE
+ ilocs = MINLOC( zalk(:,:), mask = tmask(:,:,1) == 1. )
+ ii = ilocs(1) + nimpp - 1
+ ij = ilocs(2) + njmpp - 1
+ ENDIF
+ !
+ IF(lwp) THEN
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
+ WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity <= 0 '
+ WRITE(numout,9900) kt, zmin, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
+ ENDIF
+ ENDIF
+
+
+9600 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max DIC: ',f16.10,', i j: ',2i5)
+9700 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min DIC: ',f16.10,', i j: ',2i5)
+9800 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max ALK: ',f16.10,', i j: ',2i5)
+9900 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min ALK: ',f16.10,', i j: ',2i5)
+
+ END SUBROUTINE trc_bio_check
+
+
#else
!!=====================================================================