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 12377 for NEMO/trunk/src/TOP/TRP/trcsink.F90 – NEMO

Ignore:
Timestamp:
2020-02-12T15:39:06+01:00 (4 years ago)
Author:
acc
Message:

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

Location:
NEMO/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk

    • Property svn:externals
      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r11615_ENHANCE-04_namelists_as_internalfiles_agrif@HEAD      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
  • NEMO/trunk/src/TOP/TRP/trcsink.F90

    r11536 r12377  
    2424   INTEGER, PUBLIC :: nitermax      !: Maximum number of iterations for sinking 
    2525 
     26   !! * Substitutions 
     27#  include "do_loop_substitute.h90" 
    2628   !!---------------------------------------------------------------------- 
    2729   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
     
    3537   !!---------------------------------------------------------------------- 
    3638 
    37    SUBROUTINE trc_sink ( kt, pwsink, psinkflx, jp_tra, rsfact ) 
     39   SUBROUTINE trc_sink ( kt, Kbb, Kmm, pwsink, psinkflx, jp_tra, rsfact ) 
    3840      !!--------------------------------------------------------------------- 
    3941      !!                     ***  ROUTINE trc_sink  *** 
     
    4547      !!--------------------------------------------------------------------- 
    4648      INTEGER , INTENT(in)  :: kt 
     49      INTEGER , INTENT(in)  :: Kbb, Kmm 
    4750      INTEGER , INTENT(in)  :: jp_tra    ! tracer index index       
    4851      REAL(wp), INTENT(in)  :: rsfact    ! time step duration 
     
    7073         iiter(:,:) = 1 
    7174      ELSE 
    72          DO jj = 1, jpj 
    73             DO ji = 1, jpi 
    74                iiter(ji,jj) = 1 
    75                DO jk = 1, jpkm1 
    76                   IF( tmask(ji,jj,jk) == 1.0 ) THEN 
    77                       zwsmax =  0.5 * e3t_n(ji,jj,jk) * rday / rsfact 
    78                       iiter(ji,jj) =  MAX( iiter(ji,jj), INT( pwsink(ji,jj,jk) / zwsmax ) ) 
    79                   ENDIF 
    80                END DO 
    81             END DO 
    82          END DO 
     75         DO_2D_11_11 
     76            iiter(ji,jj) = 1 
     77            DO jk = 1, jpkm1 
     78               IF( tmask(ji,jj,jk) == 1.0 ) THEN 
     79                   zwsmax =  0.5 * e3t(ji,jj,jk,Kmm) * rday / rsfact 
     80                   iiter(ji,jj) =  MAX( iiter(ji,jj), INT( pwsink(ji,jj,jk) / zwsmax ) ) 
     81               ENDIF 
     82            END DO 
     83         END_2D 
    8384         iiter(:,:) = MIN( iiter(:,:), nitermax ) 
    8485      ENDIF 
    8586 
    86       DO jk = 1,jpkm1 
    87          DO jj = 1, jpj 
    88             DO ji = 1, jpi 
    89                IF( tmask(ji,jj,jk) == 1.0 ) THEN 
    90                  zwsmax = 0.5 * e3t_n(ji,jj,jk) * rday / rsfact 
    91                  zwsink(ji,jj,jk) = MIN( pwsink(ji,jj,jk), zwsmax * REAL( iiter(ji,jj), wp ) ) 
    92                ELSE 
    93                  ! provide a default value so there is no use of undefinite value in trc_sink2 for zwsink2 initialization 
    94                  zwsink(ji,jj,jk) = 0. 
    95                ENDIF 
    96             END DO 
    97          END DO 
    98       END DO 
     87      DO_3D_11_11( 1,jpkm1 ) 
     88         IF( tmask(ji,jj,jk) == 1.0 ) THEN 
     89           zwsmax = 0.5 * e3t(ji,jj,jk,Kmm) * rday / rsfact 
     90           zwsink(ji,jj,jk) = MIN( pwsink(ji,jj,jk), zwsmax * REAL( iiter(ji,jj), wp ) ) 
     91         ELSE 
     92           ! provide a default value so there is no use of undefinite value in trc_sink2 for zwsink2 initialization 
     93           zwsink(ji,jj,jk) = 0. 
     94         ENDIF 
     95      END_3D 
    9996 
    10097      !  Initializa to zero all the sinking arrays  
     
    104101      !   Compute the sedimentation term using trc_sink2 for the considered sinking particle 
    105102      !   ----------------------------------------------------- 
    106       CALL trc_sink2( zwsink, psinkflx, jp_tra, iiter, rsfact ) 
     103      CALL trc_sink2( Kbb, Kmm, zwsink, psinkflx, jp_tra, iiter, rsfact ) 
    107104      ! 
    108105      IF( ln_timing )   CALL timing_stop('trc_sink') 
     
    110107   END SUBROUTINE trc_sink 
    111108 
    112    SUBROUTINE trc_sink2( pwsink, psinkflx, jp_tra, kiter, rsfact ) 
     109   SUBROUTINE trc_sink2( Kbb, Kmm, pwsink, psinkflx, jp_tra, kiter, rsfact ) 
    113110      !!--------------------------------------------------------------------- 
    114111      !!                     ***  ROUTINE trc_sink2  *** 
     
    121118      !!      transport term, i.e.  div(u*tra). 
    122119      !!--------------------------------------------------------------------- 
     120      INTEGER,  INTENT(in   )                         ::   Kbb, Kmm  ! time level indices 
    123121      INTEGER,  INTENT(in   )                         ::   jp_tra    ! tracer index index       
    124122      REAL(wp), INTENT(in   )                         ::   rsfact    ! duration of time step 
     
    136134      ztraz(:,:,:) = 0.e0 
    137135      zakz (:,:,:) = 0.e0 
    138       ztrb (:,:,:) = trb(:,:,:,jp_tra) 
     136      ztrb (:,:,:) = tr(:,:,:,jp_tra,Kbb) 
    139137 
    140138      DO jk = 1, jpkm1 
     
    147145      DO jn = 1, 2 
    148146         !  first guess of the slopes interior values 
    149          DO jj = 1, jpj 
    150             DO ji = 1, jpi 
    151                ! 
    152                zstep = rsfact / REAL( kiter(ji,jj), wp ) / 2. 
    153                !               
    154                DO jk = 2, jpkm1 
    155                   ztraz(ji,jj,jk) = ( trb(ji,jj,jk-1,jp_tra) - trb(ji,jj,jk,jp_tra) ) * tmask(ji,jj,jk) 
    156                END DO 
    157                ztraz(ji,jj,1  ) = 0.0 
    158                ztraz(ji,jj,jpk) = 0.0 
    159  
    160                ! slopes 
    161                DO jk = 2, jpkm1 
    162                   zign = 0.25 + SIGN( 0.25, ztraz(ji,jj,jk) * ztraz(ji,jj,jk+1) ) 
    163                   zakz(ji,jj,jk) = ( ztraz(ji,jj,jk) + ztraz(ji,jj,jk+1) ) * zign 
    164                END DO 
    165           
    166                ! Slopes limitation 
    167                DO jk = 2, jpkm1 
    168                   zakz(ji,jj,jk) = SIGN( 1., zakz(ji,jj,jk) ) *        & 
    169                      &             MIN( ABS( zakz(ji,jj,jk) ), 2. * ABS(ztraz(ji,jj,jk+1)), 2. * ABS(ztraz(ji,jj,jk) ) ) 
    170                END DO 
    171           
    172                ! vertical advective flux 
    173                DO jk = 1, jpkm1 
    174                   zigma = zwsink2(ji,jj,jk+1) * zstep / e3w_n(ji,jj,jk+1) 
    175                   zew   = zwsink2(ji,jj,jk+1) 
    176                   psinkflx(ji,jj,jk+1) = -zew * ( trb(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep 
    177                END DO 
    178                ! 
    179                ! Boundary conditions 
    180                psinkflx(ji,jj,1  ) = 0.e0 
    181                psinkflx(ji,jj,jpk) = 0.e0 
    182           
    183                DO jk=1,jpkm1 
    184                   zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 
    185                   trb(ji,jj,jk,jp_tra) = trb(ji,jj,jk,jp_tra) + zflx 
    186                END DO 
    187             END DO 
    188          END DO 
     147         DO_2D_11_11 
     148            ! 
     149            zstep = rsfact / REAL( kiter(ji,jj), wp ) / 2. 
     150            !               
     151            DO jk = 2, jpkm1 
     152               ztraz(ji,jj,jk) = ( tr(ji,jj,jk-1,jp_tra,Kbb) - tr(ji,jj,jk,jp_tra,Kbb) ) * tmask(ji,jj,jk) 
     153            END DO 
     154            ztraz(ji,jj,1  ) = 0.0 
     155            ztraz(ji,jj,jpk) = 0.0 
     156 
     157            ! slopes 
     158            DO jk = 2, jpkm1 
     159               zign = 0.25 + SIGN( 0.25, ztraz(ji,jj,jk) * ztraz(ji,jj,jk+1) ) 
     160               zakz(ji,jj,jk) = ( ztraz(ji,jj,jk) + ztraz(ji,jj,jk+1) ) * zign 
     161            END DO 
     162       
     163            ! Slopes limitation 
     164            DO jk = 2, jpkm1 
     165               zakz(ji,jj,jk) = SIGN( 1., zakz(ji,jj,jk) ) *        & 
     166                  &             MIN( ABS( zakz(ji,jj,jk) ), 2. * ABS(ztraz(ji,jj,jk+1)), 2. * ABS(ztraz(ji,jj,jk) ) ) 
     167            END DO 
     168       
     169            ! vertical advective flux 
     170            DO jk = 1, jpkm1 
     171               zigma = zwsink2(ji,jj,jk+1) * zstep / e3w(ji,jj,jk+1,Kmm) 
     172               zew   = zwsink2(ji,jj,jk+1) 
     173               psinkflx(ji,jj,jk+1) = -zew * ( tr(ji,jj,jk,jp_tra,Kbb) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep 
     174            END DO 
     175            ! 
     176            ! Boundary conditions 
     177            psinkflx(ji,jj,1  ) = 0.e0 
     178            psinkflx(ji,jj,jpk) = 0.e0 
     179       
     180            DO jk=1,jpkm1 
     181               zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t(ji,jj,jk,Kmm) 
     182               tr(ji,jj,jk,jp_tra,Kbb) = tr(ji,jj,jk,jp_tra,Kbb) + zflx 
     183            END DO 
     184         END_2D 
    189185      END DO 
    190186 
    191       DO jk = 1,jpkm1 
    192          DO jj = 1,jpj 
    193             DO ji = 1, jpi 
    194                zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 
    195                ztrb(ji,jj,jk) = ztrb(ji,jj,jk) + 2. * zflx 
    196             END DO 
    197          END DO 
    198       END DO 
    199  
    200       trb(:,:,:,jp_tra) = ztrb(:,:,:) 
     187      DO_3D_11_11( 1,jpkm1 ) 
     188         zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t(ji,jj,jk,Kmm) 
     189         ztrb(ji,jj,jk) = ztrb(ji,jj,jk) + 2. * zflx 
     190      END_3D 
     191 
     192      tr(:,:,:,jp_tra,Kbb) = ztrb(:,:,:) 
    201193      psinkflx(:,:,:)   = 2. * psinkflx(:,:,:) 
    202194      ! 
     
    216208      !!---------------------------------------------------------------------- 
    217209      ! 
    218       REWIND( numnat_ref )              ! namtrc_rad in reference namelist  
    219210      READ  ( numnat_ref, namtrc_snk, IOSTAT = ios, ERR = 907) 
    220211907   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtrc_snk in reference namelist' ) 
    221       REWIND( numnat_cfg )              ! namtrc_rad in configuration namelist  
    222212      READ  ( numnat_cfg, namtrc_snk, IOSTAT = ios, ERR = 908 ) 
    223213908   IF( ios > 0 )   CALL ctl_nam ( ios , 'namtrc_snk in configuration namelist' ) 
Note: See TracChangeset for help on using the changeset viewer.