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 8677 – NEMO

Changeset 8677


Ignore:
Timestamp:
2017-11-08T11:32:49+01:00 (6 years ago)
Author:
frrh
Message:

Commit JP/AXY's local changes to give us a stable platform.

Location:
branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO
Files:
1 added
7 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/OPA_SRC/stpctl.F90

    r6487 r8677  
    5858      INTEGER  ::   ii, ij, ik              ! temporary integers 
    5959      REAL(wp) ::   zumax, zsmin, zssh2     ! temporary scalars 
     60      REAL(wp) ::   ztmax, ztmin            ! temporary scalars 
    6061      INTEGER, DIMENSION(3) ::   ilocu      !  
    6162      INTEGER, DIMENSION(2) ::   ilocs      !  
     
    1481499500  FORMAT (' kt=',i6,' min SSS: ',1pg11.4,', i j: ',2i5) 
    149150 
     151! ==================================================================================================== 
     152! ==================================================================================================== 
     153      !                                              !AXY (25/10/17) 
     154      !                                              !* Test max/min limits of temperature 
     155      !                                              !  ---------------------------------- 
     156      ztmax =  -5.e0  ! arbitrary  low maximum value 
     157      ztmin = 100.e0  ! arbitrary high minimum value 
     158      DO jj = 2, jpjm1 
     159         DO ji = 1, jpi 
     160            IF( tmask(ji,jj,1) == 1) THEN 
     161               ztmax = MAX(ztmax,tsn(ji,jj,1,jp_tem))     ! find local maximum 
     162               ztmin = MIN(ztmin,tsn(ji,jj,1,jp_tem))     ! find local minimum 
     163            ENDIF 
     164         END DO 
     165      END DO 
     166      IF( lk_mpp )   CALL mpp_max( ztmax )                ! max over the global domain 
     167      IF( lk_mpp )   CALL mpp_min( ztmin )                ! min over the global domain 
     168      ! 
     169      IF( ztmax > 40.) THEN  ! we've got a problem 
     170         IF (lk_mpp) THEN 
     171            CALL mpp_maxloc ( tsn(:,:,1,jp_tem),tmask(:,:,1), ztmax, ii,ij ) 
     172         ELSE 
     173            ilocs = MAXLOC( tsn(:,:,1,jp_tem), mask = tmask(:,:,1) == 1.e0 ) 
     174            ii = ilocs(1) + nimpp - 1 
     175            ij = ilocs(2) + njmpp - 1 
     176         ENDIF 
     177         ! 
     178         IF(lwp) THEN 
     179            WRITE(numout,*) '***** WARNING *****' 
     180            WRITE(numout,*) 'stp_ctl : sea surface temperature > 40C' 
     181            WRITE(numout,*) '======= ' 
     182            WRITE(numout,9600) kt, ztmax, ii, ij 
     183         ENDIF 
     184      ENDIF 
     185      ! 
     186      IF( ztmin < -3.) THEN  ! we've got a problem 
     187         IF (lk_mpp) THEN 
     188            CALL mpp_minloc ( tsn(:,:,1,jp_tem),tmask(:,:,1), ztmin, ii,ij ) 
     189         ELSE 
     190            ilocs = MINLOC( tsn(:,:,1,jp_tem), mask = tmask(:,:,1) == 1.e0 ) 
     191            ii = ilocs(1) + nimpp - 1 
     192            ij = ilocs(2) + njmpp - 1 
     193         ENDIF 
     194         ! 
     195         IF(lwp) THEN 
     196            WRITE(numout,*) '***** WARNING *****' 
     197            WRITE(numout,*) 'stp_ctl : sea surface temperature < -3C' 
     198            WRITE(numout,*) '======= ' 
     199            WRITE(numout,9700) kt, ztmin, ii, ij 
     200         ENDIF 
     201      ENDIF 
     2029600  FORMAT (' kt=',i6,' max SST: ',1pg11.4,', i j: ',2i5) 
     2039700  FORMAT (' kt=',i6,' min SST: ',1pg11.4,', i j: ',2i5) 
     204! ==================================================================================================== 
     205! ==================================================================================================== 
    150206       
    151207      IF( lk_c1d )  RETURN          ! No log file in case of 1D vertical configuration 
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90

    r8442 r8677  
    6161# endif 
    6262                                   zchd, zchn, zdin, zsil 
    63       USE dom_oce,           ONLY: e3t_0, e3t_n, gphit, tmask 
     63      USE dom_oce,           ONLY: e3t_0, e3t_n, gphit, tmask, mig, mjg 
    6464# if defined key_iomput 
    6565      USE iom,               ONLY: lk_iomput 
     
    9090      USE trcoxy_medusa,     ONLY: trc_oxy_medusa 
    9191# endif 
     92      USE lib_mpp,           ONLY: ctl_stop 
     93      USE trcstat  
    9294 
    9395   !!* Substitution 
     
    121123 
    122124# if defined key_roam 
     125      !! init 
     126      f_fco2w(:,:)       = 0.0 
     127      f_fco2atm(:,:)     = 0.0 
     128      f_schmidtco2(:,:)  = 0.0 
     129      f_kwco2(:,:)       = 0.0 
     130      f_co2starair(:,:)  = 0.0 
     131      f_dpco2(:,:)       = 0.0 
     132      f_rhosw(:,:)       = 0.0 
     133      f_K0(:,:)          = 0.0 
     134      !! air pressure (atm); ultimately this will use air  
     135      !! pressure at the base of the UKESM1 atmosphere  
     136      !!                                      
     137      f_pp0(:,:)   = 1.0 
     138 
     139 
    123140      !!----------------------------------------------------------- 
    124141      !! Air-sea gas exchange 
     
    133150         DO ji = 2,jpim1 
    134151            !! OPEN wet point IF..THEN loop 
    135             if (tmask(ji,jj,1) == 1) then 
     152            IF (tmask(ji,jj,1) == 1) then 
    136153               IF (lk_oasis) THEN 
    137154                  !! use 2D atm xCO2 from atm coupling 
    138155                  f_xco2a(ji,jj) = PCO2a_in_cpl(ji,jj) 
     156                  !!!  
     157                  !!! Jpalm test on atm xCO2 
     158                  IF ( (f_xco2a(ji,jj) > 1500 ).OR.(f_xco2a(ji,jj) < 100 ) ) THEN 
     159                    IF(lwp) THEN  
     160                       WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj),          & 
     161                                       ' -- ji =', mig(ji),' jj = ', mjg(jj) 
     162                       CALL ctl_stop( 'MEDUSA - Air-Sea :', 'unrealistic atm xCO2 ' ) 
     163                    ENDIF  
     164                 ENDIF  
    139165               ENDIF 
    140166               !! 
     
    161187               'air-sea: carb-chem kt = ', kt 
    162188               CALL flush(numout) 
     189               !! JPALM add carb print: 
     190               call trc_rst_dia_stat(f_xco2a(:,:), 'f_xco2a') 
     191               call trc_rst_dia_stat(wndm(:,:), 'wndm') 
     192               call trc_rst_dia_stat(f_kw660(:,:), 'f_kw660') 
     193               call trc_rst_dia_stat(ztmp(:,:), 'ztmp') 
     194               call trc_rst_dia_stat(zsal(:,:), 'zsal') 
     195               call trc_rst_dia_stat(zalk(:,:), 'zalk') 
     196               call trc_rst_dia_stat(zdic(:,:), 'zdic') 
     197               call trc_rst_dia_stat(zsil(:,:), 'zsil') 
     198               call trc_rst_dia_stat(zpho(:,:), 'zpho') 
    163199#   endif 
    164200      DO jj = 2,jpjm1 
    165201         DO ji = 2,jpim1 
    166202            if (tmask(ji,jj,1) == 1) then 
    167                !! air pressure (atm); ultimately this will use air  
    168                !! pressure at the base of the UKESM1 atmosphere  
    169                !!                                      
    170                f_pp0(ji,jj)   = 1.0 
    171                !! 
    172                !! IF(lwp) WRITE(numout,*) ' MEDUSA ztmp    =', ztmp(ji,jj) 
    173                !! IF(lwp) WRITE(numout,*) ' MEDUSA wndm    =', wndm(ji,jj) 
    174                !! IF(lwp) WRITE(numout,*) ' MEDUSA fr_i    =', fr_i(ji,jj) 
    175203               !! 
    176204#  if defined key_axy_carbchem 
    177205#   if defined key_mocsy 
     206               !! Jpalm -- 12-09-2017 -- add extra check after reccurent 
     207               !!          carbonate failure in the coupled run. 
     208               !!          must be associated to air-sea flux or air xCO2...i 
     209               !!          Check MOCSY inputs 
     210               IF ( (zsal(ji,jj) > 75.0 ).OR.(zsal(ji,jj) < 0.0 ) .OR.        & 
     211                    (ztmp(ji,jj) > 50.0 ).OR.(ztmp(ji,jj) < -20.0 ) .OR.      & 
     212                    (zalk(ji,jj) > 35.0E2 ).OR.(zalk(ji,jj) <= 0.0 ) .OR.     & 
     213                    (zdic(ji,jj) > 35.0E2 ).OR.(zdic(ji,jj) <= 0.0 ) .OR.     & 
     214                    (f_kw660(ji,jj) > 1.0E-2 ).OR.(f_kw660(ji,jj) < 0.0 ) ) THEN 
     215                  IF(lwp) THEN  
     216                      WRITE(numout,*) ' surface T = ',ztmp(ji,jj) 
     217                      WRITE(numout,*) ' surface S = ',zsal(ji,jj) 
     218                      WRITE(numout,*) ' surface ALK = ',zalk(ji,jj) 
     219                      WRITE(numout,*) ' surface DIC = ',zdic(ji,jj) 
     220                      WRITE(numout,*) ' KW660 = ',f_kw660(ji,jj) 
     221                      WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj)    
     222                      WRITE(numout,*) ' surface pco2w  = ',f_pco2w(ji,jj) 
     223                      WRITE(numout,*) ' surface fco2w  = ',f_fco2w(ji,jj) 
     224                      WRITE(numout,*) ' surface fco2a  = ',f_fco2atm(ji,jj) 
     225                      WRITE(numout,*) ' surface co2flx = ',f_co2flux(ji,jj) 
     226                      WRITE(numout,*) ' surface dpco2  = ',f_dpco2(ji,jj) 
     227                      WRITE(numout,*) ' MOCSY input: ji =', mig(ji),' jj = ', mjg(jj),  & 
     228                                       ' kt = ', kt  
     229                      WRITE(numout,*) 'MEDUSA - Air-Sea INPUT: unrealistic surface Carb. Chemistry' 
     230                      CALL ctl_stop( 'MEDUSA - Air-Sea INPUT: ',             & 
     231                                     'unrealistic surface Carb. Chemistry -- INPUTS' ) 
     232                  ENDIF      
     233               ENDIF      
    178234               !! 
    179235               !! AXY (22/06/15): use Orr & Epitalon (2015) MOCSY-2 carbonate 
     
    200256               f_TALK(ji,jj) = (zalk(ji,jj) / f_rhosw(ji,jj)) * 1000. 
    201257               f_dcf(ji,jj)  = f_rhosw(ji,jj) 
     258               !! Jpalm -- 12-09-2017 -- add extra check after reccurent 
     259               !!          carbonate failure in the coupled run. 
     260               !!          must be associated to air-sea flux or air xCO2...i 
     261               !!          Check MOCSY inputs 
     262              IF ( (f_pco2w(ji,jj) > 1.E4 ).OR.(f_pco2w(ji,jj) < 0.0 ) .OR.     & 
     263                   (f_fco2w(ji,jj) > 1.E4 ).OR.(f_fco2w(ji,jj) < 0.0 ) .OR.     &    
     264                   (f_fco2atm(ji,jj) > 1.E4 ).OR.(f_fco2atm(ji,jj) < 0.0 ) .OR.     & 
     265                   (f_co2flux(ji,jj) > 1.E-2 ).OR.(f_co2flux(ji,jj) < -1.E-2 ) .OR.     & 
     266                   (f_dpco2(ji,jj) > 1.E4 ).OR.(f_dpco2(ji,jj) < -1.E4 ) ) THEN 
     267                 IF(lwp) THEN  
     268                     WRITE(numout,*) ' surface T = ',ztmp(ji,jj) 
     269                     WRITE(numout,*) ' surface S = ',zsal(ji,jj) 
     270                     WRITE(numout,*) ' surface ALK = ',zalk(ji,jj) 
     271                     WRITE(numout,*) ' surface DIC = ',zdic(ji,jj) 
     272                     WRITE(numout,*) ' KW660 = ',f_kw660(ji,jj) 
     273                     WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj)    
     274                     WRITE(numout,*) ' surface pco2w  = ',f_pco2w(ji,jj) 
     275                     WRITE(numout,*) ' surface fco2w  = ',f_fco2w(ji,jj) 
     276                     WRITE(numout,*) ' surface fco2a  = ',f_fco2atm(ji,jj) 
     277                     WRITE(numout,*) ' surface co2flx = ',f_co2flux(ji,jj) 
     278                     WRITE(numout,*) ' surface dpco2  = ',f_dpco2(ji,jj) 
     279                     WRITE(numout,*) ' MOCSY output: ji =', mig(ji),' jj = ', mjg(jj),  & 
     280                                       ' kt = ', kt      
     281                     WRITE(numout,*) 'MEDUSA - Air-Sea OUTPUT: unrealistic surface Carb. Chemistry' 
     282                     CALL ctl_stop( 'MEDUSA - Air-Sea OUTPUT: ',            & 
     283                                    'unrealistic surface Carb. Chemistry -- OUTPUTS' ) 
     284                 ENDIF      
     285              ENDIF      
    202286            ENDIF 
    203287         ENDDO 
    204288      ENDDO 
    205289 
     290#   if defined key_debug_medusa 
     291               !! JPALM add carb print: 
     292               call trc_rst_dia_stat(f_pco2w(:,:), 'f_pco2w') 
     293               call trc_rst_dia_stat(f_fco2w(:,:), 'f_fco2w') 
     294               call trc_rst_dia_stat(f_fco2atm(:,:), 'f_fco2atm') 
     295               call trc_rst_dia_stat(f_schmidtco2(:,:), 'f_schmidtco2') 
     296               call trc_rst_dia_stat(f_kwco2(:,:), 'f_kwco2') 
     297               call trc_rst_dia_stat(f_co2starair(:,:), 'f_co2starair') 
     298               call trc_rst_dia_stat(f_co2flux(:,:), 'f_co2flux') 
     299               call trc_rst_dia_stat(f_dpco2(:,:), 'f_dpco2') 
     300#   endif 
    206301#   else    
    207302 
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90

    r8442 r8677  
    164164      fslownflux(:,:) = 0.0  
    165165      fslowcflux(:,:) = 0.0  
     166      !! 
     167      !! JPALM -- 21-09-2017 -- needed to debug air-sea carb 
     168      f_xco2a(:,:)  = 0.0 
     169      f_pco2w(:,:)  = 0.0 
     170      f_ph(:,:)     = 0.0 
     171      f_kw660(:,:)  = 0.0 
     172      ztmp(:,:)  = 0.0 
     173      zsal(:,:)  = 0.0 
     174      zalk(:,:)  = 0.0 
     175      zdic(:,:)  = 0.0 
     176      zsil(:,:)  = 0.0 
     177      zpho(:,:)  = 0.0 
     178      f_co2flux(:,:)  = 0.0  
     179      f_pco2atm(:,:)  = 0.0 
     180      f_h2co3(:,:)    = 0.0 
     181      f_hco3(:,:)     = 0.0 
     182      f_co3(:,:)      = 0.0 
     183      f_omarg(:,:)    = 0.0 
     184      f_omcal(:,:)    = 0.0 
     185 
    166186      !! 
    167187      !! allocate and initiate 2D diag 
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90

    r8442 r8677  
    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 
    8485      USE in_out_manager,             ONLY: lwp, numout, nn_date0 
    8586# if defined key_iomput 
     
    8788# endif 
    8889      USE lbclnk,                     ONLY: lbc_lnk 
    89       USE lib_mpp,                    ONLY: ctl_stop 
     90      USE lib_mpp,                    ONLY: ctl_stop, ctl_warn 
    9091      USE oce,                        ONLY: tsb, tsn 
    9192      USE par_kind,                   ONLY: wp 
     
    115116       
    116117      PUBLIC   trc_bio_medusa    ! called in trcsms_medusa.F90 
     118      PUBLIC   trc_bio_exceptional_fix ! here  
    117119 
    118120   !!* Substitution 
     
    177179      !! temporary variables 
    178180      REAL(wp) ::    fq0,fq1,fq2,fq3,fq4 
     181      !! 
     182      !! T and S check temporary variable 
     183      REAL(wp) ::    sumtsn, tsnavg 
     184      INTEGER  ::    summask 
     185      CHARACTER(40) :: charout, charout2, charout3, charout4, charout5 
    179186      !! 
    180187      !!------------------------------------------------------------------ 
     
    450457 
    451458# if defined key_roam 
     459         !! extra MEDUSA-2 tracers 
    452460         DO jj = 2,jpjm1 
    453461            DO ji = 2,jpim1 
     
    456464                  zdtc(ji,jj) = max(0.,trn(ji,jj,jk,jpdtc)) 
    457465                  !! dissolved inorganic carbon 
    458                   zdic(ji,jj) = max(0.,trn(ji,jj,jk,jpdic)) 
     466                  zdic(ji,jj) = trn(ji,jj,jk,jpdic) 
    459467                  !! alkalinity 
    460                   zalk(ji,jj) = max(0.,trn(ji,jj,jk,jpalk)) 
     468                  zalk(ji,jj) = trn(ji,jj,jk,jpalk) 
    461469                  !! oxygen 
    462470                  zoxy(ji,jj) = max(0.,trn(ji,jj,jk,jpoxy)) 
     
    470478                  ztmp(ji,jj) = tsn(ji,jj,jk,jp_tem) 
    471479                  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 
    491480               ENDIF 
    492481            ENDDO 
    493482         ENDDO 
    494483# else 
     484         !! diagnostic MEDUSA-1 detrital carbon tracer 
    495485         DO jj = 2,jpjm1 
    496486            DO ji = 2,jpim1 
    497                if (tmask(ji,jj,jk) == 1) then 
     487               IF (tmask(ji,jj,jk) == 1) THEN 
    498488                  !! implicit detrital carbon 
    499489                  zdtc(ji,jj) = zdet(ji,jj) * xthetad 
     
    502492         ENDDO 
    503493# endif 
     494 
     495# if defined key_roam 
     496         !!================================================================ 
     497    !! AXY (03/11/17): check input fields 
     498         !! tracer values that exceed thresholds can cause carbonate system 
     499         !! failures when passed to MOCSY; temporary temperature excursions 
     500         !! in recent UKESM0.8 runs trigger such failures but are too short 
     501         !! to have physical consequences; this section checks for such 
     502         !! values and amends them using neighbouring values 
     503         !!  
     504         !! the check on temperature here is also carried out at the end of 
     505         !! each model time step and anomalies are reported in the master 
     506         !! ocean.output file; the error reporting below is strictly local 
     507         !! to the relevant ocean.output_XXXX file so will not be visible 
     508         !! unless all nodes are reporting output 
     509         !!================================================================ 
     510         !! 
     511         DO jj = 2,jpjm1 
     512            DO ji = 2,jpim1 
     513               if (tmask(ji,jj,jk) == 1) then 
     514                  !! the thresholds for the four tracers are ... 
     515                  IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT. 40.0  ) .OR.   & 
     516                       (zsal(ji,jj) .LE.  0.0) .OR. (zsal(ji,jj) .GT. 50.0  ) .OR.   & 
     517                       (zdic(ji,jj) .LE.  0.0) .OR. (zdic(ji,jj) .GT. 4.0E3 ) .OR.   & 
     518                       (zalk(ji,jj) .LE.  0.0) .OR. (zalk(ji,jj) .GT. 4.0E3 ) ) THEN 
     519                     !! 
     520                     !! all tracer values are reported in the event of any excursion 
     521                     write(charout,*)  ' Tmp = ', ztmp(ji,jj) 
     522                     write(charout2,*) ' Sal = ', zsal(ji,jj) 
     523                     write(charout3,*) ' DIC = ', zdic(ji,jj) 
     524                     write(charout4,*) ' Alk = ', zalk(ji,jj) 
     525                     write(charout5,*) mig(ji), mjg(jj), jk, kt  
     526                     IF (lwp) CALL ctl_warn( 'trc_bio_medusa: carbonate chemistry WARNING:',  & 
     527                        TRIM(charout),TRIM(charout2),TRIM(charout3),TRIM(charout4),           &  
     528                        'at i, j, k, kt:', TRIM(charout5) ) 
     529                     !! 
     530                     !! Detect, report and correct tracer excursions 
     531                     IF ( (ztmp(ji,jj) .LT. -3.0) .OR. (ztmp(ji,jj) .GT.  40.0) )             & 
     532                        CALL trc_bio_exceptional_fix(                                         & 
     533                        tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_tem), tmask(ji-1:ji+1,jj-1:jj+1,jk),    & 
     534                        'Tmp', -3.0, 40.0, ztmp(ji,jj) ) 
     535                     !! 
     536                     IF ( (zsal(ji,jj) .LE. 0.0) .OR. (zsal(ji,jj) .GT.  50.0)  )             & 
     537                        CALL trc_bio_exceptional_fix(                                         & 
     538                        tsn(ji-1:ji+1,jj-1:jj+1,jk,jp_sal), tmask(ji-1:ji+1,jj-1:jj+1,jk),    & 
     539                        'Sal', 1.0, 50.0, zsal(ji,jj) ) 
     540                     !! 
     541                     IF ( (zdic(ji,jj) .LE. 0.0) .OR. (zdic(ji,jj) .GT. 4.0E3)  )             & 
     542                        CALL trc_bio_exceptional_fix(                                         & 
     543                        trn(ji-1:ji+1,jj-1:jj+1,jk,jpdic), tmask(ji-1:ji+1,jj-1:jj+1,jk),     & 
     544                        'DIC', 100.0, 4.0E3, zdic(ji,jj) ) 
     545                     !! 
     546                     IF ( (zalk(ji,jj) .LE. 0.0) .OR. (zalk(ji,jj) .GT. 4.0E3)  )             & 
     547                        CALL trc_bio_exceptional_fix(                                         & 
     548                        trn(ji-1:ji+1,jj-1:jj+1,jk,jpalk), tmask(ji-1:ji+1,jj-1:jj+1,jk),     & 
     549                        'Alk', 100.0, 4.0E3, zalk(ji,jj) ) 
     550                  ENDIF 
     551               ENDIF 
     552            ENDDO 
     553         ENDDO 
     554# endif 
     555 
    504556# if defined key_debug_medusa 
    505557         DO jj = 2,jpjm1 
     
    657709   END SUBROUTINE trc_bio_medusa 
    658710 
     711 
     712 
     713   SUBROUTINE trc_bio_exceptional_fix(tiny_var, tiny_mask, var_nm, mini, maxi, varout) 
     714      !! JPALM (27/10/17): This function is called only when abnormal values that  
     715      !! could break the model's carbonate system are fed to MEDUSA 
     716      !!  
     717      !! The function calculates an average tracer value based on the 3x3 cell 
     718      !! neighbourhood around the abnormal cell, and reports this back 
     719      !!  
     720      !! Tracer array values are not modified, but MEDUSA uses "corrected" values 
     721      !! in its calculations 
     722      !! 
     723      !! temporary variables 
     724      REAL(wp), INTENT( in ), DIMENSION(3,3) :: tiny_var, tiny_mask 
     725      CHARACTER(25), INTENT( in )            :: var_nm 
     726      REAL(wp), INTENT( in )                 :: mini, maxi 
     727      REAL(wp), INTENT( out )                :: varout 
     728      REAL(wp)      :: sumtsn, tsnavg 
     729      INTEGER       :: summask 
     730      CHARACTER(25) :: charout1, charout2 
     731      CHARACTER(60) :: charout3, charout4 
     732      INTEGER       :: ii,ij 
     733     
     734      !! point to the center of the 3*3 zoom-grid, to check around 
     735      ii = 2 
     736      ij = 2 
     737      !! Print surounding values to check if isolated Crazy value or  
     738      !! General error 
     739      IF(lwp) WRITE(numout,*)                                 & 
     740         '----------------------------------------------------------------------' 
     741      IF(lwp) WRITE(numout,*)                                 & 
     742         'trc_bio_medusa: 3x3 neighbourhood surrounding abnormal ', TRIM(var_nm) 
     743      IF(lwp) WRITE(numout,9100)                              & 
     744         3, tiny_var(ii-1,ij+1), tiny_var(ii  ,ij+1), tiny_var(ii+1,ij+1) 
     745      IF(lwp) WRITE(numout,9100)                              & 
     746         2, tiny_var(ii-1,ij  ), tiny_var(ii  ,ij  ), tiny_var(ii+1,ij  ) 
     747      IF(lwp) WRITE(numout,9100)                              & 
     748         1, tiny_var(ii-1,ij-1), tiny_var(ii  ,ij-1), tiny_var(ii+1,ij-1) 
     749      IF(lwp) WRITE(numout,*)                                 & 
     750         'trc_bio_medusa: 3x3 land-sea neighbourhood, tmask' 
     751      IF(lwp) WRITE(numout,9100)                              & 
     752         3, tiny_mask(ii-1,ij+1), tiny_mask(ii  ,ij+1), tiny_mask(ii+1,ij+1) 
     753      IF(lwp) WRITE(numout,9100)                              & 
     754         2, tiny_mask(ii-1,ij  ), tiny_mask(ii  ,ij  ), tiny_mask(ii+1,ij  ) 
     755      IF(lwp) WRITE(numout,9100)                              & 
     756         1, tiny_mask(ii-1,ij-1), tiny_mask(ii  ,ij-1), tiny_mask(ii+1,ij-1) 
     757 
     758      !! Correct out of range values 
     759      sumtsn = ( tiny_mask(ii-1,ij+1) * tiny_var(ii-1,ij+1) ) +  & 
     760               ( tiny_mask(ii  ,ij+1) * tiny_var(ii  ,ij+1) ) +  & 
     761               ( tiny_mask(ii+1,ij+1) * tiny_var(ii+1,ij+1) ) +  & 
     762               ( tiny_mask(ii-1,ij  ) * tiny_var(ii-1,ij  ) ) +  & 
     763               ( tiny_mask(ii+1,ij  ) * tiny_var(ii+1,ij  ) ) +  & 
     764               ( tiny_mask(ii-1,ij-1) * tiny_var(ii-1,ij-1) ) +  & 
     765               ( tiny_mask(ii  ,ij-1) * tiny_var(ii  ,ij-1) ) +  & 
     766               ( tiny_mask(ii+1,ij-1) * tiny_var(ii+1,ij-1) ) 
     767      !! 
     768      summask = tiny_mask(ii-1,ij+1) + tiny_mask(ii  ,ij+1)   +  & 
     769                tiny_mask(ii+1,ij+1) + tiny_mask(ii-1,ij  )   +  & 
     770                tiny_mask(ii+1,ij  ) + tiny_mask(ii-1,ij-1)   +  & 
     771                tiny_mask(ii  ,ij-1) + tiny_mask(ii+1,ij-1) 
     772      !! 
     773      IF ( summask .GT. 0 ) THEN 
     774         tsnavg = ( sumtsn / summask ) 
     775         varout = MAX( MIN( maxi, tsnavg), mini ) 
     776      ELSE    
     777         IF (ztmp(ii,ij) .LT. mini )  varout = mini 
     778         IF (ztmp(ii,ij) .GT. maxi )  varout = maxi 
     779      ENDIF 
     780      !! 
     781      write(charout1,9200) tiny_var(ii,ij) 
     782      write(charout2,9200) varout 
     783      write(charout3,*) ' ', charout1, ' -> ', charout2 
     784      write(charout4,*) ' Tracer: ', trim(var_nm) 
     785      IF(lwp) WRITE(numout,*) 'trc_bio_medusa: ** EXCEPTIONAL VALUE SWITCHING **' 
     786      IF(lwp) WRITE(numout,*) charout4  
     787      IF(lwp) WRITE(numout,*) charout3 
     788      IF(lwp) WRITE(numout,*)                                 & 
     789         '----------------------------------------------------------------------' 
     790      IF(lwp) WRITE(numout,*)                                 & 
     791         ' ' 
     7929100  FORMAT('Row:', i1, '  ', e16.6, ' ', e16.6, ' ', e16.6) 
     7939200  FORMAT(e16.6) 
     794 
     795   END SUBROUTINE trc_bio_exceptional_fix  
     796 
    659797#else 
    660798   !!===================================================================== 
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/TOP_SRC/TRP/trctrp.F90

    r8442 r8677  
    2828   USE zpshde          ! partial step: hor. derivative       (zps_hde routine) 
    2929# if defined key_debug_medusa 
    30    USE trcrst 
     30   USE trcstat 
    3131# endif 
    3232 
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/TOP_SRC/trcrst.F90

    r8427 r8677  
    4343   USE sbc_oce, ONLY: lk_oasis  
    4444   USE oce,     ONLY: CO2Flux_out_cpl, DMS_out_cpl, chloro_out_cpl  !! Coupling variable 
     45   USE trcstat 
    4546 
    4647   IMPLICIT NONE 
     
    5253   PUBLIC   trc_rst_cal 
    5354   PUBLIC   trc_rst_stat 
    54    PUBLIC   trc_rst_dia_stat 
    55    PUBLIC   trc_rst_tra_stat 
    5655 
    5756   !! * Substitutions 
     
    697696 
    698697 
    699    SUBROUTINE trc_rst_tra_stat 
    700       !!---------------------------------------------------------------------- 
    701       !!                    ***  trc_rst_tra_stat  *** 
    702       !! 
    703       !! ** purpose  :   Compute tracers statistics - check where crazy values appears 
    704       !!---------------------------------------------------------------------- 
    705       INTEGER  :: jk, jn 
    706       REAL(wp) :: ztraf, zmin, zmax, zmean, zdrift, areasf 
    707       REAL(wp), DIMENSION(jpi,jpj) :: zvol 
    708       !!---------------------------------------------------------------------- 
    709  
    710       IF( lwp ) THEN 
    711          WRITE(numout,*) 
    712          WRITE(numout,*) '           ----SURFACE TRA STAT----             ' 
    713          WRITE(numout,*) 
    714       ENDIF 
    715       ! 
    716       zvol(:,:) = e1e2t(:,:) * fse3t_a(:,:,1) * tmask(:,:,1) 
    717       areasf = glob_sum(zvol(:,:)) 
    718       DO jn = 1, jptra 
    719          ztraf = glob_sum( tra(:,:,1,jn) * zvol(:,:) ) 
    720          zmin  = MINVAL( tra(:,:,1,jn), mask= ((tmask(:,:,1).NE.0.)) ) 
    721          zmax  = MAXVAL( tra(:,:,1,jn), mask= ((tmask(:,:,1).NE.0.)) ) 
    722          IF( lk_mpp ) THEN 
    723             CALL mpp_min( zmin )      ! min over the global domain 
    724             CALL mpp_max( zmax )      ! max over the global domain 
    725          END IF 
    726          zmean  = ztraf / areasf 
    727          IF(lwp) WRITE(numout,9001) jn, TRIM( ctrcnm(jn) ), zmean, zmin, zmax 
    728       END DO 
    729       IF(lwp) WRITE(numout,*) 
    730 9001  FORMAT(' tracer nb :',i2,'    name :',a10,'    mean :',e18.10,'    min :',e18.10, & 
    731       &      '    max :',e18.10) 
    732       ! 
    733    END SUBROUTINE trc_rst_tra_stat 
    734  
    735  
    736  
    737    SUBROUTINE trc_rst_dia_stat( dgtr, names) 
    738       !!---------------------------------------------------------------------- 
    739       !!                    ***  trc_rst_dia_stat  *** 
    740       !! 
    741       !! ** purpose  :   Compute tracers statistics 
    742       !!---------------------------------------------------------------------- 
    743       REAL(wp), DIMENSION(jpi,jpj) , INTENT(in) ::   dgtr      ! 2D diag var 
    744       CHARACTER(len=*)             , INTENT(in) ::   names     ! 2D diag name 
    745       !!--------------------------------------------------------------------- 
    746       INTEGER  :: jk, jn 
    747       CHARACTER (LEN=18) :: text_zmean 
    748       REAL(wp) :: ztraf, zmin, zmax, zmean, areasf 
    749       REAL(wp), DIMENSION(jpi,jpj) :: zvol 
    750       !!---------------------------------------------------------------------- 
    751  
    752       IF( lwp )  WRITE(numout,*) 'STAT- ', names 
    753        
    754       ! fse3t_a will be undefined at the start of a run, but this routine 
    755       ! may be called at any stage! Hence we MUST make sure it is  
    756       ! initialised to zero when allocated to enable us to test for  
    757       ! zero content here and avoid potentially dangerous and non-portable  
    758       ! operations (e.g. divide by zero, global sums of junk values etc.)    
    759       zvol(:,:) = e1e2t(:,:) * fse3t_a(:,:,1) * tmask(:,:,1) 
    760       ztraf = glob_sum( dgtr(:,:) * zvol(:,:) ) 
    761       !! areasf = glob_sum(e1e2t(:,:) * tmask(:,:,1) ) 
    762       areasf = glob_sum(zvol(:,:)) 
    763       zmin  = MINVAL( dgtr(:,:), mask= ((tmask(:,:,1).NE.0.)) ) 
    764       zmax  = MAXVAL( dgtr(:,:), mask= ((tmask(:,:,1).NE.0.)) ) 
    765       IF( lk_mpp ) THEN 
    766          CALL mpp_min( zmin )      ! min over the global domain 
    767          CALL mpp_max( zmax )      ! max over the global domain 
    768       END IF 
    769  
    770       text_zmean = "N/A" 
    771       ! Avoid divide by zero. areasf must be positive. 
    772       IF  (areasf > 0.0) THEN  
    773          zmean = ztraf / areasf 
    774          WRITE(text_zmean,'(e18.10)') zmean 
    775       ENDIF 
    776  
    777       IF(lwp) WRITE(numout,9002) TRIM( names ), text_zmean, zmin, zmax 
    778  
    779   9002  FORMAT(' tracer name :',A,'    mean :',A,'    min :',e18.10, & 
    780       &      '    max :',e18.10 ) 
    781       ! 
    782    END SUBROUTINE trc_rst_dia_stat 
    783  
    784  
    785698#else 
    786699   !!---------------------------------------------------------------------- 
  • branches/UKMO/dev_r5518_GO6_pkg_test_sst/NEMOGCM/NEMO/TOP_SRC/trcstp.F90

    r8442 r8677  
    1919   USE trcwri 
    2020   USE trcrst 
     21   USE trcstat 
    2122   USE trdtrc_oce 
    2223   USE trdmxl_trc 
Note: See TracChangeset for help on using the changeset viewer.