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 3632 for branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90 – NEMO

Ignore:
Timestamp:
2012-11-22T16:28:42+01:00 (11 years ago)
Author:
acc
Message:

Branch dev_NOC_2012_r3555. #1006. Step 9: Merge in trunk changes between revision 3385 and 3452

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90

    r3625 r3632  
    77   !!             8.5  !  02-06  (E. Durand, G. Madec)  F90 
    88   !!             9.0  !  06-07  (G. Madec)  add clo_rnf, clo_ups, clo_bat 
     9   !!        NEMO 3.4  !  03-12  (P.G. Fogli) sbc_clo bug fix & mpp reproducibility 
    910   !!---------------------------------------------------------------------- 
    1011 
     
    1819   USE oce             ! dynamics and tracers 
    1920   USE dom_oce         ! ocean space and time domain 
    20    USE phycst 
     21   USE phycst          ! physical constants 
    2122   USE in_out_manager  ! I/O manager 
    2223   USE sbc_oce         ! ocean surface boundary conditions 
    23    USE lib_mpp         ! distributed memory computing library 
    24    USE lbclnk          ! ??? 
     24   USE lib_fortran,    ONLY: glob_sum, DDPDD 
     25   USE lbclnk          ! lateral boundary condition - MPP exchanges 
     26   USE lib_mpp         ! MPP library 
     27   USE timing 
    2528 
    2629   IMPLICIT NONE 
     
    8689         SELECT CASE ( jp_cfg ) 
    8790         !                                           ! ======================= 
     91         CASE ( 1 )                                  ! ORCA_R1 configuration 
     92            !                                        ! ======================= 
     93            ncsnr(1)   = 1    ; ncstt(1)   = 0           ! Caspian Sea 
     94            ncsi1(1)   = 332  ; ncsj1(1)   = 203 
     95            ncsi2(1)   = 344  ; ncsj2(1)   = 235 
     96            ncsir(1,1) = 1    ; ncsjr(1,1) = 1 
     97            !                                         
     98            !                                        ! ======================= 
    8899         CASE ( 2 )                                  !  ORCA_R2 configuration 
    89100            !                                        ! ======================= 
     
    174185      !!      put as run-off in open ocean. 
    175186      !! 
    176       !! ** Action  :   emp   updated surface freshwater flux at kt 
     187      !! ** Action  :   emp updated surface freshwater fluxes and associated heat content at kt 
    177188      !!---------------------------------------------------------------------- 
    178189      INTEGER, INTENT(in) ::   kt   ! ocean model time step 
    179190      ! 
    180       INTEGER                     ::   ji, jj, jc, jn   ! dummy loop indices 
    181       REAL(wp)                    ::   zze2, zcoef, zcoef1 
    182       REAL(wp), DIMENSION (jpncs) ::   zfwf  
    183       !!---------------------------------------------------------------------- 
    184       ! 
     191      INTEGER             ::   ji, jj, jc, jn   ! dummy loop indices 
     192      REAL(wp), PARAMETER ::   rsmall = 1.e-20_wp    ! Closed sea correction epsilon 
     193      REAL(wp)            ::   zze2, ztmp, zcorr     !  
     194      REAL(wp)            ::   zcoef, zcoef1         !  
     195      COMPLEX(wp)         ::   ctmp  
     196      REAL(wp), DIMENSION(jpncs) ::   zfwf   ! 1D workspace 
     197      !!---------------------------------------------------------------------- 
     198      ! 
     199      IF( nn_timing == 1 )  CALL timing_start('sbc_clo') 
    185200      !                                                   !------------------! 
    186201      IF( kt == nit000 ) THEN                             !  Initialisation  ! 
     
    190205         IF(lwp) WRITE(numout,*)'~~~~~~~' 
    191206 
    192          ! Total surface of ocean 
    193          surf(jpncs+1) = SUM( e1t(:,:) * e2t(:,:) * tmask_i(:,:) ) 
    194  
    195          DO jc = 1, jpncs 
    196             surf(jc) =0.e0 
    197             DO jj = ncsj1(jc), ncsj2(jc) 
    198                DO ji = ncsi1(jc), ncsi2(jc) 
    199                   surf(jc) = surf(jc) + e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)      ! surface of closed seas 
     207         surf(:) = 0.e0_wp 
     208         ! 
     209         surf(jpncs+1) = glob_sum( e1e2t(:,:) )   ! surface of the global ocean 
     210         ! 
     211         !                                        ! surface of closed seas  
     212         IF( lk_mpp_rep ) THEN                         ! MPP reproductible calculation 
     213            DO jc = 1, jpncs 
     214               ctmp = CMPLX( 0.e0, 0.e0, wp ) 
     215               DO jj = ncsj1(jc), ncsj2(jc) 
     216                  DO ji = ncsi1(jc), ncsi2(jc) 
     217                     ztmp = e1e2t(ji,jj) * tmask_i(ji,jj) 
     218                     CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
     219                  END DO  
    200220               END DO  
    201             END DO  
    202          END DO  
    203          IF( lk_mpp )   CALL mpp_sum ( surf, jpncs+1 )       ! mpp: sum over all the global domain 
     221               IF( lk_mpp )   CALL mpp_sum( ctmp ) 
     222               surf(jc) = REAL(ctmp,wp) 
     223            END DO 
     224         ELSE                                          ! Standard calculation            
     225            DO jc = 1, jpncs 
     226               DO jj = ncsj1(jc), ncsj2(jc) 
     227                  DO ji = ncsi1(jc), ncsi2(jc) 
     228                     surf(jc) = surf(jc) + e1e2t(ji,jj) * tmask_i(ji,jj)      ! surface of closed seas 
     229                  END DO  
     230               END DO  
     231            END DO  
     232            IF( lk_mpp )   CALL mpp_sum ( surf, jpncs )       ! mpp: sum over all the global domain 
     233         ENDIF 
    204234 
    205235         IF(lwp) WRITE(numout,*)'     Closed sea surfaces' 
     
    216246      !                                                   !--------------------! 
    217247      !                                                   !  update emp        ! 
    218       zfwf = 0.e0                                         !--------------------! 
    219       DO jc = 1, jpncs 
    220          DO jj = ncsj1(jc), ncsj2(jc) 
    221             DO ji = ncsi1(jc), ncsi2(jc) 
    222                zfwf(jc) = zfwf(jc) + e1t(ji,jj) * e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj)  
    223             END DO   
    224          END DO  
    225       END DO 
    226       IF( lk_mpp )   CALL mpp_sum ( zfwf(:) , jpncs )       ! mpp: sum over all the global domain 
     248      zfwf = 0.e0_wp                                      !--------------------! 
     249      IF( lk_mpp_rep ) THEN                         ! MPP reproductible calculation 
     250         DO jc = 1, jpncs 
     251            ctmp = CMPLX( 0.e0, 0.e0, wp ) 
     252            DO jj = ncsj1(jc), ncsj2(jc) 
     253               DO ji = ncsi1(jc), ncsi2(jc) 
     254                  ztmp = e1e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj) 
     255                  CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
     256               END DO   
     257            END DO  
     258            IF( lk_mpp )   CALL mpp_sum( ctmp ) 
     259            zfwf(jc) = REAL(ctmp,wp) 
     260         END DO 
     261      ELSE                                          ! Standard calculation            
     262         DO jc = 1, jpncs 
     263            DO jj = ncsj1(jc), ncsj2(jc) 
     264               DO ji = ncsi1(jc), ncsi2(jc) 
     265                  zfwf(jc) = zfwf(jc) + e1e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj)  
     266               END DO   
     267            END DO  
     268         END DO 
     269         IF( lk_mpp )   CALL mpp_sum ( zfwf(:) , jpncs )       ! mpp: sum over all the global domain 
     270      ENDIF 
    227271 
    228272      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN      ! Black Sea case for ORCA_R2 configuration 
    229          zze2    = ( zfwf(3) + zfwf(4) ) / 2. 
     273         zze2    = ( zfwf(3) + zfwf(4) ) * 0.5_wp 
    230274         zfwf(3) = zze2 
    231275         zfwf(4) = zze2 
    232276      ENDIF 
    233277 
     278      zcorr = 0._wp 
     279 
    234280      DO jc = 1, jpncs 
    235281         ! 
    236          IF( ncstt(jc) == 0 ) THEN  
    237             ! water/evap excess is shared by all open ocean 
    238             zcoef  = zfwf(jc) / surf(jpncs+1) 
    239             zcoef1 = rcp * zcoef 
    240             emp(:,:) = emp(:,:) + zcoef 
    241             qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
    242          ELSEIF( ncstt(jc) == 1 ) THEN  
    243             ! Excess water in open sea, at outflow location, excess evap shared 
    244             IF ( zfwf(jc) <= 0.e0 ) THEN  
    245                 DO jn = 1, ncsnr(jc) 
     282         ! The following if avoids the redistribution of the round off 
     283         IF ( ABS(zfwf(jc) / surf(jpncs+1) ) > rsmall) THEN 
     284            ! 
     285            IF( ncstt(jc) == 0 ) THEN           ! water/evap excess is shared by all open ocean 
     286               zcoef    = zfwf(jc) / surf(jpncs+1) 
     287               zcoef1   = rcp * zcoef 
     288               emp(:,:) = emp(:,:) + zcoef 
     289               qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
     290               ! accumulate closed seas correction 
     291               zcorr    = zcorr    + zcoef 
     292               ! 
     293            ELSEIF( ncstt(jc) == 1 ) THEN       ! Excess water in open sea, at outflow location, excess evap shared 
     294               IF ( zfwf(jc) <= 0.e0_wp ) THEN  
     295                   DO jn = 1, ncsnr(jc) 
     296                     ji = mi0(ncsir(jc,jn)) 
     297                     jj = mj0(ncsjr(jc,jn)) ! Location of outflow in open ocean 
     298                     IF (      ji > 1 .AND. ji < jpi   & 
     299                         .AND. jj > 1 .AND. jj < jpj ) THEN  
     300                         zcoef      = zfwf(jc) / ( REAL(ncsnr(jc)) * e1e2t(ji,jj) ) 
     301                         zcoef1     = rcp * zcoef 
     302                         emp(ji,jj) = emp(ji,jj) + zcoef 
     303                         qns(ji,jj) = qns(ji,jj) - zcoef1 * sst_m(ji,jj) 
     304                     ENDIF  
     305                   END DO  
     306               ELSE  
     307                   zcoef    = zfwf(jc) / surf(jpncs+1) 
     308                   zcoef1   = rcp * zcoef 
     309                   emp(:,:) = emp(:,:) + zcoef 
     310                   qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
     311                   ! accumulate closed seas correction 
     312                   zcorr    = zcorr    + zcoef 
     313               ENDIF 
     314            ELSEIF( ncstt(jc) == 2 ) THEN       ! Excess e-p-r (either sign) goes to open ocean, at outflow location 
     315               DO jn = 1, ncsnr(jc) 
    246316                  ji = mi0(ncsir(jc,jn)) 
    247317                  jj = mj0(ncsjr(jc,jn)) ! Location of outflow in open ocean 
    248                   IF (      ji > 1 .AND. ji < jpi   & 
    249                       .AND. jj > 1 .AND. jj < jpj ) THEN  
    250                       zcoef  = zfwf(jc) / ( REAL(ncsnr(jc), wp) * e1t(ji,jj) * e2t(ji,jj) ) 
    251                       zcoef1 = rcp * zcoef 
    252                       emp(ji,jj) = emp(ji,jj) + zcoef 
    253                       qns(ji,jj) = qns(ji,jj) - zcoef1 * sst_m(ji,jj) 
    254                   END IF  
    255                 END DO  
    256             ELSE  
    257                 zcoef  = zfwf(jc) / surf(jpncs+1) 
    258                 zcoef1 = rcp * zcoef 
    259                 emp(:,:) = emp(:,:) + zcoef 
    260                 qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
    261             ENDIF 
    262          ELSEIF( ncstt(jc) == 2 ) THEN  
    263             ! Excess e-p+r (either sign) goes to open ocean, at outflow location 
    264             IF(      ji > 1 .AND. ji < jpi    & 
    265                .AND. jj > 1 .AND. jj < jpj ) THEN  
    266                 DO jn = 1, ncsnr(jc) 
    267                   ji = mi0(ncsir(jc,jn)) 
    268                   jj = mj0(ncsjr(jc,jn)) ! Location of outflow in open ocean 
    269                   zcoef  = zfwf(jc) / ( REAL(ncsnr(jc), wp) * e1t(ji,jj) * e2t(ji,jj) ) 
    270                   zcoef1 = rcp * zcoef 
    271                   emp(ji,jj) = emp(ji,jj) + zcoef 
    272                   qns(ji,jj) = qns(ji,jj) - zcoef1 * sst_m(ji,jj) 
    273                 END DO  
     318                  IF(      ji > 1 .AND. ji < jpi    & 
     319                     .AND. jj > 1 .AND. jj < jpj ) THEN  
     320                     zcoef      = zfwf(jc) / ( REAL(ncsnr(jc)) *  e1e2t(ji,jj) ) 
     321                     zcoef1     = rcp * zcoef 
     322                     emp(ji,jj) = emp(ji,jj) + zcoef 
     323                     qns(ji,jj) = qns(ji,jj) - zcoef1 * sst_m(ji,jj) 
     324                  ENDIF  
     325               END DO  
    274326            ENDIF  
    275          ENDIF  
    276          ! 
    277          DO jj = ncsj1(jc), ncsj2(jc) 
    278             DO ji = ncsi1(jc), ncsi2(jc) 
    279                zcoef  = zfwf(jc) / surf(jc) 
    280                zcoef1 = rcp * zcoef 
    281                emp(ji,jj) = emp(ji,jj) - zcoef 
    282                qns(ji,jj) = qns(ji,jj) + zcoef1 * sst_m(ji,jj) 
    283             END DO   
    284          END DO  
    285          ! 
     327            ! 
     328            DO jj = ncsj1(jc), ncsj2(jc) 
     329               DO ji = ncsi1(jc), ncsi2(jc) 
     330                  zcoef      = zfwf(jc) / surf(jc) 
     331                  zcoef1     = rcp * zcoef 
     332                  emp(ji,jj) = emp(ji,jj) - zcoef 
     333                  qns(ji,jj) = qns(ji,jj) + zcoef1 * sst_m(ji,jj) 
     334               END DO   
     335            END DO  
     336            ! 
     337         END IF 
    286338      END DO  
    287       ! 
    288       CALL lbc_lnk( emp , 'T', 1. ) 
     339 
     340      IF ( ABS(zcorr) > rsmall ) THEN      ! remove the global correction from the closed seas 
     341         DO jc = 1, jpncs                  ! only if it is large enough 
     342            DO jj = ncsj1(jc), ncsj2(jc) 
     343               DO ji = ncsi1(jc), ncsi2(jc) 
     344                  emp(ji,jj) = emp(ji,jj) - zcorr 
     345                  qns(ji,jj) = qns(ji,jj) + rcp * zcorr * sst_m(ji,jj) 
     346               END DO   
     347             END DO  
     348          END DO 
     349      ENDIF 
     350      ! 
     351      emp (:,:) = emp (:,:) * tmask(:,:,1) 
     352      ! 
     353      CALL lbc_lnk( emp , 'T', 1._wp ) 
     354      ! 
     355      IF( nn_timing == 1 )  CALL timing_stop('sbc_clo') 
    289356      ! 
    290357   END SUBROUTINE sbc_clo 
    291     
    292     
     358 
     359 
    293360   SUBROUTINE clo_rnf( p_rnfmsk ) 
    294361      !!--------------------------------------------------------------------- 
     
    314381               ii = mi0( ncsir(jc,jn) ) 
    315382               ij = mj0( ncsjr(jc,jn) ) 
    316                p_rnfmsk(ii,ij) = MAX( p_rnfmsk(ii,ij), 1.0 ) 
     383               p_rnfmsk(ii,ij) = MAX( p_rnfmsk(ii,ij), 1.0_wp ) 
    317384            END DO  
    318385         ENDIF  
     
    342409         DO jj = ncsj1(jc), ncsj2(jc) 
    343410            DO ji = ncsi1(jc), ncsi2(jc) 
    344                p_upsmsk(ji,jj) = 0.5            ! mixed upstream/centered scheme over closed seas 
     411               p_upsmsk(ji,jj) = 0.5_wp         ! mixed upstream/centered scheme over closed seas 
    345412            END DO  
    346413         END DO  
     
    380447   !!====================================================================== 
    381448END MODULE closea 
     449 
Note: See TracChangeset for help on using the changeset viewer.