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.
Changeset 9182 for branches/NERC/dev_r5518_GO6_Carb_Fail_from_GO6_9163/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90 – NEMO

Ignore:
Timestamp:
2018-01-04T17:24:35+01:00 (6 years ago)
Author:
jpalmier
Message:

JPALM - UKESM #475 -- merge with dev_r5518_GO6_Carb_Fail_from_GO6_8356 rev 9070-9133

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/NERC/dev_r5518_GO6_Carb_Fail_from_GO6_9163/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90

    r8442 r9182  
    2323   !!  -   !  2016-11  (A. Yool)              Updated diags for CMIP6 
    2424   !!  -   !  2017-05  (A. Yool)              Added extra DMS calculation 
     25   !!  -   !  2017-11  (J. Palm, A. Yool)     Diagnose tracer excursions 
    2526   !!---------------------------------------------------------------------- 
    2627   !! 
     
    8182                                            gdepw_0, gdepw_n,               & 
    8283                                            nday_year, nsec_day, nyear,     & 
    83                                             rdt, tmask 
     84                                            rdt, tmask, mig, mjg, nimpp,    & 
     85                                            njmpp  
    8486      USE in_out_manager,             ONLY: lwp, numout, nn_date0 
    8587# if defined key_iomput 
     
    8789# endif 
    8890      USE lbclnk,                     ONLY: lbc_lnk 
    89       USE lib_mpp,                    ONLY: ctl_stop 
     91      USE lib_mpp   
    9092      USE oce,                        ONLY: tsb, tsn 
    9193      USE par_kind,                   ONLY: wp 
     
    115117       
    116118      PUBLIC   trc_bio_medusa    ! called in trcsms_medusa.F90 
     119      PUBLIC   trc_bio_exceptional_fix ! here  
    117120 
    118121   !!* Substitution 
     
    177180      !! temporary variables 
    178181      REAL(wp) ::    fq0,fq1,fq2,fq3,fq4 
     182      !! 
     183      !! T and S check temporary variable 
     184      REAL(wp) ::    sumtsn, tsnavg 
     185      INTEGER  ::    summask 
     186      CHARACTER(40) :: charout, charout2, charout3, charout4, charout5 
    179187      !! 
    180188      !!------------------------------------------------------------------ 
     
    450458 
    451459# if defined key_roam 
     460         !! extra MEDUSA-2 tracers 
    452461         DO jj = 2,jpjm1 
    453462            DO ji = 2,jpim1 
     
    456465                  zdtc(ji,jj) = max(0.,trn(ji,jj,jk,jpdtc)) 
    457466                  !! dissolved inorganic carbon 
    458                   zdic(ji,jj) = max(0.,trn(ji,jj,jk,jpdic)) 
     467                  zdic(ji,jj) = trn(ji,jj,jk,jpdic) 
    459468                  !! alkalinity 
    460                   zalk(ji,jj) = max(0.,trn(ji,jj,jk,jpalk)) 
     469                  zalk(ji,jj) = trn(ji,jj,jk,jpalk) 
    461470                  !! oxygen 
    462471                  zoxy(ji,jj) = max(0.,trn(ji,jj,jk,jpoxy)) 
     
    470479                  ztmp(ji,jj) = tsn(ji,jj,jk,jp_tem) 
    471480                  zsal(ji,jj) = tsn(ji,jj,jk,jp_sal) 
    472                   !! 
    473              !! AXY (28/02/14): check input fields 
    474                   if (ztmp(ji,jj) .lt. -3.0 .or. ztmp(ji,jj) .gt. 40.0 ) then 
    475                      IF(lwp) WRITE(numout,*)                                 & 
    476                         ' trc_bio_medusa: T WARNING 2D, ',                   & 
    477                         tsb(ji,jj,jk,jp_tem), tsn(ji,jj,jk,jp_tem),          & 
    478                         ' at (', ji, ',', jj, ',', jk, ') at time', kt 
    479            IF(lwp) WRITE(numout,*)                                 & 
    480                         ' trc_bio_medusa: T SWITCHING 2D, ',                 & 
    481                         tsn(ji,jj,jk,jp_tem), ' -> ', tsb(ji,jj,jk,jp_tem) 
    482                      !! temperatur 
    483                      ztmp(ji,jj) = tsb(ji,jj,jk,jp_tem) 
    484                   endif 
    485                   if (zsal(ji,jj) .lt. 0.0 .or. zsal(ji,jj) .gt. 45.0 ) then 
    486                      IF(lwp) WRITE(numout,*)                                 & 
    487                         ' trc_bio_medusa: S WARNING 2D, ',                   & 
    488                         tsb(ji,jj,jk,jp_sal), tsn(ji,jj,jk,jp_sal),          & 
    489                         ' at (', ji, ',', jj, ',', jk, ') at time', kt 
    490                   endif 
    491481               ENDIF 
    492482            ENDDO 
    493483         ENDDO 
    494484# else 
     485         !! diagnostic MEDUSA-1 detrital carbon tracer 
    495486         DO jj = 2,jpjm1 
    496487            DO ji = 2,jpim1 
    497                if (tmask(ji,jj,jk) == 1) then 
     488               IF (tmask(ji,jj,jk) == 1) THEN 
    498489                  !! implicit detrital carbon 
    499490                  zdtc(ji,jj) = zdet(ji,jj) * xthetad 
     
    502493         ENDDO 
    503494# endif 
     495 
     496# if defined key_roam 
     497         !! --------------------------------------------- 
     498         !! JPALM -- 14-12-2017 -- Here, before any exeptionnal crazy value is 
     499         !!              removed, we want to tell to the Master Node that this  
     500         !!              Exceptionnal value did exist 
     501         !! 
     502         Call trc_bio_check(kt) 
     503 
     504         !!================================================================ 
     505    !! AXY (03/11/17): check input fields 
     506         !! tracer values that exceed thresholds can cause carbonate system 
     507         !! failures when passed to MOCSY; temporary temperature excursions 
     508         !! in recent UKESM0.8 runs trigger such failures but are too short 
     509         !! to have physical consequences; this section checks for such 
     510         !! values and amends them using neighbouring values 
     511         !!  
     512         !! the check on temperature here is also carried out at the end of 
     513         !! each model time step and anomalies are reported in the master 
     514         !! ocean.output file; the error reporting below is strictly local 
     515         !! to the relevant ocean.output_XXXX file so will not be visible 
     516         !! unless all nodes are reporting output 
     517         !!================================================================ 
     518         !! 
     519         DO jj = 2,jpjm1 
     520            DO ji = 2,jpim1 
     521               if (tmask(ji,jj,jk) == 1) then 
     522                  !! the thresholds for the four tracers are ... 
     523                  IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0  ) .OR.   & 
     524                       (zsal(ji,jj) .LE.  0.0) .OR. (zsal(ji,jj) .GT. 50.0  ) .OR.   & 
     525                       (zdic(ji,jj) .LE.  0.0) .OR. (zdic(ji,jj) .GT. 4.0E3 ) .OR.   & 
     526                       (zalk(ji,jj) .LE.  0.0) .OR. (zalk(ji,jj) .GT. 4.0E3 ) ) THEN 
     527                     !! 
     528                     !! all tracer values are reported in the event of any excursion 
     529                     write(charout,*)  ' Tmp = ', ztmp(ji,jj) 
     530                     write(charout2,*) ' Sal = ', zsal(ji,jj) 
     531                     write(charout3,*) ' DIC = ', zdic(ji,jj) 
     532                     write(charout4,*) ' Alk = ', zalk(ji,jj) 
     533                     write(charout5,*) mig(ji), mjg(jj), jk, kt  
     534                     IF (lwp) CALL ctl_warn( 'trc_bio_medusa: carbonate chemistry WARNING:',  & 
     535                        TRIM(charout),TRIM(charout2),TRIM(charout3),TRIM(charout4),           &  
     536                        'at i, j, k, kt:', TRIM(charout5) ) 
     537                     !! 
     538                     !! Detect, report and correct tracer excursions 
     539                     IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT.  40.0) )             & 
     540                        CALL trc_bio_exceptional_fix(                                         & 
     541                        tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_tem), tmask(ji-1:ji+1,jj-1:jj+1,jk),    & 
     542                        'Tmp', -3.0, 40.0, ztmp(ji,jj) ) 
     543                     !! 
     544                     IF ( (zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT.  50.0)  )             & 
     545                        CALL trc_bio_exceptional_fix(                                         & 
     546                        tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_sal), tmask(ji-1:ji+1,jj-1:jj+1,jk),    & 
     547                        'Sal', 1.0, 50.0, zsal(ji,jj) ) 
     548                     !! 
     549                     IF ( (zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3)  )             & 
     550                        CALL trc_bio_exceptional_fix(                                         & 
     551                        trn(ji-1:ji+1,jj-1:jj+1,jk,jpdic), tmask(ji-1:ji+1,jj-1:jj+1,jk),     & 
     552                        'DIC', 100.0, 4.0E3, zdic(ji,jj) ) 
     553                     !! 
     554                     IF ( (zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3)  )             & 
     555                        CALL trc_bio_exceptional_fix(                                         & 
     556                        trn(ji-1:ji+1,jj-1:jj+1,jk,jpalk), tmask(ji-1:ji+1,jj-1:jj+1,jk),     & 
     557                        'Alk', 100.0, 4.0E3, zalk(ji,jj) ) 
     558                  ENDIF 
     559               ENDIF 
     560            ENDDO 
     561         ENDDO 
     562# endif 
     563 
    504564# if defined key_debug_medusa 
    505565         DO jj = 2,jpjm1 
     
    657717   END SUBROUTINE trc_bio_medusa 
    658718 
     719 
     720 
     721   SUBROUTINE trc_bio_exceptional_fix(tiny_var, tiny_mask, var_nm, mini, maxi, varout) 
     722      !! JPALM (27/10/17): This function is called only when abnormal values that  
     723      !! could break the model's carbonate system are fed to MEDUSA 
     724      !!  
     725      !! The function calculates an average tracer value based on the 3x3 cell 
     726      !! neighbourhood around the abnormal cell, and reports this back 
     727      !!  
     728      !! Tracer array values are not modified, but MEDUSA uses "corrected" values 
     729      !! in its calculations 
     730      !! 
     731      !! temporary variables 
     732      REAL(wp), INTENT( in ), DIMENSION(3,3) :: tiny_var, tiny_mask 
     733      CHARACTER(25), INTENT( in )            :: var_nm 
     734      REAL(wp), INTENT( in )                 :: mini, maxi 
     735      REAL(wp), INTENT( out )                :: varout 
     736      REAL(wp)      :: sumtsn, tsnavg 
     737      INTEGER       :: summask 
     738      CHARACTER(25) :: charout1, charout2 
     739      CHARACTER(60) :: charout3, charout4 
     740      INTEGER       :: ii,ij 
     741     
     742      !! point to the center of the 3*3 zoom-grid, to check around 
     743      ii = 2 
     744      ij = 2 
     745      !! Print surounding values to check if isolated Crazy value or  
     746      !! General error 
     747      IF(lwp) WRITE(numout,*)                                 & 
     748         '----------------------------------------------------------------------' 
     749      IF(lwp) WRITE(numout,*)                                 & 
     750         'trc_bio_medusa: 3x3 neighbourhood surrounding abnormal ', TRIM(var_nm) 
     751      IF(lwp) WRITE(numout,9100)                              & 
     752         3, tiny_var(ii-1,ij+1), tiny_var(ii  ,ij+1), tiny_var(ii+1,ij+1) 
     753      IF(lwp) WRITE(numout,9100)                              & 
     754         2, tiny_var(ii-1,ij  ), tiny_var(ii  ,ij  ), tiny_var(ii+1,ij  ) 
     755      IF(lwp) WRITE(numout,9100)                              & 
     756         1, tiny_var(ii-1,ij-1), tiny_var(ii  ,ij-1), tiny_var(ii+1,ij-1) 
     757      IF(lwp) WRITE(numout,*)                                 & 
     758         'trc_bio_medusa: 3x3 land-sea neighbourhood, tmask' 
     759      IF(lwp) WRITE(numout,9100)                              & 
     760         3, tiny_mask(ii-1,ij+1), tiny_mask(ii  ,ij+1), tiny_mask(ii+1,ij+1) 
     761      IF(lwp) WRITE(numout,9100)                              & 
     762         2, tiny_mask(ii-1,ij  ), tiny_mask(ii  ,ij  ), tiny_mask(ii+1,ij  ) 
     763      IF(lwp) WRITE(numout,9100)                              & 
     764         1, tiny_mask(ii-1,ij-1), tiny_mask(ii  ,ij-1), tiny_mask(ii+1,ij-1) 
     765 
     766      !! Correct out of range values 
     767      sumtsn = ( tiny_mask(ii-1,ij+1) * tiny_var(ii-1,ij+1) ) +  & 
     768               ( tiny_mask(ii  ,ij+1) * tiny_var(ii  ,ij+1) ) +  & 
     769               ( tiny_mask(ii+1,ij+1) * tiny_var(ii+1,ij+1) ) +  & 
     770               ( tiny_mask(ii-1,ij  ) * tiny_var(ii-1,ij  ) ) +  & 
     771               ( tiny_mask(ii+1,ij  ) * tiny_var(ii+1,ij  ) ) +  & 
     772               ( tiny_mask(ii-1,ij-1) * tiny_var(ii-1,ij-1) ) +  & 
     773               ( tiny_mask(ii  ,ij-1) * tiny_var(ii  ,ij-1) ) +  & 
     774               ( tiny_mask(ii+1,ij-1) * tiny_var(ii+1,ij-1) ) 
     775      !! 
     776      summask = tiny_mask(ii-1,ij+1) + tiny_mask(ii  ,ij+1)   +  & 
     777                tiny_mask(ii+1,ij+1) + tiny_mask(ii-1,ij  )   +  & 
     778                tiny_mask(ii+1,ij  ) + tiny_mask(ii-1,ij-1)   +  & 
     779                tiny_mask(ii  ,ij-1) + tiny_mask(ii+1,ij-1) 
     780      !! 
     781      IF ( summask .GT. 0 ) THEN 
     782         tsnavg = ( sumtsn / summask ) 
     783         varout = MAX( MIN( maxi, tsnavg), mini ) 
     784      ELSE    
     785         IF (ztmp(ii,ij) .LT. mini )  varout = mini 
     786         IF (ztmp(ii,ij) .GT. maxi )  varout = maxi 
     787      ENDIF 
     788      !! 
     789      write(charout1,9200) tiny_var(ii,ij) 
     790      write(charout2,9200) varout 
     791      write(charout3,*) ' ', charout1, ' -> ', charout2 
     792      write(charout4,*) ' Tracer: ', trim(var_nm) 
     793      IF(lwp) WRITE(numout,*) 'trc_bio_medusa: ** EXCEPTIONAL VALUE SWITCHING **' 
     794      IF(lwp) WRITE(numout,*) charout4  
     795      IF(lwp) WRITE(numout,*) charout3 
     796      IF(lwp) WRITE(numout,*)                                 & 
     797         '----------------------------------------------------------------------' 
     798      IF(lwp) WRITE(numout,*)                                 & 
     799         ' ' 
     8009100  FORMAT('Row:', i1, '  ', e16.6, ' ', e16.6, ' ', e16.6) 
     8019200  FORMAT(e16.6) 
     802 
     803   END SUBROUTINE trc_bio_exceptional_fix  
     804 
     805   SUBROUTINE trc_bio_check(kt) 
     806      !!----------------------------------- 
     807      !! JPALM -- 14-12-2017 -- Still dealing with this micro-boil/carb failure 
     808      !!                     problem. The model is now able to correct a local 
     809      !!                     crazy value. but does it silently. 
     810      !!                     We need to spread the word to the master node. we 
     811      !!                     don't want the model  to correct values without telling us 
     812      !!                     This module will tell at least when crazy DIC or 
     813      !!                     ALK appears. 
     814      !!------------------------------------- 
     815      REAL(wp)              :: zmax, zmin    ! temporary scalars 
     816      INTEGER               :: ji,jj         ! dummy loop 
     817      INTEGER               :: ii,ij         ! temporary scalars  
     818      INTEGER, DIMENSION(2) :: ilocs         !  
     819      INTEGER, INTENT( in ) :: kt 
     820      !! 
     821      !!========================== 
     822      !! DIC Check 
     823      zmax =  -5.0  ! arbitrary  low maximum value 
     824      zmin =  4.0E4  ! arbitrary high minimum value 
     825      DO jj = 2, jpjm1 
     826         DO ji = 2,jpim1 
     827            IF( tmask(ji,jj,1) == 1) THEN 
     828               zmax = MAX(zmax,zdic(ji,jj))     ! find local maximum 
     829               zmin = MIN(zmin,zdic(ji,jj))     ! find local minimum 
     830            ENDIF 
     831         END DO 
     832      END DO 
     833      IF( lk_mpp )   CALL mpp_max( zmax )       ! max over the global domain 
     834      IF( lk_mpp )   CALL mpp_min( zmin )       ! min over the global domain 
     835      ! 
     836      IF( zmax .GT. 4.0E3) THEN  ! we've got a problem 
     837         IF (lk_mpp) THEN 
     838            CALL mpp_maxloc ( zdic(:,:),tmask(:,:,1), zmax, ii,ij ) 
     839         ELSE 
     840            ilocs = MAXLOC( zdic(:,:), mask = tmask(:,:,1) == 1. ) 
     841            ii = ilocs(1) + nimpp - 1 
     842            ij = ilocs(2) + njmpp - 1 
     843         ENDIF 
     844         ! 
     845         IF(lwp) THEN 
     846            WRITE(numout,*) 'trc_bio:tracer anomaly: *****    WARNING     *****' 
     847            WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC > 4000 ' 
     848            WRITE(numout,9600) kt, zmax, ii, ij 
     849            WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****' 
     850         ENDIF 
     851      ENDIF 
     852      ! 
     853      IF( zmin .LE. 0.0) THEN  ! we've got a problem 
     854         IF (lk_mpp) THEN 
     855            CALL mpp_minloc ( zdic(:,:),tmask(:,:,1), zmin, ii,ij ) 
     856         ELSE 
     857            ilocs = MINLOC( zdic(:,:), mask = tmask(:,:,1) == 1. ) 
     858            ii = ilocs(1) + nimpp - 1 
     859            ij = ilocs(2) + njmpp - 1 
     860         ENDIF 
     861         ! 
     862         IF(lwp) THEN 
     863            WRITE(numout,*) 'trc_bio:tracer anomaly: *****    WARNING     *****' 
     864            WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC <= 0 ' 
     865            WRITE(numout,9700) kt, zmin, ii, ij 
     866            WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****' 
     867         ENDIF 
     868      ENDIF 
     869      !! 
     870      !!========================== 
     871      !! ALKALINITY Check 
     872      zmax =  -5.0  ! arbitrary  low maximum value 
     873      zmin =  4.0E4  ! arbitrary high minimum value 
     874      DO jj = 2, jpjm1 
     875         DO ji = 2,jpim1 
     876            IF( tmask(ji,jj,1) == 1) THEN 
     877               zmax = MAX(zmax,zalk(ji,jj))     ! find local maximum 
     878               zmin = MIN(zmin,zalk(ji,jj))     ! find local minimum 
     879            ENDIF 
     880         END DO 
     881      END DO 
     882      IF( lk_mpp )   CALL mpp_max( zmax )       ! max over the global domain 
     883      IF( lk_mpp )   CALL mpp_min( zmin )       ! min over the global domain 
     884      ! 
     885      IF( zmax .GT. 4.0E3) THEN  ! we've got a problem 
     886         IF (lk_mpp) THEN 
     887            CALL mpp_maxloc ( zalk(:,:),tmask(:,:,1), zmax, ii,ij ) 
     888         ELSE 
     889            ilocs = MAXLOC( zalk(:,:), mask = tmask(:,:,1) == 1. ) 
     890            ii = ilocs(1) + nimpp - 1 
     891            ij = ilocs(2) + njmpp - 1 
     892         ENDIF 
     893         ! 
     894         IF(lwp) THEN 
     895            WRITE(numout,*) 'trc_bio:tracer anomaly: *****     WARNING    *****' 
     896            WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity > 4000 ' 
     897            WRITE(numout,9800) kt, zmax, ii, ij 
     898            WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****' 
     899         ENDIF 
     900      ENDIF 
     901      ! 
     902      IF( zmin .LE. 0.0) THEN  ! we've got a problem 
     903         IF (lk_mpp) THEN 
     904            CALL mpp_minloc ( zalk(:,:),tmask(:,:,1), zmin, ii,ij ) 
     905         ELSE 
     906            ilocs = MINLOC( zalk(:,:), mask = tmask(:,:,1) == 1. ) 
     907            ii = ilocs(1) + nimpp - 1 
     908            ij = ilocs(2) + njmpp - 1 
     909         ENDIF 
     910         ! 
     911         IF(lwp) THEN 
     912            WRITE(numout,*) 'trc_bio:tracer anomaly: *****    WARNING     *****' 
     913            WRITE(numout,*) 'trc_bio:tracer anomaly:  sea surface Alkalinity <= 0 ' 
     914            WRITE(numout,9900) kt, zmin, ii, ij 
     915            WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****' 
     916         ENDIF 
     917      ENDIF 
     918 
     919 
     9209600  FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max DIC: ',f16.10,', i j: ',2i5) 
     9219700  FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min DIC: ',f16.10,', i j: ',2i5) 
     9229800  FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max ALK: ',f16.10,', i j: ',2i5) 
     9239900  FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min ALK: ',f16.10,', i j: ',2i5) 
     924 
     925   END SUBROUTINE trc_bio_check 
     926 
     927 
    659928#else 
    660929   !!===================================================================== 
Note: See TracChangeset for help on using the changeset viewer.