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 8882 for branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90 – NEMO

Ignore:
Timestamp:
2017-12-01T18:44:09+01:00 (6 years ago)
Author:
flavoni
Message:

dev_CNRS_2017 branch: merged dev_r7881_ENHANCE09_RK3 with trunk r8864

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90

    r7646 r8882  
    33MODULE agrif_opa_sponge 
    44   !!====================================================================== 
    5    !!                ***  MODULE agrif_opa_update  *** 
    6    !! AGRIF :    
     5   !!                   ***  MODULE  agrif_opa_interp  *** 
     6   !! AGRIF: sponge package for the ocean dynamics (OPA) 
    77   !!====================================================================== 
    8    !! History :   
     8   !! History :  2.0  !  2002-06  (XXX)  Original cade 
     9   !!             -   !  2005-11  (XXX)  
     10   !!            3.2  !  2009-04  (R. Benshila)  
     11   !!            3.6  !  2014-09  (R. Benshila)  
    912   !!---------------------------------------------------------------------- 
    1013#if defined key_agrif 
     14   !!---------------------------------------------------------------------- 
     15   !!   'key_agrif'                                              AGRIF zoom 
     16   !!---------------------------------------------------------------------- 
    1117   USE par_oce 
    1218   USE oce 
    1319   USE dom_oce 
     20   ! 
    1421   USE in_out_manager 
    1522   USE agrif_oce 
    16    USE wrk_nemo   
    1723   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    1824 
     
    2430 
    2531   !!---------------------------------------------------------------------- 
    26    !! NEMO/NST 3.7 , NEMO Consortium (2015) 
     32   !! NEMO/NST 4.0 , NEMO Consortium (2017) 
    2733   !! $Id$ 
    2834   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     
    3137 
    3238   SUBROUTINE Agrif_Sponge_Tra 
    33       !!--------------------------------------------- 
    34       !!   *** ROUTINE Agrif_Sponge_Tra *** 
    35       !!--------------------------------------------- 
    36       REAL(wp) :: timecoeff 
    37       !!--------------------------------------------- 
     39      !!---------------------------------------------------------------------- 
     40      !!                 *** ROUTINE Agrif_Sponge_Tra *** 
     41      !!---------------------------------------------------------------------- 
     42      REAL(wp) ::   zcoef   ! local scalar 
     43      !!---------------------------------------------------------------------- 
    3844      ! 
    3945#if defined SPONGE 
    40       timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
    41  
     46      zcoef = REAL( Agrif_NbStepint(), wp ) / Agrif_rhot() 
     47      ! 
    4248      CALL Agrif_Sponge 
    43       Agrif_SpecialValue=0. 
     49      Agrif_SpecialValue    = 0._wp 
    4450      Agrif_UseSpecialValue = .TRUE. 
    45       tabspongedone_tsn = .FALSE. 
    46  
    47       CALL Agrif_Bc_Variable(tsn_sponge_id,calledweight=timecoeff,procname=interptsn_sponge) 
    48  
     51      tabspongedone_tsn     = .FALSE. 
     52      ! 
     53      CALL Agrif_Bc_Variable( tsn_sponge_id, calledweight=zcoef, procname=interptsn_sponge ) 
     54      ! 
    4955      Agrif_UseSpecialValue = .FALSE. 
    5056#endif 
     
    5460 
    5561   SUBROUTINE Agrif_Sponge_dyn 
    56       !!--------------------------------------------- 
    57       !!   *** ROUTINE Agrif_Sponge_dyn *** 
    58       !!--------------------------------------------- 
    59       REAL(wp) :: timecoeff 
    60       !!--------------------------------------------- 
    61  
     62      !!---------------------------------------------------------------------- 
     63      !!                 *** ROUTINE Agrif_Sponge_dyn *** 
     64      !!---------------------------------------------------------------------- 
     65      REAL(wp) ::   zcoef   ! local scalar 
     66      !!---------------------------------------------------------------------- 
     67      ! 
    6268#if defined SPONGE 
    63       timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
    64  
    65       Agrif_SpecialValue=0. 
     69      zcoef = REAL( Agrif_NbStepint(), wp ) / Agrif_rhot() 
     70      ! 
     71      Agrif_SpecialValue    = 0._wp 
    6672      Agrif_UseSpecialValue = ln_spc_dyn 
    67  
     73      ! 
    6874      tabspongedone_u = .FALSE. 
    6975      tabspongedone_v = .FALSE.          
    70       CALL Agrif_Bc_Variable(un_sponge_id,calledweight=timecoeff,procname=interpun_sponge) 
    71  
     76      CALL Agrif_Bc_Variable( un_sponge_id, calledweight=zcoef, procname=interpun_sponge ) 
     77      ! 
    7278      tabspongedone_u = .FALSE. 
    7379      tabspongedone_v = .FALSE. 
    74       CALL Agrif_Bc_Variable(vn_sponge_id,calledweight=timecoeff,procname=interpvn_sponge) 
    75  
     80      CALL Agrif_Bc_Variable( vn_sponge_id, calledweight=zcoef, procname=interpvn_sponge ) 
     81      ! 
    7682      Agrif_UseSpecialValue = .FALSE. 
    7783#endif 
     
    8187 
    8288   SUBROUTINE Agrif_Sponge 
    83       !!--------------------------------------------- 
    84       !!   *** ROUTINE  Agrif_Sponge *** 
    85       !!--------------------------------------------- 
    86       INTEGER  :: ji,jj,jk 
    87       INTEGER  :: ispongearea, ilci, ilcj 
    88       LOGICAL  :: ll_spdone 
    89       REAL(wp) :: z1spongearea, zramp 
    90       REAL(wp), POINTER, DIMENSION(:,:) :: ztabramp 
    91  
     89      !!---------------------------------------------------------------------- 
     90      !!                 *** ROUTINE  Agrif_Sponge *** 
     91      !!---------------------------------------------------------------------- 
     92      INTEGER  ::   ji, jj, ind1, ind2 
     93      INTEGER  ::   ispongearea 
     94      REAL(wp) ::   z1_spongearea 
     95      REAL(wp), DIMENSION(jpi,jpj) :: ztabramp 
     96      !!---------------------------------------------------------------------- 
     97      ! 
    9298#if defined SPONGE || defined SPONGE_TOP 
    93       ll_spdone=.TRUE. 
    9499      IF (( .NOT. spongedoneT ).OR.( .NOT. spongedoneU )) THEN 
    95          ! Define ramp from boundaries towards domain interior 
    96          ! at T-points 
     100         ! Define ramp from boundaries towards domain interior at T-points 
    97101         ! Store it in ztabramp 
    98          ll_spdone=.FALSE. 
    99  
    100          CALL wrk_alloc( jpi, jpj, ztabramp ) 
    101102 
    102103         ispongearea  = 2 + nn_sponge_len * Agrif_irhox() 
    103          ilci = nlci - ispongearea 
    104          ilcj = nlcj - ispongearea  
    105          z1spongearea = 1._wp / REAL( ispongearea - 2 ) 
    106  
     104         z1_spongearea = 1._wp / REAL( ispongearea - 1 ) 
     105          
    107106         ztabramp(:,:) = 0._wp 
    108107 
     108         ! --- West --- ! 
    109109         IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 
     110            ind1 = 1+nbghostcells 
     111            ind2 = 1+nbghostcells + (ispongearea-1) 
    110112            DO jj = 1, jpj 
    111                IF ( umask(2,jj,1) == 1._wp ) THEN 
    112                  DO ji = 2, ispongearea                   
    113                     ztabramp(ji,jj) = ( ispongearea-ji ) * z1spongearea 
    114                  END DO 
    115                ENDIF 
     113               DO ji = ind1, ind2                   
     114                  ztabramp(ji,jj) = REAL( ind2 - ji ) * z1_spongearea * umask(ind1,jj,1) 
     115               END DO 
    116116            ENDDO 
    117117         ENDIF 
    118118 
     119         ! --- East --- ! 
    119120         IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 
     121            ind1 = nlci - (1+nbghostcells) - (ispongearea-1) 
     122            ind2 = nlci - (1+nbghostcells) 
    120123            DO jj = 1, jpj 
    121                IF ( umask(nlci-2,jj,1) == 1._wp ) THEN 
    122                   DO ji = ilci+1,nlci-1 
    123                      zramp = (ji - (ilci+1) ) * z1spongearea 
    124                      ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 
    125                   ENDDO 
    126                ENDIF 
     124               DO ji = ind1, ind2 
     125                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ji - ind2 ) * z1_spongearea * umask(ind2-1,jj,1) ) 
     126               ENDDO 
    127127            ENDDO 
    128128         ENDIF 
    129129 
     130         ! --- South --- ! 
    130131         IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 
    131             DO ji = 1, jpi 
    132                IF ( vmask(ji,2,1) == 1._wp ) THEN 
    133                   DO jj = 2, ispongearea 
    134                      zramp = ( ispongearea-jj ) * z1spongearea 
    135                      ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 
    136                   END DO 
    137                ENDIF 
     132            ind1 = 1+nbghostcells 
     133            ind2 = 1+nbghostcells + (ispongearea-1) 
     134            DO jj = ind1, ind2 
     135               DO ji = 1, jpi 
     136                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - jj ) * z1_spongearea * vmask(ji,ind1,1) ) 
     137               END DO 
    138138            ENDDO 
    139139         ENDIF 
    140140 
     141         ! --- North --- ! 
    141142         IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 
    142             DO ji = 1, jpi 
    143                IF ( vmask(ji,nlcj-2,1) == 1._wp ) THEN 
    144                   DO jj = ilcj+1,nlcj-1 
    145                      zramp = (jj - (ilcj+1) ) * z1spongearea 
    146                      ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 
    147                   END DO 
    148                ENDIF 
     143            ind1 = nlcj - (1+nbghostcells) - (ispongearea-1) 
     144            ind2 = nlcj - (1+nbghostcells) 
     145            DO jj = ind1, ind2 
     146               DO ji = 1, jpi 
     147                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( jj - ind2 ) * z1_spongearea * vmask(ji,ind2-1,1) ) 
     148               END DO 
    149149            ENDDO 
    150150         ENDIF 
     
    158158         DO jj = 2, jpjm1 
    159159            DO ji = 2, jpim1   ! vector opt. 
    160                fsaht_spu(ji,jj) = 0.5_wp * visc_tra * (ztabramp(ji,jj) + ztabramp(ji+1,jj  )) 
    161                fsaht_spv(ji,jj) = 0.5_wp * visc_tra * (ztabramp(ji,jj) + ztabramp(ji  ,jj+1)) 
    162             END DO 
    163          END DO 
    164  
     160               fsaht_spu(ji,jj) = 0.5_wp * visc_tra * ( ztabramp(ji,jj) + ztabramp(ji+1,jj  ) ) 
     161               fsaht_spv(ji,jj) = 0.5_wp * visc_tra * ( ztabramp(ji,jj) + ztabramp(ji  ,jj+1) ) 
     162            END DO 
     163         END DO 
    165164         CALL lbc_lnk( fsaht_spu, 'U', 1. )   ! Lateral boundary conditions 
    166165         CALL lbc_lnk( fsaht_spv, 'V', 1. ) 
     166          
    167167         spongedoneT = .TRUE. 
    168168      ENDIF 
     
    179179            END DO 
    180180         END DO 
    181  
    182181         CALL lbc_lnk( fsahm_spt, 'T', 1. )   ! Lateral boundary conditions 
    183182         CALL lbc_lnk( fsahm_spf, 'F', 1. ) 
     183          
    184184         spongedoneU = .TRUE. 
    185185      ENDIF 
    186186      ! 
    187       IF (.NOT.ll_spdone) CALL wrk_dealloc( jpi, jpj, ztabramp ) 
    188       ! 
    189187#endif 
    190188      ! 
     
    192190 
    193191 
    194    SUBROUTINE interptsn_sponge(tabres,i1,i2,j1,j2,k1,k2,n1,n2,before) 
    195       !!--------------------------------------------- 
    196       !!   *** ROUTINE interptsn_sponge *** 
    197       !!--------------------------------------------- 
    198       INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2 
    199       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres 
    200       LOGICAL, INTENT(in) :: before 
     192   SUBROUTINE interptsn_sponge( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) 
     193      !!---------------------------------------------------------------------- 
     194      !!                 *** ROUTINE interptsn_sponge *** 
     195      !!---------------------------------------------------------------------- 
     196      INTEGER                                     , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2 
     197      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   tabres 
     198      LOGICAL                                     , INTENT(in   ) ::  before 
    201199      ! 
    202200      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    205203      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ztu, ztv 
    206204      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2) ::tsbdiff 
     205      !!---------------------------------------------------------------------- 
    207206      ! 
    208207      IF( before ) THEN 
     
    241240                        zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
    242241                        ! horizontal diffusive trends 
    243                         ztsa = zbtr * (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk)  ) 
     242                        ztsa = zbtr * (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  ) 
    244243                        ! add it to the general tracer trends 
    245244                        tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + ztsa 
     
    258257 
    259258 
    260    SUBROUTINE interpun_sponge(tabres,i1,i2,j1,j2,k1,k2, before) 
    261       !!--------------------------------------------- 
    262       !!   *** ROUTINE interpun_sponge *** 
    263       !!---------------------------------------------     
    264       INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 
    265       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 
    266       LOGICAL, INTENT(in) :: before 
    267  
    268       INTEGER :: ji,jj,jk 
    269  
     259   SUBROUTINE interpun_sponge( tabres, i1, i2, j1, j2, k1, k2, before ) 
     260      !!---------------------------------------------------------------------- 
     261      !!                 *** ROUTINE interpun_sponge *** 
     262      !!---------------------------------------------------------------------- 
     263      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2 
     264      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   tabres 
     265      LOGICAL                               , INTENT(in   ) ::   before 
     266      !! 
     267      INTEGER :: ji, jj, jk 
    270268      ! sponge parameters  
     269      INTEGER :: jmax 
    271270      REAL(wp) :: ze2u, ze1v, zua, zva, zbtr 
    272271      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ubdiff 
    273272      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff 
    274       INTEGER :: jmax 
    275       !!---------------------------------------------     
     273      !!---------------------------------------------------------------------- 
    276274      ! 
    277275      IF( before ) THEN 
    278276         tabres = un(i1:i2,j1:j2,:) 
    279277      ELSE 
    280          ubdiff(i1:i2,j1:j2,:) = (ub(i1:i2,j1:j2,:) - tabres(:,:,:))*umask(i1:i2,j1:j2,:) 
     278         ubdiff(i1:i2,j1:j2,:) = ( ub(i1:i2,j1:j2,:) - tabres(:,:,:) )*umask(i1:i2,j1:j2,:) 
    281279         ! 
    282280         DO jk = 1, jpkm1                                 ! Horizontal slab 
     
    297295               DO ji = i1,i2   ! vector opt. 
    298296                  zbtr = r1_e1e2f(ji,jj) * e3f_n(ji,jj,jk) * fsahm_spf(ji,jj) 
    299                   rotdiff(ji,jj,jk) = (-e1u(ji,jj+1) * ubdiff(ji,jj+1,jk) & 
    300                                        +e1u(ji,jj  ) * ubdiff(ji,jj  ,jk) &  
    301                                     & ) * fmask(ji,jj,jk) * zbtr  
     297                  rotdiff(ji,jj,jk) = ( -e1u(ji,jj+1) * ubdiff(ji,jj+1,jk)   & 
     298                                    &   +e1u(ji,jj  ) * ubdiff(ji,jj  ,jk) ) * fmask(ji,jj,jk) * zbtr  
    302299               END DO 
    303300            END DO 
     
    312309                     ze1v = hdivdiff(ji,jj,jk) 
    313310                     ! horizontal diffusive trends 
    314                      zua = - ( ze2u - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) )   & 
    315                            + ( hdivdiff(ji+1,jj,jk) - ze1v  ) / e1u(ji,jj) 
     311                     zua = - ( ze2u - rotdiff (ji,jj-1,jk) ) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) )   & 
     312                           + ( hdivdiff(ji+1,jj,jk) - ze1v ) * r1_e1u(ji,jj) 
    316313 
    317314                     ! add it to the general momentum trends 
     
    327324 
    328325         jmax = j2-1 
    329          IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-3) 
     326         IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-nbghostcells-2)   ! North 
    330327 
    331328         DO jj = j1+1, jmax 
     
    338335 
    339336                     ! horizontal diffusive trends 
    340                      zva = + ( ze2u - rotdiff (ji-1,jj,jk)) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) )   & 
    341                            + ( hdivdiff(ji,jj+1,jk) - ze1v  ) / e2v(ji,jj) 
     337                     zva = + ( ze2u - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) )   & 
     338                           + ( hdivdiff(ji,jj+1,jk) - ze1v ) * r1_e2v(ji,jj) 
    342339 
    343340                     ! add it to the general momentum trends 
     
    356353 
    357354 
    358    SUBROUTINE interpvn_sponge(tabres,i1,i2,j1,j2,k1,k2, before,nb,ndir) 
    359       !!--------------------------------------------- 
    360       !!   *** ROUTINE interpvn_sponge *** 
    361       !!---------------------------------------------  
    362       INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 
    363       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 
    364       LOGICAL, INTENT(in) :: before 
    365       INTEGER, INTENT(in) :: nb , ndir 
    366       ! 
    367       INTEGER  ::   ji, jj, jk 
    368       REAL(wp) ::   ze2u, ze1v, zua, zva, zbtr 
    369       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: vbdiff 
    370       REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff 
    371       INTEGER :: imax 
    372       !!---------------------------------------------  
     355   SUBROUTINE interpvn_sponge( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) 
     356      !!---------------------------------------------------------------------- 
     357      !!                 *** ROUTINE interpvn_sponge *** 
     358      !!---------------------------------------------------------------------- 
     359      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2 
     360      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   tabres 
     361      LOGICAL                               , INTENT(in   ) ::   before 
     362      INTEGER                               , INTENT(in   ) ::   nb , ndir 
     363      ! 
     364      INTEGER ::   ji, jj, jk 
     365      INTEGER ::   imax 
     366      REAL(wp)::   ze2u, ze1v, zua, zva, zbtr 
     367      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) ::   vbdiff, rotdiff, hdivdiff 
     368      !!---------------------------------------------------------------------- 
    373369 
    374370      IF( before ) THEN  
     
    376372      ELSE 
    377373         ! 
    378          vbdiff(i1:i2,j1:j2,:) = (vb(i1:i2,j1:j2,:) - tabres(:,:,:))*vmask(i1:i2,j1:j2,:) 
     374         vbdiff(i1:i2,j1:j2,:) = ( vb(i1:i2,j1:j2,:) - tabres(:,:,:) ) * vmask(i1:i2,j1:j2,:) 
    379375         ! 
    380376         DO jk = 1, jpkm1                                 ! Horizontal slab 
     
    403399         !                                                 
    404400 
    405          imax = i2-1 
    406          IF ((nbondi == 1).OR.(nbondi == 2))   imax = MIN(imax,nlci-3) 
     401         imax = i2 - 1 
     402         IF ((nbondi == 1).OR.(nbondi == 2))   imax = MIN(imax,nlci-nbghostcells-2)   ! East 
    407403 
    408404         DO jj = j1+1, j2 
     
    437433 
    438434#else 
     435   !!---------------------------------------------------------------------- 
     436   !!   Empty module                                          no AGRIF zoom 
     437   !!---------------------------------------------------------------------- 
    439438CONTAINS 
    440439   SUBROUTINE agrif_opa_sponge_empty 
    441       !!--------------------------------------------- 
    442       !!   *** ROUTINE agrif_OPA_sponge_empty *** 
    443       !!--------------------------------------------- 
    444440      WRITE(*,*)  'agrif_opa_sponge : You should not have seen this print! error?' 
    445441   END SUBROUTINE agrif_opa_sponge_empty 
Note: See TracChangeset for help on using the changeset viewer.